国内外对于准确获取运动物体与地形间相对位置这一问题，主要集中在导弹、固定翼飞机、潜艇等领域。文献[3]采用地形轮廓匹配算法准确获取导弹与地形的相对位置，但该算法不允许飞行器进行机动飞行；文献[4]采用桑迪亚惯性地形辅助导航算法修正战斗机与地形的相对位置误差，但该算法在地形梯度变化剧烈区域易造成滤波发散。对于直升机而言，经常需要在复杂的地形上空进行机动低空飞行，因此上述2种方法都不适用于直升机低空飞行告警。而文献[5]采用的互功率谱算法，可准确获取潜艇与地形的相对位置，该算法作为一种性能优良的图像匹配算法，不仅对地形适应性强，而且在匹配时允许潜艇作机动航行[6]。
由于直升机对地形匹配的技术需求与潜艇具有较大的相似性，因此，本文基于互功率谱算法，开展基于地形匹配的直升机低空飞行前视告警方法研究。首先，建立直升机逃逸轨迹模型并确定直升机前视告警边界；然后，基于互功率谱算法进行地形匹配，以准确获取直升机与地形的相对位置，并利用地形熵选择匹配区域提高匹配精度；最后，构建测试算例对本文提出的方法进行仿真验证。
1 直升机低空飞行前视告警方法 直升机前视告警方法的核心是告警边界，目前告警边界的设计方法主要有美国Honeywell公司的告警边界设计方法[7]和基于逃逸轨迹的告警边界设计方法[8]。其中，Honeywell的告警边界为考虑通用性，其边界设计趋于保守，低空飞行时虚警率较高，不适用于直升机低空飞行。基于逃逸轨迹的告警边界设计方法则基于直升机的飞行性能，能够在保证安全的基础上最大程度的缩小告警边界以减小虚警率，因此本文采用该方法。
基于逃逸轨迹的告警边界设计方法，首先需要确定直升机的逃逸轨迹，然后在此基础上设计告警边界。
1.1 直升机逃逸轨迹 直升机的逃逸过程通常可分为3个阶段：飞行员反应、拉起改出和稳定爬升。其中，飞行员反应阶段保持直升机当前状态飞行，稳定爬升阶段以拉起改出阶段最终时刻的状态为基准进行匀速飞行，这两部分均为直线运动轨迹；而拉起改出为机动飞行状态，其轨迹为复杂的空间曲线，是逃逸轨迹中最为关键的部分，直接决定了逃逸轨迹的整体形态。图 1示出了直升机逃逸轨迹的3个主要组成部分。
 |
| 图 1 直升机逃逸轨迹各部分组成示意图 Fig. 1 Schematic diagram of each part of helicopter escape trajectory |
| 图选项 |
由文献[9]可知，直升机逃逸轨迹的拉起改出部分与抛物线方程拟合度较高，该部分逃逸轨迹高度H可表示为
 | (1) |
式中：x为轨迹的横向坐标；c0、c1和c2为待定拟合系数。拟合系数的大小与直升机的飞行性能相关，通常需要基于飞行试验数据拟合求解，亦可利用高精度直升机飞行动力学模型获得，本文采用第2种方法。
在实际应用时，对于给定机型，首先离线计算各种不同组合飞行状态下逃逸轨迹拉起改出部分的拟合系数，在实际飞行过程中，通过多维插值即可获得该飞行状态下的系数值，从而得到相应逃逸轨迹的拉起改出部分。
1.2 直升机前视告警边界 参考文献[10]，将逃逸轨迹中飞行员反应阶段的时间设定为2s。此外，由于传感器误差以及地形数据库误差等因素的存在，直升机在飞行过程中应在其下方留有一段安全距离，即最小安全高度。在现有直升机近地告警产品中，如美国Honeywell公司的MXKKII，将最小安全高度的最大值设为500feet(对应高度152.4m)，对于直升机低空飞行而言，显然不能直接采用。因此，本文根据传感器误差、地形数据库误差以及地形平均斜率等因素，将最小安全高度设为50m。以UH-60直升机为例，设飞行初始状态为：前飞速度为40m/s，下降率为4m/s，则得到告警边界如图 2所示。图中：x、z为轨迹坐标。
 |
