A Non-Stationary Geometry-Based Cooperative Scattering Channel Model for MIMO Vehicle-to-Vehicle Communication Systems

Abstract

Traditional channel models for vehicle-to-vehicle (V2V) communication usually assume fixed velocity in static scattering environment. In the realistic scenarios, however, time-variant velocity for V2V results in non-stationary statistical properties of wireless channels. Dynamic scatterers with random velocities and directions have been always utilized to depict the non-stationary statistical properties of the channel. In this paper, a non-stationary geometry-based cooperative scattering channel model is proposed for multiple-input multiple-output (MIMO) V2V communication systems, where a birth-death process is used to capture the appearance and disappearance dynamic properties of moving scatterers that reflect the time-variant time correlation and Doppler spectrum characteristics. Moreover, our model has more straight and concise to study the impact of the vehicular traffic density on channel characteristics and thus avoid complicated procedure in deriving the analytical expressions of the channel parameters and functions. The numerical results validate our analysis and demonstrate that setting important parameters of our model can appropriately build up more purposeful measurement campaigns in the future.

1. Introduction

 Vehicle-to-vehicle (V2V) communication has been attracted significant research attention mainly due to be able to improve traffic safety and traffic efficiency in the context of intelligent transportation systems (ITS) over the past several years [1-2]. However, the V2V transmission links are particularly vulnerable to multi-path fading effects due to reflection, scattering, or diffraction by objects in realistic scenarios. For the design and performance evaluation of V2V communication systems, a deep understanding the underlying propagation channel characteristics is required. Channel modeling is an effective solution that provides insights into the real physical attenuation of V2V transmission links [3-5].

 In recent years, many geometry-based channel models have been proposed for multiple-input multiple-output (MIMO) V2V communication systems, where MIMO technology is used to enhance link reliability and improve spectrum efficiency by efficient exploiting multi-path propagation [6]. In [7], an adaptive geometry-based stochastic model was firstly presented for MIMO V2V communication systems, where the space-time-frequency correlation function were derived for a non-isotropic scattering environment and the corresponding three-dimensional (3D) model was proposed in [8]. Subsequently, a geometry-based stochastic two-ring scattering fading model for MIMO mobile-to-mobile (M2M) isotropic channels was proposed in [9], However, road scattering environment cannot be adequately characterized by the two-ring model. Recently, a 3D semi-ellipsoid scattering channel model to depict stationary roadside scenarios was presented in [10], where the channel model relied on the assumption of stationary scatterers, and dynamic scatterers with random velocities and directions are not taken into account on the statistical properties of the channel. In practice, it is not practical for V2V channel modeling that coexist with moving scatterers (e.g., passing vehicles, and walking pedestrians) in vehicular communication scenarios [11], which lead to a larger Doppler shift resulting in aggravating multi-path fading effects [12]. Moreover, the spacing of the transceiver MIMO antenna arrays is not large enough and the scattering of the propagation environment is insufficient, which will result in space-time correlation between the antenna signals [13-14]. Consequently, it is hard to depict the statistical properties of the MIMO V2V channel. Note that cooperative relay diversity technique has been regarded as a powerful approach to combat wireless channel fading and improve system performance [15-19]. For example, the performance of the secondary users in underlay cognitive networks with opportunistic relay selection was studied in [17]. A 3D two-cylinder model for MIMO M2M via stratospheric-relay narrowband channels was proposed in [18]. In [19], a geometry-based scattering model for MIMO V2V relay-based cooperative channels is presented. Unfortunately, the impact of time-variant velocity on channel properties was not carefully considered, which leads to erroneous simulation results in MIMO V2V system evaluation.

 Beyond the impact of dynamic moving scatterers’ considerations, realistic V2V fading channels also exhibit the general non-stationary characteristics, especially in high mobility scenarios [20]. However, in most of the existing works, a common assumption is that the channel fulfills the wide-sense stationary (WSS), and the non-stationarity has been largely neglected [21-23]. Non-stationary is a result caused by time-variant velocity of the transceiver. So far, only few studies have taken the effect of velocity variations into account [24-25]. In this paper, we present a non-stationary three-ring geometry-based cooperative scattering channel model for MIMO V2V communication systems, where the birth-death process is incorporated into the proposed model to capture the appearance and disappearance dynamic properties of moving scattering vehicles on time axes. Moreover, statistical properties of the proposed model are investigated under non-isotropic scattering conditions considering the impact of time-variant velocity.

 The remainder of the paper is organized as follows. Section 2 describes the system model. Section 3 presents the statistical properties of the proposed channel model. Section 4 provides the numerical simulation results and discussion. Finally, the conclusions are drawn in Section 5.

2. SYSTEM MODEL

 A general cooperative MIMO V2V communication model is presented in two-dimension multi-path scattering scenario as shown in Fig. 1, where there are abundant scatterers represented either moving vehicles or stationary scattering objects (e.g., parked cars, trees, buildings, etc.) between the transmitter (Tx) and the receiver (Rx) [26], and the scattering objects are treated as relay cooperative nodes. To simplify the analysis, we focus on the line-of-sight (LOS) and assume that only the single bounced scattering components are considered in non-line-of-sight (NLOS) propagation scenario. This is because that the amplitudes of scatterer’s waves are small, so the practical impact of vehicular scattering is confined to single bounced components [27-28]. In the proposed multi-path propagation channel, the channel is the relay cooperative communication network. Information is delivered from the source node to the destination node via multiple virtual relay nodes. These nodes serve as the information delivery media, which parallel scatterers in physical multi-path channels as signal power reflectors [29-30]. Fig. 2 provides the bijective mapping between multi-path scattering and relay cooperative MIMO V2V communication systems.

Fig. 1. The typical scattering V2V channel exhibits structural including moving and static scatter clusters with LOS and single bounced components.

 

Fig. 2. Analogy between scattering propagation V2V channel model and relay cooperative communication.

 

