一种基于谱反演法的离散雨滴介质的相位屏建模方法

摘要

本发明涉及一种基于谱反演法的离散雨滴介质的相位屏建模方法，属于离散雨滴介质相位屏建模领域。为了解决现有的离散雨滴介质的相位屏模型对光场信息的描述不完全，不符合虚拟试验的要求，未考虑雨滴对于光场的相位扰动的缺点，而提出一种基于谱反演法的离散雨滴介质的相位屏建模方法，包括：获取雨滴粒子的散射势能；得到经傅立叶变换后的散射势能；通过经傅立叶变换后的散射势能得到散射势能的近似值；通过公式表示雨滴粒子介电常数的二阶相关矩三维功率谱函数；通过公式表示雨滴的介电常数功率谱；通过公式表示雨滴介质折射率二阶相关矩功率谱；根据谱反演法获得传输路径上混合介质的相位屏模型；本发明适用于激光在雨中传输的仿真模型。

说明书

技术领域

本发明涉及一种基于谱反演法的离散雨滴介质的相位屏建模方法，属于离散雨滴介质相位屏建模领域。

背景技术

由于大气介质介电常数的变化，激光束在大气中传输的光波波前会受到随机相位扰动的影响。这些相位扰动主要由湍流和浑浊介质所造成，同时单个粒子的散射也会伴随一定程度的衰减。因此在虚拟试验中，针对由低散射介质所构建的大气空间，激光传输模型需要给出目标平面包括振幅和相位在内的全部光场信息。在低散射特性的介质中，降雨是一种最主要的造成激光功率衰减并引起波前相位扰动的大气介质。

目前，对于降雨影响激光束传播的研究工作大多集中在讨论降雨的衰减效应上。前期研究仅考虑遮挡效应得到激光在雨中传输的衰减模型，之后有文献结合Mie理论和M-P雨滴谱分布，分析了雨滴粒子尺度和激光波长与消光系数之间的关系，推导得到了不同波长激光的衰减系数与降雨强度之间的内在关系。

然而上述现有技术的方法均着眼于考虑降雨对激光传输的衰减作用，对光场信息的描述不完全，不符合虚拟试验的要求。而且即使一些研究中得到了目标平面的光场分布，但并未考虑雨滴对于光场的相位扰动。

发明内容

本发明的目的是为了解决现有的离散雨滴介质的相位屏模型对光场信息的描述不完全，不符合虚拟试验的要求，未考虑雨滴对于光场的相位扰动的缺点，而提出一种基于谱反演法的离散雨滴介质的相位屏建模方法包括如下步骤：

步骤一：获取在位置r处的雨滴粒子相对于在位置R_j处雨滴粒子的散射势能V(r)；

步骤二：使用傅立叶变换表示雨滴粒子的散射势能V(r)，得到经傅立叶变换后的散射势能V₀(r)；

步骤三：通过经傅立叶变换后的散射势能V₀(r)以及单个散射粒子的散射幅度函数f₀(K)，并结合Born近似方法，得到：

其中f₀(K)表示平面波经过一个半径为a的粒子后其散射振幅的空间频谱，为散射势能的近似值；

步骤四：通过公式表示雨滴粒子介电常数的二阶相关矩三维功率谱函数

其中，ε_p(r)为雨滴粒子的介电常数，为雨滴粒子介电常数的伴随矩阵，r₁和r₂为任意的两个径向距离；

步骤五：通过公式表示雨滴的介电常数功率谱：

其中，p(a,J)是雨滴尺度谱分布，与降雨强度参数J有关，降雨强度的单位为mm/h，表示每小时在单位面积上的降雨厚度，为f₀(K)的伴随矩阵；

步骤六：根据介电常数与折射率的关系结合降雨的尺度谱，通过公式表示雨滴介质折射率二阶相关矩功率谱：

其中a_m是雨滴的平均粒子半径；

步骤七：使用von Kármán模型表示湍流的折射率二阶相关矩功率谱：

其中，是湍流折射率的结构常数，K₀＝2π/L₀，K_m＝5.92/l₀，L₀和l₀分别为外尺度和内尺度；

步骤八：根据步骤六以及步骤七获得折射率起伏谱：

步骤九：根据谱反演法以及步骤八获得传输路径上混合介质的相位屏模型。

本发明的有益效果为：本发明考虑到雨滴对于光场的相位扰动，并针对以二阶矩为例的归一化闪烁系数，仿真结果与一阶微扰近似的理论结果基本一致，实时性可以满足虚拟试验的需求。

