一种海洋环境电波传导特性获取方法与流程

本发明涉及电波传播领域，具体涉及一种海洋环境电波传导特性获取方法。

背景技术：

海洋环境主要包括海洋地理、地质、水文、气象和海洋空间电磁环境等。受海洋环境各要素的影响，电波传播需涉及折射、反射、绕射和透射等不同的传播方式。研究海洋环境下电波传播问题时，需涉及电磁波理论、大气科学、传输理论等诸多领域，雷达、通信、数据链、定位系统、电子对抗和遥感遥测等无线电系统应用时需对无线电波的传播特性进行分析。研究海洋环境的变化特点及其对无线电波传播的影响机理，需建立以海洋环境气象要素、海洋环境大气折射要素和海洋环境动态海面及地形要素为主要影响因子的海洋环境电波传播预测模型和雷达探测性能评估模型。

基于民用领域的需求，现有研究人员致力于高效电波传播特性预测算法的研究，建立了许多电波传播预测模型，主要为经统计分析实测数据后归纳出的经验模型、在严格电磁理论基础上推导出的确定性模型和半经验半确定性模型三类。

dockerygd和kuttlerjr在文章animprovedimpedanceboundaryalgorithmforfouriersplit-stepsolutionsoftheparabolicwaveequation中提出了中心差分离散混合傅里叶变换算法，利用二阶中心差分方程拟合阻抗边界条件的一阶微分方程，而后通过构造辅助函数将ssft变换转换为单边离散正弦变换求解，该方法较好的解决了复杂边界条件下的电波传播问题，但是在解决动态阻抗边界条件和考虑海面粗糙度及不规则地形等混合边界条件对电波传播的影响时比较困难。

针对动态海洋表面和不规则地形组成的下表面边界条件下的电波传播问题，dockery、donohue和kuttler等人对下表面边界阻抗条件问题作了大量研究，在dockerygd,kuttlerjr.animprovedimpedanceboundaryalgorithmforfouriersplit-stepsolutionsoftheparabolicwaveequation[j].ieeetransactionsonantennasandpropagation,1996,44(12):1592-1599中，提出了基于分步傅立叶变换法的离散混合傅立叶变换(discretemixedfouriertransform，dmft)抛物方程法，可处理动态海面和不规则地表边界条件下的电波传播问题。然而dmft模型在粗糙海面等动态阻抗边界条件下，随掠射角、电波频率、海面风速变化，存在数值振荡的问题。此时，radioscience,2002,37(2):1-11.janaswamy提出了一种处理阻抗边界条件的改进算法，kuttlerjr,janaswamyr.improvedfouriertransformmethodsforsolvingtheparabolicwaveequation[j].在该文章中，算法更具稳定性，但是算法较为繁琐。

此外，现有技术中对宽角dmft抛物方程和不同地表电磁特性对电波传播影响进行技术研究。例如，胡绘斌,陈建忠.一种基于globe地图的雷达探测范围计算方法[j].电波科学学报,2010,25(1):156-160和王晶,李智,来嘉哲.雷达脉冲电波环境构建方法研究[j].系统仿真学报,2010,22(7):1566-1571。在上述两篇文章中，胡绘斌、王晶等人研究了基于数字地图的不规则地形和城市小区内的电波传播算法，但目前抛物方程法处理动态海洋表面和不规则地形的混合下表面边界条件电波传播问题还不够成熟。

技术实现要素：

本发明的目的是解决现有海洋环境电波传导算法较为复杂及数值稳定性较差的问题，提供一种海洋环境电波传导特性获取方法。

为实现以上发明目的，本发明的技术方案为：

步骤一、通过海洋表面数据获取表面边界阻抗特性α和高度步长δz；

所述海洋表面数据包括最大传播仰角θmax、边界表面的电导率δs、电磁波长λ和折射指数n；

在水平极化和垂直极化条件时，表面边界阻抗特性α计算如下：

其中，εc为边界表面相对复介电常数，由εc＝εs+60δsλi计算得到；i是复数中的虚数；k表示自由空间波数；δs是边界表面的电导率；λ表示电磁波长；θ表示传播仰角；

同时，根据判决条件选取所需的差分方法计算高度步长δz；

δz＝λ/(2sin(θmax))

其中，θmax为最大传播仰角；

