Better than Rician: modelling millimetre wave channels as two-wave with diffuse power

Abstract

This contribution provides experimental evidence for the two-wave with diffuse power (TWDP) fading model. We have conducted two indoor millimetre wave measurement campaigns with directive horn antennas at both link ends. One horn antenna is mounted in a corner of our laboratory, while the other is steerable and scans azimuth and elevation. Our first measurement campaign is based on scalar network analysis with 7 GHz of bandwidth. Our second measurement campaign obtains magnitude and phase information; it is additionally sampled directionally at several positions in space. We apply Akaike’s information criterion to decide whether Rician fading sufficiently explains the data or the generalised TWDP fading model is necessary. Our results indicate that the TWDP fading hypothesis is favoured over Rician fading in situations where the steerable antenna is pointing towards reflecting objects or is slightly misaligned at line-of-sight. We demonstrate TWDP fading in several different domains, namely, frequency, space, and time.

1 Introduction

Accurate modelling of wireless propagation effects is a fundamental prerequisite for a proper communication system design. After the introduction of the double-directional radio channel model [1], wireless propagation research (<6 GHz) started to model the wireless channel agnostic to the antennas used. More than a decade later, propagation research focusses now on millimetre wave bands to unlock the large bandwidths available in this regime [2–5]. At millimetre waves (mmWaves), omnidirectional antennas have small effective antenna areas, resulting in a high path loss [6–10]. To overcome this high path loss, researchers have proposed to apply highly directive antennas on both link ends [11–14]. Most researchers aim to achieve high directivity with antenna arrays [15–20] and a few with dielectric lenses [21–23]. When the link quality depends so much on the achieved beam-forming gain, antennas must be considered as part of the wireless channel again. Small-scale fading is then influenced by the antenna.

According to Durgin [24, p. 137], “The use of directive antennas or arrays at a receiver, for example, amplifies several of the strongest multipath waves that arrive in one particular direction while attenuating the remaining waves. This effectively increases the ratio of specular to nonspecular received power, turning a Rayleigh or Rician fading channel into a TWDP fading channel.” The mentioned two-wave with diffuse power (TWDP) fading channel describes this spatial filtering effect by two non-fluctuating receive signals together with many smaller diffuse components.

1.1 Related work

The authors of [25] investigated a simple wall scattering scenario and analysed how fading scales with various antenna directivities and different bandwidths. Increasing directivity [25], as well as increasing bandwidth [25, 26], results in an increased Rician K-factor. The authors of [27] analysed fading at 28 GHz with high gain horn antennas on both link ends. They observe high Rician K-factors even at non-line-of-sight (NLOS). This effect is explained by spatial filtering of directive antennas, as they suppress many multipath components [25]. Outdoor measurements in [28, 29], show a graphical agreement with the Rice fit, but especially Fig. 10 in [29] might be better explained as TWDP fading.

TWDP fading has already successfully been applied to describe 60 GHz near body shadowing [30]. Furthermore, as quoted above, TWDP must be considered for arrays, as they act as spatial filters [24, 31]. While theoretical work on TWDP fading is already advanced [32–37], experimental evidence, especially at millimetre waves, is still limited. For enclosed structures, such as aircraft cabins and buses, the applicability of the TWDP model is demonstrated by Frolik [38–42]. A deterministic two ray behaviour in ray tracing data of mmWave train-to-infrastructure communications is shown in [43]. A further extension of the TWDP-fading model, the so-called fluctuating two-ray fading model, was also successfully applied to fit mmWave measurement data [44–46]. This model brings in another degree of freedom and allows for common shadowing of both specular waves. The wireless channels in this present study are unblocked; thus, this model is not considered here.

Our group has conducted two measurement campaigns [47, 48] to directionally analyse receive power and small-scale fading parameters for mmWaves. This contribution is based on the measurement data gathered in [47, 48].

1.2 Outline and contributions

With this contribution, we aim to bring scientific rigour to the small-scale fading analysis of millimetre wave indoor channels. We show in Section 2 —by means of an information-theoretic approach [49] and null hypothesis testing [50]—that the TWDP model has evidence in mmWave communications.

We have conducted two measurement campaigns within the same laboratory with different channel sounding concepts. Our measurements are carried out in the V-band; the applied centre frequency is 60 GHz. For both measurement campaigns, 20 dBi horn antennas are used at the transmitter and at the receiver. The first measurement campaign (MC1) samples the channel in azimuth (φ) and elevation (θ), keeping the antenna’s (apparent) phase centre ([51, pp. 799]) at a fixed (x,y) coordinate. The transmitter is mounted in a corner of our laboratory. The sounded environment as well as the mechanical set-ups are explained in Section 3. For MC1, we sounded the channel in the frequency domain by aid of scalar network analysis, described in Section 4. These channel measurements span over 7 GHz bandwidth, supporting us to analyse fading in the frequency domain.

