一种基于椭球分层模型的电波折射修正算法

唐学海刘伯阳

①(中国人民解放军63610部队 库尔勒 841001)

②(中国人民解放军63612部队 敦煌 736200)

一种基于椭球分层模型的电波折射修正算法

唐学海*①刘伯阳②

①(中国人民解放军63610部队 库尔勒 841001)

②(中国人民解放军63612部队 敦煌 736200)

传统电波折射修正算法普遍采用大气球面分层假设，这类算法在修正高仰角测量目标时具有较好的修正精度，然而对于低仰角、远距离目标，修正精度还不高。该文提出一种折射修正算法，采用更精确的椭球面分层模型描述大气分布，利用迭代递推的方法计算修正后的目标位置，相比传统折射修正算法，计算量有所增加，但是提高了低仰角、远距离目标测量数据折射修正精度，可用于事后数据处理。

电波传播；折射修正；椭球分层模型；迭代递推

1 引言

由于地球周围充满大气，而空间各处大气成分、密度、温度、湿度和电离度、电子浓度均不相同，介质特性也相当复杂，由此导致电波或光波通过大气传播时，其传播过程不是匀速直线运动，这就造成了外测测量的折射误差。历次试验任务数据的计算结果表明，当目标到测站的距离为1000 km左右、俯仰角在2°～4°时，大气折射导致的视在位置与实际位置间的偏差能达到5000～3000 m。由此可见，远距离、低仰角的测量数据受大气折射的影响较大，折射误差修正对提高外测数据测量精度起着至关重要的作用。

折射修正首先要知道大气层的结构特性，由于地球周围大气层的结构极其复杂，地表60 km高度以上大气层受时间、天气、季节、地形地貌等因素影响，很难精确描述。然而在一定范围内其分布还是有规律可循的，传统折射修正模型[1−9]就应用了大气球面分层假设进行近似描述，此模型假定地球为正球体，大气结构在水平方向是均匀的，这样的近似描述方法在用于高仰角目标的折射修正时具有较高的精度，但对于低仰角目标，尤其是对低仰角、远距离目标的折射修正精度并不高，不能满足精密测量的需求[10−12]。其中一个重要的原因是由于地球自转及自身曲率的存在，使得大气层的整体分布在水平方向存在差异，不同纬度下的高程-折射率变化梯度明显不同，纬度越大梯度越大，纬度越小梯度越小。所以对于低仰角测量目标来说，电磁波在大气层内以较小角度远距离传播，这种水平方向的差异对电磁波折射的影响就尤为明显[13,14]。为此，本文采用大气椭球面分层假设，基于此假设提出了新的电波折射修正算法，结果表明此算法在低仰角、远距离目标的折射修正中效果优于传统球面分层模型的折射修正算法。由于此算法采用了迭代递推的方法，计算量相对于传统算法有所增加，因此可用于事后数据处理。

2 大气椭球面分层模型建立

为了提高对大气层分布的描述精度，现提出大气椭球面分层假设，即假设：

(1)地球为椭球体，地球周围的大气层也成椭球体分布；

(2)将大气层由数层同几何中心同曲率的椭球面分层，椭球面曲率与地球曲率相同(如图1所示)。

图1 大气椭球面分层模型

根据上述假设可以得到如下结论：

(1)与地球同心的任意两相邻椭球面所夹薄层内的大气物理参数，如气压、气温、湿度和电离度相同，此时在薄层内的电(光)波波速相同，折射率也相同。也就是说，电(光)波的波速v及折射率n只与所处薄层有关。

(2)相邻薄层的入射波线、折射界面法线与折射波线处于同一平面。

目前主要采用气球探空仪探测空间大气参数，计算定点上空的折射率，通常在试验任务前、后半小时到两小时释放两次探空气球，获得的大气参数可分别用于任务实时折射修正和事后折射修正。为将定点上空的大气参数转换到较大空域内，现采用下述方法：

根据大气椭球面分层假设，大气中任意点(x,y,z)所处的椭球面方程为

其中a,b,c分别为椭球面的3个半长轴，有如下关系：

其中f为地球扁率常数。

