A Mobility Model for a 3D Non-Stationary Geometry Cluster-Based Channel Model for High Speed Trains in MIMO Wireless Channels

During channel modeling for high-mobility channels, such as high-speed train (HST) channels, the velocity of the mobile radio station is assumed to be constant. However, this might not be realistic due to the dynamic movement of the train along the track. Therefore, in this paper, an enhanced Gauss–Markov mobility model with a 3D non-stationary geometry based stochastic model (GBSM) for HST in MIMO Wireless Channels is proposed. The non-isotropic scatterers within a cluster are assumed to be around the sphere in which the mobile relay station (MRS) is located. The multi-path components (MPCs) are modeled with varying velocities, whereas the mobility model is a function of time. The MPCs are represented in a death–birth cluster using the Markov process. Furthermore, the channel statistics, i.e., the space-time correlation function, the root-mean-square Doppler shift, and the quasi-stationary interval, are derived from the non-stationary model. The model shows how the quasi-stationary time increases from 0.21 to 0.451 s with a decreasing acceleration of 0.6 to 0.2 m/s2 of the HST. In addition, the impact of the distribution of the angles on the channel statistics is presented. Finally, the simulated results are compared with the measured results. Therefore, there is a close relationship between the proposed model and the measured results, and the model can be used to characterize the channel’s properties.


Introduction
The space-air-ground-sea integrated network (SAGSIN) will focus on more diversified and dynamic communication scenarios, including vehicle-to-vehicle (V2V), high-speed train (HST), unmanned aerial vehicle (UAV), satellite, and maritime communications [1]. This is the wireless communication network's vision for the B5G/6G (Beyond 5G/6G) era. Therefore, accurate and user-friendly channel models that can accurately mimic the underlying characteristics of the B5G/6G channels are essential for the successful design of this communication system. Channel modeling in these dynamic channel scenarios is essential for the evaluation and performance of a wireless communication system before and during implementation. The channel model represents how MPCs in a non-stationary wireless channel propagate in actual scattering situations [2,3]. This is crucial for assessing the effectiveness of communication systems. Some 2D geometry-based stochastic models (GBSM) with non-stationary characteristics have been presented in [4][5][6][7]. However, in 5G and beyond, elevation angles need to be considered. In [8][9][10][11], a 3D channel model for a high-speed train was proposed, and its channel statistics were derived. However, in [8,11], some of these channel models were based on ray tracing, and they only considered distinct dimensions for particular HST environments. In [10], the MPC's non-stationarity was based on tapped delay line, and yet a clustered channel model provides better insights about the channel characteristics. A cluster in this context is a collection of MPCs with identical delay, power, and angle characteristics. The ability to extract intra-cluster and inter-cluster statistics is typically provided by a cluster-based channel architecture, since statistical models with a few parameters, such as Laplacian or Gaussian distributions, can frequently capture the intra-cluster features. This provides more intuitive insights, and the model can further be compared with existing standard channel models.
In [12,13], the term "BS-Visibility Region (VR)" was used to describe the partially visible nature of clusters. The visibility nature of clusters has been mimicked using the birth-death cluster, which is also known as the Markov process. This has been documented in several studies, including [4,14]. In the B5G/6G channel model development, it is generally necessary to find a technique to effectively represent the exact channel characteristics. The non-stationary channel characteristics have to be captured in both domains of space and time [14]. Several channel features in 5G systems are intended to be represented by a GBSM, which is also termed a generic 5G channel model (MG5GCM) [15]. The model is capable of supporting numerous communication scenarios based on the general model structure, but the direction angles (azimuth and elevation), the travel times between the TX/RX, and the scattering components were determined separately. In [4], the model disregarded the non-stationary features in the frequency domain and could only change in the time and array axes. It is challenging for the model to attain spatial consistency because the angles and delays indirectly determine the positions of scatterers.
Some non-stationary GBSMs which have been proposed for HST channels and V2V channels have some characteristics in common, such as high Doppler shifts and temporal non-stationarity [16]. However, these cannot be generalized to all channels with the Doppler effect. For HST conditions where the train velocity could reach 350 km/h, widely used standard channel models are WINNER I [17], WINNER II [18], and IMT-Advanced channel models [19]. The models discussed above are only two-dimensional (2D) and can only be used in situations when the transceiver and scatterers are sufficiently apart. In addition, some of these models were generated using the temporal wide-sense stationary (WSS) assumption, and the cluster dynamics in the time domain were overlooked. By considering the time-varying angles and cluster dynamics, the HST channel model in [2] was proposed based on the IMT-Advanced channel model; however, this is a 2D model, yet elevation angles are important especially when there are ground reflection rays. Although the channel characteristics in B5G/6G systems are intended to be captured by the GBSM [3], the suggested models are generalized for all the dynamic communication scenarios. This can be seen in some of the recent models presented in [12,14,17,20,21].
The elliptical GBSM for the HST channel model which considered varying movement speed and direction was proposed in [4,22]. However, the model was two-dimensional (2D) and can only be used in situations when the transceiver and scatterers are significantly separate from each other. A non-stationary 3D deterministic model based on ray tracing for the HST channel model was proposed in [23]. However, the model could only be utilized for tunnel environments. A 3D GBSM based on a tapped delay line wireless channel model for several HST environments was proposed in [17]. The model's validity depends on the parameters of the models developed from the viaduct and cutting environments acquired during ray tracing. However, the whole procedure results in a high level of computational complexity.
The HST undergoes different trajectories of acceleration and deceleration. This happens during the taking off, change in velocity as it goes through different environments, and stopping at the station. On the other hand, beam misalignment due to dynamic vehicle traffic tends to lower quality-of-service (QoS) performance [24]. Given this trajectory, the HST dynamic velocity and motion direction are usually ignored during channel modeling. For channel models situated in 5G and beyond, mobility modeling specifications in environments with the Doppler effect have to be considered. The d hoc network simulations frequently employ the Gauss-Markov mobility model [25]. This model has lately been used to describe the UAV channel due to its accessibility and efficiency. However, since accelerating or decelerating operations, such as the starting and stopping of HSTs, cannot be modeled using the Gauss-Markov mobility model, an enhanced Gauss-Markov mobility model from [26] was used in this research work to feature the acceleration and deceleration characteristics. The enhanced Gauss-Markov mobility model is simple and available. In this paper, a GBSM for a 5G HST dynamic communication system is proposed. The model was applied in the open-case scenario of the HST. The effect of scatters on the environment was also studied. The primary contributions of this research paper are:

1.
A 3D, mobile and non-stationary cluster-based GBSM with scatterers located around the moving MRS is proposed.

2.
The HST's mobility is described by the enhanced Gauss-Makorv mobility model incorporating acceleration. 3.
The death-birth Markov model is used to model the cluster MPCs. 4.
The channel statistics, i.e., the local space-time correlation function (ST-CF), the rootmean-square Doppler shift spread, and the quasi-stationary intervals, are derived.

5.
The simulated results of the proposed model are compared with the measured results.

A Mobility Model for a 3D, Non-Stationary Cluster and Geometry-Based Channel Model
In this study, the channel impulse response (CIR) is derived from the proposed 3D non-stationary GBSM shown in Figure 1. The sum of sinusoid (SoS) simulation method in [27] corresponding to the model is derived. The MIMO channel with non-isotropic scatterers is considered. The CIR of the complex fading envelope of the MIMO channel  These parameters from Table 1 are considered for the transmission link at any time instant t from the transmitting antenna T p to R q . By considering independent components of Los and NLoS, the following formula can be used to get the time-variant channel impulse response (CIR). In this work, the CIR of the n th cluster is considered to be a Gaussian process that follows the Lindeberg-Levy theorem. By considering this, the CIR equation is deduced to have a power variance of P n (t)/2, with the expectation of 0. The condition is that the cluster number (N) tends to infinity. P n (t) denotes the normalized power which is received on the n th cluster, whereas the carrier wavelength is denoted by λ. The following equations are used to describe the two components ie., where h Los