For the second measurement campaign (MC2), described in Section 5, we improved the set-up mechanically and radio frequency (RF) wise. By adding another linear guide along the z-axis, we keep the antenna’s phase centre constant in (x,y,z) coordinate, irrespective of the antenna’s elevation. Furthermore, we changed the sounding concept to time-domain channel sounding. This approach allows us to utilise the time domain and to show channel impulse responses in Section 7. Additionally, by adjusting (x,y,z,φ,θ), we sample the channel in the spatial domain at all directions (φ,θ). These improvements enable us to show spatial correlations in Section 6, a further analysis tool to support the claims from MC1.

In summary, we demonstrate TWDP fading for directional mmWave indoor channels in the frequency-domain, in the spatial-domain, and in the time-domain.

2 Methodology—fading model identification

TWDP fading captures the effect of interference of two non-fluctuating radio signals and many smaller so-called diffuse signals [31]. The TWDP distribution degenerates to Rice if one of the two non-fluctuating radio signals vanishes. This is analogous to the well-known Rice degeneration to the Rayleigh distribution with decreasing K-factor. In the framework of model selection, TWDP fading, Rician fading, and Rayleigh fading are hence nested hypotheses [49]. Therefore, it is also obvious that among these alternatives, TWDP always allows the best possible fit of measurement data. Occam’s razor [52] asks to select, among competing hypothetical distributions, the hypothesis that makes the fewest assumptions. Different distribution functions are often compared via a goodness-of-fit test [53]. Nevertheless, the authors of [54] argue that Akaike’s information criterion (AIC) [49, 55–57] is better suited for the purpose of choosing among fading distributions. Later on, the AIC was also used in [58–62]. The AIC can be seen as a form of Occam’s razor as it penalises the number of estimable parameters in the approximating model [49] and hence aims for parsimony.

2.1 Mathematical description of TWDP fading

An early form of TWDP fading was analysed in [32]. Durgin et al. [31] introduced a random phase superposition formalism. Later, Rao et al. [35] achieved a major breakthrough and found a description of TWDP fading as conditional Rician fading. For the benefit of the reader, we will briefly repeat some important steps of [35].

The TWDP fading model in the complex-valued baseband is given as:

$$ r_{\text{complex}} = V_{1} e^{j\phi_{1}}+V_{2} e^{j\phi_{2}}+ X + jY~, $$

(1)

where V₁≥0 and V₂≥0 are the deterministic amplitudes of the non-fluctuating specular components. The phases ϕ₁ and ϕ₂ are independent and uniformly distributed in (0,2π). The diffuse components are modelled via the law of large numbers as X+jY, where $X,Y\sim \mathcal {N}\left (0,\sigma ^{2}\right)$. The K-factor is the power ratio of the specular components to the diffuse components:

$$ K = \frac{V_{1}^{2}+V_{2}^{2}}{2\sigma^{2}}~. $$

(2)

The parameter Δ describes the amplitude relationship among the specular components:

$$ \Delta=\frac{2V_{1} V_{2}}{V_{1}^{2} + V_{2}^{2}}~. $$

(3)

The Δ-parameter is bounded between 0 and 1 and equals 1 iff and only if both amplitudes are equal. The second moment of the envelope r=|r_complex| of TWDP fading is given as:

$$ \mathbb{E}\big[r^{2} \big]=\Omega=V_{1}^{2}+V_{2}^{2}+2\sigma^{2}~. $$

(4)

Expectation is denoted by $\mathbb {E}$. For bounded amplitudes V₁ and V₂, a clever choice of σ² normalises Ω, that is Ω≡1. Starting from (4), by using (2) we arrive at:

$$ \sigma^{2} = \frac{1}{2(1+K)}~. $$

(5)

Given (K,Δ,Ω≡1), the authors of [63] provide a formula for the amplitudes of both specular components:

$$ V_{1,2} = \frac{1}{2}\sqrt{\frac{K }{K+1}}\left(\sqrt{1+\Delta} \pm \sqrt{1-\Delta} \right). $$

(6)

Real-world measurement data have Ω≠1. To work with the formalism introduced above, we normalise the measurement data through estimating $\hat {\Omega }$ by the method of moments. The second moment Ω of Rician fading and TWDP fading is merely a scale factor [64, 65]. Notably, we are more concerned with a proper fit of K and Δ. Generally, estimation errors on Ω propagate to K and Δ estimates. However, Lopez-Fernandez et al. [64] achieved an almost asymptotically efficient estimator with a moment-based estimation of Ω.

Our envelope measurements are partitioned into 2 sets. We take the first set (r₁,…,r_n,…r_N) for parameter estimation of the tuple (K,Δ) as described in Section 2.2, and hypothesis testing as described in Section 2.3. The first set is carefully selected to obtain envelopes samples that are as independent as possible. The second set (r₁,…,r_m,…r_M) is the complement of the first set. We use the elements of the second set to estimate the second moment via:

$$ \hat{\Omega}=\frac{1}{M}\sum\limits_{m=1}^{M} r_{m}^{2}~, $$

(7)

where m is the sample index and M is the size of the second set. Partitioning is necessary to avoid biases through noise correlations of $\hat {\Omega }$ and $\left (\hat {K},\hat {\Delta }\right)$ [66].

By considering the estimate (7) as true parameter Ω, all distributions are parametrised by the tuple (K,Δ), solely. Example distributions are shown in Fig. 1. The cumulative distribution function (CDF) of the envelope of (1) is given in [35] as:

$$\begin{array}{*{20}l} F_{\text{TWDP}}(r ; K, \Delta) =& \\ 1 - \frac{1}{2\pi} \int\limits_{0}^{2\pi} Q_{1} \!\! & \left(\sqrt{2 K \left[ 1+ \Delta \cos(\alpha) \right]}, \frac{r}{\sigma} \right) \notag \mathrm{d}\alpha~. \end{array} $$

(8)

Fig. 1

Comparison of Rayleigh, Rician, and TWDP fading. The TWDP distribution with Δ=1 deviates from the Rice distribution. TWDP fading’s probability for deep fades is higher than for a Rayleigh distribution

The Marcum Q-function is denoted by Q₁(·,·). For Δ→0, Eq. (8) reduces to the well-known Rice CDF:

$$ F_{\text{Rice}}(r ; K) = 1 - Q_{1} \left(\sqrt{2 K }, \frac{r}{\sigma} \right)~. $$

(9)

It might sound tempting to have a second strong radio signal present; in fact, however, two waves can either superpose constructively or destructively and eventually lead to fading that is more severe than Rayleigh [38–42]. We observe the highest probability for deep fades for TWDP fading in Fig. 1.

2.2 Parameter estimation and model selection

Note that our model of TWDP fading (1) does not contain noise. Over our wide frequency range (in MC1, we have 7 GHz bandwidth), the receive noise power spectral density is not equal. A statistical noise description that is valid over our wide frequency range is frequency-dependent. To avoid the burden of frequency-dependent noise modelling, we only take measurement samples which lie at least 10 dB above the noise power and ignore noise in our estimation.

Having the envelope measurement data set (r₁,…,r_n,…r_N) at hand, we are seeking a distribution of which the observed realisations r_n appear most likely. To do so, we estimate the parameter tuple (K,Δ) via the maximum likelihood procedure:

$$\begin{array}{*{20}l} \left(\hat{K},\hat{\Delta}\right) &=\text{argmax}_{K,\Delta} \sum\limits_{n=1}^{N} \operatorname{ln} \frac{\partial F_{\text{TWDP}}(r_{n} ; K, \Delta)}{\partial r} \\ &=\text{argmax}_{K,\Delta} \sum\limits_{n=1}^{N}\operatorname{ln} f_{\text{TWDP}}(r_{n}; K, \Delta) \\ &=\text{argmax}_{K,\Delta} \sum\limits_{n=1}^{N} \operatorname{ln} \mathcal{L}(K, \Delta | r_{n})~. \end{array} $$

(10)

We denote the probability density function (PDF) by f(·), n denotes the sample index, and N the size of the set. To solve (10), we first discretise K and Δ in steps of 0.05. Next, we calculate $\frac {\partial F_{\text {TWDP}}(r ; K, \Delta)}{\partial r}$ for all parameters via numerical differentiation. Within this family of distributions, we search for the parameter vector maximising the log-likelihood function (10). For the optimal Rice fit, the maximum is searched within the parameter slice (K,Δ≡0). An exemplary fit of Rician and TWDP fading is shown in Fig. 2. As a reference, Rayleigh fading (K≡0, Δ≡0) is shown as well.

Fig. 2

CDFs: Distribution fitting for exemplary frequency domain measurement data. Illustration of the maximum likelihood fitted Rice distribution and the maximum likelihood fitted TWDP fading distribution. The Rician K-factor and the TWDP K-factor deviate significantly. Rayleigh fading is plotted as reference

To select between Rician fading and TWDP fading, we employ Akaike’s information criterion (AIC). The AIC is a rigourous way to estimate the Kullback-Leibler divergence, that is, the relative entropy based on the maximum-likelihood estimate [49]. Given the maximum-likelihood fitted parameter tuple ($\hat {K}, \hat {\Delta }$) of TWDP fading and Rician fading, we calculate the sample size-corrected AIC [49, p. 66] for Rician fading (AIC_R) or TWDP fading (AIC_T):

$$ \begin{aligned} \text{AIC}_{\text{R/T}} &= - 2 \sum_{n=1}^{N}\operatorname{ln} \mathcal{L}_{\text{R/T}}\left(\hat{K}, \hat{\Delta}|r\right) + 2U_{\text{R/T}}\\ &\quad + \frac{2U_{\text{R/T}} (U_{\text{R/T}}+1)}{N-U_{\text{R/T}}-1}, \end{aligned} $$