附图说明

图1为本发明的基于谱反演法的离散雨滴介质的相位屏建模方法的流程图；

图2为低频叠加相位屏生成示意图；

图3为对于平面波闪烁系数与已有的闪烁系数数据的对比；实线表示一阶微扰近似的闪烁系数仿真结果；离散点表示闪烁系数的理论值结果；

图4为在传输路径为4km时，不同降雨条件下，平面波的光场相干系数变化曲线；其中曲线由上至下分别代表湍流条件下、降雨量为5mm/h条件下、降雨量为10mm/h条件下以及降雨量为15mm/h条件下的相干系数变化曲线；

图5为对于高斯波闪烁系数与已有的闪烁系数数据的对比；实线表示一阶微扰近似的闪烁系数仿真结果；离散点表示闪烁系数的理论值结果。

具体实施方式

具体实施方式一：本实施方式的基于谱反演法的离散雨滴介质的相位屏建模方法包括如下步骤：

步骤一：获取在位置r处的雨滴粒子相对于在位置R_j处雨滴粒子的散射势能V(r)；

步骤二：使用傅立叶变换表示雨滴粒子的散射势能V(r)，得到经傅立叶变换后的散射势能V₀(r)；

步骤三：通过经傅立叶变换后的散射势能V₀(r)以及单个散射粒子的散射幅度函数f₀(K)，并结合Born近似方法，得到：

其中f₀(K)表示平面波经过一个半径为a的粒子后其散射振幅的空间频谱，为散射势能的近似值；

步骤四：通过公式表示雨滴粒子介电常数的二阶相关矩三维功率谱函数

其中，ε_p(r)为雨滴粒子的介电常数，为雨滴粒子介电常数的伴随矩阵，r₁和r₂为任意的两个径向距离；

步骤五：通过公式表示雨滴的介电常数功率谱：

其中，p(a,J)是雨滴尺度谱分布，与降雨强度参数J有关，降雨强度的单位为mm/h，表示每小时在单位面积上的降雨厚度，为f₀(K)的伴随矩阵；

步骤六：根据介电常数与折射率的关系结合降雨的尺度谱，通过公式表示雨滴介质折射率二阶相关矩功率谱：

其中a_m是雨滴的平均粒子半径；

步骤七：使用von Kármán模型表示湍流的折射率二阶相关矩功率谱：

其中，是湍流折射率的结构常数，K₀＝2π/L₀，K_m＝5.92/l₀，L₀和l₀分别为外尺度和内尺度；

步骤八：根据步骤六以及步骤七获得折射率起伏谱：

步骤九：根据谱反演法以及步骤八获得传输路径上混合介质的相位屏模型。

具体实施方式二：本实施方式与具体实施方式一不同的是：f₀(K)具体为：

其中k为波数，J₁表示一阶贝塞尔函数。

其它步骤及参数与具体实施方式一相同。

具体实施方式三：本实施方式与具体实施方式一或二不同的是：步骤二中，所述经傅立叶变换后的散射散射势能V₀(r)具体为：

其它步骤及参数与具体实施方式一或二相同。

具体实施方式四：本实施方式与具体实施方式一至三之一不同的是：步骤一具体为：

建立雨滴粒子的介电常数与其散射势能V(r)的关系：其中，k为波数；

将V(r)表示为：其中V₀(r-R_j)是位于r处的雨滴粒子相对于在R_j处雨滴粒子的散射势能。

其它步骤及参数与具体实施方式一至三之一相同。

具体实施方式五：本实施方式与具体实施方式一至四之一不同的是：步骤五中，所述雨滴尺度谱分布p(a,J)使用泊松分布描述，具体为：

p(a,J)＝m(J)exp[-n(J)a]；

其中，m(J)＝4.382×10⁶J^0.112m^-4；n(J)＝5932J^-0.182m^-1；

其它步骤及参数与具体实施方式一至四之一相同。

具体实施方式六：本实施方式与具体实施方式一至五之一不同的是：步骤六中，雨滴的平均粒子半径a_m具体为：

其它步骤及参数与具体实施方式一至五之一相同。

具体实施方式七：本实施方式与具体实施方式一至六之一不同的是：步骤七中，具体为

其中，Δz为两个相位屏间的距离，Δz_i为第i个相位屏与第i+1个相位屏的间距，公式中的等号表示赋值运算。

其它步骤及参数与具体实施方式一至六之一相同。