2.1 Descriptions of the Proposed Model

 Fig. 3 shows the detailed geometry-based three-ring geometric channel model based on relay cooperative scattering for MIMO V2V communication systems. Different from the existing two-ring scattering fading model [9], where static scatterers and moving scatterers cannot be distinguished, the proposed three-ring model is the combination of a two-ring around the Tx and Rx to describe moving scatterers and an ellipse ring to depict stationary roadside environments. Such a model can be easily adapted to a wide variety of scenarios by changing the shape of the scattering region (e.g., ellipse, one-ring, or two-ring) and adjusting model parameters [8]. In the proposed model, transceivers equipped with omnidirectional uniform linear arrays with L_tTx and L_rRx antenna elements are moving with varying velocities. The antenna element spacings at the Tx and Rx are denoted by \(\delta_{T}\) and \(\delta_{R}\), respectively. The angles \(\psi_{T}\) and \(\psi_{R}\) describe the tilt angle of the Tx and Rx antenna arrays relative to the x-axis, respectively. The p-th \(\left(p=1,2, \dots, L_{t}\right)\) transmit and q-th \(\left(q=1,2, \dots, L_{r}\right)\) receive antenna elements are denoted by \(A_{T}^{p}\) and \(A_{R}^{q}\), respectively. A large number of scatterers locate on three-ring, and the scatterers are distributed in a small number of clusters. For the ellipse ring, it is assumed that L effective static scatter clusters are located on an ellipse with the Tx and the Rx at its focal points. The l-th \((l=1,2, \dots, L)\) cluster is denoted by \(S_{E, l}\). The half lengths of the major and minor axes of the ellipse are denoted as a and b, respectively. The distance between the focal points is \(D=2 f\). For the two-ring, supposing there are M effective moving scatter clusters lying on a ring of radius R_t around the Tx side and the m-th \((m=1,2, \dots, M)\) cluster is denoted by \(S_{T, m}\). Similarly, N effective scatter clusters lie on a ring of radius R_r around the Rx side and the n-th \((n=1,2, \dots, N)\) cluster is denoted by \(S_{R, n}\). Furthermore, we can assume that each cluster is composed of several scattering rays which have the same delays and angle parameters. The total number of scattering rays for each cluster is N₁ for the Tx ring, N₂ for the Rx ring, N₃ and for the ellipse ring, respectively. The n₁-th scatterer for the m-th Tx cluster is denoted by \(​​S_{R, mn_{1}}\), the n₂-th scatterer for the n-th Rx cluster is denoted by \(S_{R, nn_{2}}\) and the n3-th scatterer for the l-th ellipse cluster is denoted by \(S_{E, ln_{3}}\). The angle of departure (AoD) and the angle of arrival (AoA) for the LOS component are defined as \(a_{T}^{L o S}\) and \(\beta_{R}^{\operatorname{los}}\), respectively. Parameters \(a_{T, m}\left(\beta_{T, m}\right), \quad a_{R, n}\left(\beta_{R, n}\right)\) and \(a_{E, l}\left(\beta_{E, l}\right)\) are the mean AoDs (AoAs) for the waves that impinge on the clusters \(S_{T, m}, S_{R, n}\) and \(S_{E, l}\), respectively. Parameters \(\Delta a_{T, m n_{1}}\left(\Delta \beta_{T, m n_{1}}\right)\), \(\Delta a_{R, n n_{2}}\left(\Delta \beta_{R, n n_{2}}\right)\) and \(\Delta a_{E, \ln _{3}}\left(\Delta \beta_{E, \ln _{3}}\right)\) are the relative azimuth AoDs (AoAs) offsets caused by scatterer \(S_{T, m n_{1}}, S_{R, n n_{2}}\) and \(S_{E, l n_{3}}\), respectively. Meanwhile, the reasonable assumptions \(\max \left\{R_{t}, R_{r}\right\} \ll D\) and \(\max \left\{\delta_{T}, \delta_{R}\right\} \ll \min \left\{R_{t}, R_{r}, a-f\right\}\) are utilized in this proposed model [12]. Although there exists situations when Tx, Rx, and scatter clusters are all moving, it should be noticed that those situations can be equivalent to a static Tx situation using the principles of relative motion [31]. We can assume that all moving scatter clusters are allowed to move in a relative direction γ_s with equal velocities v_s at a certain probability p_s, and the γ_s follows a uniform distribution. Besides, the Rx moves in a relative direction γ_R with a time-variant velocity of v_R(t) at time instant t.

 

Fig. 3. Proposed geometry-based cooperative scattering channel model for MIMO V2Vcommunication systems

2.2 Time-Variant Channel Impulse Response (CIR)

 As shown in Fig. 3, note that some waves from the Tx antennas may directly travel to the Rx antennas (i.e., LOS mode), whereas other waves are single bounced scattering from the clusters located around Tx before arriving at Rx (i.e., SBT mode), single bounced scattering from the clusters located around Rx before arriving at Rx (i.e., SBR mode), and single bounced scattering from the clusters located around the ellipse ring before arriving at Rx (i.e., SBE mode), respectively. Therefore, the proposed V2V channel model can be described by matrix \(H(t)=\left[h_{p q}(t)\right]_{L_{t} \times L_{r}}\) of size\(L_{t} \times L_{r}\). The CIR \(h_{p q}(t)\) between \(A_{T}^{p}\) and \(A_{R}^{q}\) can be calculated as

\(\begin{equation}\begin{aligned} h_{p q}(t) &=h_{p q}^{{LoS}}(t)+h_{p q}^{N L o S}(t) \\ &=h_{p q}^{L o S}(t)+h_{p q}^{S B T}(t)+h_{p q}^{S B R}(t)+h_{p q}^{S B E}(t) \end{aligned}\end{equation}\)       (1)

where the respective time-variant channel complex gains can be modeled as

\(\begin{equation}h_{p q}^{L o S}(t)=\sqrt{\frac{K_{p q} P_{p q}}{K_{p q}+1}} \exp \left(j\left(2 \pi \varphi_{L o S}(t)-k_{0} \varepsilon_{p q}\right)\right)\end{equation}\)       (2)

\(\begin{equation}h_{p q}^{S B T}(t)=\sqrt{\frac{\eta_{_{T}} P_{p q}}{M\left(K_{p q}+1\right)}} \sum_{m=1}^{M} \lim _{N_{1} \rightarrow \infty} \sqrt{\frac{1}{N_{1}}} \sum_{n_{1}=1}^{N_{1}} \exp \left(j 2 \pi \varphi_{m n_{1}}(t)+j\left(\vartheta_{mn_{1}}-k_{0} \varepsilon_{p q, m n_{1}}\right)\right)\end{equation}\)       (3)

\(\begin{equation}h_{p q}^{S B T}(t)=\sqrt{\frac{\eta_{_{R}} P_{p q}}{N\left(K_{p q}+1\right)}} \sum_{m=1}^{N} \lim _{N_{2} \rightarrow \infty} \sqrt{\frac{1}{N_{2}}} \sum_{n_{2}=1}^{N_{2}} \exp \left(j 2 \pi \varphi_{n n_{2}}(t)+j\left(\vartheta_{nn_{2}}-k_{0} \varepsilon_{p q, n n_{2}}\right)\right)\end{equation}\)       (4)

\(\begin{equation}h_{p q}^{S B E}(t)=\sqrt{\frac{\eta_{_E} P_{p q}}{L\left(K_{p q}+1\right)}} \sum_{l=1}^{L} \lim _{N_{3} \rightarrow \infty} \sqrt{\frac{1}{N_{3}}} \sum_{n_{3}=1}^{N_{3}} \exp \left(j 2 \pi \varphi_{ln_{3}}(t)+j\left(\vartheta_{ln_{3}}-k_{0} \varepsilon_{p q, ln_{3}}\right)\right)\end{equation}\)       (5)

where P_pq is the total power of \(A_{T}^{p}-A_{R}^{q}\) link and K_pq is the Rician factor. Here, energy-related parameters η_T ,η_R and η_E specify the numbers of NLOS components of moving scatter clusters and static scatter clusters, respectively, which can be normalized to satisfy 1 η_T+η_R+η_E=1. Note that the value of K_pq can depict vehicular traffic density to a certain extent. For a high vehicular traffic density, the value of K_pq is smaller than that in the low vehicular traffic density scenario due to dense moving vehicles, also the NLOS components bear a significant amount of power [7]. The initial phases \(\vartheta_{m n_{1}}, \vartheta_{n n_{2}}\) and \(\vartheta_{l n_{3}}\) are independent and identically distributed (i.i.d.) random variables with uniform distributions over \([0,2 \pi) . k_{0}=2 \pi / \lambda\) is the free-space wave number, where λ is carrier wavelength，which can be expressed in terms of the speed of light c and the carrier frequency f_c as \(\lambda=c / f_{c} \cdot \varepsilon_{p q}, \varepsilon_{p q, m n_{1}}, \varepsilon_{p q, n n_{2}}\) and \(\varepsilon_{p q, ln _{3}}\) are travel distances of the waves through the link \(A_{T}^{p}-A_{R}^{q}, A_{T}^{p}-S_{T, m n_{1}}-A_{R}^{q}, A_{T}^{p}-S_{R, n n_{2}}-A_{R}^{q}\) and \(A_{T}^{p}-S_{E, l n_{3}}-A_{R}^{q}\), respectively. Note that \(\max \left\{\delta_{T}, \delta_{R}\right\} \ll \min \left\{R_{t}, R_{r}\right\}\) and \(\sqrt{1+x} \approx 1+x / 2(x \ll 1)\). It is assumed that the lengths of the RX trajectory are small compared to the radii of the two-ring over the observation time. According to the assumption in [18], [32], the travel distances can be calculated as