(11)

where U is the model order. For Rician fading, the model order is U_R=1, since we estimate the K-factor, only. For TWDP fading U_T=2, as Δ is estimated additionally. The second moment Ω (estimated already with a different data set before the parameter estimation) is not part of the ML estimation (10) and therefore not accounted in the model order U. We choose between Rician fading and TWDP fading based on the lower AIC.

2.3 Validation of the chosen model

Based on (11), one of the two distributions, Rice or TWDP, will always yield a better fit. To validate whether the chosen distributions really explains the data, we state the following statistical hypothesis testing problem:

$$\begin{array}{*{20}l} \mathcal{H}_{0}: &\left\{ \begin{array}{ll} F_{\text{Rice}}\left(r; \hat{K}\right), & \text{if }\text{AIC}_{\mathrm{R}} \le \text{AIC}_{\mathrm{T}} \\ F_{\text{TWDP}}\left(r; \hat{K}, \hat{\Delta}\right), &\text{else} \end{array} \right. \\ \mathcal{H}_{1}: & \left\{ \begin{array}{ll} \neg F_{\text{Rice}}\left(r; \hat{K}\right), &\text{if }\text{AIC}_{\mathrm{R}} \le \text{AIC}_{\mathrm{T}} \\ \neg F_{\text{TWDP}}\left(r; \hat{K},\hat{\Delta}\right), &\text{else} \end{array}\right. \end{array} $$

(12)

The Boolean negation is denoted by ¬. Our statistical tool is the g-test [67, 68]¹. At a significance level α, a null hypothesis is rejected if:

$$ G=2\sum\limits_{i=1}^{m} O_{i} \ln \left(\frac{O_{i}}{E_{i}} \right) \overset{ \mathbf{?}}{>} \chi_{(1-\alpha, m-e)}^{2}~, $$

(13)

where O_i is the observed bin count in cell i and E_i is the expected bin count in cell i under the null hypothesis $\mathcal {H}_{0}$. The cell edges are illustrated with vertical lines in Fig. 2. The cell edges are chosen, such that 10 observed bin counts fall into one cell. The estimated parameters of the model are denoted by e. For Rician fading, we estimate e=2 (Ω,K) parameters, and for TWDP fading, we estimate e=3 (Ω,K,Δ) parameters in total. The (1−α) quantile of the chi-square distribution with m−e degrees of freedom is denoted by $\chi _{(1-\alpha, m-e)}^{2}$. The prescribed confidence level is 1−α=0.01.

3 Floor plan and set-ups for MC1 and MC2

Our measured environment is a mixed office and laboratory room. There are office desks in the middle of the room, and at the window side, there are laboratory desks (see Fig 3.). The main interacting objects in our channel are office desks, a metallic fridge, a wall, and the surface of the laboratory desk. These objects are all marked in Fig. 3.

Fig. 3

Floor plan of the measured environment. The floorplan indicates the multipath components that are visible in the measurement results. TX and RX switch roles for MC2. TX/RX in the right upper corner of the room is always static. RX/TX in the middle of the room is steerable, indicated by the spider’s web

Our directional measurements are carried out by using the traditional approach of mechanically steered directional antennas [69, 70]. As directional antennas, 20 dBi conical horn antennas with an 18° 3 dB opening angle are used. Our polarisation is determined by the LOS polarisation. When TX and RX are facing each other at LOS, the polarisation is co-polarised and the E-field is orthogonal to the floor. In MC1, the essential mechanical adaptation to the state-of-the-art directional channel sounding set-up [71, 72] is that the elevation-over-azimuth positioner is mounted on an xy-positioning stage. Thereby, we compensate for all linear translations caused by rotations and keep the phase centre of the horn antenna always at the same (x,y) coordinate, see Fig. 4. The z coordinate is roughly 70 cm above ground but varies 13 cm for different elevation angles.

Fig. 4

Photograph of the mechanical set-up for MC1 from the receiver point of view. The receive antenna, a conical 20 dBi horn, is mounted on a multi-axis positioning and rotating system. The azimuthal and elevation angle are controlled to scan the whole upper hemisphere. The multi-axis system moves and rotates the horn antenna such that its phase centre stays at the same (x,y) coordinate during the directional scan

For MC2, we add another linear guide along the z-axis to compensate for all introduced offsets. The horn antenna’s phase centre is thereby lifted upwards by 1 m. Now, we are able to fix the phase centre of the horn antenna at a specific (x,y,z) coordinate in space. The whole mechanical set-up and the fixed phase centre is illustrated in Fig. 5.

Fig. 5

Photograph of the improved mechanical set-up for MC2 from the receiver point of view. Our mechanical set-up consists now of five independent axes to fully compensate all offsets introduced by rotation. A schematic sketch is superimposed. All five axes are necessary to rotate the horn antenna around the phase centre at a fixed (x,y,z) coordinate. Notice that TX and RX switch roles as compared to Fig. 4

4 MC1: Scalar-valued wideband measurements

A wireless channel is said to be small-scale fading, if the receiver (RX) cannot distinguish between different multipath components (MPCs). Depending on the position of the transmitter (TX), the position of the RX and the position of the interacting components, the MPCs interfere constructively or destructively [73, pp. 27]. The fading concept only asks for a single carrier frequency, whose MPCs arrive with different phases at the RX. By spatial sampling, a statistical description of the fading process is found.

In MC1, the spatial (x,y) coordinate (of TX and RX) is kept constant. Different phases of the impinging MPCs are realised by changing the TX frequency over a bandwidth of 7 GHz. Thereby, we implicitly rely on frequency translations to estimate the parameters of the spatial fading process.

4.1 Measurement set-up

We measure the forward transfer function with an Rohde and Schwarz R&S ZVA24 vector network analyser (VNA). The VNA can measure directly up to 24 GHz. For mmWave up-conversion and down-conversion, we employ modules from Pasternack [74]. They are based on radio frequency integrated circuits described in [75]. The up-converter module and the down-converter module are operating built-in synthesiser phase-locked loops (PLLs), where the local oscillator (LO) frequency is calculated as:

$$ f_{\text{LO}}=\nicefrac{7}{4} \cdot s_{\text{PLL}} \cdot 285.714\,\text{MHz}\approx s_{\text{PLL}} \cdot 500\,\text{MHz} ~. $$

(14)

The scaling factor of the synthesiser PLL counters is denoted by s_PLL. For f_LO≈60 GHz, the scaling factor is s_PLL=120. To avoid crosstalk, we measure the transfer function via the conversion gain (mixer) measurement option of our VNA and operate the transmitter and receiver at different baseband frequencies: 601 to 1100 MHz and 101 to 600 MHz. The set-up is shown in Fig. 6.

Fig. 6

RF set-up for MC1. The combination of different PLL scaling factors allows for a measurement bandwidth of 7 GHz. The reference clock for the up-converter and the down-converter is shared. The power splitter has an isolation of 30 dB. To avoid possible leakage on the clock distribution network, attenuators additionally decouple both converters. The transfer function is measured applying the conversion gain (mixer) measurement option of the R&S VNA

4.2 Receive power and fading distributions

In Fig. 7, we show the estimated received mean power of 7 GHz bandwidth, normalised to the maximum RX power, that is

$$ P_{\text{RX,norm.}}(\varphi,\theta)=\frac{\hat{\Omega}(\varphi,\theta)}{\max_{\varphi',\theta'}\left(\hat{\Omega}(\varphi',\theta')\right)}~. $$

(15)

Fig. 7

Estimated directional receive power of MC1. There are four main interacting objects leading to stronger receive power (marked in the figure). TWDP fading occurs whenever the LOS-link is not perfectly aligned or the reflecting structure is not perfectly plain. Red diamonds mark TWDP fading, and white circles mark Rician fading. Black markers show points where both distributions are rejected by the hypothesis test. Directions less than 10 dB above the noise level are not evaluated

As already mentioned in Section 2, we partition the frequency measurements into two sets. The normalised receive power is calculated according to (7), with frequency samples spaced by 2.5 MHz. Every tenth sample is left out as these samples are used for fitting of (K,Δ) and hypothesis testing. We display the results via a stereographic projection from the south pole and use tan(θ/2) as azimuthal projection. All sampling points, lying at least 10 dB above the noise level, are subject of our study. They are displayed with red, white, or black markers. Sampling points where we decided for TWDP fading, following the procedure described in Section 2, are marked with red diamonds. White circles mark points for which AIC favours Rician fading. Four points are marked black. These points failed the null hypothesis test, and we neither argue for Rician fading nor for TWDP fading. TWDP fading occurs whenever the line-of-sight (LOS)-link is not perfectly aligned or if the interacting object cannot be described by a pure reflection.

In Fig. 8, the K-parameter of the selected hypothesis is illustrated. Figure 8 shows either the Rician K-factor or the TWDP K-factor, depending on the selected hypothesis. Note that their definitions are fully equivalent. For Rician fading, the amplitude V₂ in (1) is zero by definition. Whenever the RX power is high, the K-factor is high. Below the K-estimate, the estimate of Δ is shown. Here again, by definition, Δ≡0 whenever we decide for Rician fading. For interacting objects, the parameter Δ tends to be close to one. Note, that decisions based on AIC select TWDP fading mostly when Δ is above 0.3. Smaller Δ values do not change the distribution function sufficiently to justify a higher model order.

Fig. 8

Estimated K-factor (a) and Delta-parameter (b) of MC1. We plot the K-factor estimate of the selected hypothesis. The K-factor behaves analogously to the RX power. At LOS, the K-factor is far above 20 dB. The desk reflection has a surprisingly high K-factor of about 15 dB. Other reflections have K-factors of approximately 10 dB. The Δ-parameter for reflections tends to be close to 1. Markers have the same meaning as in Fig. 7

5 MC2: Vector-valued spatial measurements

In contrast to MC1, we no longer rely on frequency translations and are indeed sampling the channel in space. The fading results we present in Section 5.3 are evaluated at a single frequency. Fading is hence determined by the obtained spatial samples, exclusively.

5.1 Measurement set-up

At the transmitter side, a 2-GHz wide waveform is produced by an arbitrary waveform generator (AWG). A multi-tone waveform (OFDM) with Newman phases [76–78] is applied as sounding signal. The signal has 401 tones (sub-carriers) with a spacing of 5 MHz. This large spacing assures that our system is not limited by phase noise [79]. The TX sequence is repeated 2 000 times to obtain a coherent processing gain of 33 dB for i.i.d. noise. The Pasternack up-converter (the same as in MC1) shifts the baseband sequence to 60 GHz. The 20 dBi conical horn antenna, together with the up-converter is mounted on a five axis positioner to directionally steer them. As receiver, a signal analyser (SA) (R&S FSW67) with a 2-GHz analysis bandwidth is used. The received in-phase and quadrature (IQ) baseband samples are obtained from the SA. The whole system is sketched in Fig. 9.

Fig. 9

RF set-up for MC2. The VNA from Fig. 6 is replaced with an AWG and an SA. This leads to a set-up where we obtain phase information as well. An option of the SA gives us direct access to the baseband IQ samples

In MC2, for feasibility reasons, TX and RX switch places compared to MC1. The RX in form of the SA is put onto the laboratory table. The RX 20 dBi conical horn antenna is directly mounted at the RF input of the SA. The SA is located on a table close to a corner of the room; the RX antenna is not steered.

Similar to the set-ups of [80–83], proper triggering between the arbitrary waveform generator and the SA ensures a stable phase between subsequent measurements.

5.2 Receive power

For the calculation of the RX power, averaged over 2 GHz bandwidth, we perform a sweep through azimuth and elevation at a single coordinate. The LOS and wall reflection from MC1 are still visible in Fig. 10. Fading is evaluated at a single frequency in the subsection below. Nevertheless, we already indicate fading distributions by markers in Fig. 10 in order to better orient ourselves later on.

Fig. 10

Estimated directional receive power of MC2. Due to the elevated position of the steerable horn antenna, two interacting objects from Fig. 7, namely the desk and the fridge, are no longer visible. LOS and the wall reflection are still present. These regions are the only ones which are spatially sampled. Markers have the same meaning as in Fig. 7

As the steerable horn antenna is above the office desks and the fridge level, these interacting objects do not become apparent. In case the steerable TX does not hit the RX at LOS accurately, the table surface acts as reflector and a TWDP model explains the data. For wall reflections, with non-ideal alignment, TWDP also explains the data best.

5.3 Fading distributions

To obtain different spatial realisations, with the horn antenna pointing into the same direction, the coordinate of the apparent phase centre is moved to (x,y,z) positions uniformly distributed within a cube of side length 2.8λ (see Fig. 11). We realise a set of 9×9×9=729 directional measurements. This results in a spacing between spatial samples of 0.35λ in each direction. Although λ/2 sampling is quite common [25, 27], we choose the sampling frequency to be co-prime with the wavelength, to circumvent periodic effects [84]. We restrict our spatial extend to avoid changes in large-scale fading. Only at directions with strong reception levels spatial sampling is performed². Similar as in the previous section, we partition the measurements into 2 sets. The partitioning is made according to a 3D chequerboard pattern. The first set is used for the estimation of the second moment $\hat {\Omega }$, and the second set is used for the parameter tuple (K,Δ).

The best fitting K-factors, in both regions with strong reception, are illustrated in Fig. 12, top part. Below the Δ-parameters are provided. Remember, the RX in form of an SA is put on the laboratory table. In case the TX is not perfectly aligned, a reflection from the table surface yields a fading statistic captured by the TWDP model. The interaction with the wall, similar to Fig. 8, has again regions best modelled via TWDP fading.

Fig. 12

Estimated K-factor (a) and Delta-parameter (b) of MC2. Due to the elevated position, the wall reflection has a 6 dB increased K-factor as compared to MC1; see top part of Fig. 8. If the beam is not perfectly aligned, Rician fading turns again into TWDP fading. Wall reflections described by TWDP fading have a Δ-parameter of close to one. The table surface reflection leads to a significantly smaller reflected component (Δ≪1). The encircled sampling points are subject of further study in Sections 6 and 7. Markers have the same meaning as in Fig. 7

6 MC2: Efficient computation of the spatial correlation

The wall reflection from the previous section is now subject to a more detailed study. Our spatial samples are used to show spatial correlations among the drawn samples.

Our three-dimensional sampling problem, see again Fig. 11, is treated via two-dimensional slicing. For the calculation of the spatial (2D) autocorrelation function, we apply the Wiener–Khintchine–Einstein theorem that relates the autocorrelation function of a wide-sense-stationary random process to its power spectrum [85]. In two dimensions, this theorem reads [86, 87]:

$$ \mathcal{F}_{2D}\left\{\boldsymbol{C}(x,y)\right\} = \boldsymbol{S}\left(x',y'\right)\,, $$

(16)

where C is the 2D-autocorrelation and S is the power spectral density of a 2D signal. The operator $\mathcal {F}_{2D}$ denotes the 2D Fourier transform. We calculate all 2D autocorrelation functions C^(z,f) of one x−y slice at height z at a single frequency f through

$$\begin{array}{*{20}l} &\mathcal{F}_{2D}\left\{ \boldsymbol{C}^{(z,f)}(x,y) \right\} = \mathcal{F}_{2D} \left\{ \Re \left\{ \boldsymbol{H}^{(z,f)}\left(x',y'\right) \right\} \right\} \\ & \odot\text{conj} \left\{ \mathcal{F}_{2D} \left\{ \Re \left\{ \boldsymbol{H}^{(z,f)}\left(x',y'\right)\right\} \right\} \right\}~. \end{array} $$

(17)

The symbols ⊙ denotes the Hadamard multiplication. The operator conj{·} denotes complex conjugation. To ensure a real-valued autocorrelation matrix (instead of a generally complex representation [87]), from the complex-valued channel samples, only the real parts ℜ{·} are taken. The spatial autocorrelation of the imaginary parts are identical. One could also analyse the magnitude and phase individually. While the correlation of the magnitude stays almost at 1, the phase correlation patterns are similar to those of the real part.

The 2D Fourier transform $\mathcal {F}_{2D}$ is realised via a 2D discrete Fourier transform (DFT). The 2D DFT is calculated via a multiplication with the DFT matrix D from the left and the right. To mimic a linear convolution with the DFT, zero padding is necessary. We hence take the matrix ${\widetilde {\boldsymbol {H}}}^{(z,f)}$

$$ \widetilde{\boldsymbol{H}}^{(z,f)} = \left(\begin{array}{cc} \Re \left\{\boldsymbol{H}^{(z,f)} \right\}& \mathbf{0} \\ \mathbf{0} & \mathbf{0} \end{array} \right)~. $$

(18)

Furthermore, the finite spatial extend of our samples acts as rectangular window. The rectangular window leads to a triangular envelope of the the autocorrelation function. This windowed spatial correlation is denoted by:

$$\begin{array}{*{20}l} \mathbf{C}_{\text{windowed}}^{(z,f)} = \boldsymbol{D}^{H} &\left(\left(\boldsymbol{D} {\widetilde{\boldsymbol{H}}}^{(z,f)} \boldsymbol{D} \right) \right.\\ & \left. \odot \operatorname{conj} \left\{ \boldsymbol{D} {\widetilde{\boldsymbol{H}}}^{(z,f)} \boldsymbol{D} \right\} \right) \mathbf{D}^{H}. \end{array} $$

(19)

To compensate the windowing effect, we calculate the spatial correlation of the rectangular window, constructed in accordance to (18):

$$\begin{array}{*{20}l} \boldsymbol{S} = \boldsymbol{D}^{H} & \left(\left(\boldsymbol{D} \left(\begin{array}{cc} \mathbbm{1} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} \end{array} \right) \mathbf{D} \right) \right.\\ & \left. \odot \operatorname{conj} \left\{ \boldsymbol{D} \left(\begin{array}{cc} \mathbbm{1} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} \end{array} \right) \boldsymbol{D} \right\} \right) \boldsymbol{D}^{H}~. \end{array} $$

