一种电波射线上大气折射率计算方法与流程

本发明属于电波折射误差修正技术领域，具体涉及到一种电波射线上大气折射率计算方法。

背景技术：

随着电子器件精度的进一步提高，以及信号与数据处理方法的进一步发展和优化，作为现代战争中不可或缺的雷达系统得到了广泛应用，雷达系统的探测精度主要取决于电波折射误差修正精度。电波折射误差主要是由于雷达的电波通过不均匀的大气而产生的，根据目前的电波折射误差修正研究和实际应用情况，电波折射误差计算公式本身的精度较高，影响其修正精度的关键因素是电波射线上大气折射率的精确度，因此要提高电波折射误差修正精度的关键是要提高电波射线上的大气折射率的准确度。

对流层大气中的折射率主要由直接探测法、统计法、公式法和反演法等四种方法获得。直接探测法是利用“五九”型探空仪、电子探空仪、adas探测系统等探测仪进行不同高度处的大气温度、气压和湿度等参数探测，然后利用大气折射率与温度、气压和湿度的函数关系计算出对应高度处的大气折射率值，进而得到雷达所在地的大气折射率剖面，它是一种由高度和折射率组成的离散剖面。统计法是通过大量的雷达所在地的大气折射率剖面的统计，给出一种大气折射率的数学模型，目前该模型精度最高的是分段模型。公式法是利用国际上给出的hopfild、saastamoinen等公式计算大气折射率剖面。反演法是利用微波辐射计或gps等在雷达所占地进行测量，进而反演出该地的大气折射率剖面。就获得大气折射率剖面的精度而言，这四类方法中直接探测法的精度最高，其次是统计法。直接探测法获得的大气折射率剖面是实时测量，基本上不存在大气的时变误差；统计法是根据本地的大气折射率实测剖面进行统计，尽管它具有一定的误差，但建立的模型更有利于电波折射误差计算。而公式法中的公式参数是利用国外的气象参数统计得到，国外的气象参数与我国的参数存在误差，因此该方法在我国应用会存在较大的误差，利用微波辐射计或gps的反演法由于存在反演的精度和仪器本身的误差，其精度也较低。

在下垫面平坦地区，针对均匀的大气结构地区，大气折射率在水平方向上的变化很小，因此可将大气折射率视为只随高度变化的一维剖面，即将大气视为球面分层大气，这种情况下，只要利用直接探测法得到雷达所在地的大气剖面，电波射线任意位置的折射率可通过相同高度实测折射率的平移得到。针对不均匀的大气结构地区，大气折射率在水平方向上会发生变化，大气折射率的获取不能单一的视为只随高度的变化，应考虑水平方向上的变化，因此大气折射率应采用三维大气剖面。在下垫面复杂地区，由于大气折射率在水平方向上变化较大，因此大气折射率更应采用三维大气剖面。

事实上，大气在水平方向上的变化很难获得，如果采用在电波射线经过区域进行有限数量的大气折射率直接探测，不仅需要花费大量的人力和物力，而且也不可能得到精确的水平变化特征。目前在所有的电波折射误差修正中几乎都采用了大气球面分层法获得大气折射率的一维剖面，这对于下垫面大气均匀地区是可行的，但是对于下垫面大气不均匀地区，尤其是下垫面复杂的地区，会出现较大的折射率误差，进而影响电波折射误差修正精度。

技术实现要素：

本发明提供了以一种电波射线上大气折射率计算方法，能够提高不同地区任意大气结构的电波折射率的计算精度。

为实现上述目的，本发明采用如下技术方案：

一种电波射线上大气折射率计算方法，包括：

(1)根据我国地理位置范围和无线电气象环境变化的特点，采用大气环境栅格技术将全国划分成1840个栅格，建立全国大气折射率剖面数据库；

(2)对于均匀的大气结构地区，只考虑大气折射率随高度的变化，将大气视为球面分层大气，利用直接探测法得到雷达所在地的大气剖面，并通过相同高度实测折射率的平移来获取电波射线任意位置的折射率；