| 图 2 UH-60直升机逃逸轨迹及告警边界示意图 Fig. 2 Schematic diagram of escape trajectory and alert envelope for UH-60 helicopter |
| 图选项 |
2 直升机/地形相对位置的准确获取 在确定了直升机前视告警边界后，须进一步准确获取直升机与地形的相对位置，才能进行有效告警。此外，为避免匹配区域选择不当导致的匹配失败，在进行地形匹配之前，利用地形熵选择匹配区域，以提高匹配的准确性。
2.1 直升机低空飞行地形匹配方法 互功率谱算法的基本原理是利用快速傅里叶变换将2个相关矩阵由空域变换到频域，再根据傅里叶变换的位移性质，得到一个反映两者相对误差的脉冲函数[11-12]。将该方法应用于直升机与地形匹配时，须将直升机雷达探测得到的地形高程矩阵和机载地形数据库的地形高程矩阵进行快速傅里叶变换，得到一个互功率谱，由归一化的互功率谱计算出一个脉冲函数，由该函数即可得到直升机与地形的相对位置误差。
直升机雷达探测的地形高程矩阵和机载地形数据库地形高程矩阵的互功率谱G(u, v)可表示为
 | (2) |
式中：F1(u, v)和F2(u, v)为2个地形高程矩阵对应的频谱，u、v为x、y对应的频域坐标；F1*(u, v)为F1(u, v)的复共轭；Δx、Δy为2个地形高程矩阵间的水平偏差值。对式(2)进行傅里叶逆变换后，可以得到一个水平面上的脉冲函数，表示为
 | (3) |
当直升机与地形的相对位置误差为0时，该脉冲函数会产生一个明显的冲激峰，其值为1，而其他区域则为0，如图 3(a)所示。当相对位置误差较小时，2个地形高程矩阵存在的水平偏移量引起的脉冲函数将产生一个明显的冲激峰。相对位置误差越小，2个地形高程矩阵重叠区域越大，则冲激峰值越接近1，其他区域的值越接近于0，如图 3(b)所示。此时冲激峰的位置坐标即为直升机与地形的相对位置误差值。而当相对位置误差很大时，该函数的能量将从单一峰值分散到许多小峰值，如图 3(c)所示。图中：X、Y为脉冲函数坐标。
 |
| 图 3 不同相对位置误差下的脉冲函数 Fig. 3 Pulse functions with different relative position errors |
| 图选项 |
2.2 不同地形对匹配方法的影响分析 此处通过仿真测试，考察不同地形起伏程度对互功率谱算法的有效性影响。首先利用马尔可夫算法随机生成一个起伏程度较大的地形[13]，大小为220×220矩阵网格，地形高程均值为487m，极差为55m，方差为81m2，如图 4所示。
 |
| 图 4 起伏较大的随机地形 Fig. 4 Undulating random terrain |
| 图选项 |
在该地形中取出2个200×200的矩阵作为待匹配地形，两矩阵X、Y方向的平移偏差值均为20。对2个待匹配地形矩阵进行快速傅里叶变换，再计算得到脉冲函数，如图 5所示。该脉冲函数中存在一个明显的冲激峰，其坐标为(20，20)，峰值为0.44，由此可得到2个待匹配地形矩阵X、Y方向的平移偏差值均为20。
 |
| 图 5 匹配成功时的脉冲函数 Fig. 5 Pulse function on successful matching |
| 图选项 |
进一步生成一个起伏程度较小的随机地形，地形高程均值为497m，极差为5m，方差为1.6m2，如图 6所示。
 |
| 图 6 起伏较小的随机地形 Fig. 6 Small undulating random terrain |
| 图选项 |
同样在该地形中取出2个200×200的矩阵作为待匹配地形，两矩阵X、Y方向的平移偏差值均为20。计算得到脉冲函数如图 7所示。可见，该脉冲函数中不存在明显的冲激峰，此时匹配失败。
 |
