A Novel 3D Non-Stationary UAV-MIMO Channel Model and Its Statistical Properties

Qiuming Zhu＊, Kaili Jiang Xiaomin Chen Weizhi Zhong, Ying Yang

1 College of Electronic and Information Engineering, Nanjing University of Aeronautics and Astronautics, Nanjing 211106, China

2 College of Astronautics, Nanjing University of Aeronautics and Astronautics, Nanjing 211106, China

Abstract: The wireless communication systems based on Unmanned Aerial Vehicles(UAVs) have found a wide range of applications recently. In this paper, we propose a new three-dimensional (3D) non-stationary multiple-input multiple-output (MIMO) channel model for the communication links between the UAV and mobile terminal (MT). The new model originates the traditional geometry-based stochastic models (GBSMs) but considers the non-stationary propagation environment due to the rapid movements of the UAV,MT, and clusters. Meanwhile, the upgrade time evolving algorithms of time-variant channel parameters, i.e., the path number based on birth-death processes of clusters, path delays,path powers, and angles of arrival and departure, are developed and optimized. In addition,the statistical properties of proposed GBSM including autocorrelation function (ACF),cross-correlation function (CCF), and Doppler power spectrum density (DPSD) are investigated and analyzed. Simulation results demonstrate that our proposed model provides a good agreement on the statistical properties with the corresponding derived theoretical ones, which indicates its usefulness for the performance evaluation and validation of the UAV based communication systems.

Keywords: unmanned aerial vehicles (UAVs);non-stationary channel models; geometry-based stochastic models (GBSMs); statistical properties

I. INTRODUCTION

Due to the cost reduction and device miniaturization, unmanned aerial vehicles (UAVs)have found a wide range of applications such as environmental monitoring, traffic control, cargo transportation and emergency assistance [1]-[4]. Also, the UAV is an important part of the cognitive Internet of things and cognitive radio network [5]-[7]. Among these applications,the rapid deployment of emergency communication systems combing with multiple-input multiple-output (MIMO) technologies and UAV-based relay systems have gotten increasing attention [8]-[10]. To efficiently design and optimize the communication systems of UAVs with ground base stations (BSs), it is vital to have an accurate and reliable channel model for the underlying UAV-MIMO channels.

There are limited channel models for UAV communications in literatures [11-20], and they can be classified as deterministic models[11], [12], non-geometrical stochastic models (NGSMs) [13], [14] and geometry-based stochastic models (GBSMs) [15]-[20]. By using the ray-tracing method [11] or the finite difference time domain (FDTD) method[12], the deterministic models need a detailed digital map with specific trajectories, which leads to high complexity and limited generality. Based on the measurement data rather than the detailed geometrical information, the NGSMs [13], [14] characterize the channel in a stochastic way, which makes them unable to capture the geometrical-based features of propagation channels.

A novel 3D non-stationary GBSM for UAV-MIMO channels has been proposed in this paper. The proposed model includes the LoS path and NLoS paths and allows the movements of the UAV, MT and clusters.

To guarantee a good tradeoff between the complexity and accuracy, several GBSMs for high altitude platform (HAP) and UAV channels have been proposed [15]-[20]. For example, the authors in [15] presented a 2D stationary GBSM for HAP channels, but they did not consider the influence of elevation angle. The three-dimensional (3D) stationary GBSMs for HAP and UAV channels can be addressed in [16], [17] and [18], [19], respectively. The authors assumed that the scatterers were located on the surface of several cylinders with different distributions such as the uniform distribution [18] or Von Mises (VM) distribution[16], [17], [19]. It should be mentioned that the aforementioned GBSMs [15]-[19] were only designed for the stationary environments,in which all the scatterers are fixed. However,measurement results [20]-[23] have shown that the survival periods of scatterers for both UAV channels [20], [21] and cellular channels[22], [23] are finite and the scatterers could appear (birth) or disappear (death) due to the movements of UAVs and mobile stations.In [24]-[26], the Markov process has been adopted to mimic such a birth and death process for modeling the non-stationary cellular channels. Recently, a 3D non-stationary model for HAP-MIMO channels was proposed in[27]. The authors introduce an M-step 2-state Markov process to model the dynamic evolution of scatters, but lacked of details on the computation methods of channel parameters.This paper aims to fill the above research gap.Overall, the major contributions and novelties of this paper are summarized as follows:

1) This paper develops a new 3D non-stationary GBSM for UAV-MIMO channels. In order to make the new model more universal,we model the angles of arrival (AoAs) by the flexible Von Mises Fisher (VMF) distribution,and take the movements of UAV, ground mobile terminal (MT) and scatterers into account.

2) The computation and optimization methods of time-variant channel parameters, i.e.,the path number, path delay, path power and angle parameters, are proposed and analyzed.Especially, the new method can achieve a smooth transition of path power when the NLoS paths appear or disappear.

3) Statistical properties of the proposed GBSM such as the time-variant autocorrelation function(ACF), cross-correlation function (CCF) and Doppler power spectral density (DPSD) are derived and verified by simulations.

The remainder of this paper is organized as follows. In Section II, a new 3D non-stationary GBSM for UAV-MIMO channels is proposed.The computation and optimization methods of the time evolving channel parameters are presented and analyzed in Section III. The statistical properties of the proposed GBSM are derived in details in Section IV. In Section V,simulation results are given and compared with the corresponding derivation results. Finally,conclusions are given in Section VI.

II. PROPOSED 3D MIMO GBSM FOR THE UAV-GROUND LINKS

Let us consider a 3D non-stationary MIMO channel between the UAV equipped with P antennas and the MT equipped with Q antennas as shown in figure 1. The antenna arrays of both terminals are uniform linear arrays(ULA) and equipped with omni-directional antenna elements. The 3D locations of the pth transmit antenna and the qth receive antenna are denoted by rpand rq, respectively. In figure 1, the line-of-sight (LoS) component and nonline-of-sight (NLoS) components are included and the effective scatterers around the MT are distributed on the surface of several hemispheroids with different radiuses. For the nth NLoS component, the angle parameters such as the azimuth angle of departure (AAoD)αT,n,m(t), elevation angle of departure (EAoD)βT,n,m(t), azimuth angle of arrival (AAoA)αR,n,m(t), and elevation angle of arrival(EAoA) βR,n,m(t) are determined by the scattering cluster Sn. The angle parameters of LoS component such as the AAoD, EAoD, AAoA,and EAoD are denoted by(t) and(t), respectively. It is assumed that the velocities of the UAV, MT and clusters are constant during the short analytical time interval and denoted by vT, vRand vSn.

The small scale-fading channel between the UAV and MT can be represented by a complex matrix. In this paper, the matrix element hqp(t ,τ) denotes the channel impulse response(CIR) between the pth UAV antenna and the qth MT antenna, and it can be expressed by the summation of the LoS component and several NLoS components as

where N(t) is the number of NLoS components, K(t) denotes the Ricean factor,mean the complex channel coefficient, the power and delay of the nth NLoS component, respectively,andrepresent the complex channel coefficient and the delay of the LoS component. It should be mentioned that this paper only considers the power normalized small scale-fading channel model. The total power equals to one, while the power of the LoS component and all NLoS components are K(t) /(K(t)+1)and 1/(K(t)+ 1), respectively.

Fig. 1. Parameters of the proposed 3D GBSM for the UAV-ground links.

For the LoS component, the complex channel coefficient can be further expressed as (2) shown in the bottom at this page,where k=2πf0c0denotes the wave number, f0, c0represent the carrier frequency and the speed of light, and φLoSis the random initial phase.are the spherical unit vector(s) of departure and arrival signal.

For the nth NLoS component, the complex channel coefficient can be further expressed as the superposition of finite rays, i.e., as (3) shown in the bottom at this page, where M denotes the number of rays, and φn,mis the random initial phase.the spherical unit vector of the mth ray within the nth NLoS departure and arrival signal.