\(\begin{equation}\varepsilon_{p q} \approx D-k_{p} \delta_{T} \cos \left(\alpha_{T}^{L o S}-\psi_{T}\right)-k_{q} \delta_{R} \cos \left(\beta_{R}^{L o S}-\psi_{R}\right)\end{equation}\)       (6)

\(\begin{equation}\varepsilon_{p q, m n_{1}} \approx R_{r}-k_{p} \delta_{T} \cos \left(\alpha_{T, m}+\Delta a_{T, m n_{1}}-\psi_{T}\right)+D-k_{q} \delta_{R}\left(\frac{R_{t}}{D} \sin \psi_{R} \sin \left(\alpha_{T, m}+\Delta a_{T, m n_{1}}\right)-\cos \psi_{R}\right)\end{equation}\)       (7)

\(\begin{equation}\varepsilon_{p q, nn_{2}} \approx D-k_{q} \delta_{T}\left(\frac{R_{r}}{D} \sin \psi_{T} \sin \left(\beta_{R, n}+\Delta \beta_{R, nn_{2}}\right)+\cos \psi_{T}\right)+R_{r}-k_{q} \delta_{R} \cos \left(\beta_{R, n}+\Delta \beta_{R, nn_{2}}-\psi_{R}\right)\end{equation}\)       (8)

\(\begin{equation}\varepsilon_{p q, ln_{3}} \approx 2 a-k_{p} \delta_{T} \cos \left(a_{E, l}+\Delta a_{E, ln_{3}}-\psi_{T}\right)-k_{q} \delta_{R} \cos \left(\beta_{E, l}+\beta_{E, l n_{3}}-\psi_{R}\right)\end{equation}\)       (9)

where \(k_{p}=\left(L_{t}+1-2 p\right) / 2, k_{q}=\left(L_{r}+1-2 q\right) / 2, \alpha_{T}^{L o S} \approx 0\) and \(\beta_{R}^{\operatorname{los}} \approx \pi\)

The non-stationarity of the channel model is represented by the time-variant Doppler shifts parameters and the scatter clusters time evolution. According to the result in [21], the time-varying Doppler shifts for this model can be written as

\(\begin{equation}\varphi_{L o S}(t)=\int_{0}^{t} \frac{f_{c}}{c} v_{R}(t) \cos \left(\beta_{R}^{L o S}-\gamma_{R}\right) d t\end{equation}\)       (10)

\(\begin{equation}\begin{aligned} \varphi_{m n_{1}}(t)=& \int_{0}^{t} \frac{f_{c}}{c}\left[v_{R}(t)\left(\frac{R_{t}}{D} \sin \gamma_{R} \sin \left(\alpha_{T, m}+\Delta \alpha_{T, m n_{1}}\right)-\cos \gamma_{R}\right)\right) \\ &\left.-v_{s}\left(\cos \left(\alpha_{T, m}+\Delta \alpha_{T, m n_{1}}-\gamma_{s}\right)+\cos \left(\beta_{T, m}+\Delta \beta_{T, m n_{1}}-\gamma_{s}\right)\right)\right] d t \end{aligned}\end{equation}\)       (11)

\(\begin{equation}\begin{aligned} \varphi_{m n_{2}}(t)=& \int_{0}^{t} \frac{f_{c}}{c}\left[v_{R}(t) \cos \left(\beta_{R, n}+\Delta \beta_{R, nn_{2}}-\gamma_{R}\right)\right.\\ &\left.-v_{s}\left(\cos \left(\alpha_{R, n}+\Delta \alpha_{R, nn_{2}}-\gamma_{s}\right)+\cos \left(\beta_{R, n}+\Delta \beta_{R, nn_{2}}-\gamma_{s}\right)\right)\right] d t \end{aligned}\end{equation}\)       (12)

\(\begin{equation}\varphi_{l n_{3}}(t)=\int_{0}^{t} \frac{f_{c}}{c} v_{R}(t) \cos \left(\beta_{l n_{3}}+\Delta \beta_{E, l n_{3}}-\gamma_{R}\right) d t\end{equation}\)       (13)

To model the time-varying velocity of \(v_{R}(t)\), we use the following expression:

\(\begin{equation}v_{R}(t)=v_{R}\left(t_{0}\right)+a_{0} t\end{equation}\)       (14)

where the constant quantity \(v_{R}\left(t_{0}\right)\) is represented the speed of Rx motion at time t = 0, and the constant quantity a₀ is called the speed acceleration.

 In our proposed model, using the results in [27] and [29], we can describe the \(a_{E, l n_{3}}\) in terms of the \(\beta_{E, l n_{3}}\), as follows