The matrix 𝟙 denotes the all-ones matrix. Matrix S compensates the truncation effect of the autocorrelation through element-wise (Hadamard) division, denoted by ⊘. Finally, the efficient computation of the spatial correlation (17) reads:

$$ \boldsymbol{C}^{(z,f)} = \boldsymbol{C}_{\text{windowed}}^{(z,f)} \oslash \boldsymbol{S}~. $$

(20)

At a distance of 0.35λ, the measurement data is still correlated; therefore, we are able to view our correlation results on the finer, interpolated grid. The interpolation factor is 20. That means that our spatial correlations are shown on a grid of 0.35λ/20=0.0175λ distance. The very efficient implementation of (20) is applied to all (parallel) 2D slices and to all frequencies. All realisations in z and f are averaged:

$$ \bar{\boldsymbol{C}}= \frac{1}{9} \frac{1}{401} \sum\limits_{z=1}^{9} \sum\limits_{f=1}^{401} \boldsymbol{C}^{(z,f)}~. $$

(21)

Furthermore, we plot one-dimensional autocorrelation functions, evaluated along x and y, together with their two-dimensional representations. We provide two spatial correlation plots evaluated at an azimuth angle of φ=340° and φ=160° in Fig. 13, both at an elevation angle of θ=110°. The bottom part of Fig. 13 shows a correlation pattern dominated by a single wave. The spatial correlation above shows an interference pattern, which is intuitively explained by a superposition of two plane waves. The one dimensional correlations, evaluated either at the x-axis or at the y-axis, show this oscillatory behaviour as well.