具体实施方式八：本实施方式与具体实施方式一至七之一不同的是：步骤七中，具体为

其它步骤及参数与具体实施方式一至七之一相同。

采用以下实施例验证本发明的有益效果：

<实施例1>

本发明可以用于建立激光在雨中传输的仿真模型。

一)基于谱反演法的相位屏模型介绍：

本发明是基于谱反演法理论的，因此为了便于说明本发明的模型建立过程，下面先介绍基于谱反演法的相位屏模型理论。

基于谱反演法的相位屏建模方法由于其理论背景清楚，使用方便，因此应用较为广泛。其基本原理对于构造雨滴介质的相位屏模型较为重要。

在研究随机起伏介质中的光传播规律时，在传播路径上设置若干个相位屏，光束传播到一个相位屏时即把对应的相位扰动叠加到光场中，然后真空传播至下一个相位屏处，如此重复直至终点为止，这便是谱反演法。

对于随机函数f，其相关函数和三维空间谱密度的关系为Fourier变换：

其中K表示三维的空间频率向量。B_f为随机函数的相关函数，r₁和r₂为径向距离，当r₁和r₂位于垂直于发射轴的同一平面内时，一般采用ρ₁和ρ₂表示，则在一平面内，二维空间谱密度可表示为：

K_⊥表示垂直平面上的空间频率。

一个相位屏内相位相干函数与该相位屏所代表的平板内折射率相关函数之间的关系可表示为：

B_s(ρ₁-ρ₂)＝k²∫∫<n₁(ρ₁,z)n₁(ρ₂,z′)>dzdz′>

下标s表示相位的相关函数，由相位与介质折射率之间的关系及式(2)所表示的相关函数与二维空间谱密度的关系可得：

k为光波的波数，Δz为传输路径的长度。同一平面内，相位相关函数和二维空间谱密度的关系可表示为：

对于平面波，相位S的二维频谱F_s和折射率的空间谱密度Φ_n的关系为：

F_s(K_x,K_y,z)＝2π²k²ΔzΦ_n(K_x,K_y,K_z＝0,z)>

K_x,K_y,K_z分别表示在x轴、y轴、z轴上的空间频率。

对于高斯波：

根据相位起伏S的二维频谱可在相位空间构造一个二维复变函数随机场。首先得到复随机变量Φ_nm＝a_nm+ib_nm，其中a_nm，b_nm是互不相关且统计平均为0的高斯型随机变量，方差为：

<a²(n,m)>＝<b²(n,m)>＝Δk_xΔk_yΦ_n(nΔk_x,mΔk_y,z)>

对此复变函数随机场进行傅里叶变换，可以得到平面内的相位分布函数：

其中Δx×Δy表示在坐标空间上单元格的大小，N_x，N_y分别表示坐标轴上的分割数，L_x×L_y＝N_xΔx×N_yΔy表示坐标空间相位屏大小。式(9)得到的相位屏模型其二阶矩与式(5)所表示的相位起伏S的二阶矩是相同的。此方法称为谱反演法。

由谱反演法所生成的相位屏模型中，主要存在的问题是低频成分不足，并且由于傅里叶变换的周期性，所产生的随机场会给光学波前带来一部分虚假的信息。因此，有一些学者提出了基于谐波分量的低频补偿相位屏生成方法，该方法的基本思路是：

(1)对平面内包含原点的单元格进行两次细分。重复以上过程，直到达到模拟的精度要求为止。

(2)对平面内最低一级的单元格进行傅里叶变换，而其他高等级的单元格则根据它们对相位的贡献乘以相应的比例系数。比例系数与单元格面积成比例，具体的形式如式(10)所示：

其中m、n分别代表x、y坐标轴方向的级数。j、l是分别与m、n对应的比例系数p代表第几块模板。

图2低频叠加相位屏生成示意图

图2为基于谐波分量的低频补偿相位屏生成方法的示意图。最后将初级相位屏与次谐波相位屏相互叠加就可以得到经低频补偿后的相位屏模型：

θ(jΔx,lΔy)＝θ(jΔx,lΔy)+θ_sh(jΔx,lΔy)>

二)混合大气介质相位屏生成

这部分是本发明的主要部分，主要以公式推导过程来体现本发明的思路已经创新之处。

由谱反演法的基本原理，为得到雨滴粒子的相位屏模型，首先应根据雨滴粒子的尺度分布得到折射率功率谱，在推导过程中假设每一个雨滴粒子的散射特性是相同的，并且雨滴在光波传输路径上满足单散射特性。