如图2所示，在A点释放探空气球，可测得A点上空的大气参数，计算可得不同高度的折射率。也就是说，通过测量数据计算可获得A点上空任意高度的折射率n与高度h间的函数，记为n=n(h)；根据大气椭球面分层假设，不处于A点上空的B点处的折射率等于与B点处于同一椭球薄层并位于A点上空的B'点处的折射率。将A点坐标(xA,yA,zA)代入式(1)，并联立式(2)可解出(a,b,c)的值(aA,bA,cA)。对式(1)左端分别沿x,y,z方向求偏导并代入点A坐标可得点A处的切面法向量，由此可得过点A的椭球体法线方程为

图2 折射率空间变化分布

同理，可求B点所在的椭球面方程为

联立式(3)，式(4)可求出过A点的法线与B点所处的椭球面的交点B'(xB′,yB′,zB′)，向量AB'的模即为B'点距地面的高度hB′，将H(B)称为B点处的气象等效高程。

根据大气椭球面分层假设，B点与B'点处折射率相等，可求出B点处的折射率为

根据大气椭球面分层假设，利用定点的探空数据可求出周围较大范围内的大气分布情况，准确度较球面分层假设有一定提高，但通过计算获取的大气结构毕竟不能反映真实电波传播路径上的大气结构，尤其在地形复杂的试验场区差异更大，希望将来可通过更先进的大气探测方法获取更加精确的3维大气结构模型。

3 迭代递推法修正折射误差

电磁波从雷达发出后，其在大气层中的传播是在无数层介质间被不断折射后抵达探测目标(如图3所示)，由于波路可逆，回波波路与其完全相同。如图3所示，假设有一束电磁波自分层界面Ni−1传播至界面Ni，经过折射后传播至分层界面Ni+1，由第2节中大气椭球面分层假设可知，每一薄层中的折射率处处相等，即在薄层内电磁波沿直线传播。记Si为分层界面Ni与波线的交点，ni−1为分层界面Ni−1与Ni之间的折射率，vi为Si点处切平面的法向量，波线向量Si−1Si记为si−1，在分层界面Ni两边的入射角和折射角分别为θi1和θi2，在大气椭球面分层假设下，电(光)波在大气层内的传播规律满足斯涅尔定律，从而有

图3 电波在大气层中折射传播示意图

其中

在实际应用中，ni−1为分层界面Ni−1与Ni之间的折射率均值，计算方法如式(10)：

其中n(h),H(Si−1)和H(Si)的计算方法第2节已给出。

根据雷达测量数据，已知波线的起始点S0，其坐标(x0,y0,z0)为雷达阵面中心的地心系坐标(雷达站址坐标)；s0=S0S1为起始向量，其方向为波线在雷达接收面处的切线方向，即为雷达实测方向，长度为dr(选定的递推步长)；雷达实测距离为R0，波束在空间中单向传递时间为TR，由雷达测距原理知

其中c为真空中的光速。下面由上述已知量通过迭代递推法计算经折射误差修正后目标的地心系坐标。

在选定递推步长dr的前提下，将大气按椭球面分层，在每一薄层中，入射波线(折射波线)的长度均为递推步长dr，即|si|=dr 。由S0(x0,y0,z0)为起点，通过递推依次计算S1,S2,…，直到目标点S的地心系坐标。

首先计算S1点的坐标。记地心为O，则OS0+ S0S1=OS1，将OS0写作向量形式为[x0y0z0], OS1的向量形式为[x1y1z1]，从而有[x1y1z1]=[x0y0z0] +s0。

下面计算S2点的坐标：与上面的计算类似，有

[x1y1z1]已由上一步计算得出，现通过迭代法计算s

1。由于|s1|=dr，可设已知S1点坐标，由第2节的推导可计算出S1点处的椭球面法向量为

由上述推导，联立方程组

其中sin(θ11)和n0可分别由式(9)和式(10)计算得到，将s0作为s1的初值，代入式(12)，得S2的初值，再将S2代入式(10)，得n1的初值，最后将n1代入式(14)可计算得到s1的值，记为。将再代入式(12)，重复上述步骤，又可得到新的s1的值，记为。依次迭代下去，可得到,…，直到某个满足(ε为设定的迭代控制常数)，此时认为满足精度要求，可作为s1的值，由此得到S2点的坐标。