\(\begin{equation}a_{E, l n_{3}}=\left\{\begin{array}{cc} f\left(\beta_{E, ln_{3}}\right) & { if }\ 0<\beta_{E, ln _{3}} \leq \beta_{0} \\ f\left(\beta_{E, ln _{3}}\right)+\pi & { if }\ \beta_{0}<\beta_{E, ln _{3}} \leq 2 \pi-\beta_{0} \\ f\left(\beta_{E, ln _{3}}\right)+2 \pi & { if } \ 2 \pi-\beta_{0}<\beta_{E, ln _{3}} \leq 2 \pi \end{array}\right.\end{equation}\)       (16)

where

\(\begin{equation}f\left(\beta_{E, l n_{3}}\right)=\arctan \left[\frac{\left(\kappa^{2}-1\right) \sin \left(\beta_{E, ln_{13}}\right)}{2 \kappa+\left(\kappa^{2}+1\right) \cos \left(\beta_{E, l n_{3}}\right)}\right], \quad \beta_{E, ln_{3}}=\beta_{E, l}+\Delta \beta_{E, l n_{3}}\end{equation}\)       (16)

\(\begin{equation}\beta_{0}=\pi-\arctan \left(\frac{\kappa^{2}-1}{2 \kappa}\right), \kappa=\frac{a}{f}\end{equation}\)       (17)

Note that, since the direction of movement and location of clusters are not known a priori, when the number of scatterers in each cluster is approaching infinity \(\left(N_{1}, N_{2}, N_{3} \rightarrow \infty\right)\), these discrete independent random variables \(\left(\Delta \alpha_{T, m n_{1}}, \Delta \beta_{R, n n_{2}}, \text { and } \Delta \beta_{E, \ln _{3}}\right)\) can be characterized as continuous random variables \(\Delta \alpha_{T, m}, \Delta \beta_{R, n}\), and \(\Delta \beta_{E, l}\), with given probability density functions (pdfs) \(p\left(\Delta \alpha_{T, m}\right), p\left(\Delta \beta_{R, n}\right)\) and \(p\left(\Delta \beta_{E, l}\right)\), respectively. In this paper, we use the von Mises PDFs to characterize scatterers distributions because it approximates uniform, Gaussian, Laplacian distributions and admits closed-form solutions for many useful situations. The von Mises pdfs can be defined by [32]

\(\begin{equation}p\left(\Delta \alpha_{T, m}\right)=\frac{\exp \left[k_{T} \cos \left(\Delta \alpha_{T, m}\right)\right]}{2 \pi I_{0}\left(k_{T}\right)}\end{equation}\)       (18)

\(\begin{equation}p\left(\Delta \beta_{R, n}\right)=\frac{\exp \left[k_{R} \cos \left(\Delta \beta_{R, n}\right)\right]}{2 \pi I_{0}\left(k_{R}\right)}\end{equation}\)       (19)

\(\begin{equation}p\left(\Delta \beta_{E, l}\right)=\frac{\exp \left[k_{e l} \cos \left(\Delta \beta_{E, l}\right)\right]}{2 \pi I_{0}\left(k_{e l}\right)}\end{equation}\)       (20)

where \(I_{0}(\cdot)\) is the zero-order modified Bessel function of the first kind, \(\Delta \alpha_{T, m} \in[-\pi, \pi), \Delta \beta_{R, n} \in[-\pi, \pi), \Delta \beta_{E, l} \in[-\pi, \pi)\), and k_T, k_R and k_el are real-valued concentration parameters that indicate the spreads of the scatterers in each cluster. Especially when \(k_{T / \mathrm{R} / \mathrm{el}}=0\), the PDFs are approaching to the uniform distribution, and as \(k_{T / \mathrm{R}/{\mathrm{el}}}\) increases, the PDFs become increasingly non-isotropic.

 In the non-stationary scenario, the scatter clusters only exist over a certain time period. A birth-death process [34-36] is adopted to model the clusters’ appearance and disappearance. The time dependent channel fluctuations caused by the movement of Rx, i.e., \(\delta_{R}(t, \tau)\) and the movement of clusters, i.e., \(\delta_{s}(t, \tau)\) in the time span between t and t +τ can be expressed as

\(\begin{equation}\delta_{h}(t, \tau)=\delta_{R}(t, \tau)+\delta_{s}(t, \tau)\end{equation}\)       (21)

 with

\(\begin{equation}\delta_{R}(t, \tau)=\int_{t}^{t+\tau} v_{R}(t) d t=\left(v_{R}\left(t_{0}\right)+a_{0} t\right) \tau\end{equation}\)       (22)

 and

\(\begin{equation}\delta_{s}(t, \tau)=\int_{t}^{t+\tau} p_{s} v_{s} d t=p_{s} v_{s} \tau\end{equation}\)       (23)

where p_s is the probability of scatter clusters movement.

 The scatter clusters in time instant t +τ can be assumed as the sum of the surviving scatter clusters from time instant t and the newly born scatter clusters during the time interval τ . Hence, the total number of scatter clusters (M , N , and L) are time-variant revealing the non-stationarity of the V2V communication channels. The birth-death process is controlled by the generation rate of clusters λ_G and the recombination rate of clusters λ_R .The probabilities of clusters at t +τ survived from t can be modeled as

\(\begin{equation}P_{{sunvivl}}(t, \tau)=e^{-\lambda_{R} \cdot \frac{\delta_{h}(t, \tau)}{D_{c}}}\end{equation}\)       (24)

where D_c is the scenario dependent correlation factor. According to the birth-death process, the duration between clusters appearance and disappearance are exponentially distributed [35], the expectation of the number of newly generated clusters \(E\left[N_{n e w}(t+\tau)\right]\) based on the Poisson process can be calculated as [36]

\(\begin{equation}E\left[N_{n e w}(t+\tau)\right]= \frac{\lambda_{G}}{\lambda_{R}} ({1-e^{-\lambda_{R} \cdot \frac{\delta_{h}(t, \tau)}{D_{c}}}}) \end{equation}\)       (25)

3. Statistical Properties of the Theoretical Channel Model

 In this section, key statistical properties of the model, e.g., the spatial cross-correlation function (CCF), temporal autocorrelation function (ACF), Wigner-Ville spectrum, and ergodic capacity under non-isotropic scattering conditions considering the impact of velocity variations will be derived.

 The normalized spatial-temporal correlation function (ST-CF) between two complex channel impulse response \(h_{p q}(t)\) and \(h_{p^{\prime} q^{\prime}}(t+\tau)\) is defined as

\(\begin{equation}\rho_{p q, p^{\prime} q^{\prime}}\left(\delta_{T}, \delta_{R}, t, \tau\right)=\frac{E\left\{h_{p q}(t) h_{p^{\prime} q^{\prime}}^{*}(t+\tau)\right\}}{\sqrt{E\left\{\left|h_{p q}(t)\right|^{2}\right\}} \sqrt{E\left\{\left|h^*_{p^{\prime} q^{\prime}}(t+\tau)\right|^{2}\right\}}}\end{equation}\)       (26)

where E[·] is the statistical expectation operation , (·)^∗ denotes the complex conjugate operation, \(p, p^{\prime} \in\left\{1,2, \cdots L_{t}\right\}\) and \(q, q^{\prime} \in\left\{1,2, \cdots L_{r}\right\}\). Since the LOS, single bounced components are independent of each other, (26) can be separated into four terms as

\(\begin{aligned} \rho_{p q, p^{\prime} q^{\prime}}\left(\delta_{T}, \delta_{R}, t, \tau\right) &=\rho_{p q, p^{\prime} q^{\prime}}^{L o S}\left(\delta_{T}, \delta_{R}, t, \tau\right)+\rho_{p q, p^{\prime} q^{\prime}}^{S B T}\left(\delta_{T}, \delta_{R}, t, \tau\right) \\ &+\rho_{p q, p^{\prime} q^{\prime}}^{S B R}\left(\delta_{T}, \delta_{R}, t, \tau\right)+\rho_{p q, p^{\prime} q^{\prime}}^{S B E}\left(\delta_{T}, \delta_{R}, t, \tau\right) \end{aligned}\)       (27)

where \(\rho_{p q, p^{\prime} q^{\prime}}^{L o S}\left(\delta_{T}, \delta_{R}, t, \tau\right), \rho_{p q, p^{\prime} q^{\prime}}^{S B T}\left(\delta_{T}, \delta_{R}, t, \tau\right), \rho_{p q, p^{\prime} q^{\prime}}^{S B R}\left(\delta_{T}, \delta_{R}, t, \tau\right)\), and \(\rho_{p q, p^{\prime} q^{\prime}}^{S B E}\left(\delta_{T}, \delta_{R}, t, \tau\right)\) denote the normalized ST-CFs of the LOS, single bounced caused by moving and static scatter clusters components at time t, respectively.

 According to the result in [8] and [12]，and by applying the corresponding von Mises pdfs, the normalized ST-CFs can be calculated as

\(\rho_{p q, p^{\prime} q^{\prime}}^{L o S}\left(\delta_{T}, \delta_{R}, t, \tau\right)=K c_{p p^{\prime}}^{L o S}\left(\delta_{T}\right) c_{q q^{\prime}}^{L o S}\left(\delta_{R}\right) e^{j 2 \pi\left(\varphi_{L o s}(t)-\varphi_{L o s}(t+\tau)\right)}\)       (28)

\(\rho_{p q, p^{\prime} q^{\prime}}^{S B T}\left(\delta_{T}, \delta_{R}, t, \tau\right)=\frac{\eta_{_T}}{M} \sum_{m=1}^{M} \int_{-\pi}^{\pi} \int_{-\pi}^{\pi} e^{-j k_{0}\left(\varepsilon_{pq, m n1}-\varepsilon_{p^{\prime}q^{\prime} ,m n1}\right)} e^{j 2 \pi j\left(\varphi_{m n1}(t)-\varphi_{m n1}(t+\tau)\right)} p\left(\gamma_{s}\right) p\left(\Delta \alpha_{T, m}\right) d \gamma_{s} d \Delta \alpha_{T, m}\)       (29)

\(\rho_{p q, p^{\prime} q^{\prime}}^{S B R}\left(\delta_{T}, \delta_{R}, t, \tau\right)=\frac{\eta_{_R}}{N} \sum_{n=1}^{N} \int_{-\pi}^{\pi} \int_{-\pi}^{\pi} e^{-j t_{0}\left(\varepsilon_{pq, nn_{2}}-\varepsilon_{p' q',nn_{2}}\right)} e^{j 2 \pi j\left(\varphi_{nn_{2}}(t)-\varphi_{nn_{2}}(t+\tau)\right)} p\left(\gamma_{s}\right) p\left(\Delta \beta_{R, n}\right) d \gamma_{s} d \Delta \beta_{R, n}\)       (30)

\(\rho_{p q, p^{\prime} q^{\prime}}^{S B E}\left(\delta_{T}, \delta_{R}, t, \tau\right)=\frac{\eta_{_E}}{L} \sum_{l=1}^{L} \int_{-\pi}^{\pi} c_{p p^{\prime}}^{S B E}\left(\delta_{T}, \Delta \alpha_{E, k}\right) c_{q q^{\prime}}^{S B E}\left(\delta_{R}, \Delta \beta_{E, l}\right) e^{j 2 \pi\left(\varphi_{ln_3}(t)-\varphi_{ln_3}(t+\tau)\right)} p\left(\Delta \beta_{E, l}\right) d \Delta \beta_{E, l}\)       (31)

where \(c_{p p^{\prime}}^{L o S}\left(\delta_{T}\right)=e^{j k_{0} \delta_{T}\left(p-p^{\prime}\right) \cos \psi_{T}}=e^{-j k_{0} \delta_{R}\left(q-q^{\prime}\right) \cos \psi_{R}} c_{p p^{\prime}}^{S B E}\left(\delta_{T}, \Delta \alpha_{E, l}\right)=e^{j k_{0} \delta_{T}\left(p-p^{\prime}\right) \cos \left(\alpha_{E, l}+\Delta \alpha_{E, l}-\psi_{T}\right)},c_{q q^{\prime}}^{S B E}\left(\delta_{R}, \Delta \beta_{E, k}\right)\)\(=e^{j k_{0} \delta_{R}\left(q-q^{\prime}\right) \cos \left(\beta_{E, k}+\Delta \beta_{E, k}-\psi_{R}\right)}, p\left(\gamma_{s}\right)=1 / 2 \pi\).

3.1 Spatial CCF

 The spatial CCF \(\rho_{p q, p^{\prime} q^{\prime}}\left(\delta_{T}, \delta_{R}, t\right)\) is equal to the ST-CFs \(\rho_{p q, p^{\prime} q^{\prime}}\left(\delta_{T}, \delta_{R}, \tau, t\right)\) in (27) by setting τ = 0 . Therefore, the spatial CCF can be obtained as

\(\rho_{p q, p^{\prime} q^{\prime}}^{L o S}\left(\delta_{T}, \delta_{R}, t\right)=K c_{p p^{\prime}}^{L o S}\left(\delta_{T}\right) c_{q q^{\prime}}^{L o S}\left(\delta_{R}\right)\)       (32)

\(\rho_{p q, p^{\prime} q^{\prime}}^{S B T}\left(\delta_{T}, \delta_{R}, t\right)=\frac{\eta_{_T}}{M} \sum_{m=1}^{M} \int_{-\pi}^{\pi} e^{-j k_{0}\left(\varepsilon_{p q, m n1}-\varepsilon_{p' q' ,m n1}\right)} p\left(\Delta \alpha_{T, m}\right) d \Delta \alpha_{T, m}\)       (33)

\(\rho_{p q, p^{\prime} q^{\prime}}^{S B R}\left(\delta_{T}, \delta_{R}, t, \tau\right)=\frac{\eta_{R}}{N} \sum_{n=1}^{N} \int_{-\pi}^{\pi} e^{-j k_{0}\left(\varepsilon_{p q, n n2}-\varepsilon_{p^{\prime} q^{\prime}, nn_{2}}\right)} p\left(\Delta \beta_{R, n}\right) d \Delta \beta_{R, n}\)       (34)

\(\rho_{p q, p^{\prime} q^{\prime}}^{S B E}\left(\delta_{T}, \delta_{R}, t\right)=\frac{\eta_{E}}{L} \sum_{l=1}^{L} \int_{-\pi}^{\pi} c_{p p^{\prime}}^{S B E}\left(\delta_{T}, \Delta \alpha_{E, l}\right) c_{q q^{\prime}}^{S B E}\left(\delta_{R}, \Delta \beta_{E, l}\right) p\left(\Delta \beta_{E, l}\right) d \Delta \beta_{E, l}\)       (35)

 The integral in (35) has to be solved numerically. However, by assuming that \(R_{r} \sin \psi_{R} \sin \left(\alpha_{T, m}+\Delta \alpha_{T, m}\right) / D\) and \(R_{t} \sin \psi_{T} \sin \left(\beta_{R, n}+\Delta \beta_{R, n}\right) / D\) are small values, by using the equality \(\int_{-\pi}^{\pi} \exp (a \sin c+b \cos c) d c=2 \pi I_{0}(\sqrt{a^{2}+b^{2}})\) and trigonometric transformations, the integral in (33) and (34) can be obtained the closed-form solutions for (36) and (37) , respectively.

\(\rho_{p q, p^{\prime} q^{\prime}}^{S B T}\left(\delta_{T}, \delta_{R}, t\right) \approx \frac{\eta_{T}}{I_{0}\left(k_{T}\right) M} \sum_{m=1}^{M} e^{-j k_{0}\left(q^{\prime}-q\right) \delta_{R} \cos \psi_{R}} I_{0}(\sqrt{x_{m_{1}}^{2}+y_{m_{1}}^{2}})\)       (36)

\(\rho_{p q, p^{\prime} q^{\prime}}^{S B R}\left(\delta_{T}, \delta_{R}, t\right) \approx \frac{\eta_{R}}{I_{0}\left(k_{R}\right) N} \sum_{n=1}^{N} e^{-j k_{0}\left(p^{\prime}-p\right) \delta_{T} \cos \psi_{T}} I_{0}(\sqrt{x_{m_{2}}^{2}+y_{m_{2}}^{2}})\)       (37)

where

\(x_{m_{1}}=j k_{0}\left(p^{\prime}-p\right) \delta_{T} \cos \left(\psi_{T}-\alpha_{T, m}\right)+j k_{0}\left(q^{\prime}-q\right) \delta_{R} R_{t} / D \sin \psi_{R} \sin \alpha_{T, m}+k_{T}\)       (38)

\(y_{m_{1}}=j k_{0}\left(p^{\prime}-p\right) \delta_{T} \sin \left(\psi_{T}-\alpha_{T, m}\right)+j k_{0}\left(q^{\prime}-q\right) \delta_{R} R_{t} / D \sin \psi_{R} \cos \alpha_{T, m}\)       (39)

\(x_{m_{2}}=j k_{0}\left(q^{\prime}-q\right) \delta_{R} \cos \left(\psi_{R}-\beta_{R, n}\right)+j k_{0}\left(p^{\prime}-p\right) \delta_{T} R_{r} / D \sin \psi_{T} \sin \beta_{R, n}+k_{R}\)       (40)

\(y_{m_{1}}=j k_{0}\left(q^{\prime}-q\right) \delta_{R} \sin \left(\psi_{R}-\beta_{R, n}\right)+j k_{0}\left(p^{\prime}-p\right) \delta_{T} R_{r} / D \sin \psi_{T} \cos \beta_{R, n}\)       (41)

3.2 Temporal ACF

The temporal ACF \(r_{p q}(t, \tau)\) is defined as \(r_{p q}(t, \tau)=E\left\{h(t) h^{*}(t+\tau)\right\}\), which is equal to the \(\rho_{p q, p^{\prime} q^{\prime}}\left(\delta_{T}, \delta_{R}, t, \tau\right)\) in (27) by setting \(\delta_{T}=0\) and \(\delta_{R}=0\). It has to be mentioned that the temporal ACF should consider the channel fluctuations characteristics caused by the appearance and disappearance dynamic properties of moving scatterers during the time interval τ .Therefore, the temporal ACF can be obtained as

\(r_{p q}(t, \tau)=P_{ {survival}}(t, \tau) \cdot\left(r_{p q}^{{Los}}(t, \tau)+r_{p q}^{S B T}(t, \tau)+r_{p q}^{S B R}(t, \tau)+r_{p q}^{S B E}(t, \tau)\right)\)       (42)

where

\(r_{p q}^{\operatorname{Los}}(\tau, t)=K e^{j 2 \pi\left(\varphi_{{LoS}}(t)-\varphi_{LoS}(t+\tau)\right)}=K e^{\left.j 2 \pi \cos \gamma_{R} \tau\left(f_{\max_0 }+\frac{1}{2} \tau_{\max_{a}}+tf_{\max _{a}}\right)\right)}\)       (43)