雨滴粒子的介电常数与其散射势能V(r)的关系为：

其中，k为波数。

在位置r处雨滴粒子的散射势能可表示为：

V₀(r-R_j)是位于r处的雨滴粒子相对于在R_j处雨滴粒子的散射势能。

使用傅里叶变换表示雨滴粒子的散射势能：

K表示空间频率。通过雨滴粒子散射势能的傅里叶变换与单个散射粒子的散射幅度函数f₀(K)，并结合Born近似，可得：

f₀(K)表示平面波经过一个半径为a的粒子后其散射振幅的空间频谱。由于在一般的雨滴谱分布中，绝大部分的雨滴尺度相对于激光波长要大很多，因此f₀(K)可以用它的衍射部分来进行表示：

其中，k表示波数，J₁表示一阶贝塞尔函数。

雨滴粒子介电常数的二阶相关矩三维功率谱函数可表示为：

由以上的推导过程雨滴的介电常数功率谱可表示为：

其中p(a,J)是雨滴尺度谱分布，其与降雨强度参数J有关，降雨强度的单位为mm/h，表示每小时在单位面积上的降雨厚度，本文中降雨的尺度谱分布采用泊松分布来进行描述：

p(a,J)＝m(J)exp[-n(J)a] (19)

其中m(J)＝4.382×10⁶J^0.112m^-4；n(J)＝5932J^-0.182m^-1。根据介电常数与折射率的关系结合降雨的尺度谱分布可得雨滴介质折射率二阶相关矩功率谱可表示为：

根据如上的雨滴粒子泊松尺度分布函数可得的各个参数。其中，a_m是雨滴的平均粒子半径，由雨滴的尺度谱分布可表示为：

除降雨外，还需要考虑湍流对激光传输的影响，湍流的折射率二阶相关矩功率谱采用von Kármán模型表示：

是湍流折射率的结构常数，K₀＝2π/L₀，K_m＝5.92/l₀，L₀和l₀分别为外尺度和内尺度。由式(20)，(22)可得混合介质的折射率起伏谱为：

根据谱反演法和式(23)可得到在传输路径上混合介质的相位屏模型。

三)大气介质相位屏模型参数选择

这部分是对混合大气介质相位屏生成中的进一步解释，主要研究如下参数：湍流谱内尺度l₀，外尺度L₀，折射率结构常数屏间距离。

传输路径上混合介质的相位屏模型与湍流谱内尺度l₀，外尺度L₀，折射率结构常数屏间距离和降雨强度等参数有关。

其中相位屏间的距离应大于屏间大气介质所有非均匀元的空间尺度，由于一般湍流的非均匀元的尺度要远远大于雨滴，则：

Δz>L₀>

L₀为屏间介质的外尺度。

同时，相位屏之间的传输应满足几何光学近似的条件，要小于湍流的Fresnel尺度，则有：

激光每传输到一个相屏时即把相位屏间湍流大气引起的相位扰动叠加到光波波前上,然后在真空中传输至下一个相位屏处,再叠加下一个相位屏对应的相位扰动,再在真空中传输直至终点为止。在数值求解两相位屏间湍流大气引起的相位

扰动时一般均假设相位屏间湍流强度均匀,即折射率结构常数相等。对于均匀湍流路径的光传输该假设不会给模拟结果带来误差。而对于非均匀湍流路径的光传输，该假设将会造成一定的模拟误差。减小这种误差的途径之一是在传输路径上设置尽可能多的相位屏。然而一方面相位屏间距的限制条件使得相位屏数目有限以致不可能对所有的湍流起伏区进行采样；另一方面过多的相位屏增加了计算负担，降低了计算效率。因此在上述假设和现有的计算条件下，为了减小模拟误差应着眼如何合理选择相位屏间的问题上。在以往的数值模拟过程中，人们常把相位屏间的近似为前一个相位屏处或后一个相位屏处的折射率结构常数，这都将带来较大的模拟误差，为了减小误差，常用的选择方法如下：

取为两相位屏处折射率结构常数的几何平均值：

取为两相位屏间湍流大气折射率结构常数的路径平均值：

四)实验结果与分析

光场分布的一阶矩为平均振幅，表明了光场经浑浊介质传播后未经破坏的部分；单点的二阶矩为场的平均强度；两点的二阶矩可用来表示光场的相干特性。此处主要根据目标平面的二阶矩所表示的光场相干特性和归一化闪烁指数对仿真结果进行分析和验证。