继续同样的递推，可依次得到点S3,S4,…,Sk的坐标，k的值由如下方法确定：

已知TR是雷达实测出的脉冲单向传递时间，设脉冲波束从雷达发出后在每个薄层中的传播时间为ti，则

其中k为雷达至目标间大气被分层后的层数，也为总的递推次数。根据式(15)求得薄层内波线向量si的模长为

在递推过程中，将从雷达发射机到Si点处电磁波传播时间记为iT，则

持续递推，当TR−Ti≤ti时，递推总次数k=i+1。因此，由上述递推方法计算出的点Sk的坐标即为目标点S的坐标。

4 模型实例修正效果及分析

为检验本文提出模型的实际修正效果，本节选取某次试验任务中的雷达实测数据作为目标测试数据，分别使用本文介绍的折射修正算法与传统球面分层折射修正算法对原始数据进行修正。此次任务中被测目标相对于雷达距离远、仰角低，雷达测量精度受折射影响较大，且本次任务中被测目标携带GPS定位装置，可将GPS遥测航迹作为标准航迹，使用标准航迹与经过两种不同折射修正算法修正的位置数据进行比较，可以更加直观地体现不同修正算法之间的差异。修正过程中使用的大气折射率数据来自同一数据源，均由分析计算此次任务中的气象探空数据获得。

将雷达测得的目标原始地心系坐标记为(xp,yp,zp)，即视在位置；将目标GPS遥测地心系坐标记为(xg,yg,zg)，即标准位置；将经过椭球面分层折射修正的目标地心系坐标记为(xe,ye,ze)；将使用传统球面折射修正算法修正过的目标地心系坐标记为(xr,yr,zr)。则视在位置与标准位置之间的位置差为

图4所示为使用目标数据计算得到的δp-g(t)随时间的变化关系图，可以看出雷达测量数据精度受折射影响极大，未经过折射修正前，同一时刻目标视在位置与实际位置间最大距离大于4000 m。

经椭球面分层折射修正算法修正过的坐标与标准位置之间的位置差为

经传统球面分层折射修正算法修正过的坐标与标准位置之间的位置差为

图5所示为对目标数据分别使用两种修正算法修正得到的坐标与标准位置之间的位置差，即δe-g,δr-g随时间变化关系图。对比图4可知，经过折射修正后的位置精度较修正之前得到了很大提高；使用椭球面分层折射修正算法获得的目标位置更接近标准位置。

下面分别将经过两种算法修正后的目标地心系坐标转换到测站系下，称为修正后测站系坐标，与修正前的测站系坐标(即雷达原始测量数据)进行比较，分析两种算法的区别。图6所示为不同算法折射修正结果与原始雷达测量数据在测站系下的坐标差值随时间的变化关系，其中图6(a)，图6(b)和图6(c)分别为两种不同修正算法在斜距、方位角和俯仰角上的修正差异。

通过在测站系下的比较，可直观地体现出两种不同修正算法的差别，即修正过程中对电磁波传播波路恢复的差异。具体来说，由于折射率取值方式不同，使得两种修正模型在同一空间位置使用的折射率不同，主要造成了两种修正模型在斜距和方位角上的修正差异；由于不同模型假设的分层形状存在差异，传统球面分层折射算法不会改变目标的方位角[15]，也就是说认为电磁波全程都在同一平面内传播，而本文算法在修正的过程中会使波线在水平方向上发生偏转；通过上文的比较可知这样的修正结果与电磁波实际传播路径更为吻合。

另外，由于本文所介绍的折射修正算法采用迭代递推的计算方法，所以其修正误差与步长有直接关系。由于减小步长会增加运算量，所以需要一个较为合适的递推步长dr，使得在满足精度要求的前提下尽量提高运算速度。将步长dr=1 m时计算出的修正后坐标作为基准位置，经过实验得出，当步长dr≤5000 m 时，修正后的位置与基准位置间距离小于10 m；当步长dr≤1000 m ，修正后的位置与基准位置间距离小于1.5 m。