\(r_{p q}^{S B T}(t, \tau)=\frac{\eta_{T}}{M} \sum_{m=1}^{M} \int_{-\pi}^{\pi} \int_{-\pi}^{\pi} e^{j 2 \pi j\left(\varphi_{m n_1}(t)-\varphi_{m n_1}(t+\tau)\right)} p\left(\gamma_{s}\right) p\left(\Delta \alpha_{T, m}\right) d \gamma_{s} d \Delta \alpha_{T, m}\)       (44)

\(r_{p q}^{S B R}(t, \tau)=\frac{\eta_{R}}{N} \sum_{n=1}^{N} \int_{-\pi}^{\pi} \int_{-\pi}^{\pi} e^{j 2 \pi j\left(\varphi_{nn_{2}}(t)-\varphi_{nn_{2}}(t+\tau)\right)} p\left(\gamma_{s}\right) p\left(\Delta \beta_{R, n}\right) d \gamma_{s} d \Delta \beta_{R, n}\)       (45)

\(r_{p q}^{S B E}(t, \tau)=\frac{\eta_{E}}{L} \sum_{l=1}^{L} \int_{-\pi}^{\pi} e^{j 2 \pi j\left(\varphi_{l n_3}(t)-\varphi_{l n_3}(t+\tau)\right)} p\left(\Delta \beta_{E, l}\right) d \Delta \beta_{E, l}\)       (46)