(3)对于不均匀的大气结构地区以及下垫面复杂地区，采用三维大气剖面，当电波射线的开始位置是在雷达所在的栅格内，利用直接探测法得到雷达所在地的大气折射率剖面；当电波射线落到其他栅格内时，先计算出射线点的位置，然后再利用已建立的全国大气剖面模型数据库得到该射线点的大气折射率剖面，利用此折射率剖面进行计算，得到电波射线上的大气折射率。

优选的，全国大气折射率剖面数据库的建立是根据我国地理位置范围和无线电气象环境变化的特点，采用大气环境栅格技术将全国划分成1840个栅格，然后利用分段模型对大气折射率进行计算，其模型为：

其中，h0为地面海拔高度，km；n0为地面大气折射率；g为地面到1km高度的折射率梯度，1/km；n1为地面1km高度上的大气折射率；n9为海拔9km高度上的大气折射率；c为地面1km至海拔9km大气层的衰减系数，1/km；c9为海拔9km以上高度大气层的衰减系数，1/km，下标0为区分作用，h表示目标指定的海拔高度，单位为km；

利用近20年全国大气探空数据对每个栅格按月份进行统计，各个探空站每天都有2次大气剖面实际探测数据，因此每个栅格内每月共有1200组探空数据，在这些数据中除去特殊天气外的数据后按式(1)进行统计，可得到1200组参数n0i、gi、ci和c9，再将所有参数进行回归统计，就可建立各个栅格内每个月的g与n0、c与n0、c9与n0之间的关系式，最后建立由栅格、地面高度h0、月份、g与n0、c与n0、c9与n0的回归系数组成的全国大气折射率剖面数据库。

优选的，当栅格内没有探空数据时，选用相邻有探空站的四个栅格内的同时间数据，利用拉格朗日插值公式求得该栅格内的同时间大气剖面数据，这样就使得各个栅格内都具有大气剖面数据。

优选的，当雷达电波射线的开始位置是在雷达所在的栅格内时，大气折射率剖面由直接探测方法得到，由于电波在大气传播时会发生折射现象，所以通常用射线来表示电波传播的轨迹，大气折射率剖面是描述大气折射率沿高度变化的数学公式或相关数据，因此计算大气折射率剖面时，需要首先计算电波传播的射线轨道，其具体计算过程如下：

已知射线初始位置上的大气折射指数n0、矢径r0和初始仰角θ0，矢径r0表示当地地球半径，则由斯奈尔定律的变形式可得出：

θ表示目标在海拔h处电波射线上的仰角，n表示大气折射指数，r表示当地地球半径；

根据三角关系，可推导出：

假设射线上某点的位置用水平距离d和高度h表示，当该点沿射线移动到另外一点时，则该点的半径和高度都会有所变化，即：

其中，是地心张角的微变量，根据数学几何关系有：

其中，dd是水平距离d的微变量，结合式(5)和(6)，且r＝r0+h，则有：

将等式两边同时求积分，则当海拔高度h取0至ht时所对应地面距离d可表示为：

其中，ht表示实测大气折射率剖面最高高度，由上式得知，d是ht的函数，当高度h连续做微小的变化，即h依次取0、δh、2δh、…、nδh时，根据上述函数关系式可得到高度对应的地面距离d分别为0、d1、d21、…、dn，当δh→h时，相邻的离散点可视为连续，进而形成了射线轨道，δh为随机选取的数据，表示高度h的变化量。

优选的，在射线点的位置的计算过程中，采用空间大地直角坐标系(x，y，z)和大地直角坐标系(l，b，h)来描述电波射线上的任一点，在空间大地直角坐标系中，地球中心为坐标原点o，z轴是地球的旋转轴，指向地球自转轴方向，x轴指向起始大地子午面与赤道的交点，y轴与x轴、z轴构成右手坐标系，x轴和y轴构成赤道面，空间的任一点p在大地坐标系中的位置用大地经度l、大地维度b和大地高度h表示，其中，过p点大地子午面与起始大地子午面的夹角为大地经度l，p点对椭球面的法线与赤道面的夹角为大地维度b，p点沿法线至椭球面的距离为该点的大地高度h；