图4 视在位置与标准位置间的位置差

图5 使用两种修正算法得到的坐标与标准位置的差异

图6 不同算法折射修正结果与原始雷达测量数据在测站系下的坐标差

传统折射修正算法采用球体模型，球面上每一点的法线方向均指向球心，因此可以先在地心矩方向对电磁波的传播路径进行整体积分，然后对运算结果进行整体迭代运算，修正折射方向估计偏差；而本文方法基于椭球体模型，球面上每一点的法线不都指向地心，且传播路径并不位于同一空间平面，因此无法在地心矩方向整体积分，所以采用了先分步迭代再积分的方法，在每一步积分前先进行迭代，修正折射方向估计偏差，导致运算效率较传统算法更低。表1为使用不同步长、不同修正算法对某次航天任务的雷达全程测量数据进行折射修正的效率和精度比较。

5 结束语

为了提高雷达外测数据的处理精度，本文在对大气层分布结构细致分析的基础上，充分利用了任务时球探空获取的高空大气数据，提出了大气椭球面分层模型，并利用迭代递推的方法推导了电波折射修正算法。实验数据分析表明，本文算法修正精度较传统修正算法有了较大程度的提高，然而由于运算量较大，本文算法的实时性较差，更适合任务结束后对雷达数据的精确分析，目前已初步应用于部分雷达数据事后处理软件。后续工作将进一步对算法进行优化改进，在不降低修正精度的基础上提高运算效率，争取逐步将此算法应用于雷达数据实时处理软件中。

表1 不同修正算法对相同测量数据进行修正的效率和精度比较

[1] 黄家贵, 赵华, 刘元. 电波折射误差实时修正方法研究[J]. 飞行器测控学报, 2011, 30(5): 6-10.

Huang Jia-gui, Zhao Hua, and Liu Yuan. A method for correction of radio wave refraction error in real-time data processing[J]. Journal of Spacecraft TT & C Technology, 2011, 30(5): 6-10.

[2] 王红光, 吴振森, 朱庆林. 大气折射对雷达低仰角跟踪误差的影响分析[J]. 电子与信息学报, 2012, 34(8): 1893-1896.

Wang Hong-guang, Wu Zhen-sen, and Zhu Qing-lin. Influence analysis of atmospheric refraction on low-angle radar tracking errors[J]. Journal of Electronics & Information Technology, 2012, 34(8): 1893-1896.

[3] 武征, 潘佳梁, 胡梦中. 测量雷达对流层大气折射误差修正方法研究[J]. 无线电工程, 2014, 44(2): 73-76.

Wu Zheng, Pan Jia-liang, and Hu Meng-zhong. Research on troposphere refraction error correction of measurement radar[J]. Radio Engineering, 2014, 44(2): 73-76.

[4] 韩先平. 光电经纬仪大气折射误差修正方案[J]. 光电技术应用, 2013, 28(3): 11-13.

Han Xian-ping. Error correction methods for atmospheric refraction of electro-optical theodolite[J]. Electro-optic Technology Application, 2013, 28(3): 11-13.

[5] 王东鹤, 王喜军, 何昕, 等. 光学测量中光波大气折射的精密修正[J]. 光电工程, 2013, 40(4): 45-51.

Wang Dong-he, Wang Xi-jun, He Xin, et al.. High-precision correction of atmospheric refraction in optical measurement applications[J]. Opto-electronic Engineering, 2013, 40(4): 45-51.

[6] 邹士迁, 刘鑫. 火控雷达距离测量误差分析及修正方法[J]. 火力与指挥控制, 2013, 38(5): 148-151.

Zou Shi-qian and Liu Xin. Error analysis and correction of range error for fire control radar[J]. Fire Control & Command Control, 2013, 38(5): 148-151.

[7] 张捍卫, 李爱国, 丁安民. 适应于大气折射研究的压高公式[J].导航定位学报, 2013, 1(2): 11-14.

Zhang Han-wei, Li Ai-guo, and Ding An-min. Restudy on the barometric height formula in atmospheric refraction theory[J]. Journal of Navigation and Positioning, 2013, 1(2): 11-14.