Fig. 13

Spatial correlation plot at φ=160°,340° and θ=110°. For the (a) wall reflection at φ=160°, the pattern shows an interference of two plane waves, supporting the TWDP fading assumption. For (b) LOS at φ=340°, we observe a spatial correlation pattern dominated by one wave. The white dashed lines illustrate plane wave phase fronts

7 MC2: Time-gated fading results

To confirm that our observations are not artefacts of our measurement set-up, for example back-lobes of the horn antenna, we now study the wireless channel in the time domain. Our 2-GHz wide measurements from MC2 allow for a time resolution of approximately 0.5 ns. This corresponds to a spatial resolution of 15 cm. We plot the channel impulse responses (CIRs) as a function of distance, namely the LOS excess length Δs, that is:

$$ h(\Delta s) = h((\tau - \tau_{\text{LOS}})c_{0})~. $$

(22)

The scatter-plot of the CIRs for φ=160° is shown in Fig. 14. The LOS CIR at φ=340° is displayed as reference as well. The steerable TX is positioned more than a metre apart from the wall. This amounts in an excess distance of approximately two to three metres. At this excess distance, a cluster of multipath components is present. Note, if the horn antenna points towards the wall, the wave emitted by the back-lobe of the horn antenna is received at zero excess distance. Still, the receive power of the back-lobe is far below the components arriving from the wall reflection. Fading is hence determined by the wall scattering behaviour.