It should be noticed that our proposed model allows all the channel parameters such as N(t), K(t), hand τqp,n(t) to be time-variant, which can capture the non-stationary characteristic of the real propagation channels due to the movements of UAVs, MTs and clusters.

III. COMPUTATION AND OPTIMIZATION OF TIME-INVARIANT CHANNEL PARAMETERS

3.1 Path number of valid NLoS paths

Due to the movements of the UAV, MT and scatterers, the propagation environment between the UAV and MT is time dependent.Consequently, the number of valid NLoS paths is time-varying, which means some propagation paths can disappear and some new paths can appear. We adopt the Markov process to model such a birth and death process [24] and the process is determined by a birth rate λBand a death rate λD. Thus, the survival probability of each path is related with the time dependent channel fluctuation ε(t,∆t) and can be calculated by

Fig. 2. Survival probabilities with different Dc.

where Dcdenotes the scenario dependent coefficient, which is related with the space correlation. Typical values of Dccan be chosen with the same order of correlation distances such as 10m, 50m, 100m and 500m. Since ε(t,∆t)is mainly caused by the movement of the UAV εT(t,∆t), the MT εR(t,∆t) and clusters εS(t,∆t) in the time span from t to t+∆t, it can be expressed as

For simplicity, we use the mean velocity of clusters Sn,n =1,...,N denoted by vS(t) and the probability of movement denoted by Pcto represent the fluctuation factor of clusters. Thus,ε(t,∆t) can be further simplified as

The mean number of newly generated clusters at time instant t+∆t can be computed by

In order to demonstrate the birth-death process of NLoS paths in the proposed model, we adopt the measured parameters in [21], where the MT is above the roofs of the buildings.The simulation parameters are as follows,λB=28m, λD=14m, ∆t=1ms, and Pc=0.3.figure 2 depicts the survival probability as a function of the velocity of the UAV for correlation factors Dc=10,50,100,500m. As we can see that the velocity of the UAV increases,the survival probability decreases. Meanwhile,as the correlation factor increases, the survival probability increase as well. Note that increasing Dcbeyond 50m has a negligible effect on the survival probability when the velocity in the range of 150 km/h. figure 3 shows the appearance and disappearance of each NLoS path for different velocities with Dc=100m.It is clearly shown that when the velocity increases, the survival interval shorten as the survival probability decreases, which results in the birth-death process more frequently.

3.2 Path delay and power

During the survival interval, the delay and power of each path are related with the geometrical parameters. The 3D location vectors of the UAV, MT and cluster Snat the time instant t can be calculated by

where LT(t), LR(t) and LSn(t) denote the time-variant 3D location vectors of the UAV,MT and cluster Sn, respectively. Thus, the time-variant distances can be calculated by

where DT,n(t )，DR,n(t) and DLoS(t) denote the time-variant distances between the UAV and cluster Sn, the MT and cluster Sn, and the UAV and MT.

By combing the delays between the UAV and cluster Snand between Snand the MT, the total delay of the nth NLoS path at time instant t can be calculated as

while the delay of LoS path can be calculated as

The power of the nth NLoS path can be calculated by the corresponding path delay as [25]

where rτ,τσ mean the delay scalar and delay spread, and Ynfollows the Gaussian distribution. However, this method would lead to the discontinue path powers when the NLoS path appear or disappear. In order to achieve a smooth transition which is more realistic, in this paper we adopt a filter to smooth the output path power. The upgrade method can be expressed as

where ∗ denotes the convolution symbol and W(t) represents a cosine filter

Fig. 3. Birth and death process of NLoS paths.

where αwmeans the roll off factor and fsis the sampling frequency. Finally, the total power of all NLoS paths should be 1/(K(t)+1), so the path power should be normalized as

