高精度预测甚低频电波场强随时间变化的方法

1.本发明属于计算电磁学技术领域，涉及一种高精度预测甚低频电波场强随时间变化的方法。


背景技术：

2.按照国际电信联盟(international telecommunication union，iut)对于无线电波频段的划分，3khz～30khz频段范围内的无线电波称为甚低频(very low frequency，vlf)电磁波，其对应的波长为10km～100km，是无线电频谱中极具特点的频段，在应用过程中的传播具有十分明显的特点。在此频段，电波波长较长的无线电波沿地球表面传输时，信号的传播衰减小、传播特性稳定，可在全球范围内进行传播，同时能够渗透一定深度的海水和土壤为水下或者地下设备所接收，因而甚低频电磁波被广泛应用于超远程导航、授时、通信、定位，特别是对潜通信和导航等方面。同时在地震预测、电离层扰动监测等领域也有着较为广泛的应用。
3.由于地面和电离层都具有良好的反射特性，并且这两个反射壁之间的距离可与甚低频电波波长可比拟，因此电波在地面和电离层之间被来回多次反射，且被两个反射壁引导向前传播，此传播过程的机理被称为波导传播。对于在地面和低电离层所构成的地-电离层波导中传播的甚低频电磁波而言，随地球经纬度呈现出复杂空间变化的低电离层、地球磁场、空气层、地层等成为影响其传播特性的主要因素；尤其是低电离层和地球磁场，其除了随昼夜、季节、太阳黑子周期等相对规则的变化外，还存在偶发e层电离层扰动、d层的消失和出现、电离层暴、磁暴等随机的突发情况。这些复杂的电磁环境的变化使得甚低频电磁波呈现出异常复杂的空间变化特性和时间变化特性。受低电离层反射的影响，甚低频电波传播过程中会发生干涉，再加上电波传播本身的衰减，以及地层媒质复杂的空时变化等，使得甚低频场强计算变的十分困难。因此构建更为精确的电离层参数模型，掌握甚低频电波在地-电离层波导结构中的传播特性对于提高甚低频定位导航授时系统的精度等具有重要意义。
4.研究甚低频电磁波在地-电离层波导中的传播特性，现有的解析理论包括：1、奥斯汀经验公式：此公式的应用范围一般在2000～16000km，应用条件是波在白天，海水上传播，因此奥斯汀经验公式在近距离上的计算误差较大；没有近距离上海陆交界的过渡情况；且甚低频的传播始终与接收点电导率有关，而奥斯汀经验公式只考虑了海水上传播，因此要严格分析场强大小就需要探索新的计算方法。2、波跳理论：是基于射线理论原理推导而来，总场有天波和地波两部分组成，其中，天波会受到电离层复杂变化的影响，导致传播不稳定，且波跳理论只适用于近距离的vlf电波传播计算。3、波导模方法：波导模方法是基于甚低频波在地-电离层波导中的传播机理推导而来，地-电离层波导中的vlf总场可以表示为几种模式场的叠加，且每种模式都必须满足亥姆霍兹方程和上下阻抗边界。与奥斯汀经验公式和波跳理论相比，波导模方法具有计算精度高，计算距离远，计算效率高等特点。其中可通过方法3掌握电波传播特点，利用现有的电离层参数模型，进一步验证算法精确度并提
高模型精度。
5.ccir建议低电离层电子密度随高度变化采用指数模型以及国际参考电离层模型(international ionospheric reference model，iri)，指数模型中电离层的参考高度h以及梯度系数β的确定是由实验和探测结果得到的，在特定的时间以及季节，其值是一致的，因此只能反映简单的电离层变化；iri模型只能提供65-2000km范围内的随时间，空间变化的电离层电子密度，电子温度等值，对于65km以下的电离层参数视为零，因此两者都存在一定的缺陷，不能更精确的反映实际电离层的变化。


技术实现要素：