步骤二、获取场衰减函数初始场条件；

利用天线辐射方向图f(θ)和口径场分布的傅里叶变换计算dmft算法的场衰减函数初始场条件u(x0,mδz)；

当x＝x0时，场衰减函数初始场条件u(x0,mδz)的计算公式为式(1)；

其中，n表示自然数；γχ表示水平极化或垂直极化fresnel反射系数；δp表示傅立叶变换所需的步进频率，即分辨频率；z0表示发射天线高度；z表示传播区域高度；j表示累加数，从0到n-1；

步骤三、构建辅助函数和离散辅助函数；

3.1)l＝0，即xl＝x0时，根据场衰减函数初始场条件构建初始场下的等效辅助函数w(x0,mδz)和初始离散辅助函数u(x0,jδp)；

其中，m为高度步进数；

3.2)l＞0即xl≠x0时，根据上一距离步进的离散辅助函数结果u(xl-1,jδp)，计算当前离散辅助函数结果u(xl,jδp)，并进行辅助函数w(xl,mδz)的计算；

其中，δx为距离步长；p表示傅里叶变换频率；

步骤四、获取特解辅助函数；

将离散辅助函数w(xl,mδz)的结果代入递推公式(6)计算特解辅助函数up(xl,mδz)；

up(xl,mδz)＝w(xl,mδz)+r·u(xl,(m-1)δz)m＝1，2，…，n(6)

其中，r表示当前位置距发射天线的直线距离；

步骤五、获取中间因子c(xl)；

5.1)计算距离因数a：

5.2)l＝0，即xl＝x0时，根据场衰减函数初始场条件计算中间因子c(xl)，即c(xl)＝c(x0)；

l＞0即xl≠x0时，根据上一步进中间因子c(xl-1)结果迭代计算下一步进中间因子c(xl)；

其中，r^m表示当前位置距发射天线的直线距离的m次方；

步骤六、获取场衰减函数u(xl,mδz)；

通过式(9)获取计算因子b(xl)，代入式(10)计算场衰减函数u(xl,mδz)；

u(xl,mδz)＝up(xl,mδz)+b(xl)r^mm＝0,1,2,…,n(10)

其中，b(xl)表示计算因子；∑'表示m＝0和m＝n项加和时各乘以0.5；

步骤七、令m＝m+1，l＝l+1，返回步骤三，获得不同高度步进数和距离步进数的衰减函数值；

步骤八、通过传播场分布函数获得场强；

通过步骤七获取的不同距离步进数l和高度步进数m的衰减函数值，获取不同高度和距离下的场强，直至口径场内所有高度和距离下的场强值计算完成；场强通过传播场分布函数ψ(xl,mδz)获得，该传播场分布函数的输出即为场强；

fr＝e^{ik[n(z)-1]δx}

ψ(xl,mδz)＝fr·u(xl,mδz)

其中，n表示折射指数，fr表示折射指数。

与现有技术相比，本发明方法的优点是：

1.为了便于数值计算，本发明方法在传统的pe算法基础上采用中心差分dmft算法，利用二阶中心差分方程拟合阻抗边界条件的一阶微分方程，而后通过构造辅助函数将ssft变换转换为单边离散正弦变换求解，该方法较好的解决了复杂边界条件下的电波传播问题。

2.本发明方法采用前向—后向差分代替传统的中心差分方法，克服alpha盲点的问题，提高了算法的稳定性和精度，同时以一阶差分方法拟合阻抗边界条件，对减少算法复杂度提高计算效率有重要意义。

附图说明

图1为电磁波主能量在轴向锥形区域内传播示意图；

图2为传统抛物方程法的badalpha现象示意图；

图3为后向差分的dmft算法badalpha现象示意图；

图4为本发明改进dmft方法结果示意图；

图5为水平极化电波传播因子分布示意图；

图6为水平极化idmft算法与two-ray算法对比示意图；

图7为垂直极化电波传播因子分布示意图；

图8为垂直极化idmft算法与two-ray算法对比示意图；

图9为传播距离为10km出两种算法计算结果示意图；

图10为传播距离为50km出两种算法计算结果示意图；

图11为传播距离为100km出两种算法计算结果示意图；

图12为传播距离为150km出两种算法计算结果示意图。

具体实施方式