假设雷达探测到的视在距离、仰角和方位角分别为re、θe和下标e为区分作用，并非变量，由于大气的水平方向变化很小，方位角产生的折射误差很小，假设雷达电波射线的起始点的大地坐标为(l0，b0，h0)，空间大地直角坐标为(x0，y0，z0)，设定在电波射线的高度步长为δh，则电波射线上任意一点i离地面的高度hi＝hi-1+δh，下标i、i-1和0均为区分作用，并非变量；

由雷达测量参数可以得到电波射线上任意一点上的空间大地直角坐标为：

这里，i＝1，2，3，…，i-1；

根据电波传播理论可以得到

式中，n0和ni分别为地面和海拔i处的大气折射率，r0为雷达到地心的距离；

将空间大地直角坐标转换成大地坐标系对应参数(li，bi)，即

其中，a为地球长半轴的长度，a＝6378140m；e²为地球第一偏心率的平方量值，e²＝0.0066943848，

利用式(12)即可计算出射线点的位置(li，bi)。

优选的，计算电波射线上大气折射率时，

①在雷达电波射线开始点处采用当地直接探测到的地面大气折射率n0＝1+n0×10^-6；

②其他位置的大气折射率计算时，首先由步骤5得到电波射线任意点的大地坐标(li，bi)，然后根据(li，bi)判别该点所处的栅格；

当该点处于雷达所在地的栅格时，采用该地直接探测的大气折射率剖面；当该点处于其他栅格时，首先利用水平插值方法获得该栅格的地面大气折射率n0，然后利用全国大气剖面模型数据库得到gi、ci和c9，最后利用公式(1)得到该栅格的大气折射率模型，对该模型进行高度差值计算即可获得电波射线处的大气折射率值。

本发明针对不同的大气结构地区，对电波射线上大气折射率进行计算，能够有效地提高平坦地区大气结构不均匀情况以及下垫面复杂地区的大气折射率获取和电波折射误差修正精度；针对均匀的大气结构地区，由于大气折射率在水平方向上的变化很小，只考虑大气折射率随高度的变化，将大气视为球面分层大气，利用直接探测法得到雷达所在地的大气剖面，电波射线任意位置的折射率可通过相同高度实测折射率的平移获取；针对不均匀的大气结构地区以及下垫面复杂地区，大气折射率在水平方向上会发生变化，大气折射率的获取不能单一的视为只随高度的变化，应考虑水平方向上的变化，因此大气折射率采用三维大气剖面，从而较为精确地获得电波射线上大气折射率的组合大气折射率测量方法。本发明在工程应用上容易实现，其计算误差小，进而能够有效地提高电波折射率的计算精度。

附图说明

图1为本发明的流程图；

图2为本发明所述射线描迹法中电波传播射线的路径图；

图3为本发明所述的大地坐标系与空间大地直角坐标系的示意图；

图4为本发明所述实施例中的目标飞行航迹图。

具体实施方式

以下结合附图对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本实发明的部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的其他所有实施例，都属于本发明的保护范围。

如图1所示，本发明所述的一种电波射线上大气折射率计算方法，包括：

(1)根据我国地理位置范围和无线电气象环境变化的特点，采用大气环境栅格技术将全国划分成1840个栅格，建立全国大气折射率剖面数据库；

全国大气折射率剖面数据库的建立是根据我国地理位置范围和无线电气象环境变化的特点，采用大气环境栅格技术将全国划分成1840个栅格，目前我国能够常年进行大气参数探测的气象站共有一百多个，绝大多数栅格内没有探空数据，只有少量的栅格内有实测探空数据，当栅格内没有探空数据时，选用相邻有探空站的四个栅格内的同时间数据，利用拉格朗日插值公式求得该栅格内的同时间大气剖面，这样就使得各个栅格内都具有大气剖面数据。利用拉格朗日插值公式求得该栅格内的同时间大气剖面的过程为现有技术，不再赘述。

利用分段模型对大气折射率进行计算，分段模型具有较高精度，其模型为：

其中，h0为地面海拔高度，km；n0为地面大气折射率；g为地面到1km高度的折射率梯度，1/km；n1为地面1km高度上的大气折射率；n9为海拔9km高度上的大气折射率；c为地面1km至海拔9km大气层的衰减系数，1/km；c9为海拔9km以上高度大气层的衰减系数，1/km，下标0为区分作用，h表示目标指定的海拔高度，单位为km；