6.本发明的目的是提供一种高精度预测甚低频电波场强随时间变化的方法，解决了现有技术中存在的不能精确预测甚低频电波场强随时间变化的问题。
7.本发明所采用的技术方案是高精度预测甚低频电波场强随时间变化的方法，具体按照以下步骤实施：
8.步骤1、结合iri模型获取随高度变化的电离层电子密度ne；
9.步骤2、将步骤1得到的电离层电子密度ne代入指数模型电子密度的计算公式中，计算各个时刻指数模型电离层参考高度h和梯度系数β；
10.步骤3、基于波导模方法，推导甚低频电磁波的电场分量er；
11.步骤4、通过步骤2的各个时刻指数模型电离层参考高度h和梯度系数β求解电离层的复介电常数，进一步求解电离层的反射系统和表面阻抗；
12.步骤5、通过步骤4的表面阻抗和步骤3的电场分量er，获得不同时间段对应的甚低频电波的场强，既而可以精确预测甚低频电波场强随时间的变化情况。
13.本发明的特点还在于：
14.步骤1结合iri2016模型获取随高度变化的电子密度ne，iri2016模型在线连接地址为http://irimodel.org/。
15.iri2016模型的具体参数设置如下：
16.其中输入参数：
17.电离层高度z：60～150km；时间间隔、分层厚度、接收点位置、时间选取依据实际情况来确定；
18.输出参数为电子密度ne；
19.将所得数据保存为.txt或者.dat格式，方便在程序中调用。
20.步骤2中，结合指数模型电子密度的计算公式：
21.n(z)＝1.43
×
10
13
e-0.15he(β-0.15)(z-h)
(1/m3)
ꢀꢀ
(1)
22.对上式两边化简并取对数得到：
[0023][0024]
上式为求解修正指数模型电离层参考高度和梯度系数的计算公式；
[0025]
将步骤1中iri模型底层的电离层电子密度代入公式2中，求解各个时间点的参考高度h和梯度系数β。
[0026]
步骤3中推导甚低频电磁波的电场分量的公式如下：
[0027]
垂直电偶极子在地-电离层波导中产生电场的表达式具体为：
[0028][0029]
计算公式中的具体参数为：
[0030][0031]
其中idl为垂直电偶极子的电矩，是空气中的波阻抗，k为真空中的波数，θ表示的是接收点与发射点之间的大圆角弧度，p
kw
为发射点处的辐射功率(kw)，a为地球半径，z0为垂直电偶极子距离地面高度，λ为波长；
[0032]
其中x为：
[0033][0034]
λ
′n为波导中垂直电偶极子n阶模产生的激励因子；可表示为：
[0035][0036]
其中tn表示模方程的第n个根，模方程为：
[0037]
a(tn)b(tn)＝1
[0038]
其中：
[0039][0040][0041][0042][0043][0044]
其中，δg为地层归一化表面阻抗，δi为电离层归一化表面阻抗，w1(t)和w2(t)表示第一类和第二类airy函数；
[0045][0046]
式中的积分路径γ
1,2
分别为u的复平面上从-∞e
±
i2π/3
到原点，然后沿实轴延伸至+∞；
[0047]
式中zn(z)表示为高度函数，其表达式为：
[0048]
[0049]
其中
[0050][0051]
其中z为观察点距离地面的高度。
[0052]
步骤4的求解电离层的反射系数和表面阻抗是基于传播矩阵的分层思想进行的。
[0053]
步骤4具体为：
[0054]
在离地高度为h+t处向上看(即z-h≥t处，第m+1层)的电离层表面阻抗z
m+1
为：
[0055][0056]
根据各层之间电磁场的切向分量连续的特性，通过求解可以得到以下递推公式：
[0057][0058][0059][0060]
当第一层的表面阻抗z1求出后，则垂直极化波对电离层的反射系数容易求出，若以z-h＝-z0为参考面，则反射系数可以表示为：
[0061][0062]
由于z-h＝-z0处的电子密度非常小，且通常以z-h＝0平面为电离层的参考面，因此，根据反射系数的定义，可得到z-h＝0参考面的等效反射系数：
[0063][0064]
所以电离层的表面阻抗为：
[0065][0066]
本发明基于iri模型反演指数模型参数，得到随时间变化的电离层参考高度和梯度系数，在波导模理论的基础上，对甚低频电波场强随时间变化的传播特性进行分析预测。基于甚低频波导模理论，利用iri模型反演指数模型参数，获取更为精确的电离层参数模型，计算甚低频电波随时间变化的电场强度并与实测数据进行对比，为研究甚低频电波的传播特性提供参考。
附图说明
[0067]
图1是本发明的流程示意图；
[0068]
图2是本发明实施例中本发明方法、结合指数模型的传统解析方法和实测值的甚低频波信号幅值的日变化对比图；
[0069]
图3为电离层下边界高度对甚低频波信号幅值的影响图。
具体实施方式
[0070]
下面结合附图和具体实施方式对本发明进行详细说明。
[0071]
本发明高精度预测甚低频电波场强随时间变化的方法，如图1所示。其主要推导过程如下：
[0072]
步骤1结合iri2016模型获取随高度变化的电子密度，具体为：
[0073]
iri2016在线模型连接为：http://irimodel.org/；
[0074]
iri模型的具体参数设置为：
[0075]
输入参数：
[0076]
电离层高度z：60～150km，分层厚度为：100m
[0077]
接收点位置(gwalior)：8
°
23
′
n,77
°
45
′e[0078]
时间选取：2008年12月25日00:00-24:00(ut)，时间间隔：0.5h。
[0079]
输出参数：
[0080]
电子密度：ne[0081]
将数据保存为.txt格式，方便在程序中调用。
[0082]
其中，分层厚度为电离层上边界高度和分层厚度，可以根据实际情况来确定；接收点位置也可以选取不同的地理位置，时间间隔可根据实际需要来确定，如0.5-1h均可。
[0083]
步骤2结合指数模型中低电离层电子密度的计算公式：
[0084]
n(z)＝1.43
×
10
13
e-0.15he(β-0.15)(z-h)
(1/m3)
ꢀꢀ
(13)
[0085]
对上式两边化简并取对数得到:
[0086][0087]
上式为求解指数模型电离层参考高度的计算公式。
[0088]
将步骤1得到的电离层电子密度的值代入步骤2得到的计算公式中，得到两个含有两个未知数的方程，求解二元一次方程，即可得到电离层的参考高度和梯度系数，具体步骤如下：
[0089]
白天的电离层参考高度和梯度系数：
[0090]
将iri模型底层65.1km和65.5km处的电子密度代入步骤2得到的计算公式中，求解01:00～18:00各个时间点的h,β。
[0091]
夜晚的电离层参考高度和梯度系数：
[0092]
将iri模型地层80.1km和80.5km处的电子密度值代入步骤2得到的计算公式中，求解19:00-24:00各个时间点的h,β。
[0093]
步骤4基于波导模方法，推导甚低频电磁波的电场分量er：
[0094]
1)求解高度增益函数
[0095]
垂直电偶极子idl所产生的各场的分量可以用位函数u来表示，即
[0096][0097]
位函数u必须满足亥姆霍兹方程
[0098][0099]
利用分离变量法，将位函数表示为
[0100][0101]
其中z(r)为r的高度增益函数；φ(θ)为传递函数，根据分离变量法，结合亥姆霍兹方程，上式可以转换为：
[0102][0103][0104]
其中传递函数φ(θ)的解为勒让德函数
[0105][0106]
高度增益函数的求解即就是求解满足边界条件的微分方程，利用airy函数近似法求解高度增益函数，具体为：
[0107]
对高度增益函数进行简单的变量置换后即能化为标准的stokes方程，令
[0108][0109]
则高度增益函数可以转换为
[0110][0111]
上式为标准的stokes方程，其两种独立解为第一类和第二类airy函数，方程的解可以表示为：
[0112]
z(z)＝aw1(t-y)+bw2(t-y)
ꢀꢀ
(19)
[0113]
考虑到高度增益函数必须满足边界条件：
[0114][0115]
[0116]
将边界条件代入到高度增益函数中，可得到模方程：
[0117]
a(tn)b(tn)＝1
ꢀꢀ
(20)
[0118]
其中tn表示模方程的第n个根，a(tn)和b(tn)的计算公式分别如下：
[0119][0120][0121]
其中
[0122][0123]
为了归一化起见，规定在地面上的高度函数为1，故取
[0124][0125]
2)求解位函数
[0126]
波导中每一个可能存在的波型的位函数可以表示为
[0127][0128]
当观察点不接近于场源与对极点时有
[0129][0130]
其中λn为第n阶模式的激励因子，其表达式为
[0131][0132]
其中
[0133][0134]
3)基于airy函数近似下模方程根的求解
[0135]
理想情况下模方程根的求解方法：
[0136]
模方程公式为：
[0137][0138]
由于在理想情况下，vlf电波传播时没有损耗，即衰减率为0，因此模方程的根tn为实数，此时
[0139]
w1(t)＝u(t)-iv(t),w2(t)＝u(t)+iv(t)
ꢀꢀ
(25)
[0140]
进一步化简模方程可得：
[0141][0142]
令：
[0143][0144]
求解以上方程即可得到模方程根tn。
[0145]
非理想情况下求解方法：
[0146]
一般情况下的模方程为：
[0147]
a(tn)b(tn)＝1
ꢀꢀ
(27)
[0148]
取p＝1/qi，则模方程的根可看成q和p的函数tn(q,p)。
[0149]
当q＝0，p＝0时，tn(0,0)可由理想条件下的方程求出，由a(tn)b(tn)＝1可得，tn(q,p)应满足
[0150]
(w1(t
n-y0)+pw1′
(t
n-y0))(w2′
(tn)-qw2(tn))
[0151]
＝(w2(t
n-y0)+pw2′
(t
n-y0))(w1′
(tn)-qw1(tn))
[0152]
对上式两端q求导，并根据airy函数两阶导的性质和伏龙斯基行列式进行简化后可得
[0153][0154]
同理，对上式两端对p求导，简化后可得
[0155][0156]
将以上两个方程联立，利用初始条件tn(0,0)＝t进行求解，即可得出一般阻抗边界条件下模方程的根。
[0157]
步骤4电离层表面阻抗和反射系数的求解：
[0158]
当z-h较小或者为负值时，不能简单的求解处电磁场高度变化的规律，因此反射系数和整个电离层的等效表面阻抗就不能求解处其解析表达式。因此基于传播矩阵的分层思想，认为将电离层进行分层，将每一层看成是均匀介质。若分层足够多，间距相当小，则可认为和实际的连续变化的模型差别不大。
[0159]
将电离层分为m层，每一层的厚度记为hm，h为电离层等效高度，β为梯度系数，z为观测点距离地面的高度。当β(z-h)≥2时，即深入电离层的内部时，根据麦克斯韦方程可以求得，在离地面高度为z处向上看的电离层表面阻抗zi(z)为：
[0160]
[0161]
其中η0为真空中的波阻抗，k0为真空中的波数，k1为
[0162][0163]
在离地高度为h+t处向上看(即z-h≥t处，第m+1层)的电离层表面阻抗z
m+1
为：
[0164][0165]
根据各层之间电磁场的切向分量连续的特性，通过求解可以得到以下递推公式：
[0166][0167][0168][0169]
其中θ为电离层入射角，可取为第m层中的中值，即
[0170][0171][0172]
其中zm表示第m层离地面高度。
[0173]
当第一层的表面阻抗z1求出后，则垂直极化波对电离层的反射系数容易求出，若以z-h＝-z0为参考面，则反射系数可以表示为：
[0174][0175]
由于z-h＝-z0处的电子密度非常小，且通常以z-h＝0平面为电离层的参考面，因此，根据反射系数的定义，可得到z-h＝0参考面的等效反射系数：
[0176][0177]
所以电离层的归一化表面阻抗为：
[0178][0179]
低电离层碰撞频率的计算公式为：
[0180]
v(z)＝1.83
×
10
11
exp(-0.15z)
ꢀꢀ
(35)
[0181]
电离层的等效复介电常数为：
[0182][0183]
电离层中的传播波数ki为：
[0184][0185]
步骤5将步骤4更新得到的电离层表面阻抗的值代入电场强度的计算公式中，进一步求解不同时间段对应的甚低频电波的电场强度。
[0186]
实施例
[0187]
为了检验本发明的正确性和高效性，基于波导模理论，利用iri模型反演指数模型参数，计算了vtx-gwalior传播路径上接收点处的场强的日变化规律并与实测数据进行了对比。
[0188]
计算参数的设置如下：电离层采用修正的指数模型，电离层取0～150km水平分层各向同性电离层，电离层的分层厚度取为100米。发射点为印度甚低频发射站vtx，其地理坐标为北纬08
°
26
′
，东经77
°
44
′
，接收的位于印度的gwalior，其地理坐标为北纬26
°
14
′
，东经78
°
10
′
。发射频率为f＝18.2khz，接收点与发射点之间的大圆距离为1976km。时间选取为2008年12月25日，电离层参数选取00:00-24:00的电离层数据。空气层取为真空；地面厚度为60km的平均陆地，地面参数选取为平均陆地其介电常数为ε＝22，电导率为σ＝3e-3s/m；辐射功率p＝1000w，空间网格δr＝1km，其电离层参考高度和梯度系数有本发明方法计算而来。
[0189]
采用本发明方法、采用传统的电离层指数模型的波导模方法，并采取实测结果进行对比，结果如图2所示。从图中可见，本发明的计算结果与实测数据变化趋势吻合较好。同传统的指数模型相比，提高了预测精度，并且分析了甚低频波场强随时间的变化特性，验证了本发明的正确性和高效性。
[0190]
如图3所示，比较了不同电离层下边界高度的设置对甚低频电波场强的影响，发现在65km以下，电离层的下边界高度设置对于甚低频电波的场强大小，以及在传播路径上的场强分布都有一定的影响。
[0191]
本发明与只用iri模型进行甚低频波信号的场强预测方法相比，修正之后的指数模型考虑到了iri模型在65km和80km以下电子密度为零的情况，同时修正了传统指数模型对于在冬季白天和夜晚(夏季白天和夜晚)的电离层参考高度和梯度系数的值是确定的，不随时间及地理位置的变化而变化。而修正后的指数模型在一天内的不同时间的电离层参考高度和梯度系数不同，更为接近实际的电离层。