Fig. 14

Scatter-plot of the CIRs. We plot the CIRs as a function of spatial distance, where Δs=0 corresponds to the LOS distance. Our spatial resolution (a channel tap) is 15 cm. The spatial extend of our sampled cube (729 samples) is 2.8λ=1.4 cm, a magnitude smaller than the spatial resolution. The scatter-plot is evaluated at a (a) wall reflection (φ=160°) and at (b) LOS (φ=340°). The mean power is plotted with a continous red line. We observe that the arrival cluster centred at 2.5 m fades very deeply. The gray highlighted region around 2.5 m is further analysed in Fig. 15

The gray-highlighted region of Fig. 14 (top part) shows a reflection cluster that corresponds to the excess distance of the wall reflection. The distributions of each channel tap are represented by a violin plot in Fig. 15. A violin plot illustrates the distribution estimated via Gaussian kernels [88]. Figure 15 clearly demonstrates that the TWDP-decided distributions have multiple modes. The AIC decisions are plotted as markers at the mean power levels.

Fig. 15

Violin-plot of the CIR time-gated for the wall reflection. This figure shows a zoom-in of the gray-highlighted region in Fig. 14. In contrast to Fig. 14, the y-axis is in linear scale. Thereby, the violin plot indicates the distribution at each tap. The marker shows the mean value. The marker style codes best fitting distributions