| 图 7 匹配失败时的脉冲函数 Fig. 7 Pulse function on failed matching |
| 图选项 |
上述仿真测试表明，互功率谱算法的匹配效果受地形起伏程度影响很大，地形起伏程度越大，地形越独特，匹配准确度越高，因此山地等地形较适合地形匹配，而平原、沙漠等地形易导致匹配失败。
2.3 地形熵选择匹配区域 针对平原、沙漠等包含信息量较少的区域，本文参考文献[14]，采用地形熵选择匹配区域方法来解决此问题。
地形熵能够用来描述地形起伏的剧烈程度，当某区域的局部地形高度值变化较为剧烈时，地形所包含的信息量比较丰富，此时计算得到的地形熵较小，而越小的地形熵越有利于进行匹配运算。地形熵的定义式为
 | (4) |
式中：p为无量纲的地形高度值；h为地形高度值；i为地形采样点的序号; M为地形采样点个数的总量；Hterrain为地形熵。
在进行地形匹配前，首先计算直升机当前位置附近区域的地形熵，若地形熵小于预先设定阈值，则进行地形匹配。这样，避免了在平原等区域进行匹配导致的匹配失败，有效提高了地形匹配的准确性。
以图 4地形为例，可取出共441个200×200的高程矩阵，参考文献[15]将地形熵阈值设为9，计算得到有25个高程矩阵的地形熵大于阈值，而其中有19个匹配失败。将该25个高程矩阵筛除后，剩余416个矩阵均能够匹配成功。
3 算法仿真测试 本文利用MATLAB开展对基于地形匹配的直升机低空飞行前视告警方法的仿真测试。首先给出算法性能的评价指标，然后给出仿真测试的参数设置，最后给出仿真结果及分析。
3.1 算法性能评价指标 根据文献[16]，此处采用虚警率和告警成功率作为评价告警方法的性能指标。
虚警是指告警系统在直升机安全飞行时给出了不必要的告警，虚警率定义为
 | (5) |
式中：Nu为虚警算例个数；N为算例总数。
告警成功率则是相对告警失败而言的。告警失败是指告警提示过晚，飞行员得到告警信号时已经无法及时改出而导致撞地事故。告警成功率定义为
 | (6) |
式中：Nf告警失败算例个数。
表 1给出了前视告警方法的性能指标判断标准。
表 1 前视告警方法性能指标判断标准 Table 1 Judgement criteria for performance index of forward looking alert method
| 告警边界与地形是否相交 | 正常轨迹与地形是否相交 | 逃逸轨迹与地形是否相交 | 告警情况 |
| 是 | 是 | 否 | 告警成功 |
| 是 | 是 | 是 | 告警失败 |
| 是 | 否 | 否 | 虚警 |
表选项
需要说明的是，计算虚警率和告警成功率的算例总数仅包括产生告警的算例。
3.2 仿真参数设置 为对比验证本文提出的基于地形匹配的直升机低空飞行前视告警方法的优势，仿真中设置了2个测试组，一个为采用本文方法的匹配组，另一个为采用常规告警方法的对照组，2个测试组的告警边界相同，均采用基于逃逸轨迹的告警边界，最小安全高度设为50m，但对照组未结合地形匹配算法，无法修正直升机与地形的相对位置误差。
1) 地形参数设置：类似图 4，仿真地形采用马尔可夫算法随机生成的500组地形，地形高程均值约为400m，地形矩阵大小均为220×220，设每个矩阵网格大小为90m，总大小即为19800m×19800m，2个测试组均在500组地形进行仿真飞行测试。
2) 直升机飞行参数设置：仿真基于UH-60直升机飞行动力学模型。直升机由地形矩阵坐标(20，100)处出发，在500m海拔高度以40m/s的速度平飞，无线电高度均值约为100m，终点为矩阵坐标(160，100)处，每0.5s主动探测直升机周围10×10矩阵网格的地形数据，并计算地形熵，当地形熵小于9时，则进行地形匹配。飞行期间发生告警或到达终点则停止该次仿真测试。
3) 误差参数设置：测试开始时，2个测试组直升机与地形的位置误差均为3个矩阵网格，即270m，匹配组与对照组中各有一半测试组的误差方向沿直升机飞行方向，即直升机导航位置位于实际位置前方，另一半误差方向相反，即直升机导航位置位于实际位置后方。
3.3 仿真结果与分析 由3.2节可知，在仿真测试中，对照组存在导航位置误差，告警边界以导航位置为准，而正常飞行轨迹和逃逸轨迹则以实际位置为准。取仿真测试中的典型算例为例，如图 8和图 9所示。图 8(a)、图 8(b)为测试中同一地形下的2个测试组，其中匹配组告警成功，而对照组产生虚警。图 8(a)中，告警边界、正常飞行轨迹与地形相交，而逃逸轨迹与地形不相交，此时告警成功。图 8(b)中导航位置在实际位置前方，告警边界与地形相交，而正常飞行轨迹、逃逸轨迹与地形不相交，此时为虚警。图 8(c)为该地形下，2个测试组直升机与地形相对位置误差随时间变化的曲线。从图中可以看出，仿真开始时2组的相对位置误差均为3个矩阵网格，但匹配组通过地形匹配修正了这一误差。
 |