以下结合附图和具体实施例对本发明的内容作进一步详细描述。

本发明提供一种海洋环境电波传导特性获取方法，该方法针对海洋环境的电波传播，主要是指在海面上方非均匀分布的大气结构及气象因素条件下，电磁波在动态海洋表面及不规则地形上的传播。本发明方法在分析传统傅里叶(ssft)和中心差分离散混合傅里叶变换(dmft)算法的基础上，提出前向-后向混合差分dmft算法及其具体实现步骤，解决了在计算动态阻抗边界和混合介质边界条件下的数值震荡问题，且修正传统模型在宽角抛物方程模式下求解待定系数和递推公式中存在的错误，在提高计算效率的同时解决传统算法在计算动态混合复杂边界条件下的电波传播问题时由于差分法求解过程中收敛域限制导致数值震荡而出现bad-alpha现象的缺陷。

如图1所示，当电磁波主能量以很小的角度沿x轴正向传播时，可认为其能量主要集中在一个锥形的区域内，传播区域的高度用z表示，定义沿x轴传播的场衰减函数为u(x,z)，设定距离步长δx和高度步长δz，空间坐标(x,z)通过距离步进数l和高度步进数m表示为(xl,mδz)(即xl＝lδx，z＝mδz)，其中l为大于等于0的自然数，m为大于0的自然数，其中x0＝0，表示发射天线位置。

本发明提供的海洋环境电波传导特性获取方法具体为计算复杂边界条件下的电波传播预测模型中的前-后向混合差分dmft算法，其包括以下步骤：

步骤一、通过海洋表面数据获取表面边界阻抗特性α和高度步长δz；

海洋表面数据具体包括最大传播仰角θmax、边界表面的电导率δs、电磁波长λ、折射指数n等；

利用电磁波长及阻抗边界的电磁参数计算表面边界阻抗特性α，在水平极化和垂直极化条件时，表面边界阻抗特性α计算如下：

上述，εc为边界表面相对复介电常数，由εc＝εs+60δsλi计算得到；

其中，i是复数中的虚数；k表示自由空间波数；δs是边界表面的电导率；λ表示电磁波长；θ表示传播仰角；

同时，根据判决条件选取所需的差分方法计算高度步长δz；

δz＝λ/(2sin(θmax))

其中，θmax为最大传播仰角；

步骤二、获取场衰减函数初始场条件；

利用天线辐射方向图f(θ)和口径场分布的傅里叶变换计算dmft算法的场衰减函数初始场条件u(x0,mδz)；

当x＝x0时，场衰减函数初始场条件u(x0,z)的计算公式为式(1)；

其中，n表示自然数；γχ表示水平极化或垂直极化fresnel反射系数；δp表示傅立叶变换所需的步进频率，即分辨频率；z0表示发射天线高度；z表示传播区域高度；j表示累加数，从0到n-1；

步骤三、构建辅助函数和离散辅助函数

3.1)l＝0(即xl＝x0)时，根据场衰减函数初始场条件构建初始场下的等效辅助函数w(x0,mδz)和初始离散辅助函数u(x0,jδp)；

其中，m为高度步进数；

3.2)l＞0即(xl≠x0)时，根据上一距离步进(即xl-1)的离散辅助函数结果u(xl-1,jδp)，计算当前离散辅助函数结果u(xl,jδp)，并进行辅助函数w(xl,mδz)的计算；

其中，δx为距离步长，p表示傅里叶变换频率；

步骤四、获取特解辅助函数；

将离散辅助函数w(xl,mδz)的结果代入递推公式(6)计算特解辅助函数up(xl,mδz)；

up(xl,mδz)＝w(xl,mδz)+r·u(xl,(m-1)δz)m＝1，2，…，n(6)

其中，r表示当前位置距发射天线的直线距离；

步骤五、获取中间因子c(xl)；

5.1)计算距离因数a：

5.2)l＝0(即xl＝x0)时，根据场衰减函数初始场条件计算中间因子c(xl)(即c(xl)＝c(x0))；

l＞0即(xl≠x0)时，根据上一步进中间因子c(xl-1)结果迭代计算下一步进中间因子c(xl)；

其中，r^m表示当前位置距发射天线的直线距离的m次方；

步骤六、获取场衰减函数u(xl,mδz)；