by using the equality \(\int_{-\pi}^{\pi} \exp (a \sin c+b \cos c) d c=2 \pi I_{0}(\sqrt{a^{2}+b^{2}})\), the integral in (44) , (45) and (46) can be obtained the closed-form solutions for (47) ,(48), and (49) , respectively.

\(r_{p q}^{S B T}(t, \tau)=\frac{\eta_{T}}{I_{0}\left(k_{T}\right) M} J_{0}\left(4 \pi \tau f_{\max }^{S}\right) \sum_{m=1}^{M} I_{0}(\sqrt{z_{m_{1}}^{2}+w_{m_{1}}^{2}}) e^{-j 2 \pi \tau\left(f_{\max _{0}}+f_{\max _{0}} t\right) \cos \Psi_{R}}\)       (47)

\(r_{p q}^{S B R}(t, \tau)=\frac{\eta_{R}}{I_{0}\left(k_{R}\right) N} J_{0}\left(4 \pi \tau f_{\max }^{S}\right) \sum_{n=1}^{N} I_{0}(\sqrt{z_{m_{2}}^{2}+w_{m_{2}}^{2}})\)       (48)

\(r_{p q}^{S B E}(t, \tau)=\frac{\eta_{E}}{I_{0}\left(k_{R}\right) L} \sum_{l=1}^{L} I_{0}(\sqrt{\left.\left.{k_{e l}}^{2}-4 \pi^{2} \tau^{2}\left(f_{\max _{0}}+f_{\max _{a}} t\right)\right)^{2}+4 \pi j k_{el} \tau\left(f_{\max _{0}}+f_{\max _{0}} t\right) \cos \left(\beta_{E, l}-\gamma_{R}\right)\right)}\)       (49)

where

\(f_{\max _{0}}^{R}=f_{c} \mathrm{v}_{R}\left(t_{0}\right) / c, \quad f_{\max _{a}}^{R}=f_{c} a_{0} / c, \quad f_{\max }^{S}=f_{c} \mathrm{v}_{s} / c\)       (50)

\(z_{m_{1}}=-j 2 \pi \tau\left(f_{\max _{0}}^{R}+f_{\max _{a}}^{R} t\right) \sin \gamma_{R} \cos \alpha_{T, m} R_{r} / D\)       (51)

\(w_{m_{1}}=-j 2 \pi \tau\left(f_{\mathrm{max}_{0}}^{R}+f_{\mathrm{max}_{a}}^{R} t\right) \sin \gamma_{R} \sin \alpha_{T, m} R_{r} / D+k_{T}\)       (52)

\(z_{m_{2}}=j 2 \pi \tau\left(f_{\max _{0}}^{R}+f_{\max _{a}}^{R} t\right) \sin \left(\beta_{R, n}-\gamma_{R}\right)\)       (53)

\(w_{m_{2}}=-j 2 \pi \tau\left(f_{\max _{0}}^{R}+f_{\max _{a}}^{R} t\right) \cos \left(\beta_{R, n}-\gamma_{R}\right)+k_{R}\)       (54)

3.3 Wigner-Ville Spectrum

 The Wigner-Ville spectrum S(f,t) of the proposed channel model is obtained by the Fourier transform of the temporal ACF \(r_{p q}(t, \tau)\) with respect to the time interval τ , which is given by

\(S(f, t)=\int_{-\infty}^{\infty} r_{p q}(t, \tau) \exp (-j 2 \pi f \tau) d \tau\)       (55)

By substituting (49) into (55), the expression for the Wigner-Ville spectrum of single bounced components for static scatter clusters can be written as

