Spatially correlated gamma-gamma scintillation in atmospheric optical channels Jos´e Mar´ıa Garrido-Balsells,∗ Antonio Jurado-Navas, Jos´e Francisco Paris, Miguel Castillo-V´azquez, and Antonio Puerta-Notario Dpt. Communications Engineering, University of M´alaga Campus de Teatinos s/n, E-29071 M´alaga, Spain ∗ [email protected]
Abstract: In this paper, novel analytical closed-form expressions are derived for the probability density function of the sum of identically distributed correlated gamma-gamma random variables that models an optical atmospheric channel communication with receiver spatial diversity. The mathematical expressions here proposed provide a general procedure to obtain information about the scintillation effects induced by turbulence over a diversity reception scheme implementing equal-gain combining method. Both, validity and accuracy of the obtained statistical distribution are corroborated by comparing the analytical results to numerical results obtained by Monte-Carlo simulations. These simulations are particularized for constant, exponential and circular correlation models, corresponding to three different receivers spatial configurations. In addition, the extreme situations of no correlation and fully correlated received signals are also studied. The presented expressions lead to a simple and easy-to-compute analytical procedure of analyzing atmospheric optical communications systems with correlated spatial diversity. © 2014 Optical Society of America OCIS codes: (010.1300) Atmospheric propagation; (010.1330) Atmospheric turbulence; (060.2605) Free-space optical communication; (290.5930) Scintillation.
References and links 1. L. C. Andrews and R. L. Phillips, Laser Beam Propagation through Random Media (SPIE, 2005). 2. M.A. Al-Habash, L.C. Andrews and R.L. Phillips, “Mathematical model for the irradiance probability density function of a laser beam propagating through turbulent media,” Opt. Engineering 40, 1554–1562 (2001). 3. A. Jurado-Navas, J. M. Garrido-Balsells, J. F. Paris, M. Castillo-V´azquez, and A. Puerta-Notario, “General analytical expressions for the bit error rate of atmospheric optical communication systems,” Opt. Lett. 36, 4095–4097 (2011). 4. N. Wang and J. Cheng, “Moment-based estimation for the shape parameters of the gamma-gamma atmospheric turbulence model,” Opt. Express 18, 12824–12831 (2010). 5. A. Dang, “A closed-form solution of the bit-error rate for optical wireless communication systems over atmospheric turbulence channels,” Opt. Express 19, 3494–3502 (2011). 6. M. Toyoshima, H. Takenaka, and Y. Takayama, “Atmospheric turbulence-induced fading channel model for space-to-ground laser communications links,” Opt. Express 19, 15965–15975 (2011). 7. M. M. Ibrahim and A. M. Ibrahim, “Performance analysis of optical receivers with space diversity reception,” Proc. IEEE Commun. 143, 369–372 (1996). 8. Z. Chen, S. Yu, T. Wang, G. Wu, S. Wang, and W. Gu, “Channel correlation in aperture receiver diversity systems for free-space optical communication,” J. Opt. 14, 125710 (2012).
#214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21820
9. M. Razavi and J. H. Shapiro, “Wireless optical communications via diversity reception and optical preamplification,” IEEE Trans. Wireless Commun. 4, 975–983 (2005). 10. H. Samimi, “Distribution of the sum of K-distributed random variables and applications in free-space optical communications,” IET Optoelectronics 6, 1–6 (2012). 11. T.A. Tsiftsis, H. G. Sandalidis, G. K. Karagiannidis, M. Uysal, “Optical wireless links with spatial diversity over strong atmospheric turbulence channels,” IEEE Trans. Wireless Commun. 8, 951–957 (2009). 12. H. Moradi, H.H. Refai, P.G. LoPresti, “Switch-and-stay and switch-and-examine dual diversity for high-speed free-space optics links,” IET. Optoelectronics 6, 34–42 (2012). 13. H. Moradi, H.H. Refai, P.G. LoPresti, “Thresholding-based optimal detection of wireless optical signals,” J. Opt. Commun. Net. 2, 689–700 (2010). 14. S. Kotz, J. Adams, “Distribution of sum of identically distributed exponentially correlated gamma variables,” Ann. Math. Statist. 35, 277–283 (1964). 15. V. A. Aalo, “Performance of maximal-ratio diversity systems in a correlated Nakagami-fading environment,” IEEE Trans. Commun. 43, 2360–2369 (1995). 16. M. S. Alouini, A. Abdi, M. Kaveh, “Sum of gamma variates and performance of wireless communication systems over Nakagami-fading channels,” IEEE Trans. Vehic. Tech. 50, 1471–1480 (2001). 17. P.G. Moschopoulos, “The distribution of the sum of independent gamma random variables,” Ann. Inst. Statist. Math. (Part A) 37, 541–544 (1985). 18. P. Lombardo, G. Fedele, M. M. Rao, “MRC performance for binary signals in Nakagami fading with general branch correlation,” IEEE Trans. Commun. 47, 44–52 (1999). 19. A. Erdelyi, Tables of Integrals Transforms, vol. I (McGraw Hill, 1954). 20. H. Exton, Multiple Hypergeometric Functions and Applications (John Wiley & Sons, 1976). 21. M. Abramowitz, I.A. Stegun, Handbook of Mathematical Functions (Dover Publications, 1972). 22. A. Oppenheim, R. W. Schafer, Discrete-time Signal Processing (Prentice Hall, 1999). 23. I. S. Gradshteyn and I.M. Ryzhik, Table of Integrals, Series and Products (Elsevier, 2007). 24. M.C. Jeruchim, P. Balaban, K.S. Shanmugan, Simulation of Communication Systems (Plenum, 1992). 25. K. Zhang, Z. Song, Y. L. Guan, “Simulation of Nakagami fading channels with arbitrary cross-correlation and fading parameters,” IEEE Trans. Wireless Commun. 3, 1463–1468 (2004).
1.
Introduction and problem statement
In the past few decades, scientific community has given considerable attention to investigations in atmospheric optical communication (AOC) systems. Due to the complexity associated with phase or frequency modulation, current AOC systems typically use intensity modulation with direct detection (IM/DD). In these systems, atmospheric turbulence leads to fading of the received signal intensity (scintillation), degrading the average bit error rate of such links. This performance degradation is deduced from the probability density function (pdf) of the irradiance. Thus, most widely accepted irradiance pdf models have led to consider optical scintillations as a conditional random process [1–3]. Among them, the gamma-gamma (GG) distribution [1, 2] has gained wide acceptance [4–6]. A promising solution to mitigate the degrading effects of atmospheric turbulence is the use of spatial diversity reception [7, 8]. In this technique, by employing multiple apertures at the receiver side, the inherent redundancy of spatial diversity has the potential to significantly enhance the performance of AOC systems. To evaluate this performance, it is essential to provide an analytical expression for the pdf of the sum of correlated irradiances seen at the receiver side. In [7], the fading induced by turbulence was assumed to be uncorrelated at each of the optical receivers. However, in order for this assumption to hold true, the spacing between receivers should be greater than the fading correlation length, which may be difficult to achieve in practice due to the limited available physical space or because the receiver spacing required for uncorrelated fading may exceed the beam diameter in power-limited links with well-collimated beams. As in [7], many other works have assumed that scintillation is uncorrelated at each of the optical receivers [9–11]. Alternatively, in a further step, correlation was incorporated to different statistical models and simulations [8]. In this sense, a plausible approximation consists in considering identically distributed variables, as assumed in [12–14]. To the best of the author’s knowledge, in the bibliography there are no closed-form math-
#214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21821
Fig. 1. Single-input multiple-output (SIMO) system model with a laser transmitter and N optical receivers.
ematical expressions for the pdf of correlated gamma-gamma distributed random variables (RVs). Thus, the main goal in this paper is to develop an analytical pdf model of the received optical irradiance, consistent with the scintillation theory [1], considering correlated spatial diversity reception with N receivers, and employing equal-gain combining (EGC) diversity technique. Figure 1 shows the system model with spatial diversity reception considered in this work. In this model, assuming on-off keying (OOK), IM/DD modulation and detection technique and that the N photodetectors are identical, the received optical irradiance, Irx , can be written as a sum of the individual i-th received irradiances as N
Irx = ∑ Irx,i .
(1)
i=1
Additionally, Irx,i can be expressed as Irx,i = I0 Ii , where I0 the irradiance in absence of atmospheric turbulence and Ii is the normalized irradiance associated to the atmospheric scintillation produced in the i-th channel. To characterize only the effect of the turbulence, I0 = 1 is assumed in this analysis. Thus, the normalized global irradiance I can be expressed as N
I = ∑ Ii .
(2)
i=1
As explained before, the scintillation induced by the atmospheric turbulence on the optical signal is modeled here by the GG statistical distribution. Thus, Ii a GG random variable (RV) which models the normalized irradiance observed by the i-th receiver. Furthermore, as described in [1, 2], this normalized irradiance is the result of a modulation process and can be written as the product Ii = XiYi , where Xi and Yi are RVs that arise from large-scale and smallscale effects over the received signal in the i-th photodetector. Both, Xi and Yi , are governed by a Gamma statistical distribution, denoted as G (a, b), and have a pdf given by f γ (γ ) =
γ γ a−1 exp − , ba Γ(a) b
(3)
where a and b are the shape parameters of the distribution, and Γ(·) is the Gamma function. The first and second moments of these Gamma RVs are, E[γ ] = ab and var[γ ] = ab2 , respectively. On the one hand, as described in [1], large-scale fluctuations in the irradiance are generated by turbulent eddies with sizes ranging from the scattering disk to the outer scale L0 . Since the #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21822
scattering disk is given by L/ρ0 k, where L is the path length, ρ0 is the spatial coherent radius and k is the wave number, it increases when L grows, ρ0 decreases (i.e. the strength of turbulence increases), or both. When this occurs, the average size of the large-scale eddies grows and, as a consequence, the probability that signals received in each photodetector pass through the same eddies (common eddies) also rises, resulting in an increase in the channel correlation [8]. Taking into account this channel behavior and in order to achieve compact expressions, in this work, we assume that signals of the N receivers are deflected by the same eddies and, therefore, largescale effects are common for all of them, i.e. we assume that Xi = X, ∀i = 1...N. Note that, with this assumption, we consider the most unfavorable scenario for large-scale correlation, which is also useful as a benchmark in optical systems design. On the other hand, small-scale contributions to scintillation are primarily diffractive in nature and are associated with turbulent eddies smaller than ρ0 . In this case, unlike Xi , small-scale effects of Yi are dependent on each receiver. However, even in this case, it can be assumed that Yi , ∀i = 1...N, will exhibit the same statistical properties and, therefore, all Yi will be identically distributed. In this situation, Eq. (2) can be rewritten as N
I = X ∑ Yi = XV.
(4)
i=1
Here, X is governed by a Gamma distribution, with parameters αx and βx , denoted as X ∼ G (αx , βx ). Likewise, the N small-scale identically distributed variates Yi also follow a Gamma distribution, but with the i-th channel parameters αi = α and βi = β , i.e. Yi ∼ G (α , β ). Finally, V = ∑Ni=1 Yi is the RV that describes the small scale atmospheric effect over the combined received optical irradiance. Note that since I is a normalized optical irradiance, the large and small scale parameters must fulfill that αx = 1/βx and α = 1/(N β ), respectively, and therefore, E[I] = E[X]E[V ] = 1. Note that, both αx and α are related to the effective number of large-scale and small-scale eddies, respectively, and they range from values close to 1, for strong turbulence conditions, and values >> 1, for weak turbulence intensities, as detailed in [2]. When channel correlation is considered, the pdf of V is obtained by solving the pdf of the sum of correlated Gamma RVs. It is assumed that the correlation between the variates Yi and Y j , with i, j = 1..N is given by the correlation coefficient ρi j defined as cov (Yi ,Y j ) , ρi j = var (Yi ) var (Y j )
(5)
where 0 ≤ ρi j ≤ 1. Thus, the correlation information between the small scale in the N optical signals received can be expressed by means of a NxN correlation matrix Cy , given as ⎡ ⎤ ρ12 · · · ρ1N 1 ⎢ ρ21 1 · · · ρ2N ⎥ ⎢ ⎥ Cy = ⎢ . , (6) . .. ⎥ .. .. ⎣ .. . . ⎦ ρN1 ρN2 · · · 1 N×N where the value of the correlation coefficient ρi j depends on the distance between i-th and j-th receivers. Although the development presented here is valid for any system whose correlation matrix is expressed as Eq. (6), in order to analyze the accuracy of the proposed analytical expressions three different correlation models, which are related to three typical receivers geometries, are considered in this work. In particular, constant, exponential and circular models are applied. The first one, the constant correlation model is described in [15, 16], and it corresponds to the case of evenly spaced receivers, with ρi j = ρ ∀i, j = 1..N, being ρ the correlation coefficient #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21823
between the Gamma distributed signals received at adjacent photodetectors. The second one, the exponential correlation model is applied to a lineal equidistant photodetectors array, being ρi j = ρ |i− j| [15]. Finally, the circular correlation model corresponds to the case of N receivers situated on a circumference, presenting the correlation matrix a N-th circular symmetry [16]. 2.
Sum of correlated Gamma variates
As explained in the previous section, the first goal is to obtain a closed-form expression for the pdf of the sum of correlated Gamma RVs Yi , with i = 1...N. To this end, the Moschopoulos theorem [17] is extended in the way proposed by Alouini et al. [16], but with the restriction of being every Yi in Eq. (4) identically distributed, as previously supposed. Assuming that for an arbitrary function φ (x) its Laplace transform is denoted as L {φ (x); s}, then, from [16, Eq(9)], the moment generating function (MGF) of V can be expressed by using the statistical averaging operator E[·], as N
MV (s) = E[esV ] = L { fV (v); −s} = ∏(1 − λi s)−α ,
(7)
i=1
where {λi }Ni=1 are the eigenvalues of the matrix A = DC, being D a N × N diagonal matrix with the entries β for i = 1...N, and C a N × N positive definite correlation matrix given by ⎡ ⎤ √ √ ρ12 · · · ρ1N 1 √ √ ⎢ ρ21 1 ··· ρ2N ⎥ ⎢ ⎥ C=⎢ , (8) ⎥ .. .. .. .. ⎣ ⎦ . . . . √ √ ρN1 ρN2 · · · 1 N×N whose elements are the correlation coefficients of the underlying Gaussian processes that lead to the Gamma distributed fading, as described in [15, 16, 18]. Now, the cumulative density function (cdf) of V , FV (v), can be written in the following form 1 1 1 N L {FV (v); s} = L { fV (v); s} = MV (−s) = ∏(1 + λi s)−α , s s s i=1 which can be algebraically modified to obtain the identity N α ⎞−α ⎛ 1 1 ∏ N N − λi λi Γ (N α ) ∏ (1 + λi s)−α = i=1Γ (N α ) · sN α · ∏ ⎝1 − s ⎠ i=1 i=1
(9)
(10)
and by comparing it to [19, Eq. (5)], this new expression is obtained N
1
∏ (1 + λi s)−α = det (A)α Γ (N α ) · L
(N)
vN α −1 Φ2
α , . . . , α ; N α ; −vλ1−1 , . . . , −vλi−1
,
i=1
(11) (N) where the determinant det(A) = ∏Ni=1 λi and Φ2 is the confluent Lauricella function [19, 20]. Finally, from Eqs. (9) and (11), the pdf and cdf of V , can be obtained respectively by using the inverse Laplace transform vN α −1 −1 × fV (v) = L {MV (−s); v} = [det(A)]α Γ(N α ) (12) v v (N) , × Φ2 α , ..., α ; N α ; − , ..., − λ1 λN #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21824
and −1
FV (v) = L
3.
1 vN α MV (−s); v = × s [det(A)]α Γ(N α + 1) v v (N) . × Φ2 α , ..., α ; N α + 1; − , ..., − λ1 λN
(13)
Statistical distribution of the received irradiance I
Now, starting from the normalized received optical irradiance I = XV , as indicated in Eq. (4), the unconditional pdf for the normalized combined received irradiance, fI (I), is obtained by calculating the mixture of a Gamma distribution for X defined by Eq. (3) as G (αx , 1/αx ), fX (x) =
αx (αx x)αx −1 exp (−αx x), Γ(αx )
(14)
and the function fV (v) presented above in Eq. (12). Therefore, by first fixing x and writing v = I/x, the conditional pdf, fV (I|x), is given by fV (I|x) 1 = fV (I|x), |∂ I/∂ V | x
(15)
I N α −1 × xN α [det(A)]α Γ(N α ) I I (N) . α , ..., α ; N α ; − , ..., − × Φ2 xλ1 xλN
(16)
fI (I|x) = and then, fI (I|x) =
To obtain the unconditional irradiance distribution, fI (I), the average of Eq. (16) over the Gamma distribution given in Eq. (14) must be calculated, which leads to fI (I) =
∞ 0
fI (I|x) fX (x)dx =
I N α −1 αxαx I, α [det(A)] Γ(N α ) Γ(αx )
(17)
where I=
∞ 0
(N)
xαx −N α −1 · Φ2 exp (−αx x)dx
(18)
(N) (N) with Φ2 = Φ2 α , ..., α ; N α ; − xλI , ..., − xλIN . 1 At this point, a clever way to solve fI (I) is obtained by applying the Gauss-Laguerre quadrature [21, Eq. (25.4.45)] to approximate the value of the integral I. Following this approach, the irradiance pdf can be expressed as fI (I)
xiαx −N α I N α −1 αxαx n × ∑ [det(A)]α Γ(N α ) Γ(αx ) i=1 (n + 1)2 [Ln+1 (xi )]2 I αx I αx (N) × Φ2 α , ..., α ; N α ; − , ..., − xi λ1 xi λN
(19)
where xi is the i-th root of the Laguerre polynomial with degree n, denoted as Ln (x).
#214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21825
4.
Closed-form expression for fV (v) and fI (I) when α ∈ N
Due to the fact that the pdf equation derived in previous section can only be solved by numerical analysis and in order to achieve closed-form expressions for fV (v) and fI (I), in this work we assume that α ∈ N. It is worth noting that, from a physical viewpoint, this is a reasonable assumption since, as previously mentioned, in a GG distribution this parameter represents the effective number of small-scale turbulent eddies for each of the N receiving channels. Note that, for real values of α , the pdf functions will be bounded by the pdfs corresponding to the nearest natural values of α . Under this premise (α ∈ N), the integral I can be solved in a closed-form expression. In this way, by using the variable change t = 1/x, the Eq. (18) can be expressed as ∞ α x dt, f (t)t −αx exp − I= (20) t 0 (N) where f (t) = t N α −1 · Φ2 α , ..., α ; N α ; − λIt , ..., − λItN . To solve this expression, the function 1 f (t) can be rewritten by applying the Laplace transform. Thus, from [19, Eq. (4.24-5)], the Laplace transform of f (t) is given by F(s) = L { f (t); s} =
I −α Γ(N α ) N 1+ , sN α ∏ λi s i=1
which, after some manipulations, can be also expressed as −α N −I F(s) = Γ(N α ) ∏ s − . λi i=1
(21)
(22)
At this point, by using the partial fraction expansion method described in [22], the following conversion from product to summation of fractions can be applied N
N αi 1 cmi ∏ (s − di )αi = ∑ ∑ (s − di )m , i=1 i=1 m=1
(23)
where N is the number of different eigenvalues of matrix A, di = −I/λi . In addition, to consider the multiplicity of some eigenvalues, αi is the product of the algebraic multiplicity of the eigenvalue, denoted as μA (λi ), with the parameter α , i.e. αi = μA (λi )α . Finally cmi is the corresponding fraction coefficient, whose detailed calculation is described in Appendix A. Then, the application of this expansion technique over the Eq. (22) results in the following expression N αi cmi . (24) F(s) = Γ(N α ) ∑ ∑ (s − di )m i=1 m=1 Now, the inverse Laplace transform of F(s) can be calculated. For this purpose, knowing that 1 tn −1 exp (−at)u(t), (25) ;t = L (s + a)n+1 n! with u(t) being the Heaviside step function, the inverse Laplace transform of the double sum argument in Eq. (24) is derived as 1 t m−1 −1 exp (dit)u(t), (26) L m ;t = (s − di ) (m − 1)! #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21826
and, consequently, f (t) = L−1 {F(s);t} is finally expressed as follows, noting that di was previously defined as di = −I/λi , N αi I t m−1 exp − t . f (t) = Γ(N α ) ∑ ∑ cmi (27) (m − 1)! λi i=1 m=1 Now, including f (t) in Eq. (20), and by employing [23, Eq. (3.471-9)], the Eq. (20) can be solved as N αi m−αx m−αx cmi αx I − m−2αx 2 2 . (28) λi αx I Km−αx 2 I = 2Γ(N α ) ∑ ∑ λi i=1 m=1 (m − 1)! And finally, the pdf of the combined received irradiance is expressed, from Eq. (17), as
fI (I) =
N αi 2 cmi m−2αx m+2αx λ αx × ∑ ∑ α [det(A)] Γ(αx ) i=1 m=1 Γ(m) i
×I
N α −1− m−2αx
αx I , Km−αx 2 λi
(29)
providing an easy-to-compute closed-form analytical method to obtain the statistical properties of the combined received irradiance I from the turbulence parameters and the channel correlation information. It must be noted that, from Eq. (50), the parameter cmi is also dependent on the normalized irradiance I. Additionally, from the assumption of α ∈ N and returning to the definition of the pdf of the sum of correlated Gamma RVs, given in Eq. (12), the mathematical procedure shown in Eqs. (23) to (27) can be applied to obtain an equivalent closed-form expression of fV (v). Then, from Eq. (7) and Eq. (23) N αi cmi 1 −1 (30) fV (v) = L ∑ ∑ (s − di )m ;t , [det(A)]α i=1 m=1 where, in this case, di = −λi−1 and the coefficients cmi are conveniently calculated with the procedure described in Appendix A. By using again the inverse Laplace transform given in Eq. (26), the pdf for the sum of Gamma RVs can be calculated as N αi 1 v cmi m−1 (31) fV (v) = ∑ ∑ Γ(m) v exp − λi . [det(A)]α i=1 m=1 An interesting result is reached when previous equation is manipulated to obtain the following equivalent expression N αi cmi fV (v) = ∑ ∑ λ m gm,i (v), (32) α i V [det(A)] i=1 m=1 where gVm,i (v) ∼ G (m, λi ). Then, fV (v) is said to be a mixture of Gamma RVs with shaping parameters m and λi . The summation of the corresponding coefficients verifies that N
αi
cmi λ m = 1, α i [det(A)] i=1 m=1
∑∑
(33)
since, in the particular case of Eq. (31), the coefficients cmi are independent from I. By contrast, regarding to fI (I), as previously mentioned, these coefficients are dependent on I and the equivalent mixture definition cannot be applied to Eq. (29). #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21827
5.
Numerical results
In this section, the expressions developed in the previous section are validated. To this end, analytical results obtained from these expressions are compared to numerical results obtained by simulations. In particular, several cases corresponding to different levels of correlation among the N channels have been considered. First, the extreme cases of no correlated and total correlated channels are analyzed. Latter, three cases of partially correlated channels are also studied. In all these cases, EGC method is assumed for combining the received irradiance in each photodetector. It is worth noting that, in order to obtain the numerical simulation pdf results, the histogram pdf estimator described in [24] is employed. Cases of no correlated and total correlated channels: When there is no correlation among the small scale variates, the Yi signals are independent and Cy is then a diagonal correlation matrix. Consequently, the eigenvalues of A are λi = λ = β , ∀i = 1..N, and the det(A) = ∏Ni=1 λi = λ N α . This case corresponds to a single-input singleoutput (SISO) channel where the small scale effects are modeled by a Gamma distribution whose shape parameter is α = N α , i.e. V ∼ G (N α , 1/(N α )). Thus, by applying the calculus procedure described in previous sections, the pdf of the received normalized irradiance is expressed as fI (I) =
N α + αx N α + αx 2 αx 2 λ − 2 × Γ(N α )Γ(αx )
×I
N α + αx −1 2
KN α −αx
2
αx I λ
(34) .
and, by substituting here λ = β = 1/N α , leads to a GG pdf with large and small scale parameters αx and N α , respectively, given by N α + αx
fI (I) =
√ N α + αx 2 (αx N α ) 2 I 2 −1 KN α −αx 2 αx N α I . Γ(N α )Γ(αx )
(35)
Likewise, when the optical irradiance at the N receivers is completely correlated, then, every element of the correlation matrix ρi j = 1 and, thus, all the random sequences Yi are exactly the same, i.e. Yi = Y ∀i = 1..N, and V = NY , where Y ∼ G (α , 1/(N α )). Then, by applying the scaling of RVs, fV (v) = fY (v/N)/N, the pdf of V results in a Gamma distribution such V ∼ G (α , 1/α ) and the pdf of I follows the corresponding gamma-gamma distribution with large and small scale parameters αx and α , respectively, expressed as α + αx
√ 2 (αx α ) 2 α +αx −1 I 2 fI (I) = Kα −αx 2 αx α I . Γ(α )Γ(αx )
(36)
Figure 2 shows the obtained results for the pdf of the normalized combined optical irradiance I. In this figure, analytical results calculated from Eq.(35) and Eq. (36) for ρ = 0 and ρ = 1 are represented with solid line and dashed line, respectively, whereas numerical results, obtained by generating the corresponding Gamma RVs are depicted with circles. Typical moderate turbulence intensity (αx = 10 and αy = 5) and N = 4 optical receivers are assumed in the comparison. From this figure, it is clearly observed that analytical and numerical results curves fit perfectly, corroborating the validity of the mathematical expressions proposed for these two limit cases.
#214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21828
Numerical results Analytical pdf ρ=1 Analytical pdf ρ=0
I
Probability density function, f (I)
1
0.8
0.6
0.4
0.2
0
0.5
1
1.5
2
2.5
3
Normalized total received irradiance, I
Fig. 2. Analytical and simulation results for the pdf of the normalized irradiance I for the extreme correlation cases of ρ = 0 and ρ = 1. N = 4 optical receivers and typical moderate turbulence intensity with αx = 10 and αy = 5 are assumed.
Cases of partially correlated channels: In contrast to the previous cases, to analyze realistic ρi j values, the generation of correlated Gamma sequences must be considered. For this purpose, a method to generate correlated Gamma RVs from linear combinations of independent Gamma variates, proposed in [25], is employed. The method is based on the Cholesky decomposition and on matching the first and second order moments of both independent and correlated Gamma set. Firstly, by applying the Cholesky decomposition, the correlation matrix Cy can be expressed as the product of a lower triangular matrix L by its transposed LT , i.e. Cy = LLT . If y = [y1 , . . . , yN ]T is the set of N correlated Gamma RVs with correlation matrix Cy , and g = [g1 , . . . , gN ]T is a set of N independent Gamma variates with parameters (αg,i , βg,i ), being its correlation matrix a N × N identity matrix, then y = Lg or yi =
i
∑ lik gk ,
(37)
k=1
where lik denotes the (i, k)-th position of L. Secondly, the moment matching procedure is applied, which is based on the knowledge of the moments values for the correlated Gamma variables yi , i.e. the first order moment of yi is E[yi ] = αβ , ∀i = 1..N, whereas the second order moment is var[yi ] = αβ 2 , ∀i = 1..N, being α and β the parameters defined in section 1. In addition, from Eq. (37) and taking into account that the variance of the sum of independent RVs is the sum of each RV variance, the first and second order moments of yi can be respectively expressed as E[yi ] =
i
∑ lik E[gk ],
(38)
k=1
and var[yi ] =
i
∑ lik2 var[gk ].
(39)
k=1
Then, the characteristics parameters αg,i and βg,i , are iteratively obtained from Eq. (38) and Eq. (39) as αβ − ∑i−1 k=1 lik E[gk ] = αg,i βg,i , (40) E[gi ] = lii #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21829
1.4
ρ=0.2 ρ=0.4 ρ=0.8
1.2
0.8
fI(I)
0.8
V
f (v)
1
ρ=0.2 ρ=0.4 ρ=0.8
1
0.6
0.6 0.4
0.4 0.2
0.2
0
0 0.2
0.4
0.6
0.8
1
1.2
1.4
1.6
1.6
1.2
0.5
1
1.5
2
2.5
Normalized irradiance, I (b) 1.2
σ 2=0.17 I σ 2=0.32 I σ 2=1
1
I
I
0.8
fI(I)
1 V
2
σI2=0.17 σ 2=0.32 I σ 2=1
1.4
f (v)
1.8
Sum of correlated Gamma RVs, V (a)
0.8
0.6
0.6 0.4 0.4 0.2
0.2 0.5
1
1.5
Sum of correlated Gamma RVs, V (c)
2
2.5
0
0.5
1
1.5
2
2.5
Normalized irradiance, I (d)
Fig. 3. Analytical and simulation results for the pdf of V and I RVs obtained with the constant correlation model and N = 3 equidistant optical receivers. In (a) and (b), the results are shown for different correlation levels with a typical turbulence intensity σI2 = 0.32. In (c) and (d) the results are shown as a function of turbulence conditions (weak, moderate and strong) with a correlation level given by ρ = 0.6.
and var[gi ] =
2 αβ 2 − ∑i−1 2 k=1 lik var[gk ] = αg,i βg,i , 2 lii
(41)
and thus, the parameters are calculated respectively as βg,i = var[gi ]/E[gi ] and αg,i = E[gi ]/βg,i . Finally, once the correlated RVs are generated and combined with the EGC technique, the numerical pdf is obtained in each case, by applying a histogram based pdf estimator, as mentioned before. The obtained results are compared to the analytical pdf calculated with the closed-form expressions proposed in this paper. Note that in the following figures, the turbulence intensity is given by the scintillation index, defined in [1] and particularized for the GG statistical distribution here studied as σI2 = αx−1 + α −1 + (αx α )−1 . Analytical and numerical results obtained for the constant correlation model are shown in Fig. 3. In this figure, dotted, dashed and solid lines represent analytical results, whereas circles depict numerical results. For this correlation model, N = 3 equidistant optical receivers are assumed in all cases. The simple correlation matrix is given by ⎤ ⎡ 1 ρ ρ (42) Cy = ⎣ ρ 1 ρ ⎦ . ρ ρ 1 In Fig. 3(a) and Fig. 3(b), the analytical and the simulation results for the pdf of V and I RVs are presented for different levels of correlation among the receivers: from low correlation (ρ = 0.2), which corresponds to separated photodetectors, to high correlation (ρ = 0.8) which corresponds to closer optical receivers. Likewise, in Fig. 3(c) and Fig. 3(d), the same pdfs are plotted but now for different levels of turbulence intensity, in particular, from weak (σI2 = 0.17, with {α = 15, αy = 10}), to moderate (σI2 = 0.32, with {αx = 10, αy = 5}) and #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21830
1.2
ρ=0.2 ρ=0.4 ρ=0.8
1.6 1.4
0.8
1
I
f (I)
V
f (v)
1.2
ρ=0.2 ρ=0.4 ρ=0.8
1
0.8
0.6 0.4
0.6 0.4
0.2
0.2 0 0.2
0.4
0.6
0.8
1
1.2
1.4
1.6
1.8
2
0.5
Sum of correlated Gamma RVs, V (a)
2
2
σI =0.17 2 σ =0.32 I2 σ =1 I
1 0.8
f (I)
I
1
I
V
f (v)
1.5
1.2
2
σI =0.17 2 σ =0.32 I2 σ =1
1.5
1
Normalized irradiance, I (b)
0.6 0.4
0.5
0.2 0
0.5
1
1.5
2
Sum of correlated Gamma RVs, V (c)
2.5
0
0.5
1
1.5
2
2.5
Normalized irradiance, I (d)
Fig. 4. Analytical and simulation results for the pdf of V and I RVs obtained with the exponential correlation model and N = 5 optical receivers. In (a) and (b), the results are shown for different correlation levels with a typical turbulence intensity σI2 = 0.32. In (c) and (d) the results are shown as a function of turbulence conditions (weak, moderate and strong) with a correlation level given by ρ = 0.6.
strong (σI2 = 1, with {α = 15, αy = 10}) turbulence intensities, respectively, considering a fixed correlation level with ρ = 0.6. From Fig. 3, it follows that, regardless of the correlation levels or the turbulence conditions, there is a very clear agreement between obtained analytical results and numerical simulations. Thus, these results corroborate the high accuracy of proposed closed-form expressions. Now, the same comparison is extended for the exponential and circular correlation models. In Fig. 4, analytical an numerical results for the exponential model are depicted. Here, a linear photodetector array with N = 5 receivers is assumed, with the following correlation matrix, ⎡ ⎤ 1 ρ ρ2 ρ3 ρ4 2 ρ3 ⎥ ⎢ ρ ⎢ 2 1 ρ ρ ⎥ ρ 1 ρ ρ2 ⎥ (43) Cy = ⎢ ⎢ ρ ⎥. ⎣ ρ3 ρ2 ρ 1 ρ ⎦ 1 ρ4 ρ3 ρ2 ρ Again, from Figs. 4(a)–4(d) the high level of accuracy obtained with the proposed closedform expressions is observed, regardless of the correlation levels or the turbulence conditions. Finally, in Fig. 5 the results for the circular correlation model are presented. In this case, the same numerical simulations carried out in [16] are performed. Thus, a N = 5 circular receivers array is considered with a correlation matrix defined as ⎡ ⎤ 1 ρ2 ρ3 ρ4 ρ5 ⎢ ρ5 1 ρ2 ρ3 ρ4 ⎥ ⎢ ⎥ ⎥ (44) Cy = ⎢ ⎢ ρ4 ρ5 1 ρ2 ρ3 ⎥ , ⎣ ρ3 ρ4 ρ5 1 ρ2 ⎦ ρ2 ρ3 ρ4 ρ5 1 #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21831
0.8
ρ =ρ =0.64 12 15 ρ =ρ =0.36
0.6
0.6
13
ρ12=ρ15=0.8 ρ =ρ =0.6 13
14
14
ρ =ρ =0.64 12 15 ρ =ρ =0.36
I
0.8
f (I)
14
13
1
fV(v)
1
ρ =ρ =0.8 12 15 ρ =ρ =0.6
1.2
13
14
0.4 0.4 0.2
0.2
Numerical results
Numerical results 0
0.5
1
1.5
0
2
Sum of correlated Gamma RVs, V (a) 1.6
fV(v)
2
2.5
2
σI =0.17 2 σ =0.32 I2 σ =1 I
0.8
I
Numerical results
1.5
1
f (I)
1 0.8
1.2
1
Normalized irradiance, I (b)
1.2
2
σ =0.17 I2 σ =0.32 I2 σ =1 I
1.4
0.5
0.6
Numerical results
0.6 0.4
0.4 0.2
0.2 0 0.5
1
1.5
2
Sum of correlated Gamma RVs, V (c)
2.5
0
0.5
1
1.5
2
2.5
3
Normalized irradiance, I (d)
Fig. 5. Analytical and simulation results for the pdf of V and I RVs obtained with the circular correlation model and N = 5 optical receivers. In (a) and (b), the results are shown for different correlation levels with a typical turbulence intensity σI2 = 0.32. In (c) and (d) the results are shown as a function of turbulence conditions (weak, moderate and strong) with a correlation level given by ρ2 = ρ12 = ρ15 = 0.8 and ρ3 = ρ13 = ρ14 = 0.6.
where, for symmetry, ρ2 = ρ5 and ρ3 = ρ4 . In particular, Fig. 5(a) and Fig. 5(b), show the pdf of V and I RVs for two correlation cases. In the first one, ρ2 = ρ12 = ρ15 = 0.8 and ρ3 = ρ13 = ρ14 = 0.6; in the second case, ρ2 = ρ12 = ρ15 = 0.64 and ρ3 = ρ13 = ρ14 = 0.36. Furthermore, Fig. 5(c) and 5(d) depict the effect of turbulence conditions on the accuracy of the proposed expressions. Again, as expected, a very high agreement between the numerical simulations and the analytical results is clearly observed. 6.
Concluding remarks
In this paper, the analysis of the probability density function of the received optical irradiance for a spatially diversity free space optical communication system has been performed. The gamma-gamma statistical distribution has been assumed to model the turbulence induced scintillation that produce the fading on the optical received signal, whereas the combination procedure used is the equal-gain combining scheme. In this context, novel analytical closedform expressions have been derived for the pdf of the sum of identically distributed spatially correlated Gamma random variables and then, for the corresponding sum of correlated gammagamma variates. The mathematical expressions here proposed have been corroborated with Monte-Carlo numerical simulations, resulting in an excellent agreement between numerical and analytical results. This comparison have been performed for the limiting cases of no correlation (ρ = 0) and fully correlated received signals (ρ = 1), as well as with three realistic correlations models (constant, exponential and circular), corresponding to three different receivers configurations. The study here presented corroborates the high accuracy of the proposed closed-form expressions, leading to a mathematically tractable and easily computable method to characterize the statistical behavior of the received optical irradiance in a spatial diversity context.
#214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21832
Appendix A In this Appendix, the calculation procedure of the partial fraction expansion coefficients, cmi , is described. From the method presented in [22], these coefficients are obtained from ⎡ ⎤
cmi =
dαi −m ⎢ N 1 1 ⎥ · α −m ⎣ ∏ αj ⎦ i (αi − m)! dw (w − d ) j j=1 j =i
,
(45)
w=di
where, from Eq. (23), N is the number of different eigenvalues of matrix A, d j = −I/λ j , and αi = μA (λi )α . Note that μA (λi ) was defined as the algebraic multiplicity of the eigenvalue i. Here, for the computing of each coefficient cmi , an extension of the Leibniz product rule based on the multinomial theorem can be applied. In this sense, from the generalized Leibniz rule defined as d n (u · v) = (u + v)(n) (46) dxn being u = u(x) and v = v(x), and the multinomial theorem given by (u1 + . . . + um ) =
∑
n
k1 +...+km =n
n k1 , . . . , km
N
∏ utkt ,
(47)
j=1
the n-th derivative of a product operator can be expanded in the following form N N dn n (k ) (n) = (u u + . . . + u ) = j N 1 ∑ ∏ uj j . k dxn ∏ , . . . , k N 1 j=1 j=1 k +...+kN =n
(48)
1
In this case, identifying Eq. (48) in Eq. (45), u j = (w − d j )−α j and its n-th derivative is (n) easily computed as u j = (−1)n (α j )n (w − d j )−α j −n , where (α j )n = Γ(α j + n)/Γ(α j ) is the Pochhammer symbol. Therefore, N
(k j )
∏ uj
=
j=1 j =i
N
∏
j=1 j =i
(−1)k j (α j )k j (w − d j )−α j −k j ,
(49)
and finally, the coefficient cmi can be expressed as 1 cmi = (αi − m)! k ×
N
∏
j=1 j =i
∑
1 +... +kN =αi −m iˆ
αi − m × k1 , . . ., kN iˆ
−α j −k j
(−1) (α j )k j (di − d j ) kj
(50)
,
where iˆ indicates that ki is omitted in the previous corresponding sequences, and w has been substituted for di , as indicated in Eq. (45). Acknowledgment This work was supported by the Spanish Ministerio de Econom´ıa y Competitividad, Project TEC2012-36737.
#214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21833
Abstract: In this paper, novel analytical closed-form expressions are derived for the probability density function of the sum of identically distributed correlated gamma-gamma random variables that models an optical atmospheric channel communication with receiver spatial diversity. The mathematical expressions here proposed provide a general procedure to obtain information about the scintillation effects induced by turbulence over a diversity reception scheme implementing equal-gain combining method. Both, validity and accuracy of the obtained statistical distribution are corroborated by comparing the analytical results to numerical results obtained by Monte-Carlo simulations. These simulations are particularized for constant, exponential and circular correlation models, corresponding to three different receivers spatial configurations. In addition, the extreme situations of no correlation and fully correlated received signals are also studied. The presented expressions lead to a simple and easy-to-compute analytical procedure of analyzing atmospheric optical communications systems with correlated spatial diversity. © 2014 Optical Society of America OCIS codes: (010.1300) Atmospheric propagation; (010.1330) Atmospheric turbulence; (060.2605) Free-space optical communication; (290.5930) Scintillation.
References and links 1. L. C. Andrews and R. L. Phillips, Laser Beam Propagation through Random Media (SPIE, 2005). 2. M.A. Al-Habash, L.C. Andrews and R.L. Phillips, “Mathematical model for the irradiance probability density function of a laser beam propagating through turbulent media,” Opt. Engineering 40, 1554–1562 (2001). 3. A. Jurado-Navas, J. M. Garrido-Balsells, J. F. Paris, M. Castillo-V´azquez, and A. Puerta-Notario, “General analytical expressions for the bit error rate of atmospheric optical communication systems,” Opt. Lett. 36, 4095–4097 (2011). 4. N. Wang and J. Cheng, “Moment-based estimation for the shape parameters of the gamma-gamma atmospheric turbulence model,” Opt. Express 18, 12824–12831 (2010). 5. A. Dang, “A closed-form solution of the bit-error rate for optical wireless communication systems over atmospheric turbulence channels,” Opt. Express 19, 3494–3502 (2011). 6. M. Toyoshima, H. Takenaka, and Y. Takayama, “Atmospheric turbulence-induced fading channel model for space-to-ground laser communications links,” Opt. Express 19, 15965–15975 (2011). 7. M. M. Ibrahim and A. M. Ibrahim, “Performance analysis of optical receivers with space diversity reception,” Proc. IEEE Commun. 143, 369–372 (1996). 8. Z. Chen, S. Yu, T. Wang, G. Wu, S. Wang, and W. Gu, “Channel correlation in aperture receiver diversity systems for free-space optical communication,” J. Opt. 14, 125710 (2012).
#214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21820
9. M. Razavi and J. H. Shapiro, “Wireless optical communications via diversity reception and optical preamplification,” IEEE Trans. Wireless Commun. 4, 975–983 (2005). 10. H. Samimi, “Distribution of the sum of K-distributed random variables and applications in free-space optical communications,” IET Optoelectronics 6, 1–6 (2012). 11. T.A. Tsiftsis, H. G. Sandalidis, G. K. Karagiannidis, M. Uysal, “Optical wireless links with spatial diversity over strong atmospheric turbulence channels,” IEEE Trans. Wireless Commun. 8, 951–957 (2009). 12. H. Moradi, H.H. Refai, P.G. LoPresti, “Switch-and-stay and switch-and-examine dual diversity for high-speed free-space optics links,” IET. Optoelectronics 6, 34–42 (2012). 13. H. Moradi, H.H. Refai, P.G. LoPresti, “Thresholding-based optimal detection of wireless optical signals,” J. Opt. Commun. Net. 2, 689–700 (2010). 14. S. Kotz, J. Adams, “Distribution of sum of identically distributed exponentially correlated gamma variables,” Ann. Math. Statist. 35, 277–283 (1964). 15. V. A. Aalo, “Performance of maximal-ratio diversity systems in a correlated Nakagami-fading environment,” IEEE Trans. Commun. 43, 2360–2369 (1995). 16. M. S. Alouini, A. Abdi, M. Kaveh, “Sum of gamma variates and performance of wireless communication systems over Nakagami-fading channels,” IEEE Trans. Vehic. Tech. 50, 1471–1480 (2001). 17. P.G. Moschopoulos, “The distribution of the sum of independent gamma random variables,” Ann. Inst. Statist. Math. (Part A) 37, 541–544 (1985). 18. P. Lombardo, G. Fedele, M. M. Rao, “MRC performance for binary signals in Nakagami fading with general branch correlation,” IEEE Trans. Commun. 47, 44–52 (1999). 19. A. Erdelyi, Tables of Integrals Transforms, vol. I (McGraw Hill, 1954). 20. H. Exton, Multiple Hypergeometric Functions and Applications (John Wiley & Sons, 1976). 21. M. Abramowitz, I.A. Stegun, Handbook of Mathematical Functions (Dover Publications, 1972). 22. A. Oppenheim, R. W. Schafer, Discrete-time Signal Processing (Prentice Hall, 1999). 23. I. S. Gradshteyn and I.M. Ryzhik, Table of Integrals, Series and Products (Elsevier, 2007). 24. M.C. Jeruchim, P. Balaban, K.S. Shanmugan, Simulation of Communication Systems (Plenum, 1992). 25. K. Zhang, Z. Song, Y. L. Guan, “Simulation of Nakagami fading channels with arbitrary cross-correlation and fading parameters,” IEEE Trans. Wireless Commun. 3, 1463–1468 (2004).
1.
Introduction and problem statement
In the past few decades, scientific community has given considerable attention to investigations in atmospheric optical communication (AOC) systems. Due to the complexity associated with phase or frequency modulation, current AOC systems typically use intensity modulation with direct detection (IM/DD). In these systems, atmospheric turbulence leads to fading of the received signal intensity (scintillation), degrading the average bit error rate of such links. This performance degradation is deduced from the probability density function (pdf) of the irradiance. Thus, most widely accepted irradiance pdf models have led to consider optical scintillations as a conditional random process [1–3]. Among them, the gamma-gamma (GG) distribution [1, 2] has gained wide acceptance [4–6]. A promising solution to mitigate the degrading effects of atmospheric turbulence is the use of spatial diversity reception [7, 8]. In this technique, by employing multiple apertures at the receiver side, the inherent redundancy of spatial diversity has the potential to significantly enhance the performance of AOC systems. To evaluate this performance, it is essential to provide an analytical expression for the pdf of the sum of correlated irradiances seen at the receiver side. In [7], the fading induced by turbulence was assumed to be uncorrelated at each of the optical receivers. However, in order for this assumption to hold true, the spacing between receivers should be greater than the fading correlation length, which may be difficult to achieve in practice due to the limited available physical space or because the receiver spacing required for uncorrelated fading may exceed the beam diameter in power-limited links with well-collimated beams. As in [7], many other works have assumed that scintillation is uncorrelated at each of the optical receivers [9–11]. Alternatively, in a further step, correlation was incorporated to different statistical models and simulations [8]. In this sense, a plausible approximation consists in considering identically distributed variables, as assumed in [12–14]. To the best of the author’s knowledge, in the bibliography there are no closed-form math-
#214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21821
Fig. 1. Single-input multiple-output (SIMO) system model with a laser transmitter and N optical receivers.
ematical expressions for the pdf of correlated gamma-gamma distributed random variables (RVs). Thus, the main goal in this paper is to develop an analytical pdf model of the received optical irradiance, consistent with the scintillation theory [1], considering correlated spatial diversity reception with N receivers, and employing equal-gain combining (EGC) diversity technique. Figure 1 shows the system model with spatial diversity reception considered in this work. In this model, assuming on-off keying (OOK), IM/DD modulation and detection technique and that the N photodetectors are identical, the received optical irradiance, Irx , can be written as a sum of the individual i-th received irradiances as N
Irx = ∑ Irx,i .
(1)
i=1
Additionally, Irx,i can be expressed as Irx,i = I0 Ii , where I0 the irradiance in absence of atmospheric turbulence and Ii is the normalized irradiance associated to the atmospheric scintillation produced in the i-th channel. To characterize only the effect of the turbulence, I0 = 1 is assumed in this analysis. Thus, the normalized global irradiance I can be expressed as N
I = ∑ Ii .
(2)
i=1
As explained before, the scintillation induced by the atmospheric turbulence on the optical signal is modeled here by the GG statistical distribution. Thus, Ii a GG random variable (RV) which models the normalized irradiance observed by the i-th receiver. Furthermore, as described in [1, 2], this normalized irradiance is the result of a modulation process and can be written as the product Ii = XiYi , where Xi and Yi are RVs that arise from large-scale and smallscale effects over the received signal in the i-th photodetector. Both, Xi and Yi , are governed by a Gamma statistical distribution, denoted as G (a, b), and have a pdf given by f γ (γ ) =
γ γ a−1 exp − , ba Γ(a) b
(3)
where a and b are the shape parameters of the distribution, and Γ(·) is the Gamma function. The first and second moments of these Gamma RVs are, E[γ ] = ab and var[γ ] = ab2 , respectively. On the one hand, as described in [1], large-scale fluctuations in the irradiance are generated by turbulent eddies with sizes ranging from the scattering disk to the outer scale L0 . Since the #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21822
scattering disk is given by L/ρ0 k, where L is the path length, ρ0 is the spatial coherent radius and k is the wave number, it increases when L grows, ρ0 decreases (i.e. the strength of turbulence increases), or both. When this occurs, the average size of the large-scale eddies grows and, as a consequence, the probability that signals received in each photodetector pass through the same eddies (common eddies) also rises, resulting in an increase in the channel correlation [8]. Taking into account this channel behavior and in order to achieve compact expressions, in this work, we assume that signals of the N receivers are deflected by the same eddies and, therefore, largescale effects are common for all of them, i.e. we assume that Xi = X, ∀i = 1...N. Note that, with this assumption, we consider the most unfavorable scenario for large-scale correlation, which is also useful as a benchmark in optical systems design. On the other hand, small-scale contributions to scintillation are primarily diffractive in nature and are associated with turbulent eddies smaller than ρ0 . In this case, unlike Xi , small-scale effects of Yi are dependent on each receiver. However, even in this case, it can be assumed that Yi , ∀i = 1...N, will exhibit the same statistical properties and, therefore, all Yi will be identically distributed. In this situation, Eq. (2) can be rewritten as N
I = X ∑ Yi = XV.
(4)
i=1
Here, X is governed by a Gamma distribution, with parameters αx and βx , denoted as X ∼ G (αx , βx ). Likewise, the N small-scale identically distributed variates Yi also follow a Gamma distribution, but with the i-th channel parameters αi = α and βi = β , i.e. Yi ∼ G (α , β ). Finally, V = ∑Ni=1 Yi is the RV that describes the small scale atmospheric effect over the combined received optical irradiance. Note that since I is a normalized optical irradiance, the large and small scale parameters must fulfill that αx = 1/βx and α = 1/(N β ), respectively, and therefore, E[I] = E[X]E[V ] = 1. Note that, both αx and α are related to the effective number of large-scale and small-scale eddies, respectively, and they range from values close to 1, for strong turbulence conditions, and values >> 1, for weak turbulence intensities, as detailed in [2]. When channel correlation is considered, the pdf of V is obtained by solving the pdf of the sum of correlated Gamma RVs. It is assumed that the correlation between the variates Yi and Y j , with i, j = 1..N is given by the correlation coefficient ρi j defined as cov (Yi ,Y j ) , ρi j = var (Yi ) var (Y j )
(5)
where 0 ≤ ρi j ≤ 1. Thus, the correlation information between the small scale in the N optical signals received can be expressed by means of a NxN correlation matrix Cy , given as ⎡ ⎤ ρ12 · · · ρ1N 1 ⎢ ρ21 1 · · · ρ2N ⎥ ⎢ ⎥ Cy = ⎢ . , (6) . .. ⎥ .. .. ⎣ .. . . ⎦ ρN1 ρN2 · · · 1 N×N where the value of the correlation coefficient ρi j depends on the distance between i-th and j-th receivers. Although the development presented here is valid for any system whose correlation matrix is expressed as Eq. (6), in order to analyze the accuracy of the proposed analytical expressions three different correlation models, which are related to three typical receivers geometries, are considered in this work. In particular, constant, exponential and circular models are applied. The first one, the constant correlation model is described in [15, 16], and it corresponds to the case of evenly spaced receivers, with ρi j = ρ ∀i, j = 1..N, being ρ the correlation coefficient #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21823
between the Gamma distributed signals received at adjacent photodetectors. The second one, the exponential correlation model is applied to a lineal equidistant photodetectors array, being ρi j = ρ |i− j| [15]. Finally, the circular correlation model corresponds to the case of N receivers situated on a circumference, presenting the correlation matrix a N-th circular symmetry [16]. 2.
Sum of correlated Gamma variates
As explained in the previous section, the first goal is to obtain a closed-form expression for the pdf of the sum of correlated Gamma RVs Yi , with i = 1...N. To this end, the Moschopoulos theorem [17] is extended in the way proposed by Alouini et al. [16], but with the restriction of being every Yi in Eq. (4) identically distributed, as previously supposed. Assuming that for an arbitrary function φ (x) its Laplace transform is denoted as L {φ (x); s}, then, from [16, Eq(9)], the moment generating function (MGF) of V can be expressed by using the statistical averaging operator E[·], as N
MV (s) = E[esV ] = L { fV (v); −s} = ∏(1 − λi s)−α ,
(7)
i=1
where {λi }Ni=1 are the eigenvalues of the matrix A = DC, being D a N × N diagonal matrix with the entries β for i = 1...N, and C a N × N positive definite correlation matrix given by ⎡ ⎤ √ √ ρ12 · · · ρ1N 1 √ √ ⎢ ρ21 1 ··· ρ2N ⎥ ⎢ ⎥ C=⎢ , (8) ⎥ .. .. .. .. ⎣ ⎦ . . . . √ √ ρN1 ρN2 · · · 1 N×N whose elements are the correlation coefficients of the underlying Gaussian processes that lead to the Gamma distributed fading, as described in [15, 16, 18]. Now, the cumulative density function (cdf) of V , FV (v), can be written in the following form 1 1 1 N L {FV (v); s} = L { fV (v); s} = MV (−s) = ∏(1 + λi s)−α , s s s i=1 which can be algebraically modified to obtain the identity N α ⎞−α ⎛ 1 1 ∏ N N − λi λi Γ (N α ) ∏ (1 + λi s)−α = i=1Γ (N α ) · sN α · ∏ ⎝1 − s ⎠ i=1 i=1
(9)
(10)
and by comparing it to [19, Eq. (5)], this new expression is obtained N
1
∏ (1 + λi s)−α = det (A)α Γ (N α ) · L
(N)
vN α −1 Φ2
α , . . . , α ; N α ; −vλ1−1 , . . . , −vλi−1
,
i=1
(11) (N) where the determinant det(A) = ∏Ni=1 λi and Φ2 is the confluent Lauricella function [19, 20]. Finally, from Eqs. (9) and (11), the pdf and cdf of V , can be obtained respectively by using the inverse Laplace transform vN α −1 −1 × fV (v) = L {MV (−s); v} = [det(A)]α Γ(N α ) (12) v v (N) , × Φ2 α , ..., α ; N α ; − , ..., − λ1 λN #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21824
and −1
FV (v) = L
3.
1 vN α MV (−s); v = × s [det(A)]α Γ(N α + 1) v v (N) . × Φ2 α , ..., α ; N α + 1; − , ..., − λ1 λN
(13)
Statistical distribution of the received irradiance I
Now, starting from the normalized received optical irradiance I = XV , as indicated in Eq. (4), the unconditional pdf for the normalized combined received irradiance, fI (I), is obtained by calculating the mixture of a Gamma distribution for X defined by Eq. (3) as G (αx , 1/αx ), fX (x) =
αx (αx x)αx −1 exp (−αx x), Γ(αx )
(14)
and the function fV (v) presented above in Eq. (12). Therefore, by first fixing x and writing v = I/x, the conditional pdf, fV (I|x), is given by fV (I|x) 1 = fV (I|x), |∂ I/∂ V | x
(15)
I N α −1 × xN α [det(A)]α Γ(N α ) I I (N) . α , ..., α ; N α ; − , ..., − × Φ2 xλ1 xλN
(16)
fI (I|x) = and then, fI (I|x) =
To obtain the unconditional irradiance distribution, fI (I), the average of Eq. (16) over the Gamma distribution given in Eq. (14) must be calculated, which leads to fI (I) =
∞ 0
fI (I|x) fX (x)dx =
I N α −1 αxαx I, α [det(A)] Γ(N α ) Γ(αx )
(17)
where I=
∞ 0
(N)
xαx −N α −1 · Φ2 exp (−αx x)dx
(18)
(N) (N) with Φ2 = Φ2 α , ..., α ; N α ; − xλI , ..., − xλIN . 1 At this point, a clever way to solve fI (I) is obtained by applying the Gauss-Laguerre quadrature [21, Eq. (25.4.45)] to approximate the value of the integral I. Following this approach, the irradiance pdf can be expressed as fI (I)
xiαx −N α I N α −1 αxαx n × ∑ [det(A)]α Γ(N α ) Γ(αx ) i=1 (n + 1)2 [Ln+1 (xi )]2 I αx I αx (N) × Φ2 α , ..., α ; N α ; − , ..., − xi λ1 xi λN
(19)
where xi is the i-th root of the Laguerre polynomial with degree n, denoted as Ln (x).
#214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21825
4.
Closed-form expression for fV (v) and fI (I) when α ∈ N
Due to the fact that the pdf equation derived in previous section can only be solved by numerical analysis and in order to achieve closed-form expressions for fV (v) and fI (I), in this work we assume that α ∈ N. It is worth noting that, from a physical viewpoint, this is a reasonable assumption since, as previously mentioned, in a GG distribution this parameter represents the effective number of small-scale turbulent eddies for each of the N receiving channels. Note that, for real values of α , the pdf functions will be bounded by the pdfs corresponding to the nearest natural values of α . Under this premise (α ∈ N), the integral I can be solved in a closed-form expression. In this way, by using the variable change t = 1/x, the Eq. (18) can be expressed as ∞ α x dt, f (t)t −αx exp − I= (20) t 0 (N) where f (t) = t N α −1 · Φ2 α , ..., α ; N α ; − λIt , ..., − λItN . To solve this expression, the function 1 f (t) can be rewritten by applying the Laplace transform. Thus, from [19, Eq. (4.24-5)], the Laplace transform of f (t) is given by F(s) = L { f (t); s} =
I −α Γ(N α ) N 1+ , sN α ∏ λi s i=1
which, after some manipulations, can be also expressed as −α N −I F(s) = Γ(N α ) ∏ s − . λi i=1
(21)
(22)
At this point, by using the partial fraction expansion method described in [22], the following conversion from product to summation of fractions can be applied N
N αi 1 cmi ∏ (s − di )αi = ∑ ∑ (s − di )m , i=1 i=1 m=1
(23)
where N is the number of different eigenvalues of matrix A, di = −I/λi . In addition, to consider the multiplicity of some eigenvalues, αi is the product of the algebraic multiplicity of the eigenvalue, denoted as μA (λi ), with the parameter α , i.e. αi = μA (λi )α . Finally cmi is the corresponding fraction coefficient, whose detailed calculation is described in Appendix A. Then, the application of this expansion technique over the Eq. (22) results in the following expression N αi cmi . (24) F(s) = Γ(N α ) ∑ ∑ (s − di )m i=1 m=1 Now, the inverse Laplace transform of F(s) can be calculated. For this purpose, knowing that 1 tn −1 exp (−at)u(t), (25) ;t = L (s + a)n+1 n! with u(t) being the Heaviside step function, the inverse Laplace transform of the double sum argument in Eq. (24) is derived as 1 t m−1 −1 exp (dit)u(t), (26) L m ;t = (s − di ) (m − 1)! #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21826
and, consequently, f (t) = L−1 {F(s);t} is finally expressed as follows, noting that di was previously defined as di = −I/λi , N αi I t m−1 exp − t . f (t) = Γ(N α ) ∑ ∑ cmi (27) (m − 1)! λi i=1 m=1 Now, including f (t) in Eq. (20), and by employing [23, Eq. (3.471-9)], the Eq. (20) can be solved as N αi m−αx m−αx cmi αx I − m−2αx 2 2 . (28) λi αx I Km−αx 2 I = 2Γ(N α ) ∑ ∑ λi i=1 m=1 (m − 1)! And finally, the pdf of the combined received irradiance is expressed, from Eq. (17), as
fI (I) =
N αi 2 cmi m−2αx m+2αx λ αx × ∑ ∑ α [det(A)] Γ(αx ) i=1 m=1 Γ(m) i
×I
N α −1− m−2αx
αx I , Km−αx 2 λi
(29)
providing an easy-to-compute closed-form analytical method to obtain the statistical properties of the combined received irradiance I from the turbulence parameters and the channel correlation information. It must be noted that, from Eq. (50), the parameter cmi is also dependent on the normalized irradiance I. Additionally, from the assumption of α ∈ N and returning to the definition of the pdf of the sum of correlated Gamma RVs, given in Eq. (12), the mathematical procedure shown in Eqs. (23) to (27) can be applied to obtain an equivalent closed-form expression of fV (v). Then, from Eq. (7) and Eq. (23) N αi cmi 1 −1 (30) fV (v) = L ∑ ∑ (s − di )m ;t , [det(A)]α i=1 m=1 where, in this case, di = −λi−1 and the coefficients cmi are conveniently calculated with the procedure described in Appendix A. By using again the inverse Laplace transform given in Eq. (26), the pdf for the sum of Gamma RVs can be calculated as N αi 1 v cmi m−1 (31) fV (v) = ∑ ∑ Γ(m) v exp − λi . [det(A)]α i=1 m=1 An interesting result is reached when previous equation is manipulated to obtain the following equivalent expression N αi cmi fV (v) = ∑ ∑ λ m gm,i (v), (32) α i V [det(A)] i=1 m=1 where gVm,i (v) ∼ G (m, λi ). Then, fV (v) is said to be a mixture of Gamma RVs with shaping parameters m and λi . The summation of the corresponding coefficients verifies that N
αi
cmi λ m = 1, α i [det(A)] i=1 m=1
∑∑
(33)
since, in the particular case of Eq. (31), the coefficients cmi are independent from I. By contrast, regarding to fI (I), as previously mentioned, these coefficients are dependent on I and the equivalent mixture definition cannot be applied to Eq. (29). #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21827
5.
Numerical results
In this section, the expressions developed in the previous section are validated. To this end, analytical results obtained from these expressions are compared to numerical results obtained by simulations. In particular, several cases corresponding to different levels of correlation among the N channels have been considered. First, the extreme cases of no correlated and total correlated channels are analyzed. Latter, three cases of partially correlated channels are also studied. In all these cases, EGC method is assumed for combining the received irradiance in each photodetector. It is worth noting that, in order to obtain the numerical simulation pdf results, the histogram pdf estimator described in [24] is employed. Cases of no correlated and total correlated channels: When there is no correlation among the small scale variates, the Yi signals are independent and Cy is then a diagonal correlation matrix. Consequently, the eigenvalues of A are λi = λ = β , ∀i = 1..N, and the det(A) = ∏Ni=1 λi = λ N α . This case corresponds to a single-input singleoutput (SISO) channel where the small scale effects are modeled by a Gamma distribution whose shape parameter is α = N α , i.e. V ∼ G (N α , 1/(N α )). Thus, by applying the calculus procedure described in previous sections, the pdf of the received normalized irradiance is expressed as fI (I) =
N α + αx N α + αx 2 αx 2 λ − 2 × Γ(N α )Γ(αx )
×I
N α + αx −1 2
KN α −αx
2
αx I λ
(34) .
and, by substituting here λ = β = 1/N α , leads to a GG pdf with large and small scale parameters αx and N α , respectively, given by N α + αx
fI (I) =
√ N α + αx 2 (αx N α ) 2 I 2 −1 KN α −αx 2 αx N α I . Γ(N α )Γ(αx )
(35)
Likewise, when the optical irradiance at the N receivers is completely correlated, then, every element of the correlation matrix ρi j = 1 and, thus, all the random sequences Yi are exactly the same, i.e. Yi = Y ∀i = 1..N, and V = NY , where Y ∼ G (α , 1/(N α )). Then, by applying the scaling of RVs, fV (v) = fY (v/N)/N, the pdf of V results in a Gamma distribution such V ∼ G (α , 1/α ) and the pdf of I follows the corresponding gamma-gamma distribution with large and small scale parameters αx and α , respectively, expressed as α + αx
√ 2 (αx α ) 2 α +αx −1 I 2 fI (I) = Kα −αx 2 αx α I . Γ(α )Γ(αx )
(36)
Figure 2 shows the obtained results for the pdf of the normalized combined optical irradiance I. In this figure, analytical results calculated from Eq.(35) and Eq. (36) for ρ = 0 and ρ = 1 are represented with solid line and dashed line, respectively, whereas numerical results, obtained by generating the corresponding Gamma RVs are depicted with circles. Typical moderate turbulence intensity (αx = 10 and αy = 5) and N = 4 optical receivers are assumed in the comparison. From this figure, it is clearly observed that analytical and numerical results curves fit perfectly, corroborating the validity of the mathematical expressions proposed for these two limit cases.
#214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21828
Numerical results Analytical pdf ρ=1 Analytical pdf ρ=0
I
Probability density function, f (I)
1
0.8
0.6
0.4
0.2
0
0.5
1
1.5
2
2.5
3
Normalized total received irradiance, I
Fig. 2. Analytical and simulation results for the pdf of the normalized irradiance I for the extreme correlation cases of ρ = 0 and ρ = 1. N = 4 optical receivers and typical moderate turbulence intensity with αx = 10 and αy = 5 are assumed.
Cases of partially correlated channels: In contrast to the previous cases, to analyze realistic ρi j values, the generation of correlated Gamma sequences must be considered. For this purpose, a method to generate correlated Gamma RVs from linear combinations of independent Gamma variates, proposed in [25], is employed. The method is based on the Cholesky decomposition and on matching the first and second order moments of both independent and correlated Gamma set. Firstly, by applying the Cholesky decomposition, the correlation matrix Cy can be expressed as the product of a lower triangular matrix L by its transposed LT , i.e. Cy = LLT . If y = [y1 , . . . , yN ]T is the set of N correlated Gamma RVs with correlation matrix Cy , and g = [g1 , . . . , gN ]T is a set of N independent Gamma variates with parameters (αg,i , βg,i ), being its correlation matrix a N × N identity matrix, then y = Lg or yi =
i
∑ lik gk ,
(37)
k=1
where lik denotes the (i, k)-th position of L. Secondly, the moment matching procedure is applied, which is based on the knowledge of the moments values for the correlated Gamma variables yi , i.e. the first order moment of yi is E[yi ] = αβ , ∀i = 1..N, whereas the second order moment is var[yi ] = αβ 2 , ∀i = 1..N, being α and β the parameters defined in section 1. In addition, from Eq. (37) and taking into account that the variance of the sum of independent RVs is the sum of each RV variance, the first and second order moments of yi can be respectively expressed as E[yi ] =
i
∑ lik E[gk ],
(38)
k=1
and var[yi ] =
i
∑ lik2 var[gk ].
(39)
k=1
Then, the characteristics parameters αg,i and βg,i , are iteratively obtained from Eq. (38) and Eq. (39) as αβ − ∑i−1 k=1 lik E[gk ] = αg,i βg,i , (40) E[gi ] = lii #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21829
1.4
ρ=0.2 ρ=0.4 ρ=0.8
1.2
0.8
fI(I)
0.8
V
f (v)
1
ρ=0.2 ρ=0.4 ρ=0.8
1
0.6
0.6 0.4
0.4 0.2
0.2
0
0 0.2
0.4
0.6
0.8
1
1.2
1.4
1.6
1.6
1.2
0.5
1
1.5
2
2.5
Normalized irradiance, I (b) 1.2
σ 2=0.17 I σ 2=0.32 I σ 2=1
1
I
I
0.8
fI(I)
1 V
2
σI2=0.17 σ 2=0.32 I σ 2=1
1.4
f (v)
1.8
Sum of correlated Gamma RVs, V (a)
0.8
0.6
0.6 0.4 0.4 0.2
0.2 0.5
1
1.5
Sum of correlated Gamma RVs, V (c)
2
2.5
0
0.5
1
1.5
2
2.5
Normalized irradiance, I (d)
Fig. 3. Analytical and simulation results for the pdf of V and I RVs obtained with the constant correlation model and N = 3 equidistant optical receivers. In (a) and (b), the results are shown for different correlation levels with a typical turbulence intensity σI2 = 0.32. In (c) and (d) the results are shown as a function of turbulence conditions (weak, moderate and strong) with a correlation level given by ρ = 0.6.
and var[gi ] =
2 αβ 2 − ∑i−1 2 k=1 lik var[gk ] = αg,i βg,i , 2 lii
(41)
and thus, the parameters are calculated respectively as βg,i = var[gi ]/E[gi ] and αg,i = E[gi ]/βg,i . Finally, once the correlated RVs are generated and combined with the EGC technique, the numerical pdf is obtained in each case, by applying a histogram based pdf estimator, as mentioned before. The obtained results are compared to the analytical pdf calculated with the closed-form expressions proposed in this paper. Note that in the following figures, the turbulence intensity is given by the scintillation index, defined in [1] and particularized for the GG statistical distribution here studied as σI2 = αx−1 + α −1 + (αx α )−1 . Analytical and numerical results obtained for the constant correlation model are shown in Fig. 3. In this figure, dotted, dashed and solid lines represent analytical results, whereas circles depict numerical results. For this correlation model, N = 3 equidistant optical receivers are assumed in all cases. The simple correlation matrix is given by ⎤ ⎡ 1 ρ ρ (42) Cy = ⎣ ρ 1 ρ ⎦ . ρ ρ 1 In Fig. 3(a) and Fig. 3(b), the analytical and the simulation results for the pdf of V and I RVs are presented for different levels of correlation among the receivers: from low correlation (ρ = 0.2), which corresponds to separated photodetectors, to high correlation (ρ = 0.8) which corresponds to closer optical receivers. Likewise, in Fig. 3(c) and Fig. 3(d), the same pdfs are plotted but now for different levels of turbulence intensity, in particular, from weak (σI2 = 0.17, with {α = 15, αy = 10}), to moderate (σI2 = 0.32, with {αx = 10, αy = 5}) and #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21830
1.2
ρ=0.2 ρ=0.4 ρ=0.8
1.6 1.4
0.8
1
I
f (I)
V
f (v)
1.2
ρ=0.2 ρ=0.4 ρ=0.8
1
0.8
0.6 0.4
0.6 0.4
0.2
0.2 0 0.2
0.4
0.6
0.8
1
1.2
1.4
1.6
1.8
2
0.5
Sum of correlated Gamma RVs, V (a)
2
2
σI =0.17 2 σ =0.32 I2 σ =1 I
1 0.8
f (I)
I
1
I
V
f (v)
1.5
1.2
2
σI =0.17 2 σ =0.32 I2 σ =1
1.5
1
Normalized irradiance, I (b)
0.6 0.4
0.5
0.2 0
0.5
1
1.5
2
Sum of correlated Gamma RVs, V (c)
2.5
0
0.5
1
1.5
2
2.5
Normalized irradiance, I (d)
Fig. 4. Analytical and simulation results for the pdf of V and I RVs obtained with the exponential correlation model and N = 5 optical receivers. In (a) and (b), the results are shown for different correlation levels with a typical turbulence intensity σI2 = 0.32. In (c) and (d) the results are shown as a function of turbulence conditions (weak, moderate and strong) with a correlation level given by ρ = 0.6.
strong (σI2 = 1, with {α = 15, αy = 10}) turbulence intensities, respectively, considering a fixed correlation level with ρ = 0.6. From Fig. 3, it follows that, regardless of the correlation levels or the turbulence conditions, there is a very clear agreement between obtained analytical results and numerical simulations. Thus, these results corroborate the high accuracy of proposed closed-form expressions. Now, the same comparison is extended for the exponential and circular correlation models. In Fig. 4, analytical an numerical results for the exponential model are depicted. Here, a linear photodetector array with N = 5 receivers is assumed, with the following correlation matrix, ⎡ ⎤ 1 ρ ρ2 ρ3 ρ4 2 ρ3 ⎥ ⎢ ρ ⎢ 2 1 ρ ρ ⎥ ρ 1 ρ ρ2 ⎥ (43) Cy = ⎢ ⎢ ρ ⎥. ⎣ ρ3 ρ2 ρ 1 ρ ⎦ 1 ρ4 ρ3 ρ2 ρ Again, from Figs. 4(a)–4(d) the high level of accuracy obtained with the proposed closedform expressions is observed, regardless of the correlation levels or the turbulence conditions. Finally, in Fig. 5 the results for the circular correlation model are presented. In this case, the same numerical simulations carried out in [16] are performed. Thus, a N = 5 circular receivers array is considered with a correlation matrix defined as ⎡ ⎤ 1 ρ2 ρ3 ρ4 ρ5 ⎢ ρ5 1 ρ2 ρ3 ρ4 ⎥ ⎢ ⎥ ⎥ (44) Cy = ⎢ ⎢ ρ4 ρ5 1 ρ2 ρ3 ⎥ , ⎣ ρ3 ρ4 ρ5 1 ρ2 ⎦ ρ2 ρ3 ρ4 ρ5 1 #214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21831
0.8
ρ =ρ =0.64 12 15 ρ =ρ =0.36
0.6
0.6
13
ρ12=ρ15=0.8 ρ =ρ =0.6 13
14
14
ρ =ρ =0.64 12 15 ρ =ρ =0.36
I
0.8
f (I)
14
13
1
fV(v)
1
ρ =ρ =0.8 12 15 ρ =ρ =0.6
1.2
13
14
0.4 0.4 0.2
0.2
Numerical results
Numerical results 0
0.5
1
1.5
0
2
Sum of correlated Gamma RVs, V (a) 1.6
fV(v)
2
2.5
2
σI =0.17 2 σ =0.32 I2 σ =1 I
0.8
I
Numerical results
1.5
1
f (I)
1 0.8
1.2
1
Normalized irradiance, I (b)
1.2
2
σ =0.17 I2 σ =0.32 I2 σ =1 I
1.4
0.5
0.6
Numerical results
0.6 0.4
0.4 0.2
0.2 0 0.5
1
1.5
2
Sum of correlated Gamma RVs, V (c)
2.5
0
0.5
1
1.5
2
2.5
3
Normalized irradiance, I (d)
Fig. 5. Analytical and simulation results for the pdf of V and I RVs obtained with the circular correlation model and N = 5 optical receivers. In (a) and (b), the results are shown for different correlation levels with a typical turbulence intensity σI2 = 0.32. In (c) and (d) the results are shown as a function of turbulence conditions (weak, moderate and strong) with a correlation level given by ρ2 = ρ12 = ρ15 = 0.8 and ρ3 = ρ13 = ρ14 = 0.6.
where, for symmetry, ρ2 = ρ5 and ρ3 = ρ4 . In particular, Fig. 5(a) and Fig. 5(b), show the pdf of V and I RVs for two correlation cases. In the first one, ρ2 = ρ12 = ρ15 = 0.8 and ρ3 = ρ13 = ρ14 = 0.6; in the second case, ρ2 = ρ12 = ρ15 = 0.64 and ρ3 = ρ13 = ρ14 = 0.36. Furthermore, Fig. 5(c) and 5(d) depict the effect of turbulence conditions on the accuracy of the proposed expressions. Again, as expected, a very high agreement between the numerical simulations and the analytical results is clearly observed. 6.
Concluding remarks
In this paper, the analysis of the probability density function of the received optical irradiance for a spatially diversity free space optical communication system has been performed. The gamma-gamma statistical distribution has been assumed to model the turbulence induced scintillation that produce the fading on the optical received signal, whereas the combination procedure used is the equal-gain combining scheme. In this context, novel analytical closedform expressions have been derived for the pdf of the sum of identically distributed spatially correlated Gamma random variables and then, for the corresponding sum of correlated gammagamma variates. The mathematical expressions here proposed have been corroborated with Monte-Carlo numerical simulations, resulting in an excellent agreement between numerical and analytical results. This comparison have been performed for the limiting cases of no correlation (ρ = 0) and fully correlated received signals (ρ = 1), as well as with three realistic correlations models (constant, exponential and circular), corresponding to three different receivers configurations. The study here presented corroborates the high accuracy of the proposed closed-form expressions, leading to a mathematically tractable and easily computable method to characterize the statistical behavior of the received optical irradiance in a spatial diversity context.
#214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21832
Appendix A In this Appendix, the calculation procedure of the partial fraction expansion coefficients, cmi , is described. From the method presented in [22], these coefficients are obtained from ⎡ ⎤
cmi =
dαi −m ⎢ N 1 1 ⎥ · α −m ⎣ ∏ αj ⎦ i (αi − m)! dw (w − d ) j j=1 j =i
,
(45)
w=di
where, from Eq. (23), N is the number of different eigenvalues of matrix A, d j = −I/λ j , and αi = μA (λi )α . Note that μA (λi ) was defined as the algebraic multiplicity of the eigenvalue i. Here, for the computing of each coefficient cmi , an extension of the Leibniz product rule based on the multinomial theorem can be applied. In this sense, from the generalized Leibniz rule defined as d n (u · v) = (u + v)(n) (46) dxn being u = u(x) and v = v(x), and the multinomial theorem given by (u1 + . . . + um ) =
∑
n
k1 +...+km =n
n k1 , . . . , km
N
∏ utkt ,
(47)
j=1
the n-th derivative of a product operator can be expanded in the following form N N dn n (k ) (n) = (u u + . . . + u ) = j N 1 ∑ ∏ uj j . k dxn ∏ , . . . , k N 1 j=1 j=1 k +...+kN =n
(48)
1
In this case, identifying Eq. (48) in Eq. (45), u j = (w − d j )−α j and its n-th derivative is (n) easily computed as u j = (−1)n (α j )n (w − d j )−α j −n , where (α j )n = Γ(α j + n)/Γ(α j ) is the Pochhammer symbol. Therefore, N
(k j )
∏ uj
=
j=1 j =i
N
∏
j=1 j =i
(−1)k j (α j )k j (w − d j )−α j −k j ,
(49)
and finally, the coefficient cmi can be expressed as 1 cmi = (αi − m)! k ×
N
∏
j=1 j =i
∑
1 +... +kN =αi −m iˆ
αi − m × k1 , . . ., kN iˆ
−α j −k j
(−1) (α j )k j (di − d j ) kj
(50)
,
where iˆ indicates that ki is omitted in the previous corresponding sequences, and w has been substituted for di , as indicated in Eq. (45). Acknowledgment This work was supported by the Spanish Ministerio de Econom´ıa y Competitividad, Project TEC2012-36737.
#214451 - $15.00 USD Received 4 Jul 2014; revised 14 Aug 2014; accepted 16 Aug 2014; published 2 Sep 2014 (C) 2014 OSA 8 September 2014 | Vol. 22, No. 18 | DOI:10.1364/OE.22.021820 | OPTICS EXPRESS 21833