[8] 莫剑飞, 张荣涛. 一种对流层电波折射误差修正方法[J]. 火力与指挥控制, 2010, 35(增刊): 164-165.

Mo Jian-fei and Zhang Rong-tao. A method for radio refraction error correction of troposphere[J]. Fire Control & Command Control, 2010, 35(Suppl.): 164-165.

[9] 赵军, 王西京, 张华, 等. 外测数据对流层折射误差修正及精度分析[J]. 飞行器测控学报, 2014, 33(1): 25-29.

Zhao Jun, Wang Xi-jing, Zhang Hua, et al.. Correction of tropospheric refraction error of tracking data and error analysis[J]. Journal of Spacecraft TT & C Technology, 2014, 33(1): 25-29.

[10] 张瑜, 闫雷兵. 低层大气电波折射修正公式误差研究[J]. 火力与指挥控制, 2007, 32(9): 119-121.

Zhang Yu and Yan Lei-bing. The formula error research of the radiowave rafractive correction in low-layer atmosphere[J]. Fire Control and Command Control, 2007, 32(9): 119-121.

[11] 赵永波, 张守宏. 雷达低角跟踪环境下的最大似然波达方向估计方法[J]. 电子学报, 2004, 32(9): 1520-1523.

Zhao Y B and Zhang S H. Maximum likelihood DOA estimation in radar low-angle tracking environment[J]. Acta Electronica Sinica, 2004, 32(9): 1520-1523.

[12] Boman K and Stoica P. Low angle estimation: models, methods, and bounds[J]. Digital Signal Processing, 2001, 11(1): 35-79.

[13] 杨志强, 陈详明, 赵振维. 对流层电波折射误差修正经验模型研究[J]. 电波科学学报, 2008, 23(3): 580-584.

Yang Zhi-qiang, Chen Xiang-ming, and Zhao Zhen-wei. Empirical model for radio wave refractive error correction of troposphere[J]. Chinese Journal of Radio Science, 2008, 23(3): 580-584.

[14] 韩燕, 强希文, 冯建伟, 等. 大气折射率高度分布模式及其应用[J]. 红外与激光工程, 2009, 38(2): 267-271.

Han Yan, Qiang Xi-wen, Feng Jian-wei, et al.. Height distribution profiles and its application of atmosphere refractive index[J]. Infrared and Laser Engineering, 2009, 38(2): 267-271.

[15] 张瑜, 胡笑君. 雷达方位角折射误差修正方法研究[J]. 现代雷达, 2008, 30(12): 18-23.

Zhang Yu and Hu Xiao-jun. Study on the method of correction for radar azimuth angle refraction error[J]. Modern Radar, 2008, 30(12): 18-23.

唐学海： 男，1984年生，博士，工程师，从事雷达外测数据处理、雷达成像方面的研究.

刘伯阳： 男，1990年生，学士，工程师，从事雷达数据处理方面的研究.

Radio Wave Refraction Correction Algorithm Based on Ellipsoid Stratification Atmospheric Structure Model

Tang Xue-hai①Liu Bo-yang②

①(PLA Unit 63610, Korla 841001, China)

②(PLA Unit 63612, Dunhuang 736200, China)

Traditional radio wave refraction correction algorithms commonly use spherical stratification atmospheric structure model, such kind of model has better results for the measurement data of high elevation angle target. However, it is not good enough for the long-range and low elevation one. This paper proposes a refraction correction algorithm using an ellipsoid stratification model to describe the atmospheric structure. Furthermore, the corrected results are computed by iteration and recursion method. Compared with the traditional refraction correction algorithms, the computational cost of the new model increases, but its low elevation and long-range target refraction correction accuracy is improved, and it can be used in post-flight data processing.

Radio wave propagation; Refraction correction; Ellipsoid stratification atmospheric structure model; Iteration and recursion

TN011

： A

：1009-5896(2015)06-1507-06

10.11999/JEIT141195

2014-09-12收到，2015-01-14改回

国家863计划创新基金(2011AA219205)资助课题

*通信作者：唐学海 txh0203@163.com