\(\begin{array}{l} S_{E}(f, t)=\frac{\eta_{E}}{I_{0}\left(k_{R}\right) L} \sum_{l=1}^{L} \int_{-\infty}^{\infty} \exp (-j 2 \pi f \tau) \\ \times I_{0}(\sqrt{\left.\left.{k_{e l}}^{2}-4 \pi^{2} \tau^{2}\left(f_{\max _{0}}+f_{\max _{a}} t\right)\right)^{2}+4 \pi j k_{e l} \tau\left(f_{\max _{0}}+f_{\max _{0}} t\right) \cos \left(\beta_{E, l}-\gamma_{R}\right)\right)} d \tau \end{array}\)       (56)

Let \(k_{e l}^{2}=k_{e l}^{2}\left(\cos ^{2}\left(\beta_{E, l}-\gamma_{R}\right)+\sin ^{2}\left(\beta_{E, l}-\gamma_{R}\right)\right)\), (56) can be expressed as

\(S_{E}(f, t)=\frac{\eta_{E}}{I_{0}\left(k_{R}\right) L} \sum_{l=1}^{L} \int_{-\infty}^{\infty} I_{0}\left(k_{0} v_{R}(t) j \sqrt{\left(\tau-\frac{k_{e l} \cos \left(\beta_{E, l}-\gamma_{R}\right)}{j k_{0} v_{R}(t)}\right)^{2}+\left(\frac{k_{e l} \sin \left(\beta_{E, l}-\gamma_{R}\right)}{j k_{0} v_{R}(t)}\right)^{2}}\right) \exp (-j 2 \pi f \tau) d \tau\)       (57)

Let , \(\zeta=\tau-\omega, \omega=\frac{k_{R} \cos \left(\beta_{E, l}-\gamma_{R}\right)}{j k_{0} v_{R}(t)}\), and according to \(I_{0}(x)=J_{0}(j x)\), (57) can be expressed as

\(S_{E}(f, t)=\frac{2 \eta_{E}}{I_{0}\left(k_{R}\right) L} \sum_{l=1}^{L} \exp (-j 2 \pi f \omega) \int_{0}^{\infty} J_{0}\left(k_{0} v_{R}(t) \sqrt{\zeta^{2}+\left(\omega \tan \left(\beta_{E, l}-\gamma_{R}\right)\right)^{2}}\right) \exp (-j 2 \pi f \zeta) d \zeta\)       (58)

According to \((\sqrt{\alpha^{2}-\beta^{2}}) \int_{0}^{\infty} J_{0}(\alpha \sqrt{x^{2}+y^{2}}) \cos (\beta x) d x=\cosh (y \sqrt{\alpha^{2}-\beta^{2}})\) and the integral nature of the symmetric even function, (58) can thus be rewritten as

\(\begin{array}{c} S_{E}(f, t)=\frac{2 \eta_{E}}{I_{0}\left(k_{e l}\right) L} \exp (-j 2 \pi f \omega) \sum_{l=1}^{L} \int_{0}^{\infty} J_{0}\left(k_{0} v_{R}(t) \sqrt{\zeta^{2}+\left(\omega \tan \left(\beta_{E, l}-\gamma_{R}\right)\right)^{2}}\right) \cos (2 \pi f \zeta) d \zeta \\ \quad=\frac{2 \eta_{E}}{I_{0}\left(k_{R}\right) L} \frac{1}{\sqrt{\left(k_{0} v_{R}(t)\right)^{2}-(2 \pi f)^{2}})} \exp (-j 2 \pi f \omega) \sum_{l=1}^{L} \cosh \left(\omega \tan \left(\beta_{E, l}-\gamma_{R}\right) \sqrt{\left(k_{0} v_{R}(t)\right)^{2}-(2 \pi f)^{2}}\right) \end{array}\)       (59)

3.4 Ergodic Capacity

 Channel capacity is defined as the mutual information of the channel maximized over all possible input distributions, which is usually one of the important indicators that fully describes the characteristics of the channel [37]. Since the channel capacity relates to the space-time correlation characteristics and non-uniform angular spread of scatterers. Furthermore, the ergodic capacity has been generally used as one of channel capacity measures for MIMO V2V commutation systems due to less complicated to study analytically. The ergodic capacity is expressed in b/s/Hz as [38-39]

\(C=E\left(C_{t}\right)=E\left(\log _{2} \operatorname{det}\left(\mathbf{I}_{L_{r}}+\frac{\bar{\gamma}}{L_{t}} \mathbf{H}(t) \mathbf{H}^{H}(t)\right)\right)\)       (60)

where it is assumed that \(L_{t} \geq L_{r}\), the transmitter has no channel knowledge, and the receiver has perfect channel knowledge. (·)^H denotes the transpose conjugate operation, det(·) denotes the matrix determinant, \(\mathbf{I}_{L_{r}}\) is the \(L_{r} \times L_{t}\) identity matrix, \(\bar{\gamma}=P_{T} / \sigma^{2}\) is the average signal-to-noise ratio (SNR), and \(\mathbf{H}(t)\) are generated as

\(\mathbf{H}(t)=\mathbf{R}_{r}^{1 / 2} \mathbf{H}_{w}\left(\mathbf{R}_{t}^{1 / 2}\right)^{T}\)       (61)

where \(\mathbf{H}_{w}\) is an \(L_{r} \times L_{t}\) stochastic matrix with complex Gaussian independent identically distributed entries, (·)^1/2 denotes the matrix square root operation, (·)^T denotes the transpose operation, \(\mathbf{R}_{r}\) and \(\mathbf{R}_{t}\) are the correlation matrices. For the sake simplicity, assuming that \(L_{t}=L_{r}=N_{L}\) and the normalized channel correlation coefficient between any two channels is ρ . According to the exponential correlation matrix model [40], a simple upper bound closed-form analytical expression can be expressed as

\(C=N_{L} \cdot \log _{2}\left(1+\bar{\gamma}\left(1-\rho^{2}\right) / N_{L}\right)+\log _{2}\left(\frac{1+\bar{\gamma} / N_{L}}{1+\bar{\gamma}\left(1-\rho^{2}\right) / N_{L}}\right)\)       (62)

4. Numerical Results

 In this section, we investigate the cooperative scattering channel model in detail for each statistical property. The following main parameters are used for our numerical analysis unless stated otherwise [24],[36],[41]: \(R_{t}=R_{r}=15 m, \psi_{T}=\psi_{R}=\pi / 4, L_{r}=L_{r}=N_{L}=4, \gamma_{R}=\pi / 3\),\(p=q=1, p^{\prime}=q^{\prime}=2 \quad, \quad a=180 m, \quad f=150 m, v_{R}\left(t_{0}\right)=5 \mathrm{km} / h \quad, v_{s}=40 \mathrm{km} / h, f_{c}=5.9 \mathrm{GHz}\)\(D_{c}=10, a_{0}=2.5 \mathrm{m} / \mathrm{s}^{2}, \lambda_{R}=0.04, p_{c}=0.3, M=N=L=4, K=0.156, \eta_{T}=\eta_{R}=0.2, \eta_{E}=0.6\), and \(k_{T}=k_{R}=k_{e l}=0\) . Meanwhile, considering the scatter clusters distributed on three-ring uniformly, we choose extended method of exact Doppler spread (EMEDS) [42] to derive the \(\alpha_{T, m}, \beta_{R, n}\), and \(\beta_{E, l}\), in the simulation.

 Fig. 4 illustrates the relation between Spatial CCF \(\rho_{11,22}\left(\delta_{T}, \delta_{R}\right)\) and normalized antenna element spacings of TX and RX\(\left(\delta_{T} / \lambda, \delta_{R} / \lambda\right)\). It can be observed that the normalized correlation coefficient is equal to 1 in case of \(\delta_{T} / \lambda=0\) and \(\delta_{R} / \lambda=0\), which means the fading between different antenna links is completely correlated. Meanwhile, the spatial correlation gradually decreases as normalized transceiver antenna element spacing increases. Especially, the spatial correlation is basically maintained a shock at 0.2 attachments and gradually stabilized when \(\delta_{T} / \lambda\) is greater than or equal to 2.5. As a result, the increasing in transceiver antenna element spacing does not always improve the performance of the system.