D
The horizontal distance between the center of the MRS and BS at initial time R (n,m) radius of the sphere around MRS f s (t) half spacing between the two foci of the ellipse δ T ,δ R antenna spacing at the MRS and BS θ T ,θ R orientation of the MRS and BS antenna array in the x − y plane, respectively ϕ T ,ϕ R angles of elevation of the MRS and BS antenna array relative to the x − y plane, respectively v R MRS velocity γ R motion direction of the MRS α Los R , β Los R AAoA and EAoA of the Los path, respectively, S (n,m) The m th scatterer in the n th cluster EAoD of the wave that impinges from effective scatterers S (n,m) of the n th cluster γ S , φ Horizontal and elevation moving direction of scatterers , respectively K Ricean K factors V S Velocity of moving scatterers In this case, τ Los (t) is the propagation delay of the LOS component, and τ n,m (t) designates the propagation delay of the m th resolvable subpath in the n th cluster. The phase angle φ m,n (t) is randomly distributed with a uniform distribution over [−π, π]. The direct link between T x and R x is T p → R q which is the wave traveling distance denoted by ε p,q (t) of LOS Component. The link through the S n,m (t) scatterer is ε p,S n,m (t) + ε S nm,q (t), which represents the of the path distance from the p th antenna element to the S n,m scatterer and the path distance from the S n,m scatterer to the q th antenna element. The derivations of both travel paths are as follows.
The derivation of the travel-distance terms is written as . We can deduce the following relation seen between AoDs and AoAs of the SB rays as from the . The Doppler shifts of the m th LOS sub-path within the n th cluster are represented by f LoS D (t) and f nm,D (t) in (2a) and (2b), respectively. The maximum Doppler shift resulting from the The cluster sub-path number grows towards infinity. Thus, the latency difference between a cluster's sub-paths becomes very minimal. As a result, the envelope of the CIR exhibits a Rayleigh distribution. Considering the locations of the moving vehicle R x and static BS T x , the derived CIR of the LOS link can be taken as a deterministic process. If the variables that change with time are introduced into the reference model, the suggested 3D model is capable of describing the non-stationary properties of the HST propagation channels.

Enhanced Gauss-Markov Mobility Model
In the case of the enhanced Gauss-Markov mobility model from [26],v R (t) represents the velocity of the HST at any time t. Therefore, it is possible to represent the non-stationary velocity of the HST at time instant t as In this case, v R (t 0 ) represents the initial velocity of the HST at t 0 . The time separation and the acceleration of the HST are represented as ∆ t and α T , respectively. At any time instant t − ∆t, we shall have v R (t − ∆ t ). The zero-means Gaussian distribution will represent the variable χ t−∆ t , which is also random. This will have a variance σ v . The range [0, 1] of Υ determines the timeline of the movement of the HST. This movement presents the changing time instant t and the preceding time, especially whenever Υ is zero.

Cluster Process Evolution
In a time-variant (i.e., non-stationary) environment, clusters can only be present for a short while. With time, new taps constantly appear, survive for a while, or "survive", and then finally disappear, or "die". Discrete Markov systems provide a reasonable description for this generation-recombination pattern. The mobility of the MRS in a network design supported by relays is the leading cause of the time variance of the HST communication system. When modeling the clusters of the MPCs, a genetic appearance (birth) and disappearance (death) mechanism has been popularly used to describe changing clusters. A Markov process in [28] was used for representing the time-varying cluster M(t) formation in the suggested model. The time-varying clusters can be represented as The symbols λ G and λ R represent the clusters' rates of birth/generation and death/ recombination. In addition, the fluctuations caused by the movement of the MRS can also be derived as Given that the factor δ t,∆t of the channel fluctuates, the approximation can be represented as The two types of cluster states were analyzed. These were analyzed during the process of development as time changed from t to t + ∆ t , i.e., the remaining clusters at time t and the new clusters being formed. The following formula can be used to determine the cluster's survival likelihood from t to t + ∆t.
where the type of situation factor Dc quantifies the correlation coefficient of the scenario movement for the evolving clusters, which can also be called the clustered birth-death process correlation gap. The survival probability is shown to be dependent on the channel fluctuation factor δ(t, ∆t). The average number of emerging new clusters throughout the same evolution process: where E(.) denotes the expectation.
Consider a scenario where the cluster birth-death process is significantly impacted by the movements of HSTs, and Equations (8) and (9) are functions for the MRS velocity v R . In general, the survival probability of the cluster increases as v R (t) decreases, and vice versa. Within the ∆t period, λG and λR determine the mean maximum numbers of newly developing clusters.

The HST Time-Varying Distances
The minimum distance between the base station and the truck is denoted by D min , and the D proj is the projection distance of D s (t) on the railway truck, which can be derived as in Figure 3 of [9]. By using trigonometry, D s (t) = D proj (t) 2 + h 2 , where h = H BS − H and D proj (t) = D 2 min (t) + D verticle (t). The vertical distance between the base of the base station and the projection of the mobile relay station at time t is denoted by D verticle (t). The antenna is placed on the rooftop of the train.