| 图 8 对照测试A Fig. 8 Controlled trial A |
| 图选项 |
 |
| 图 9 对照测试B Fig. 9 Controlled trial B |
| 图选项 |
图 9(a)、图 9(b)为仿真测试中同一地形下的另外2个测试组，其中匹配组告警成功而对照组告警失败。图 9(a)所述状态与图 8(a)相似，为告警成功。而图 9(b)中导航位置在实际位置后方，告警边界、正常飞行轨迹、逃逸轨迹均与地形相交，此时告警失败。图 9(c)为该地形下，2个测试组直升机与地形相对位置误差随时间变化的曲线。同样可以从图中看出，匹配组通过地形匹配修正了相对位置误差。
完成对照组、匹配组各500组仿真测试算例后，得到的测试结果如表 2所示。
表 2 告警方法仿真测试结果 Table 2 Simulation test results of alert method
| 性能指标 | 匹配组 | 对照组 |
| 告警成功率/% | 95.47 | 65.93 |
| 虚警率/% | 3.77 | 19.34 |
表选项
从表 2的仿真结果可以看出，对于采用地形匹配的本文告警方法，相比于常规告警方法，虚警率明显降低(约16%)，告警成功率显著提高(近30%)。
4 结论 本文将基于互功率谱算法的图像匹配方法应用于直升机地形匹配，有效地修正了直升机与地形的相对位置误差，并据此提出了直升机低空飞行前视告警方法，有效实现了直升机低空飞行防撞告警。主要结论如下：
1) 互功率谱算法的匹配效果受地形起伏程度影响，地形起伏程度越大，地形越独特，匹配准确度越高。
2) 对于直升机低空飞行状态，未结合地形匹配的常规告警方法由于直升机与地形的相对位置不准确，导致告警成功率低、虚警率高。结合地形匹配后的前视告警方法告警成功率显著提高，虚警率明显降低，证明了本文提出的前视告警方法能够很好地适应直升机低空飞行的防撞告警要求。
参考文献
| [1] | 刘仪. 防范可控飞行撞地事故的最新进展[J]. 国际航空, 2001(12): 41-42. LIU Y. The new development of preventing CFIT[J]. International Aviation, 2001(12): 41-42. (in Chinese) |
| [2] | 杨富民.直升机地形防撞告警理论与技术研究[D].郑州: 解放军信息工程大学, 2012: 1-3. YANG F M.Research on theory and technology of helicopter terrain collision warning[D]. Zhengzhou: Information Technology University of PLA, 2012: 1-3(in Chinese). http://cdmd.cnki.com.cn/Article/CDMD-90005-1013161203.htm |
| [3] | GOLDEN J P.Terrain contour matching(TERCOM): A cruise missile guidance aid[C]//Proceedings of 24th SPIE Annual Technical Symposium.Bellingham: SPIE, 1980, 238: 10-18. |
| [4] | BAIRD C A.Terrain-aided altitude computations on the AFTI/F-16[C]//Proceedings of IEEE Position Location & Navigation Symposium.Piscataway, NJ: IEEE Press, 1990: 474-481. |
| [5] | 徐遵义, 黄东武, 赵洪銮. 水下地形辅助导航新方法仿真[J]. 舰船科学技术, 2012, 34(6): 53-56. XU Z Y, HUANG D W, ZHAO H L. A new method for underwater terrain-aided navigation[J]. Ship Science and Technology, 2012, 34(6): 53-56. DOI:10.3404/j.issn.1672-7649.2012.06.012 (in Chinese) |
| [6] | 曹闻, 李弼程, 邓子建. 一种基于小波变换的图像配准方法[J]. 测绘通报, 2004(2): 16-19. CAO W, LI B C, DENG Z J. An image registration method based on wavelet transformation[J]. Bulletin of Surveying and Mapping, 2004(2): 16-19. DOI:10.3969/j.issn.0494-0911.2004.02.006 (in Chinese) |
| [7] | ISHIHARA Y, CONNER K.Method and computer program product for helicopter ground proximity warning system: US6583733B2[P]. 2003-12-03. |
| [8] | 周成中.基于逃逸轨迹告警原理的直升机地形感知与告警算法研究[D].南京: 南京航空航天大学, 2015: 7-9. ZHOU C Z.Study on helicopter terrain awareness and alert algorithm based on escape trajectory alert[D]. Nanjing: Nanjing University of Aeronautics and Astronautics, 2015: 7-9(in Chinese). |
| [9] | ANDERSON T, JONES W, BEAMON K.Design and implementation of TAWS for rotary wing aircraft[C]//Proceedings of Aerospace Conference Aerospace Conference.Piscataway, NJ: IEEE Press, 2011: 1-7. |
| [10] | WILLIAMS S. Towards terrain-aided navigation for underwater robotics[J]. Advanced Robotics, 2001, 15(5): 533-549. DOI:10.1163/156855301317033559 |
| [11] | 郑志彬, 叶中付. 基于相位相关的图像配准算法[J]. 数据采集与处理, 2006, 21(4): 444-449. ZHENG Z B, YE Z F. Image registration algorithm based on phase-correlation[J]. Journal of Data Acquisition & Processing, 2006, 21(4): 444-449. DOI:10.3969/j.issn.1004-9037.2006.04.014 (in Chinese) |
| [12] | 朱佳媛. 基于FFT相位相关的图像匹配算法与实现[J]. 科技风, 2013(14): 11. ZHU J Y. Image matching algorithm and implementation based on phase correlation of FFT[J]. Technology Wind, 2013(14): 11. DOI:10.3969/j.issn.1671-7341.2013.14.007 (in Chinese) |
| [13] | 张继贤, 柳健, 李德仁. 地形生成技术与方法的研究[J]. 中国图象图形学报, 1997, 2(9): 638-645. ZHANG J X, LIU J, LI D R. Technical approach to realistic terrain generation[J]. Journal of Image and Graphics, 1997, 2(9): 638-645. (in Chinese) |
| [14] | 王英钧. 地形辅助导航综述[J]. 航空电子技术, 1997, 5(1): 24-30. WANG Y J. A review of terrain-assisted navigation[J]. Avionics Technology, 1997, 5(1): 24-30. (in Chinese) |
| [15] | 席雷平, 陈自力, 李小民. 基于地形熵的地形适配区选择准则研究[J]. 弹箭与制导学报, 2010, 30(4): 247-249. XI L P, CHEN Z L, LI X M. Discussion on the navigable terrain selection algorithm based on terrain entropy[J]. Journal of Projectiles Rockets Missiles and Guidance, 2010, 30(4): 247-249. DOI:10.3969/j.issn.1673-9728.2010.04.071 (in Chinese) |
| [16] | FITZGERALD T R, BRUNNER M T. Use of high-fidelity simulation in the development of an F/A-18 active ground collision avoidance system[J]. AIAA Journal, 1995, 33(6): 68-78. |