We evaluated the fading statistic in space for φ=160° at the channel tap corresponding to approximately 2.5 m excess distance. This channel tap is mid in the cluster belonging to the wall reflection. Fig. 16 clearly shows a TWDP fading behaviour, confirmed by AIC.

Fig. 16

CDF: Distribution fitting for spatial measurement data, time gated by the channel tap at 2.5 m, φ=160°, and θ=110°. Note that, similar to the fitting result in Fig. 2, the estimated Rician K-factor is again much smaller than the TWDP K-factor

8 Conclusion

We demonstrate, by means of model selection and hypothesis testing, that TWDP fading explains observed indoor millimetre wave channels. Rician fits of reported studies must be considered with caution. As two exemplary fits, in Figs. 2 and 16, show Rician K-factors tend to be much smaller than their TWDP companions. There is more power in the specular components than is predicted by the Rician fit. The TWDP fading fit accounts for a possible cancellation of two specular waves. Our results are verified through two independent measurement campaigns. For MC1 and MC2 we even used different RF hardware. While MC1 was limited to results in the frequency domain, MC2 allowed a careful study in the spatial-domain and the time-domain.

Having this strong evidence at hand, we claim that the TWDP fading model is more accurate to describe mmWave indoor channels. The flexibility of this model allows furthermore to obtain Rician fading (Δ≡0) and Rayleigh fading (K≡0) results with the same channel model.