利用近20年全国大气探空数据对每个栅格按月份进行统计，各个探空站每天都有2次大气剖面实际探测数据，则每个栅格内每月共有1200组探空数据，在这些数据中除去特殊天气外的数据后按式(1)进行统计，可得到1200组参数n0i、gi、ci和c9，再将所有参数进行回归统计，就可建立各个栅格内每个月的g与n0、c与n0、c9与n0之间的关系式，最后建立由栅格、地面高度h0、月份、g与n0、c与n0、c9与n0的回归系数组成的全国大气折射率剖面数据库。其中，将所有参数进行回归统计，以及建立由栅格、地面高度h0、月份、g与n0、c与n0、c9与n0的回归系数组成的全国大气折射率剖面数据库的过程均为现有技术，不再赘述。在实际应用中，只要确定月份，再通过插值法得到每个栅格内的地面大气折射率n0，利用该数据库可得到gi、ci和c9，最后利用公式(1)即可得到该区域的大气折射率模型n(h)。

(2)对于均匀的大气结构地区，只考虑大气折射率随高度的变化，将大气视为球面分层大气，利用直接探测法得到雷达所在地的大气剖面，并通过相同高度实测折射率的平移来获取电波射线任意位置的折射率；

(3)对于不均匀的大气结构地区以及下垫面复杂地区，采用三维大气剖面，当电波射线的开始位置是在雷达所在的栅格内，利用直接探测法得到雷达所在地的大气折射率剖面；当电波射线落到其他栅格内时，先计算出射线点的位置，然后再利用已建立的全国大气剖面模型数据库得到该射线点的大气折射率剖面，利用此折射率剖面进行计算，得到电波射线上的大气折射率。

当雷达电波射线的开始位置是在雷达所在的栅格内时，大气折射率剖面由直接探测方法得到，由于电波在大气传播时会发生折射现象，所以通常用射线来表示电波传播的轨迹，大气折射率剖面是描述大气折射率沿高度变化的数学公式或相关数据，因此计算大气折射率剖面时，需要首先计算电波传播的射线轨道，在本实施例中采用射线描迹法进行计算，具体计算过程如下：

如图2所示，已知射线初始位置上的大气折射指数n0、矢径r0和初始仰角θ0，矢径r0表示当地地球半径，则由斯奈尔定律的变形式可得出：

θ表示目标在海拔h处电波射线上的仰角，n表示大气折射指数，r表示当地地球半径；

根据三角关系，可推导出：

在图2中，假设射线上某点的位置用水平距离d和高度h表示，当该点沿射线移动到另外一点时，则该点的半径和高度都会有所变化，即：

其中，是地心张角的微变量，根据数学几何关系有：

其中，dd是水平距离d的微变量，结合式(5)和(6)，根据图2可知r＝r0+h，则有：

将等式两边同时求积分，则当海拔高度h取0至ht时所对应地面距离d可表示为：

其中，ht表示实测大气折射率剖面最高高度，由上式得知，d是ht的函数，当高度h连续做微小的变化，即h依次取0、δh、2δh、…、nδh时，根据上述函数关系式可得到高度对应的地面距离d分别为0、d1、d21、…、dn，当δh→h时，相邻的离散点可视为连续，进而形成了射线轨道，δh为随机选取的数据，表示高度h的变化量。

如图3所示，在射线点的位置的计算过程中，采用空间大地直角坐标系(x，y，z)和大地直角坐标系(l，b，h)来描述电波射线上的任一点，在空间大地直角坐标系中，地球中心为坐标原点o，z轴是地球的旋转轴，指向地球自转轴方向，x轴指向起始大地子午面与赤道的交点，y轴与x轴、z轴构成右手坐标系，x轴和y轴构成赤道面，空间的任一点p在大地坐标系中的位置用大地经度l、大地维度b和大地高度h表示，其中，过p点大地子午面与起始大地子午面的夹角为大地经度l，p点对椭球面的法线与赤道面的夹角为大地维度b，p点沿法线至椭球面的距离为该点的大地高度h；