To visualize the time evolving process of the path power, we set the same simulation parameters as figure 3 and the rest parameters are configured as follows, M=100,K=5, αR=90°, αT=0°, DLoS(t0) = 150m,βR,n(t0) = random, αR,n(t0) = random and αw=0.8. By applying the traditional method in [28] and our proposed method, figure 4 compares the output time-variant powers for all components. In this figure, the total power is normalized and the ratio of the power of the LoS component and all NLoS components is fixed as 5. It is clearly shown that the power of each NLoS path changes during the survival time interval. For the traditional method in [28], the output power has abrupt change when the path disappear or appear. In contrast,for our proposed method, a smooth transition of output power can be observed at the appearance and disappearance of each path.

3.3 3D angle parameters

The azimuth and elevation angles of AoA and AoD are assumed to be independent and obey the VM distribution [16], [17], [19] and uniform distribution [18], respectively. These assumptions are very helpful to get the closedform expressions of channel statistical properties. However, the authors in [29] have demonstrated that the azimuth and elevation angles are dependent under some scenarios and the jointly VMF distribution fit the measurement data well. To make our channel model more universal, we model the joint distribution of the azimuth angle α and elevation angle β by the VMF distribution, which is defined as

whereα andβ indicate the mean values of azimuth angle and elevation angle, Im/2−1(⋅) denotes the zero-order modified Bessel function and m=3 denotes the VMF distribution, and κ controls the concentration of the distribution.In a scattering environment, as κ→∞ the scattering becomes a point source and when κ=0 the isotropic scattering occurs.

Fig. 4. The time evolving process of path powers with different methods.

The initial angles of all paths can be generated randomly or extracted from the environment dependent measurement data. During the time evolving stage, the mean angles of the time-variant EAoD, EAoA, AAoD and AAoA,denoted by βT,n(t), βR,n(t), αT,n(t) and αR,n(t),are real-time tracked. Let us set(t) andn(t) as the x, y, and z components of the 3D locations LT/R/Sn(t). By using the transform from the Cartesian coordination to spherical coordination. The time-variant EAoD, EAoA and AAoD, AAoA of the nth NLoS path can be calculated by

In order to demonstrate the time evolving process of angle parameters, it is assumed that the AoA and AoD of the nth NLoS path follow the VM distribution with κ=27.3 for the AoA [30] and κ=100 for the AoD. The AoAs and AoDs of the LoS path and the nth NLoS path are tracked by (18)-(21), respectively.The instantaneous normalized PDFs of AoAs and AoDs at three time instants, i.e., t=0s, 3s,5s are given in figure 5(a) and figure 5(b). As we can see that both the AoA and AoD are time-variant.

IV. STATISTICAL PROPERTIES OF THE PROPOSED GBSM

4.1 Time-variant transfer function

The time-variant transfer function (TF) of the channel is defined as the Fourier transform of the time-variant CIR hqp(t ,τ) with respect to the propagation delay can be defined as [32]

where r={δT,δR} is the space lag vector and δT,δRdenote the space lags at the UAV and MT, respectively.

4.2 Time-variant ACF

In this paper, we only consider the underlying channel being non-stationary in the time domain. The normalized 3D spatial-temporal correlation function of the proposed GBSM can be defined as

where (⋅)*denotes the complex conjugation operation.

The effect of survival probability from t to t+∆t should be taken into account and the ACF can be xpressed by setting q1= q2=q and p1= p2= p in (23) as

Fig. 5. The time evolving process of angle parameters.

In the similar way, we can derive the ex-pression of(∆t;t ) and finally obtain the ACFs of the proposed GBSM.

4.3 Time-variant CCF

By substituting ∆t=0 in (23), we can obtain the CCF between two different pairs of antennas as (28), whereshown at the bottom, wheredenote the CCFs at the UAV and MT, respectively. Takeas an example and it can be expressed as (31).

In the similar way, we can derive the expression of rand finally obtain the CCFs of the proposed GBSM.

4.4 Time-variant DPSD

The DPSD of the proposed channel model isalso time-variant and can be defined by the Fourier transform of ACF as