Fig. 4. Spatial CCF versus normalized transceiver antenna element spacings.

Fig. 5. Spatial CCF \(\rho_{11,22}\left(\delta_{T}, \delta_{R}\right)\) versus normalized RX antenna element spacing \(\delta_{R} / \lambda\) under different Rician factor K.

 Fig. 5 presents spatial CCF \(\rho_{11,22}\left(\delta_{T}, \delta_{R}\right)\) versus normalized RX antenna element spacing \(\delta_{R} / \lambda\) under different Rician factor K. We set \(\delta_{T} / \lambda=\delta_{R} / \lambda\). It is shown that the spacing of antenna element spacing corresponding to the first closely to zero point of the correlation coefficient is about 0.4, and the spatial correlation is the lowest at this time. The spatial CCF curve is in accord with the fluctuation characteristic of the Bessel function, which is consistent with Fig. 4. Meanwhile, The spatial CCF increase with the increasing of K. This is because the value of K can depict vehicular traffic density to a certain extent. For a low vehicular traffic density, the value of K is larger than that in the high vehicular traffic density scenario due to sparse moving vehicles (relay nodes), also the LOS components bear a significant amount of power. There exist a large number of single bounced relay scatting links without direct transmission link in case of K=0 (Rayleigh fading). Moreover, the existence of multi-path (multi-relay) scattering is beneficial to reduce the correlation between the antenna links in the case of urban traffic congestion. To a certain extent, the performance of the MIMO V2V communication system can be improved.

Fig. 6. Temporal ACF \(r_{p q}(t, \tau)\) versus time interval .. under various values of the moving scatter clusters velocity \(v_{s}\) for SBR mode.

 Fig. 6 plots the temporal ACF \(r_{p q}(t, \tau)\) versus time interval τ under various values of the moving scatter clusters velocity \(v_{s}\) for SBR mode. From Fig. 6, it is clear that higher value of \(v_{s}\) leads to lower time correlation property. This is because the channel fluctuations caused by the movement of moving scatter clusters increase as the moving scatter clusters velocity \(v_{s}\) increases, and thus the probabilities of clusters during the time interval τ decreases. Finally, the value of temporal correlation property decreases according to the theoretical analysis of equation (42).

Fig. 7. Temporal ACF \(r_{p q}(t, \tau)\) versus time interval τ under various speed acceleration α₀ and observation time instant t for SBE mode

Fig. 8. Doppler Wigner-Ville spectra \(S_{E}(f, t)\) ft versus observation time instant t and Doppler frequency f for SBE mode

 The non-stationary statistical characteristics of temporal ACF \(r_{p q}(t, \tau)\) and Doppler Wigner-Ville spectrum \(S_{E}(f, t)\) ft for SBE mode are presented in Fig. 7 and Fig. 8, respectively. From Fig. 7, we observe that the non-stationary of \(r_{p q}(t, \tau)\)τ is not only related to time interval τ , but also depending \(\alpha_{0}\) and t. As the τ increases, the fluctuation amplitude of the temporal ACF gradually converges. And the temporal ACF value decreases with the increasing of \(\alpha_{0}\) or t. This is due to the fact that the velocity of the mobile vehicle RX increases with t and \(\alpha_{0}\) , which leads to a larger Doppler frequency shift and channel fluctuations. Thereby the temporal autocorrelation is reduced. The trend of temporal ACF is consistent with the literature [25], which verifies the authenticity of the proposed model. Meanwhile, the non-stationary characteristic of the \(S_{E}(f, t)\) can be referred to the closed-form solution in (59), which is displayed in Fig. 8. Similar to the classic Jakes (Clarke) spectrum, but the value of \(S_{E}(f, t)\) at different t also has a non-stationary time-varying characteristic, which is characterized by an increase in motion speed due to acceleration with increasing time. The bandwidth of \(S_{E}(f, t)\) can be gradually broadened.

Fig. 9. Ergodic capacity C as a function of the average SNR for different angular spread \(k_{e l}\)

Fig. 10. Ergodic capacity C versus effective scatter clusters number L

 The effect of angular spread \(k_{e l}\) and number of clusters L on the ergodic capacity under various SNR for SBE mode are presented in Fig. 9 and Fig. 10, respectively. Fig. 9 plots ergodic capacity C as a function of the average SNR for different angular spread \(k_{e l}\) . It can be observed that ergodic capacity is dependent to the value of \(k_{e l}\) . Obviously the ergodic capacity decreases as the scattering becomes more non-isotropic. For example, When we set the average SNR equal to 20 dB , the numerical results of the ergodic capacity are 34.2 bit/s/Hz , 17.8 bit/s/Hz, 13.3 bit/s/Hz and 13.1 bit/s/Hz for \(k_{e l}=0, k_{e l}=3, k_{e l}=6\) and \(\) , respectively. The von Mises degenerates into a uniform distribution when \(k_{e l}=0\), the ergodic capacity is the largest. A large value of \(k_{e l}\) results in low ergodic capacity, the underlying physical reason is that it results in highly directional channels where multi-path cooperative scattering components are densely distributed in certain directions. However, as \(k_{e l}\) continues to increase to a certain extent, the tendency of reduction becomes slower. This is due to the value of correlation tends to stabilization. Fig. 10 shows the ergodic capacity C versus effective scatter clusters number L . It can be easily deduced from this figure that as the number of effective scatter clusters increases, the ergodic capacity increases due to the fact that multi-path (multi-relay) scattering components decrease the correlation.

 Fig. 11 demonstrates ergodic capacity versus time interval τ under various speed acceleration α₀ and observation time instant t . We observe that ergodic capacity rises rapidly at the beginning and then gradually tends to stabilization with increasing time interval τ from all considered speed acceleration and observation time instant. The reason is that temporal autocorrelation value gradually drops and tends to stabilization as the value of τ increases. At the same time, the greater the moving velocity of the transceiver at the current time instant t, the bigger the ergodic capacity and the faster the convergence, which can be seen from the analysis of Fig. 7. For example, under the same value of τ , the ergodic capacity for α₀ = 2 and t = 3 is always greater than or equal the ergodic capacity for α₀ = 0 and t = 1.

Fig. 11. Ergodic capacity C versus time interval τ under various speed acceleration α₀ and observation time instant t .

5. Conclusion

 In this paper, we have proposed a non-stationary three-ring geometric channel model based on relay cooperative scattering for MIMOV2V communication systems，wherein the scatter clusters are treated as the relay cooperative nodes. The non-stationary V2V channel temporal correlation, Wigner-Ville spectrum, and channel ergodic capacity are then investigated for non-isotropic scattering environments. The numerical results demonstrate that time-variant velocities and dynamically moving scatterers results in non-stationarity channel characteristics. Moreover, the proposed model provides a more straight and concise to study the impact of some parameters (e.g., antenna element spacing, the angular spread, effective scatter clusters number) and the vehicular traffic density on channel characteristics. For future work, the proposed model can be considered as a theoretical guidance for establishing V2V measurement campaigns.