假设雷达探测到的视在距离、仰角和方位角分别为re、θe和下标e为区分作用，并非变量，由于大气的水平方向变化很小，方位角产生的折射误差很小，可以忽略不计，假设雷达电波射线的起始点的大地坐标为(l0，b0，h0)，空间大地直角坐标为(x0，y0，z0)，设定电波射线的高度步长为δh，则电波射线上任意一点i离地面的高度hi＝hi-1+δh，下标i、i-1和0均为区分作用，并非变量；

由雷达测量参数可以得到电波射线上任意一点上的空间大地直角坐标为：

这里，i＝1，2，3，…，i-1；

根据电波传播理论可以得到：

其中，n0和ni分别为地面和海拔i处的大气折射率，r0为雷达到地心的距离；

为了从全国大气剖面模型数据库中选择电波射线上任意点所对应的位置，将空间大地直角坐标转换成大地坐标系对应参数(li，bi)，即

其中，a为地球长半轴的长度，a＝6378140m；e²为地球第一偏心率的平方量值，e²＝0.0066943848，

利用式(12)即可计算出射线点的位置(li，bi)。

在计算电波射线上大气折射率时，

①在雷达电波射线开始点处采用当地直接探测到的地面大气折射率n0＝1+n0×10^-6；

②其他位置的大气折射率计算时，首先由步骤5得到电波射线任意点的大地坐标(li，bi)，然后根据(li，bi)判别该点所处的栅格；

当该点处于雷达所在地的栅格时，采用该地直接探测的大气折射率剖面；当该点处于其他栅格时，首先利用水平插值方法获得该栅格的地面大气折射率n0，然后利用全国大气剖面模型数据库得到gi、ci和c9，最后利用公式(1)得到该栅格的大气折射率模型，对该模型进行高度差值计算即可获得电波射线处的大气折射率值。此过程为现有技术，不再赘述。

本发明的具体使用过程为：为了验证针对不同地区任意大气结构的电波射线上大气折射率获取方法，根据下垫面复杂地区大气结构的特点，利用某基地地面雷达对目标的实际测量参数和目标的真实位置、以及当时实测雷达处的大气剖面和全国大气折射率剖面模型数据库进行折射误差的计算和比对。首先在该地区用带有gps的直升机作为目标，该直升机飞过的位置作为真实位置，利用大气环境栅格技术建立全国大气折射率剖面数据库，在雷达所在地采用直接探测法，雷达电波射线的开始位置是在雷达所在的栅格内，大气折射率剖面由直接探测方法得到；在其他电波射线上的位置，先计算出射线点的位置，然后再利用所根据得我国地理位置范围(经、纬度)和无线电气象环境变化的特点，采用大气环境栅格技术建立全国大气折射率剖面数据库；选已建立的全国大气剖面模型数据库得到该位置的大气折射率，从而较为精确地获得电波射线上的大气折射率；当电波射线落到其他栅格内时，电波射线点位置的计算则通过数据库获取该栅格内的电波波射线点大气折射率剖面，再用此折射率剖面进行折射误差计算，从而减小电波射线上的折射率误差，从而得到在工程应用上容易实现、其修正精度高、又能有效地提高电波折射误差修正精度的电波射线上大气折射率计算方法。

以下结合具体实施例对本发明的技术方案进行进一步描述。

以新乡为例，用一个携带有gps的飞机作为目标，并以下垫面复杂地区为例来分析，试验中的目标是一个携带有gps的飞机，它能给出目标的真实位置，其目标飞行航线如图4所示，计算出航线点的位置；雷达所在地则采用直接探测法，然后再利用大气环境栅格技术建立全国大气折射率剖面数据库，根据已建立的全国大气剖面模型数据库得到该位置的大气折射率，从而较为精确地获得电波射线上的大气折射率。通过飞机上所在地的大气剖面直接探测法与利用栅格技术建立全国大气折射率剖面模型数据库形成的组合方法，再用某基地地面雷达对目标的实际测量参数和目标的真实位置、以及当时实测雷达处的大气剖面和全国大气折射率剖面模型数据库进行折射误差的计算和比对，从而得到电波射线上大气折射率，该发明能够提高电波射线上折射率的准确度及电波折射误差修正的精度。