模拟激光在湍流中传输的过程，仿真实验采用如下的参数：发射光束的波长为λ＝1.06μm。网格宽度选择为湍流内尺度的1/2，湍流谱内尺度为1.2mm，外尺度为80m，折射率结构常数屏间距离为100m，降雨强度为30mm/h。

由以上参数，可分别得到湍流谱相位屏和叠加雨滴粒子后的混合介质相位屏。

对于发射光波的波形，分别选择平面波和基膜高斯波束的形式，仿真结果对比已有的光强闪烁系数和相位相干系数的数据来验证模型的正确性。

(1)平面波

对于平面波，目标平面相位相关函数谱与折射率功率谱之间的关系为：

Φ_φ(K_⊥,0,J)＝2π²k²ΔzΦ_n(K_⊥,0,J)>

其中，Δz为两个相位屏之间的距离。因此由式(23)所给出的混合介质的折射率功率谱可得雨滴与湍流介质混合的相位屏模型，将其应用于数据模拟方法中，并结合真空传输模型，可到目标的光场分布，进而得到目标平面的闪烁系数。

光强闪烁系数的定义如式(29)所示：

常通过一阶微扰近似得到了平面波在浑浊介质中传输时闪烁系数类似于Rytov近似的推导结果，其表明：在光学厚度比较低的时，目标平面的光强闪烁系数基本等于光学厚度τ：

而随着光学厚度的增加，闪烁系数逐渐趋近于饱和，闪烁系数基本等于1：

闪烁系数由多次试验后所得到的光场进行统计平均得到，其表达式为：

M表示网格数，N表示样本总数。

对于平面波，本文方法所得到的闪烁系数与式(30)和式(31)进行比较，结果如图3所示，闪烁系数的仿真结果与理论值吻合的程度存在一定的偏差，当光学厚度较小时，吻合的程度较好，而随着光学厚度的增加，闪烁指数呈上升的趋势，这是由于已有的闪烁系数的相关理论结果是在一阶微扰近似的条件下得到，当进入强起伏区域时，所得到的理论结果将不再适用。

由于雨滴粒子的散射效应，降雨对接收平面的相干因子有很大的影响。随着降雨强度的增加，相干系数逐渐降低，目标平面的光场会逐渐趋于相位无关场。图4表示在传输路径为4km时，不同降雨条件下，平面波的光场相干系数变化曲线。

(2)基膜高斯波束

基膜高斯光束的定义可表示为：

E(r)＝U(r)e^ikz>

对于高斯型的波束，目标平面相位相关函数谱与折射率功率谱之间的关系可表示为：

其中，Δz为两个相位屏之间的距离。

针对于基膜高斯光束，根据一阶微扰近似所得到的闪烁系数的推导结果可表示为：

其中，k＝2π/λ为空间波数。a_m＝6.19×10^-4J^0.182为雨滴的平均半径。τ＝2.638×10^-4zJ^0.658为雨滴介质的在传输距离为z时的光学厚度。为菲涅耳数。a₀＝0.3cm为光束的初始半径。本文方法所得到的闪烁指数与式(30)和(31)的比较结果如图5所示。

光学厚度逐渐增加时，也会出现类似于输入波形已有为平面波时的偏差，原因与平面波的原因类似。

由于本模型的主要计算部分均采用FFT进行计算，因此计算效率相对较高，在硬件环境为(WindowsXP，INTEL Core(TM)2，2.66GHZ)的条件下，在每个相位屏之间的仿真时间约为0.008s。在传输距离为1千米左右的距离上，相位屏总数为10～7个，总仿真运行时间T_total为0.08～0.056s，基本可满足虚拟试验对实时仿真的需求。

结论：本发明根据数值模拟方法建立了一种描述激光在雨中传输的仿真模型。由雨滴的尺度分布谱，基于谱反演法建立了雨滴的相位屏模型。仿真结果表明，针对以二阶矩为例的归一化闪烁系数，仿真结果与一阶微扰近似的理论结果基本一致，随着降雨强度的增加，目标处光场的相干性逐渐下降，光场分布逐渐发散。实时性基本可以满足虚拟试验的需求。

本发明还可有其它多种实施例，在不背离本发明精神及其实质的情况下，本领域技术人员当可根据本发明作出各种相应的改变和变形，但这些相应的改变和变形都应属于本发明所附的权利要求的保护范围。