通过式(9)获取计算因子b(xl)，代入式(10)计算场衰减函数u(xl,mδz)；

u(xl,mδz)＝up(xl,mδz)+b(xl)r^mm＝0,1,2,…,n(10)

其中，b(xl)表示计算因子；∑'表示m＝0和m＝n项加和时各乘以0.5；

步骤七、令m＝m+1，l＝l+1，返回步骤三，获得不同高度步进数和距离步进数的衰减函数值；

步骤八、获取传播场分布函数，获得场强；

通过步骤七获取的不同距离步进数l和高度步进数m的衰减函数值，获取不同高度和距离下的场强，直至口径场内所有高度和距离下的场强值计算完成；场强通过传播场分布函数ψ(xl,mδz)获得，该传播场分布函数的输出即为场强；

fr＝e^{ik[n(z)-1]δx}

ψ(xl,mδz)＝fr·u(xl,mδz)

其中，n表示折射指数，fr表示折射指数。

以下对本发明方法的结果进行验证。

计算理想均匀大气条件下的电波传播特性，折射指数n＝1；地表边界为平静海水表面，边界表面相对复介电常数εc＝70，边界表面的电导率δs＝3.76s/m，发射天线高z0＝150m，电波频率f＝2ghz。高斯天线，天线增益30db，天线3db带宽θbw＝7°，天线仰角θ0＝0°，最大传播仰角θmax＝15°，电波水平极化，最大传播距离为50km，dmft算法的wape距离发射天线的水平距离的步长δx＝150m。

传统的中心差分dmft算法由于自身缺陷，将出现严重的数值震荡，如图2所示，算法失效。一阶后向差分dmft算法即可解决大多数条件下的电波计算需求，但是仍存在解不稳定问题。如图3为利用后向差分dmft的计算结果，随着计算距离增加，掠射角逐渐变小，100km后的计算区域出现了明显的数值振荡；图4为相同条件下的本发明方法前-后向混合差分dmft计算结果，与图2和图3相比，具有很好的收敛性。因此，在大尺度风驱海面的动态阻抗条件下，改进dmft算法较传统的dmft算法具有更好的适用性和稳定性

计算理想均匀大气条件下的电波传播特性，折射指数n＝1；地表边界为动态海洋表面边界，海面风速12m/s，边界表面相对复介电常数εc＝70，边界表面的电导率δs＝3.76s/m，发射天线高z0＝150m，电波频率f＝2ghz。高斯天线，天线增益30db，天线3db带宽θbw＝7°，天线仰角θ0＝0°，最大传播仰角θmax＝15°，电波水平极化，最大传播距离为200km，dmft算法的wape距离步长δx＝150m。

图5～图8显示了理想均匀大气及光滑海洋表面边界条件下，水平极化和垂直极化时改进dmft(idmft)算法计算的电波传播因子分布伪彩图及传播距离x＝50km处idmft与two-ray算法计算的传播因子剖面分布曲线。图5和图7中，受地表反射及电波场强干涉叠加效应的影响，远区电波场呈波瓣分布，且波瓣宽度随距离增加逐渐增大；固定距离剖面处以该区域场中值为中心呈规则震荡分布，如图6和图8所示；受地表电磁吸收效应影响，水平极化电波传播衰减较大，较垂直极化情况有更大的峰值振荡和较小的极大值。由图6和图8可以看出，idmft与two-ray算法吻合的很好，验证了idmft算法的正确性。

计算理想均匀大气条件下的电波传播特性，折射指数n＝1；地表边界为平静海水表面，边界表面相对复介电常数εc＝70，边界表面的电导率δs＝3.76s/m,，发射天线高z0＝150m，电波频率f＝2ghz。高斯天线，天线增益30db，天线3db带宽θbw＝7°，天线仰角θ0＝0°，最大传播仰角θmax＝15°，电波水平极化，最大传播距离为50km，dmft算法的wape距离步长δx＝150m。

图9～图12显示出利用改进dmft算法(标注为idmft)和miller-brown模型(图中虚线)在距离为10km、50km、100km和150km垂直剖面上的传播因子分布。两种算法在中值位置处吻合较好，在峰值位置处存在偏差，这是由于miller-brown模型是基于相干波近似得出，其精度随海面风速增大而降低，随传播距离增加而升高。