-
近20年来,研究者对大气湍流折射率的性质表现出了极大地兴趣,提出了各种理论模型来研究其特性。其中,针对Kolmogorov谱的分析最为广泛,Kolmogorov的湍流理论认为,对于局域均匀且各向同性的大气湍流,存在一个惯性区域,该区域的折射率波动的特征可以用大气湍流折射率谱密度函数${{{\varPhi }}_n}\left( {{k_r}} \right)$来描述,其中${k_r} = {{2\pi }}/{\lambda }$是空间波函数,单位为$ {{\rm{m}}}^{-1},\lambda $表示光波长。Kolmogorov谱可以适用于绝大多数近地FSO系统,其具体表达式为[15]:
$$ \varPhi_n\left(k_r, h\right)=\frac{0.033 C_n^2}{k_r^{\tfrac{11}{3}}},\; \frac{2 \pi}{L_0}< k_r <\frac{2 \pi}{l_0} $$ (1) 式中:11/3为低空湍流谱指数;${L_0}$为折射率波动的外部相关长度,称为外尺度;${l_0}$为内尺度,表示内部相关长度[16];$C_n^2$为折射率结构常数,是大气光学中的一个非常重要的参数,它描述了大气折射率起伏的强度,可以用它来表示湍流的强弱。在近地面的FSO系统中,大气湍流在传输路径上是均匀分布的,即同一时刻大气折射率结构常数$C_n^2$在整条传输链路上保持不变[3]。
而在地星链路上,激光传输贯穿整个大气层,需要采用大气折射率结构常数的海拔分布模型,即$C_n^2\left( h \right)$。Hufnagel-Valley在实测的基础上,根据实验数据,提出了大气湍流强度$C_n^2$随高度变化的模型,表示为[17]:
$$ \begin{split} C_n^2\left( h \right) =& 0.005\;94{\left( {\frac{v}{{27}}} \right)^2} \times {\left( {{{10}^{ - 5}}h} \right)^{10}}\mathop {\rm{e}}\nolimits^{\tfrac{{ - h}}{{1\;000}}} + 2.7 \times {10^{ - 16}} \times \\ & \mathop {\rm{e}}\nolimits^{\tfrac{{ - h}}{{1\;500}}} + A\mathop {\rm{e}}\nolimits^{\tfrac{{ - h}}{{100}}} \left[ {{{\rm{m}}^{ - {2}/{3}}}} \right] \end{split} $$ (2) 式中:$h$为距离地面的海拔高度;$v\left( {{\rm{m}}/{\rm{s}}} \right)$为均方根风速。当$A = 1.7 \times {10^{ - 14}}\left( {{{\rm{m}}^{ - {2}/{3}}}} \right)$,$v = 21\left( {{\text{m}}/{\rm{s}}} \right)$ 时称为Hufnagel-Valley 5/7模型,可以用来表示白天的大气湍流强度随高度的变化,该模型被广泛地应用于大气湍流仿真建模研究中[18]。因此,文中有关大气湍流强度的计算也是基于该模型。
高层大气中,Kolmogorov谱不能正确地描述大气湍流统计[12, 19-20],特别是在垂直传输路径上,湍流是不均匀的,其湍流谱的谱指数也随高度变化[19],而不是固定的11/3。也就是说,无论是在近地面受地面因素影响明显的情况下还是在高空自由的大气中,大气湍流大都与各向同性特性有一定的偏离,任何偏离局地均匀各向同性的湍流都可称为Non-Kolmogorov湍流。大气湍流的功率谱指数在不同的高度、气压、温度和湿度等因素下可能会有所不同,通常来说,随着高度的增加,大气湍流的功率谱指数会逐渐减小[21]。
Zilberman等人将卫星激光通信的大气信道看作是由地面到2 km的Kolmogorov湍流、2~10 km的谱指数为10/3的Non-Kolmogorov湍流以及10 km以上谱指数为5的Non-Kolmogorov湍流组成,同时给出了在该模型下谱指数随高度变化的表达式[22]:
$$\begin{split} \alpha \left( h \right) =& \frac{{{\alpha _1}}}{{1 + {{(h/{H_1})}^{{b_1}}}}} + \frac{{{\alpha _2} \cdot {{(h/{H_1})}^{{b_1}}}}}{{1 + {{(h/{H_1})}^{{b_1}}}}} \cdot \frac{1}{{1 + {{(h/{H_2})}^{{b_2}}}}} +\\ & \frac{{{\alpha _3} \cdot {{(h/{H_2})}^{{b_2}}}}}{{1 + {{(h/{H_2})}^{{b_2}}}}} \end{split}$$ (3) 式中:$ {\alpha _1} $、$ {\alpha _2} $和$ {\alpha _3} $分别为边界层、对流层和平流层内的湍流谱指数;${{{H}}_1}$和${{{H}}_2}$为大气层边界;$ {b_1} $和$ {b_2} $为层间平滑系数。
但是,该模型的准确性尚未得到验证,实际上,具体的谱指数随高度的变化趋势可能因不同的环境条件和研究区域而有所不同,因此文献[22]中的湍流谱模型不具有普适性,具体的研究需要依据实际情况和观测数据来确定。根据现有的测量与实验研究,自由对流层中的湍流比平流层中更加强烈,以我国青藏地区和青海湖地区为例,青藏地区边界层的谱指数在2.4~2.7之间,自由对流层的谱指数大约在2.9~3.1之间,平流层的谱指数大约在2.9~2.7之间[23]。青海湖地区在低空、中空和高空的谱指数则分别为3.67、3.47和3.20[24]。
因此,文中采用可以表示任意谱指数的Non-Kolmogorov湍流折射率波动的三维功率谱模型 ,不考虑湍流内外尺度时,具体地表达式如下[12, 25-26]:
$$ {\varPhi _n}\left( {{k_r},h,\alpha } \right) = A\left( \alpha \right)\beta (h){k_r}^{ - \alpha } $$ (4) $$ A\left( \alpha \right) = \varGamma \left( {\alpha {{ - 1}}} \right)\cos \left( {\alpha \pi {\text{/2}}} \right){{\bigg/}}\left( {{\text{4}}{\pi ^{\text{2}}}} \right) $$ (5) 式中:$\; \beta \left(h\right) $为Non-Kolmogorov谱的折射率结构常数,是随高度和湍流谱指数变化的参数,与$C_n^2\left( h \right)$的关系为:
$$ \beta (h) = \frac{{A \left({{11}}/{3}\right)}}{{A(\partial )}}C_n^2(h){( \sqrt{k / L})^{(\partial - 11/3)}} $$ (6) 式中:$L$为总传输距离。
-
分步传输算法(split-tep beam propagation method)被广泛的用于光束在随机介质中的数值仿真。主要原理是将光在随机介质中的传输分解成一段段真空与相位屏的共同作用。当传输距离$ \Delta z $大于外尺度$ {L_0} $时,可以将公式(4)空间折射率谱${\varPhi _n}\left( {{k_r},h,\alpha } \right)$写成路径依赖的形式:
$$ {\varPhi _n}\left( {{k_p}} \right) = {\varPhi _n}({k_p},{k_z} = 0)\int_{{z_i}}^{{z_{i + 1}}} {\sigma _n^2\left( z \right){\rm{d}}z} $$ (7) 式中:$ \displaystyle\int_{{z_i}}^{{z_{i + 1}}} {\sigma _n^2\left( z \right){\rm{d}}z} $为传输距离$ \left[ {{z_i},{z_{i + 1}}} \right] $内折射率波动引起的方差。
“薄”相位屏是图1中湍流层的厚度远远小于沿层传播的距离,采用“薄”相位屏时不需要对层的厚度进行积分,从而消除了上式对传输距离的依赖。因此,在数值仿真中,不需要产生全三维的折射率$\mathop \eta \nolimits_1 \left( {x,y,z} \right)$,而是采用二维的随机相位屏$\varphi \left( {{x_j},{y_j}} \right)$对光场进行调制。随机相位屏的相位功率谱密度函数与大气折射率的功率谱密度的关系为:
$$ {\varPhi _\varphi }\left( {{k_p}} \right) = 2\pi {\varPhi _n}({k_p},{k_z} = 0) $$ (8) 根据数值仿真的采样点数$N$与采样间隔$x$,可以得到随机相位屏频谱的方差${\sigma ^2}\left( {{k_x},{k_y}} \right)$:
$$ {\sigma ^2}\left( {{k_x},{k_y}} \right) = {\left( {\frac{{2\pi }}{{Nx}}} \right)^2}{\varPhi _\varphi }\left( {{k_x},{k_y}} \right) $$ (9) 再通过傅里叶反演法(FFT)[27-28],对相位屏频谱的方差进行傅里叶反变换,得到随机相位屏在x, y平面的时域形式:
$$ \varphi \left( {x,y} \right) = IFFT\left( {C\sqrt {{\sigma ^2}\left( {{k_x},{k_y}} \right)} } \right) $$ (10) 式中:$C$为一个均值为0、方差为1的$N \times N$的复随机数矩阵。
在软件仿真中,公式(10)需要写成傅里叶级数的形式:
$$ \begin{split} {{\varPhi }}\left( {jx,ly} \right) =& \mathop \sum \limits_{m = 0}^{{N_x}} \mathop \sum \limits_{n = 0}^{{N_y}} C\sqrt {{k_x}{k_y}} \times \sqrt {{{\varPhi }}\left( {m{k_x},{{n}}{k_y}} \right)} \times\\ &{\text{exp}}\left( {2\pi i\left( {\frac{{jm}}{{{N_x}}} + \frac{{ln}}{{{N_y}}}} \right)} \right) \end{split} $$ (11) 式中:${k_x} = \dfrac{{2\pi }}{{Nx}}$, ${k_y} = \dfrac{{2\pi }}{{Ny}}$为波数增量。
傅里叶变换法模拟湍流相位屏的方法虽然简单,但是生成的相位屏会产生低频不足的现象,因为该方法生成的相位屏的最小频率为${k_x}$,而不包含$\left( { - \dfrac{{{k_x}}}{2},\dfrac{{{k_x}}}{2}} \right)$和$\left( { - \dfrac{{{k_y}}}{2},\dfrac{{{k_y}}}{2}} \right)$这部分频率的贡献,从而造成相位屏低频不足。Lane等提出了一种增加次谐波的方法能够对 FFT 法生成的相位屏的低频进行补偿,该方法生成的次谐波表示为:
$$ \begin{split} {\varphi _{S H}}\left( {j{{\Delta }}x,l{{\Delta }}y} \right) =& \mathop \sum \limits_{p = 1}^{{N_p}} \mathop \sum \limits_{m = 0}^{{N_x}} \mathop \sum \limits_{n = 0}^{{N_y}} C\sqrt {{k_x}{k_y}} \times \sqrt {{{\varPhi }}\left( {m{k_x},{\text{n}}{k_y}} \right)} \times\\ &{\text{exp}}\left( {2\pi i\left( {\frac{{jm}}{{{3^p}{N_x}}} + \frac{{ln}}{{{3^p}{N_y}}}} \right)} \right) \end{split} $$ (12) 式中:$P$为次谐波级数,次谐波频率间隔分别为${k_{xp}} = {k_x}/{3^p}$, ${k_{yp}} = \dfrac{{{k_y}}}{{{3^p}}}$。值得注意的是,不同的湍流屏对低频补偿的要求不同,特别是Non-Kolmogorov谱,低频补偿的级数取决于谱指数的大小。在具体的数值仿真中,应根据仿真条件选择合适的补偿级数[28]。
因此,最终观测面的光场输出为:
$$ {U_{out}} = U\left( {{Z_n},{Z_{n - 1}}} \right)U\left( {{Z_{n - 1}},{Z_{n - 2}}} \right) \cdots U\left( {{Z_1},{Z_0}} \right){U_{in}} $$ (13) 式中:$\left( {{Z_{j + 1}},{Z_j}} \right) \cong {\text{exp}}\left( {i{\alpha _j}{{\Delta }}{z_j}A} \right) {\text{exp}}\left( {i\varphi \left( {{x_j},{y_j}} \right) } \right) {\text{exp}} ( i( {1 -} { {\alpha _j}} ) \cdot {{\Delta }}{z_j}A )$,表示光束在第j的传输,${\text{exp}} \left( {i{\alpha _j}{{\Delta }}{z_j}A} \right)$和${\text{exp}} \left( {i\left( {1 - {\alpha _j}} \right){{\Delta }}{z_j}A} \right)$为相位屏前以及相位屏后的传播算子,${{\Delta }}{z_j}$第$j$次传输距离,${\alpha _j}$为第$j$ 次传输的比例因子,$\varphi \left( {{x_j},{y_j}} \right)$表示第$j$次传输时的湍流影响;${U_{in}}$表示原始光束。
-
在地对卫星自由空间光通信链路中,大气湍流的影响主要集中在海拔20 km以下部分,海拔20 km以上可以近似为真空传输[6]。假设地星激光传输的天顶角(光线入射方向与天顶方向的夹角)为$\zeta $,则地星链路实际长度为:
$$ {L_T} = \sec \left( \xi \right) \times {{H}} $$ (14) 针对大气湍流谱指数随高度变化的特性,不再考虑整条路径上单一的大气折射率强度,而是按照大气特性对卫星与地面之间连续的湍流路径进行分层-边界层(0~2 km)、对流层(2~10 km)、平流层(10~20 km),每一层都采用一个独立的湍流强度,边界层的湍流强度最大,平流层的湍流强度最小。每一层采用谱指数不同的Non-kolmogorov谱模拟。
同时,为了降低计算机的运行时间,提高仿真效率,分步传输次数越少越好。假设每一层的大气湍流只采用一个相位屏来模拟,每一个相位屏都放在大气层的边界,如图2所示,至少也需要3层相位屏,文中在每一个大气层的边界,即2、10、20 km各放一个相位屏,分步传输三次,具体的传输模型如图2所示。
-
分步传输模型绝大多数都基于菲涅耳衍射,菲涅耳积分内的二次相位因子不受带宽限制,为了防止混叠,保证仿真精度,需要对单次传输距离以及采样间距等参数进行约束,使其满足奈奎斯特准则。
假设源平面和观测平面的采样间隔分别为$ {\delta _1} $、$ {\delta _n} $($ {\delta _n} > {\delta _1} $),采样点数$ N $,光束的波长为$ \lambda $,则相位屏的采样间隔以及分步传输的最大传输距离要满足:
$$ \Delta {z_i} \leqslant \frac{{\min {{({\delta _1},{\delta _n})}^2}N}}{\lambda } $$ (15) 观察公式(15)可以看出,分步传输的最大传输距离仅取决于$ {\delta _1} $、$ N $以及$ \lambda $。当激光波长$ \lambda $以及采样点数$ N $唯一确定时,最大传输距离仅取决于源平面的采样间隔$ {\delta _1} $。因此,只需选择合适的$ {\delta _1} $就可以保证仿真的精度。
在采用分步传输算法进行数值模拟时,是将路径中湍流大气引起的相位扰动叠加到光波波前上。对于均匀湍流路径的光传输,一般均假设相位屏间湍流强度均匀即折射率结构函数$\hat {C_n^2}$相等,该假设不会给模拟结果带来误差。而对于非均匀湍流路径的光传输,该假设将会造成一定的模拟误差,因此需要合理选择$\hat {C_n^2}$。根据文献[29]的分析,发现当相位屏间$\hat {C_n^2}$为实际湍流路径折射率结构常数的路径平均时,模拟误差最小,因此,在文中的仿真中,大气湍流折射率结构常数的选取均采用该方法。
$$ \hat{{{C}_{n}^{2}}(\Delta {z}_{i})}=\frac{1}{\Delta {z}_{i}}{\displaystyle \underset{{z}_{i-1}}{\overset{{z}_{i}}{\int }}{C}_{n}^{2}(z){\rm{d}}z} $$ (16) 为简单起见,文中只考虑从地面到海拔20 km湍流部分的光传输特性。假设天顶角$\zeta = 0$,激光在湍流大气中的垂直向上传播距离为20 km,激光波长选择550 nm,采样点数512,观测面的望远镜孔径为0.3 m,折射率结构参数模型采用HV5/7模型,初始光束的波束半径为0.05 m,每个相位屏均放置在大气层的边界位,根据仿真约束条件,仿真参数设置如表1所示。
表 1 仿真参数
Table 1. Simulation parameters
Parameters Notion Value Wavelength of Gauss beam $ \lambda $ 0.5 μm Gauss beam radius ${{w} }_{0}$ 0.05 m Zenith angle $ \theta $ 0° Outer scale ${{L} }_{0}$ 50 m Inner scale ${{l} }_{0}$ 0.001 m Turbulence model ${\rm{HV}}5/7$ ${\rm{HV}}5/7$ Transmitter grid spacing $ {\delta }_{t} $ 0.0035 m Receiver grid spacing $ {\delta }_{n} $ 0.005 m Sample points $ N $ 512 Scaling factor $ {\alpha }_{j} $ 0 Diameter of the observation aperture $ {D}_{2} $ 0.5 m 文中使用的发射光束为最低阶横向电磁(TEM)高斯束波$TM{E_{00}}$波,其发射孔径位于$z = 0$平面,在该平面内其振幅分布为高斯分布, 具有单位振幅且曲率相位前半径无穷大的高斯光束表示为 [12]:
$$ {U_0}\left( {r,0} \right) = {\text{exp}}\left( { - \frac{{{x^2} + {y^2}}}{{w_0^2}}} \right) $$ (17) 式中:${w_0}$为光束的初始腰斑半径,当${w_0} = 5 \;{\rm{cm}}$时,在腰斑处的初始辐度分布如图3所示。
-
对Kolmogorov湍流以及谱指数为10/3的Non-Kolmogorov湍流分别进行模拟,并对其进行低频补偿,仿真参数为:大气湍流强度为${{C}}_n^2 = 2.01 \times {10^{ - 17}}{{\rm{m}}^{ - 2/3}}$,相位屏尺度为G=1 m,采样点数为 N=512,仿真结果见图4和图5。
图 4 不同补偿次数下Kolmogorov湍流相位屏
Figure 4. Kolmogorov turbulent phase screens with different low frequency compensation
图 5 不同补偿次数下10/3幂律Non-Kolmogorov湍流相位屏
Figure 5. 10/3 Non-Kolmogorov turbulent phase screens with different low frequency compensation
相位结构函数(Phase Structure Fuinction, PSF)描述了波前相位的二阶统计特征,可以用来检验相位屏的准确性[30]。在数值仿真中,随机相位屏相位结构函数的统计量为[31]:
$$ {D_p}\left( r \right) = < {\left[ {\phi \left( {r + \rho } \right) - \phi \left( r \right)} \right]^2} > $$ (18) 式中: $\;\rho $为两点之间的距离;$ < \cdot > $为系统平均,表示在已生成的相位屏上对所有距离相等的点。
为了更加直观地反映低频补偿前后相位屏的变化,同时验证湍流模拟的准确性,图6给出了两种湍流相位屏未叠加次谐波、叠加不同级次谐波的相位结构函数模拟值与理论值比较。可以看出未叠加次谐波补偿的相位结构函数与理论值差异较大,经过次谐波补偿后得到明显改善。
-
光场的互相关因子(Mutual Coherence Factor, MCF)是一个二阶场矩,常用于验证传输模型的准确性[32-33]。假设接收端的光场为U(r, L),其中L是沿正z轴从发射机的发射孔径到接收机的传播距离,r是接收机平面上到传播轴的矢量。则互相关函数为:
$$ {{{\varGamma }}_2}\left( {{r_1},{r_2},L} \right) = \left\langle {U\left( {{r_1},L} \right){U^*}\left( {{r_2},L} \right)} \right\rangle $$ (19) 式中:$ < > $为统计平均;${r_1}$、${r_2}$为观测面上两个不同的点。
相干因子$ \mu \left( {{r_1},{r_2},z} \right) $表示初始光束的相干损失,是互相干函数${{\varGamma }}\left( {{r_1},{r_2},L} \right)$的模量,$\;\rho {\text{ = }}\left| {{r_1} - {r_2}} \right|$是两观测点之间的距离:
$$ \mu \left( {{r_1},{r_2},z} \right){\text{ = }}\mu \left( {\rho ,z} \right) = \frac{{\left| {{{\varGamma }}\left( {{r_1},{r_2},L} \right)} \right|}}{{{{\left[ {{{\varGamma }}\left( {{r_1},{r_{\text{1}}},L} \right){{\varGamma }}\left( {{r_{\text{2}}},{r_2},L} \right)} \right]}^{1/2}}}} $$ (20) 理论上,光场的互相干因子与波结构函数有关[3]:
$$ \mu \left( {r,z} \right) = {\text{exp}}\left( { - D\left( {r,z} \right)/2} \right) $$ (21) 有关Kolmogorov谱的波结构函数已经得到了广泛的研究[3, 15, 34-35],而关于Non-Kolmogorov谱的波结构函数还没有确切的表达式,因此本小节的验证也是基于传统的Kolmogorov谱。在Kolmogorov谱下平面波的波结构函数为:
$$ D\left( r \right) = 6.68{\left(\frac{r}{{{r_{0sw}}}}\right)^{5/3}} $$ (22) 式中:$ {r_{0,s}} $表示球面波的相干长度;$ {r_{0,s}} = {(0.423{k^2}{C_n}^2} \cdot \dfrac{3}{8}{L_T})^{ - 3/5} $。
均方误差(Mean Squared Error, MSE)是最常用的一种失真度指标,可以用来计算仿真曲线与理论曲线之间差的平方的均值。其数学表达式为[36]:
$$ {{MSE}} = \left( {\frac{1}{n}} \right)\times \sum {({y_i} - {y_t})} $$ (23) 式中:$n$为数据点的数量;${y_i}$为第$i$个数据点的真实值;$ {y_t} $为第$i$个数据点的理论预测值。MSE越小,表示仿真值与理论值之间的差异越小,仿真效果越好。
为了验证文中所提出的3层传输模型的准确性,图7数值仿真了1000次随机实现下,基于Kolmogorov谱的不同湍流层数(2层、3层、6层、11层以及21层)下观测面光束的互相干因子与理论值的比较。在此基础上,表2计算了不同湍流层数下互相干因子的仿真曲线与理论曲线的均方误差。
表 2 光场互相干因子的理论与仿真误差比较
Table 2. Comparison of theoretical and simulation errors of mutual coherence factor
Number of layers 2 3 6 11 21 Mean squared error 0.1471 6.87×10−4 3.43×10−5 3.38×10−5 3.15×10−5 仿真结果表明,2层模型与理论相差较大,对于地星链路这样的湍流环境,至少也需要3层相位屏。理论上,当传输距离一定,使用的湍流层数越多,越贴近理论,即模拟越准确[37],文中仿真也证明了这一点。但是,湍流层数量越多,分层模型越复杂,计算机运行时间非常久。对于复杂的地星传输系统,使用3层模型既保证了传输的准确性,也提高了仿真效率。
-
根据上面两小节的分析,地星链路的湍流是不均匀的,$C_n^2$随高度变化,并且近地面的湍流谱不适用于高空。本小节在3层传输模型下,对比了谱指数随高度变化的 Non-Kolmogorov谱与Kolmogorov谱对上行高斯光束的影响。每一层的大气湍流强度都是按照公式(16)取的路径平均值,假设Non-Kolmogorov谱在每一层的谱指数分别为11/3、3.5和3.3。每一层的相位屏参数如表3所示。
表 3 相位屏参数
Table 3. Parameters of each phase screen
Atmospheric layer Turbulence intensity$ {C}_{n}^{2}/{{\rm{m}}}^{-2/3} $ Non-Kolmogorov spectral index $ \alpha $ Kolmogorov
spectral index $ \alpha $Boundary $ 9.99\times {10}^{-16} $ 11/3 11/3 Troposphere $ 2.01\times {10}^{-17} $ 3.5 11/3 Stratosphere $ 7.51\times {10}^{-18} $ 3.3 11/3 图8和图9分别为采用不同的湍流谱模型(Kolmogorov谱与Non-Kolmogorov谱)模拟的3个大气层的相位屏模拟结果,图9是对图8中每一个相位屏进行低频补偿后的结果。观测图9可以发现,与图8相比,补偿后的相位屏更加光滑。对于边界层的相位屏,即图(a)和图(d)的大气扰动基本一致,均为[−6, 6] rad,这是因为Non-Kolmogorov和Kolmogorov的谱指数在这一层都为11/3。对流层图(b)和图(e)的大气扰动相差也不大,最大的区别在于平流层,在Kolmogorov谱中,大气扰动几乎为0,而Non-Kolmogorov谱生成的相位屏相位扰动比较大。
图 8 不同湍流下次谐波补偿前三个大气层对应的相位屏模拟结果。(a)~(c) Kolmogorov; (d)~(f) Non-Kolmogorov
Figure 8. Simulation results of phase screen simulations corresponding to three atmospheres without harmonic compensation. (a)-(c) Kolmogorov; (d)-(f) Non-Kolmogorov
图 9 不同湍流下次谐波补偿后3个大气层对应的相位屏模拟结果。(a)~(c) Kolmogorov;(d)~(f) Non-Kolmogorov
Figure 9. Simulation results of phase screen simulations corresponding to three atmospheres with harmonic compensation. (a)-(c) Kolmogorov; (d)-(f) Non-Kolmogorov
为了保证湍流模拟的可靠性,图10计算了图8和图9中对应的相位屏补偿前后的相位结构函数与理论值的比较曲线,仿真结果表明,经过次谐波补偿,相位屏均得到了很好的改善,补偿后的相位屏结构函数与理论值更加接近。
图11是高斯光束分别在Kolmogorov谱和Non-Kolmogorov谱下的传输仿真结果。图中,从左到右分别对应边界层2 km处,对流层10 km处以及平流层20 km处的二维光幅度轮廓分布,图(a)~(c)为Kolmogorov湍流下的传输结果,图(d)~(f)为Non-Kolmogorov湍流下的传输结果。可以发现,无论是Kolmogorov谱还是Non-Kolmogorov谱,在传输过程中,光束中心辐度都会随着传输距离的增加逐渐减小,光束的轮廓也越来越大,即光束扩展现象。在相同的仿真条件下,与谱指数固定不变的Kolmogorov谱相比,当谱指数随传输距离变化时,湍流对传输光的影响更明显,主要体现在更大的光幅度衰减以及更明显的光束扩展。不过,在文中的激光传输数值仿真中,由于参数设置的原因,所以整体高斯光束幅值的改变很小。为了进一步分析非kolmogorov湍流以及Kolmogorov湍流在非均匀的地星链路上对相干光束相干性的影响,图12根据公式(20)绘制了Kolmogorov与Non-Kolmogorov谱下输出光束的相干度在横向距离$r$上的变化。结果表明,在靠近光轴中心的位置,Non-Kolmogorov的相干度要比Kolmogorov谱大,而随着$r$值的增加,NK谱下的光束相干度要比K谱小很多,说明谱指数的动态变化对传输光束的相干性产生了影响。
Numerical simulation of three-layer transmission model for ground-to-satellite turbulent path
-
摘要: 星地激光通信具有高带宽通信的潜力,但大气湍流会显著影响卫星对地和地对卫星这类通信系统的能力。星地通信系统的室外实验费用昂贵且难以重现,现有的数值仿真大都基于水平均匀路径,不适用于星地链路长距离且不均匀的湍流路径,为了评估湍流对通信系统的影响,因此开发适用于星地链路的数值仿真是非常重要的。在数值仿真中,模拟湍流的相位屏数量太多会增加系统的复杂度,通过计算传统的Kolmogorov谱下相位屏为2层、3层、6层、11层以及21层时观测面光场的互相干因子并将其与理论值比较,发现三层模型能够在保证传输准确性的同时降低系统的复杂度。在此基础上,提出了地星大气激光传输三层模型,分析了不同湍流谱模型下(Kolmogorov谱与Non-Kolmogorov谱)的高斯光束辐度轮廓变化以及接收光束相干性的变化。仿真结果表明,在相同的大气条件下,与Kolmogorov谱相比,谱指数随高度变化的Non-Kolmogorov谱对传输光束的幅度以及相干性的影响更大。
-
关键词:
- 数值模拟 /
- 大气湍流 /
- 光传播 /
- Kolmogorov谱 /
- Non-Kolmogorov谱 /
- 分步传播算法
Abstract:Objective Satellite-ground laser communication has the potential for high-bandwidth communication, but atmospheric turbulence can significantly affect the capabilities of communication systems. Outdoor experiments of satellite-ground communication systems are expensive and difficult to reproduce. Most of the existing numerical simulations are based on horizontal uniform paths, and are not suitable for turbulence paths of non-uniform satellite-ground links. In order to evaluate the impact of turbulence on communication systems, it is very important to develop numerical simulations suitable for satellite-ground links. In numerical simulation, simulating with too many phase screens will increase the complexity of the system. Therefore, it is necessary to develop a simple and reliable numerical simulation model for ground-satellite or satellite-ground turbulence path. Methods The phase screens under the Kolmogorov spectrum and the non-Kolmogorov spectrum are simulated by Fourier inversion method, and the subharmonic compensation method is used to compensate the phase screens. By calculating the constraints of the numerical simulation, a three-layer transmission simulation model is proposed. In order to verify the model, the laser transmission under the traditional Kolmogorov turbulence spectrum model with 2, 3, 6, 11 and 21 layers is simulated on the basis of the split-step method. Results and Discussions Simulations of phase screens under the Kolmogorov spectrum and the non-Kolmogorov spectrum show that the spectral index has a great influence on the phase screen simulation (Fig.4-5). Different spectral indices have different requirements for subharmonic compensation (Fig.6). The mutual coherence factor of the optical field on the observation surface is calculated, and compared with the theoretical value (Fig.7). The mean square error was used to measure the distortion between simulation and theory for different layers of phase screens (Tab.2). It is found that the three-layer model can ensure the accuracy of transmission and reduce the complexity of the system. The changes in the amplitude profile of Gaussian beam (Fig.11) and the changes in the coherence of the received beam under different turbulence spectrum models (Kolmogorov spectrum and non-Kolmogorov spectrum) (Fig.12) are analyzed, the simulation results show that under the same atmospheric conditions, the Gaussian beam in the non-Kolmogorov spectrum will produce obvious amplitude attenuation and beam spread, and the coherence of the optical field at the observation plane will decrease faster. Conclusions In this study, a simple and reliable numerical simulation model for ground-satellite turbulence path system is designed. In this model, the turbulent path is divided into three layers, and each layer is described by different spectral indices according to the characteristics of atmospheric turbulence. The mutual coherence factors were calculated and compared with the theoretical value under the traditional Kolmogorov spectrum with 2, 3, 6, 11 and 21 phase screens. The feasibility of the three-layer transmission model was verified, the proposed model can ensure the simulation accuracy and the efficiency of computer simulation. Based on the proposed three-layer transmission model, the transmission process of the traditional Kolmogorov spectrum and the non-Kolmogorov spectrum with the spectral index varying with the height are simulated, and the changes of the optical amplitude and the observation plane coherence are compared between the two. The results show that under the same atmospheric conditions, the Kolmogorov spectrum in the near-surface region is basically consistent with the results of the non-Kolmogorov spectrum, while with the increase of transmission distance, the Gaussian beam uplink will generate obvious amplitude attenuation and beam spread in the non-Kolmogorov spectrum. The coherence of the light field on the observation surface decrease faster at the same time. The proposed simulation model can be used to evaluate the optical transmission characteristics of the ground-satellite uplink under the influence of complex atmospheric turbulence with fewer phase screens and simulation time, which provides convenience for predicting the expected wavefront behavior at the aperture of the satellite optical receiver. -
表 1 仿真参数
Table 1. Simulation parameters
Parameters Notion Value Wavelength of Gauss beam $ \lambda $ 0.5 μm Gauss beam radius ${{w} }_{0}$ 0.05 m Zenith angle $ \theta $ 0° Outer scale ${{L} }_{0}$ 50 m Inner scale ${{l} }_{0}$ 0.001 m Turbulence model ${\rm{HV}}5/7$ ${\rm{HV}}5/7$ Transmitter grid spacing $ {\delta }_{t} $ 0.0035 m Receiver grid spacing $ {\delta }_{n} $ 0.005 m Sample points $ N $ 512 Scaling factor $ {\alpha }_{j} $ 0 Diameter of the observation aperture $ {D}_{2} $ 0.5 m 表 2 光场互相干因子的理论与仿真误差比较
Table 2. Comparison of theoretical and simulation errors of mutual coherence factor
Number of layers 2 3 6 11 21 Mean squared error 0.1471 6.87×10−4 3.43×10−5 3.38×10−5 3.15×10−5 表 3 相位屏参数
Table 3. Parameters of each phase screen
Atmospheric layer Turbulence intensity$ {C}_{n}^{2}/{{\rm{m}}}^{-2/3} $ Non-Kolmogorov spectral index $ \alpha $ Kolmogorov
spectral index $ \alpha $Boundary $ 9.99\times {10}^{-16} $ 11/3 11/3 Troposphere $ 2.01\times {10}^{-17} $ 3.5 11/3 Stratosphere $ 7.51\times {10}^{-18} $ 3.3 11/3 -
[1] Wang Y, Zhu L, Feng W K. Performance study of wavelength diversity serial relay OFDM FSO system over exponentiated Weibull channels [J]. Optics Communications, 2021, 478: 126470. doi: 10.1016/j.optcom.2020.126470 [2] Long L A , Hs A , Rz A , et al. Increasing system tolerance to turbulence in a 100-Gbit/s QPSK free-space optical link using both mode and space diversity [J]. Optics Communications, 2021, 480: 126488. [3] Andrews L C, Phillips R L. Laser Beam Propagation Through Random Media[M]. US: SPIE Press Book, 2010, PM152: 212. [4] 姜义君. 星地激光通信链路中大气湍流影响的理论和实验研究 [D]. 哈尔滨: 哈尔滨工业大学, 2010. Jiang Yijun. Theoretical and experimental study on atmospheric turbulence in satellite-ground laser communication link[D]. Harbin: Harbin Institute of Technology, 2010. (in Chinese) [5] 都文和. 星地激光通信中非柯尔莫哥洛夫湍流影响研究 [D]. 哈尔滨: 哈尔滨工业大学, 2010. Du Wenhe. Study on the influence of Non-Kolmogorov turbulence on satellite to ground laser communication [D]. Harbin: Harbin Institute of Technology, 2010. (in Chinese) [6] Liu X, Zhang Q, Xin X, et al. Numerical simulation of ground-to-satellite laser transmission based on unequal spacing phase screen[C]//2019 18th International Conference on Optical Communications and Networks (ICOCN), IEEE, 2019: 1-3. [7] 钱仙妹, 朱文越, 饶瑞中. 地空激光大气斜程传输湍流效应的数值模拟分析[J]. 红外与激光工程, 2008(05): 787-792. doi: 10.3969/j.issn.1007-2276.2008.05.009 Qian Xianmei, Zhu Wenyue, Rao Ruizhong. Numerical simulation and analysis of turbulence effect of ground-to-air laser propagation in atmospheric oblique path [J]. Infrared and Laser Engineering, 2008, 37(5): 787-792. (in Chinese) doi: 10.3969/j.issn.1007-2276.2008.05.009 [8] Zilberman A, Golbraikh E, Kopeika N S. Lidar studies of aerosols and non-kolmogorov turbulence in the mediterranean troposphere[C]//Electro-Optical and Infrared Systems: Technology and Applications II, SPIE, 2005, 5987: 15-26. [9] Shchepakina E, Korotkova O. Second-order statistics of stochastic electromagnetic beams propagating through non-Kolmogorov turbulence [J]. Optics Express, 2010, 18(10): 10650-10658. doi: 10.1364/OE.18.010650 [10] Toselli I, Andrews L C, Phillips R L, et al. Free space optical system performance for laser beam propagation through non Kolmogorov turbulence for uplink and downlink paths[C]//Atmospheric Optics: Models, Measurements, and Target-in-the-Loop Propagation, SPIE, 2007, 6708: 17-28. [11] Chen M, Gao T, Hu S, et al. Simulating non-Kolmogorov turbulence phase screens based on equivalent structure constant and its influence on simulations of beam propagation [J]. Results in Physics, 2017, 7: 3596-3602. doi: 10.1016/j.rinp.2017.09.034 [12] Du W, Yao Z, Liu D, et al. Influence of non-kolmogorov turbulence on intensity fluctuations in laser satellite communication [J]. Journal of Russian Laser Research, 2012, 33: 90-97. doi: 10.1007/s10946-012-9263-1 [13] Huang Y, Gao Z, Wang F, et al. Effective radius of curvature of partially coherent Hermite–Gaussian beams propagating through non-Kolmogorov turbulence [J]. Journal of Modern Optics, 2012, 59(19): 1674-1679. doi: 10.1080/09500340.2012.735712 [14] 熊准. 星地光通信中大气湍流影响的模拟与实验研究 [D]. 武汉: 武汉大学, 2014. Xiong Zhun. Simulation and experimental study on atmospheric turbulence in satellite-ground optical communication [D]. Wuhan: Wuhan University, 2014. (in Chinese) [15] Lane R G, Glindemann A, Dainty J C. Simulation of a kolmogorov phase screen [J]. Waves in Random Media, 1992, 2(3): 209. doi: 10.1088/0959-7174/2/3/003 [16] Lyras N K, Kourogiorgas C I, Kapsis T T, et al. Ground-to-satellite optical link turbulence effects: propagation modelling & transmit diversity performance[C]//2019 13th European Conference on Antennas and Propagation (EuCAP), IEEE, 2019: 1-5. [17] 程知, 侯再红, 靖旭, 等. Hufnagel-Valley湍流强度廓线的高精度实时反演方法[J]. 红外与激光工程, 2013, 42(06): 1562-1567. doi: 10.3969/j.issn.1007-2276.2013.06.033 Cheng Zhi, Hou Zaihong, Jing Xu, et al. High-precision real-time inversion method of turbulence intensity profile in Hufnagel-valley [J]. Infrared and Laser Engineering, 2013, 42(6): 1562-1567. (in Chinese) doi: 10.3969/j.issn.1007-2276.2013.06.033 [18] Pudasaini P R, Vera M, Pokheral M. Numerical simulation of optical propagation through atmospheric turbulence[C]//APS Texas Sections Fall Meeting Abstracts, 2010. [19] Zilberman A, Golbraikh E, Kopeika N S. Propagation of electromagnetic waves in Kolmogorov and non-Kolmogorov atmospheric turbulence: three-layer altitude model [J]. Applied Optics, 2008, 47(34): 6385-6391. doi: 10.1364/AO.47.006385 [20] 马春波, 黄可, 敖珺. Kolmogorov与Non-Kolmogorov湍流相位屏的数值仿真方法研究[J]. 计算机与数字工程, 2019, 47(12): 2982-2988. Ma Chunbo, Huang Ke, Ao Jun. Numerical simulation of kolmogorov and Non-kolmogorov turbulence phase screens [J]. Computer and Digital Engineering, 2019, 47(12): 2982-2988. (in Chinese) [21] Lilly D, Lester P F. Waves and turbulence in the stratosphere [J]. Journal of Atmospheric Sciences, 1974, 31(3): 800-812. doi: 10.1175/1520-0469(1974)031<0800:WATITS>2.0.CO;2 [22] Yi X, Liu Z J, Yue P. Uplink laser satellite-communication system performance for a Gaussian beam propagating through three-layer altitude spectrum of weak-turbulence [J]. Optik, 2013, 124(17): 2916-2925. doi: 10.1016/j.ijleo.2012.08.082 [23] Hu X, Wu X, Yang Q, et al. Estimation and characterization of atmospheric turbulence in the free atmosphere above the Tibetan Plateau using the thorpe method [J]. Appl Opt, 2023, 62(4): 1115-1137. doi: 10.1364/AO.483677 [24] Chen X, Anel J A, Su Z, et al. The deep atmospheric boundary layer and its significance to the stratosphere and troposphere exchange over the Tibetan Plateau [J]. PloS One, 2013, 8(2): e56909. doi: 10.1371/journal.pone.0056909 [25] Na Z, Xin S, Yan Z, et al. Simulated experiment of the light intensity influenced by non-Kolmogorov turbulence [J]. Acta Photonica Sinica, 2018, 47(6): 0601002. (in Chinese) [26] Du W, Yuan Q, Cheng X, et al. Scintillation index of a spherical wave propagating through Kolmogorov and non-Kolmogorov turbulence along laser-satellite communication uplink at large zenith angles [J]. Journal of Russian Laser Research, 2021, 42: 198-209. doi: 10.1007/s10946-021-09951-8 [27] Xu T, Li Y, Du W, et al. Simulating atmospheric turbulence using a spatial light modulator based on fourier transform[C]//2014 Conference on Lasers and Electro-Optics (CLEO)-Laser Science to Photonic Applications, IEEE, 2014: 1-2. [28] 李玉杰, 朱文越, 饶瑞中. 非Kolmogorov大气湍流随机相位屏模拟[J]. 红外与激光工程, 2016, 45(12): 169-176. Li Yujie, Zhu Wenyue, Rao Ruizhong. Stochastic phase screen simulation of non-Kolmogorov atmospheric turbulence [J]. Infrared and Laser Engineering, 2016, 45(12): 169-176. (in Chinese) [29] 钱仙妹, 朱文越, 饶瑞中. 非均匀湍流路径光传播数值模拟中相位屏间C_n~2的选取[J]. 光学学报, 2008(10): 1856-1860. doi: 10.3321/j.issn:0253-2239.2008.10.004 Qian Xianmei, Zhu Wenyue, Rao Ruizhong. Selection of C_n~2 between phase screens in numerical simulation of light propagation in non-uniform turbulent path [J]. Acta Optica Sinica, 2008(10): 1856-1860. (in Chinese) doi: 10.3321/j.issn:0253-2239.2008.10.004 [30] Coles W A, Filice J P, Frehlich R G, et al. Simulation of wave propagation in three-dimensional random media [J]. Applied Optics, 1995, 34(12): 2089-2101. doi: 10.1364/AO.34.002089 [31] Liu T, Zhang J, Lei Y, et al. Optimal subharmonics selection for atmosphere turbulence phase screen simulation using the subharmonic method [J]. Journal of Modern Optics, 2019, 66(9): 986-991. doi: 10.1080/09500340.2019.1596323 [32] Jera E, Mohamed A. Gaussian beam propagation through turbulent atmosphere using second-order split-step algorithm[C]//2019 IEEE National Aerospace and Electronics Conference (NAECON), IEEE, 2019: 675-677. [33] Belmonte A. Feasibility study for the simulation of beam propagation: consideration of coherent lidar performance [J]. Applied Optics, 2000, 39(30): 5426-5445. doi: 10.1364/AO.39.005426 [34] Dun L, Yu N, Wuming W, et al. Numerical simulation and validation method of atmospheric turbulence of phase screen in rotation [J]. Infrared and Laser Engineering, 2017, 46(12): 1211003. (in Chinese) [35] Frehlich R. Simulation of laser propagation in a turbulent atmosphere [J]. Applied Optics, 2000, 39(3): 393-397. doi: 10.1364/AO.39.000393 [36] Hastie T, Tibshirani R, Friedman J H, et al. The Elements of Statistical Learning: Data Mining, Iinference, and Prediction[M]. New York: Springer, 2009. [37] Schmidt J D. Numerical Simulation of Optical Wave Propagation: With Examples in MATLAB[M]. US: SPIE Press Book, 2010, PM199: 808.