V. SIMULATION RESULTS AND VALIDATION

In order to verify the statistical properties of proposed GBSM, the carrier frequency is 2.4GHz and te AoA and AoD of the nth NLoS path follow the VMF distribution with κ=27.3 for the AoA [30] and κ=100 for the AoD. In the simulation, the random AoAs and AoDs can be generated by following steps.Firstly, by using (18) and (19), we calculated the mean angles of AoD and AoA. Then, by using the generating method in [29], we generated M=100 angles for each NLoS path,which obey the VMF distribution. Finally,we updated all angle values based on the time-variant geometrical relationship. In order to demonstrate the influence of angle parameters on the statistical properties, we assumed that there is only one NLoS path and it exists during the short simulation period. The rest simulation parameters are as follows,, K=10,αS=random, βS=random, βR,n(t0) = random,αR,n(t0) = random and αw=0.8.

By using (24)-(27), the absolute values of the theoretical ACF are calculated and given in figure 6(a). For comparison purpose, the theoretical ACFs at three time instants t = 0s, 1s,and 2s are extracted from figure 6(a) and the corresponding simulated ACFs are shown in figure 6(b). It shows that the simulated ACF of our model fits well to the corresponding theoretical results, which verifies the theoretical derivations as well as simulations. By using(28)-(31), the absolute values of the theoretical CCF are calculated and given by figure 7(a). In this figure, we assume antenna spaces on the UAV and MT are the same, and the axis has been normalized with respect to half wavelength. In addition, comparison between theoretical results and simulated ones at different time instants is presented in figure 7(b).A good consistent between theoretical results and simulation results ensures the correctness of our model and derivations.

Fig. 6. Time-variant absolute value of ACF.

Fig. 7. Time-variant absolute value of CCF.

Fig. 8. Time-variant DPSD.

Fig. 9. Comparison of ACFs between the proposed model and the measurement data.

By using (32), the theoretical DPSDs are calculated and given by figure 8(a). For comparison purpose, the DPSDs of the proposed channel model at three different time instants are illustrated in figure 8(b). The DPSDs drift over time due to the movements of the UAV,MT and clusters can be clearly observed.Again, the simulation results align well with the theoretical results, illustrating the correctness of our model and derivations. From figure 6 - figure 8, we can also see that the statistical properties of our proposed model i.e.,the ACF, CCF and DPSD change over time.It means that our proposed model can reproduce the non-stationary characteristic of the non-stationary channels.

Figure 9 compares the theoretical ACF of the proposed model with the measurement data[25]. The simulation parameters are as follows,shown in figure 9, the excellent agreement between the theoretical result and the measurement data confirms the effectiveness of the proposed model.

VI. CONCLUSIONS

A novel 3D non-stationary GBSM for UAV-MIMO channels has been proposed in this paper. The proposed model includes the LoS path and NLoS paths and allows the movements of the UAV, MT and clusters.Considering the smooth power transition of birth-death NLoS paths, the time evolving algorithms of time-variant channel parameters have been developed and optimized. Based on the proposed non-stationary GBSM, time-variant statistical properties such as the TF, ACF,CCF, and DPSD have also been derived and validated by simulations. Therefore, the new model as well as the theoretical results are very helpful on the performance evaluation and system optimization of the UAV-based wireless communication systems. Although the proposed model cannot be applied on the millimeter-wave (mmWave) communication systems directly, but note that some mmWave channel models [32], [33] adopt the similar frame work of the transitional models such as WINNER+ channel model [28]. In the future work, we will extend the proposed non-stationary GBSM to the mmWave band channels.

ACKNOWLEDGEMENTS

This work was supported by the National Key Scientific Instrument and Equipment Development Project (Grant No. 2013YQ200607),China NSF Grants (Grant No. 61631020),Aeronautical Science Foundation of China(Grant No.2017ZC52021) and Open Foundation for Graduate Innovation of NUAA (Grant No. kfjj20170405 and kfjj20180408).