Method of Equal Volume and the Proposed Sum of the Sinusoidal Simulation Model
The AoAs and AoDs can be integrally computed to determine the channel's statistical properties, and the analytical model assuming infinite rays in each scatterer is utilized. Unfortunately, the design of a channel simulator restricts the realization of integral computing. Therefore, the SoS simulation model is constructed by substituting an integral calculation with a summing calculation to make the channel model more applicable. The sum of sinusoidal simulation (SoS) model from [29] introduces discrete AoAs and AoDs, since the variables of integration in the integral computation are AoAs and AoDs, and the remaining parameters are the same as those in the analytical model. The von Mises distribution has been proposed to generate the AoAs and AoDs. This method's fundamental concept is to use the inverse function of integration to create a set of {α i , β i } that can satisfy the condition This method's comprehensive description is presented in [30]. Based on the MMEV and the suggested model, the corresponding channel characteristics of the SoS simulation model can be produced. However, to formulate the deterministic simulation model, a fixed number of intercluster subpaths with zero phase is assumed. By using the method described of equal areas (MMEV) in [30], the discrete set of α R n,m , β R n,m , R n,m is developed based on these distributions in (11). These new scatterers are consequently dispersed all around the MRS. Then, the other variables, particularly AAoDs and EAoDs, can be determined following the geometric relationship described above. The variables α T n,m and β T n,m can be deduced, as we can derive the following relationship between the AoDs and AoAs of the SB rays: α T n,m ≈ arccos We continue to derive α T n,m and β T n,m from (11) and (12) to understand the development of AAoDs and EAoDs. The basic characteristics of HST channels, such as the time-variant speed, AoA and AoD, power delay profile (PDP), and delay variables, are realized through this development process.

Cluster Delay Update
Once the departure and arrival angles MPC have been defined at a specific moment, The propagation distance for each propagation path from the T x to the R x can be determined. The LOS path's propagation delay is determined using the following equation: In this case, c denotes light speed. In addition, given the m th sub path for n th cluster, the propagation delay is given by τ n,m (t) = ε p,S nm (t) + ε S nm ,q (t))/c (15) The sub-paths of the propagation delay in (14) become generic for both emerging clusters and survival clusters. The cluster center has to be estimated first to determine a cluster's delay. The power means, which is denoted by the K framework, is commonly employed to evaluate the disparities between MPCs, as in [31]. For this model, the difference between AAoAs of intra-cluster and sub-paths traveling through the center is considered to be not greater than π. The total multipath component distance along a specific sub-path for i ∈ M(t) with other intra-cluster sub-paths i is based on an approximation where the multipath component distance angular terms g R n,j and g R n,i are given by g R n,j = [cos β R n,i cos α R n,j , cos β R n,j sin α R n,j , sin β R n,j ] By minimizing (15), the multipath component distances between every sub-path and the rest of the subpaths are determined from minimization. The time of arrival delay of the subpath τ n,j from the n th cluster with a delay of τ n . The cluster's power generation is dependent on the WINNER II channel model [18]. The normalized cluster delays indicated by the symbol τ n are arranged in increasing order.

Cluster Power Update
Cluster powers are computed with the assumption of a power delay profile with a single slope exponential. The parameters provided by the initial IMT-A channel model can be utilized to compute the randomly mean power for the n th cluster, whereas at time t 0 , the power delay is denoted as P n (t 0 ), as in [32]. Therefore, using the time-varying delays determined earlier, the random average powers for the n th cluster at time t may be represented as The exponential distribution is used for the cluster power, and the normal distribution is for the shadowing per cluster.
where σ τ denotes delay spread, Z n is the per cluster shadowing, and r τ is the delay distribution proportionality factor. The sum of the cluster is normalized to one, and thus, the individual cluster powers have to be normalized using The cluster power is modeled in the death-birth process, as in [27], to have a linear increase or decrease towards a threshold value. This is to prevent interruptions during the process. When the cluster power is averaged, the sub-path power is also assumed to be identical. Therefore, the sum of all subpaths is one due to normalization.

Local Space-Time Correlation Function
The normalized local complex space-time correlation function (ST-CF) of the suggested model is described as In this case, E. represents the expectation operator, and the complex conjugate operator is denoted by (.) * . The superposition of the NLOS component and LOS component is used to determine the ST-CF as follows: By applying the corresponding distribution, NLos × e j2π f nm,D(t+∆t)∆t × e j2π f (τ n,m (t)−τ n,m (t+∆t)) (23b) where P NLos = P n (t)P n (t + ∆t), A NLos = ε p,S nm (t) + ε S nm ,q (t) − ε p,S nm (t + ∆t) − ε S nm ,q (t + ∆t),

Channel Time-Variant Transfer Function
The Fourier transform of the time-variant CIR is obtained from the time-variant channel transfer function (TCTF) as where h Los pq (t, f ) and h NLos pq (t, f ) represent the two components i.e., the LOS and NLOS components of the TCTF, defined as

The RMS Delay Spread
The root-mean-square-delay spread is the measure of wideband delay dispersion which is derived as the square root of the PDP channel in [28]: Whereas the summation of the power-weighted averages of the mean Doppler shifts is derived from the mean Doppler sifts. The time-varying mean Doppler shift at time t is denoted by τ 2 , and therefore, For the m th sub-path for n th cluster, the associated power is given by P n,m (t).
While considering a non-WSS channel, at any time t, the mean Doppler shifts τ(t) are used to derive the time-variant. The root mean square-delay spread denoted by σ τ (t) is In the same manner, the azimuth/elevation angular spread of departure and the azimuth/elevation angular spread of arrival can all be calculated with the same procedure.

The Stationary Time Interval
To obtain the quasi-stationary interval, in this work the channel is often classified into the WSS and non-WSS, as in [33]. In this paper, a quasi-stationary channel is assumed when the parameter w, which is also equal to (T s ), is 10%. The mean Doppler shift τ(t) and RMS-DS σ τ (t) are both time variants in non-WSS channels. These are assumed to be constant for a WSS channel. The quasi-stationary time interval is used to evaluate the relative error of the RMS-DS across time. The equation shows that the shortest time interval, Ts, is the quasi-stationary time interval, which is described by

Numerical Simulation and Results
The statistical characteristics of the suggested theoretical and simulation models are analyzed and evaluated in this section. Through numerical simulations, we look into how several important model parameters affect the channel characteristics. In this case, some of the numerical values used are in Table 2. Table 2. Numerical values used for the simulation model.

Parameter
Numerical Value Figure 2 displays the time-variation of the total number of clusters with the Markov process. In this case, the initial cluster number is set to be N(t 0 ) = 20, which is also used in the IMT-A channel model [34]. From the NLoS UMa environment, the following parameters are considered, where the rate for cluster appearance is λ A = 0.8/m and the rate for cluster disappearance is λ D = 0.04/m. Then, P s = 0.3, v R = 60 m/s .  Figure 3 depicts the time-varying absolute value of the TCFs at different instances. It is observed that there is a decrease in time correlation with an increase in time. The suggested model is capable of describing the non-WSS channel, since the time correlation varies for different time instances. Further, another factor that affects the time correlation is changing the angle direction of the MRS and the locations of the scatterers for any given trajectory of the HST. The correlation between different time separations is shown to be time-varying. At time t = 0, it has a higher correlation than at other times due to increasing acceleration.  Figure 4 shows how the parameter of k impacts the AAoA. The absolute value of the TCFs for different k is observed. For example, at t = 0 s and t = 1 s , it is seen that a large k causes a higher channel correlation. The explanation is that as k increases, the intra-cluster scatterers gradually get more densely concentrated in particular directions.    Figure 6 shows the percentage relative error of RMS-DSS for different accelerations α R . In addition, the increasing MRS acceleration affects the quasi-stationary time interval by decreasing it. This is a result of the increasing Doppler shift. This indicates that during take-off of the HST, or increasing the acceleration at any point of the trajectory, the channel is more non-stationary. For example, by taking an increase in acceleration of the HST from 0.2 to 0.3 m/s 2 , it is observed that the quasi-stationary time can decrease by around 0.09s. From Table 3, the quasi-stationary time for varying acceleration is shown. Considering a w equal to 10%, the results are shown in Figure 6.   Figure 7 depicts how the angle parameters affect the stationary interval more than the cluster power. Figure 8 shows the stationary interval of the proposed model, which used the simulation parameters that were chosen based on the measurement configuration in [35] and the IMT-A channel model in [34]. As opposed to the measured and the IMT-A channel model, the stationary interval of the proposed model is equal to 14 ms for 80% and 20 ms for 60%. These are significantly shorter compared to the measured data: 9 ms for 80% and 20 ms for 60%. Then, the original IMT-A channel model produced 22.5 ms for 80% and 38.3 ms for 60%.

Conclusions
A 3D mobility non-stationary cluster-based model has been proposed for the HST channel. The channel statistics, i.e., the space-time correlation function, the root-meansquare Doppler shift, and the quasi-stationary interval were derived from the non-stationary model. The model shows how the quasi-stationary interval decreases with an increasing acceleration of the HST. Compared with the IMT-A advanced channel, the measured data and proposed model's stationary intervals are significantly shorter. The proposed model has a good agreement with measured results. Furthermore, this interval decreases with increasing acceleration of the HST making it important to incorporate the mobility model in the HST channel. Reason being, the trajectory of the train can be dynamic as it moves along the truck gradient during starting and stopping. For future research, the model can be further extended from MIMO to massive MIMO. This can be considered for higher frequencies in the SAGSIN wireless channel modeling application.