引言

人工增雨飞机在携带催化作业装备实施增雨作业的同时, 能搭载云降水探测仪器进行大气环境及云物理参数探测, 具有目标明晰, 机动性强的特点, 故在进行人工增雨研究试验时, 具有条件的地方多采用飞机进行人工增雨探测作业, 在空中追逐具有增雨潜力的云, 探测其状况, 根据云的具体条件, 按照增雨目的要求, 实施相应的作业手段, 并检测作业效果, 同时收集云降水的物理状况数据, 供研究分析使用^{[1, 2]}。

目前, 北京地区区域可供水量已难以满足工农业生产与城市生活用水需要, 北京市政府在积极开源节流, 引水济京; 并在以往高炮、火箭和气球作业的基础上, 积极开展飞机人工增雨作业, 以最大能力来缓解北京地区淡水资源的紧缺^{[2, 3]}。

飞机人工增雨作业中, 飞机探测作业航迹是制定和分析云降水环境状况, 讨论飞机作业效果的重要分析数据。在截获实时显示的飞机航迹图之外, 由于还需结合其它飞行探测资料, 如PMS (Particle Measurement System)粒子探测系统、温度、含水量, 以及相对应的气象雷达、气象卫星等多种观测资料, 进行更为细致的飞行探测作业资料的处理和分析, 来进一步明确云降水过程中的物理机理、制定更为科学的探测作业方案和进行更具说服力的飞机人工增雨效果评估^[4-6]。

GPS全球定位系统(Global Position System)是美国1970年代开始研制的导航定位系统, 随着其技术的发展, 其应用范围已不局限于导航定位^[7]。在气象领域, 我国利用GPS的空间定位能力, 测量高空风并发展为GPS探空设备^{[8, 9]};利用地基GPS接收机测量整层大气的水汽含量^[10];利用低轨卫星上的GPS接收机探测大气折射率廓线, 并进一步计算温度和湿度的垂直分布^[11]。这些数据可为短期天气预报、水循环、业务数值天气预报等提供准确的气象资料。

美国自1989年将GPS技术应用到WP3型环境飞机探测中^[12]。我国1996年开始也在飞机人工增雨探测作业中使用这项技术, 通过10年的应用发展, 目前我国开展的飞机人工增雨作业都在使用这项技术, 同时还对GPS数据的拓展使用进行了积极尝试^{[13, 14]}。

在以往的飞机增雨探测作业中, 对GPS飞行状态及飞机运动轨迹记录资料的分析, 使用微软的Excel进行处理。使用中发现微软的Excel虽具有良好的数据整合和图形显示功能, 但由于航迹数据为1s或5s一个记录, 对于一个架次, 飞行近2小时, 有近万行记录数据的飞机探测作业轨迹记录, 用手工方法或Excel软件判断数据断点和重复段, 工作量较大。而且在Excel三维图表中, 每个数据列最多允许有4000个数据点, 如数据点超过这一数量, 则需创建2个或更多的数据列, 才能作图。

MATLAB是矩阵实验室(Matrix Laboratory)之意, MATLAB语言是当今国际上科学界(尤其是自动控制领域)最具影响力、也是最有活力的软件, 它起源于矩阵运算, 并已经发展成一种高度集成的计算机语言。MATLAB提供了强大的科学运算、灵活的程序设计流程、高质量的图形可视化与界面设计、便捷的与其他程序和语言接口的功能。除具备卓越的数值计算能力外, 它还提供了专业水平的符号计算, 文字处理, 可视化建模仿真和实时控制等功能。MATLAB的基本数据单位是矩阵, 它的指令表达式与数学、工程中常用的形式十分相似, 故用MATLAB来解算问题要比用C、FORTRAN语言等完成相同的事情要简捷得多。MATLAB语言目前已在各国高校与研究单位起着重大的作用^{[15, 16]}。

本文利用Matlab, 分析整理了2004年北京增蓄型飞机人工增雨探测作业试验中夏延-ⅢA飞机21架次飞行航迹, 从中选择问题较多的2004年8月12日10:09 —12: 56近2个小时的飞机航迹数据, 利用Matlab讨论GPS航迹资料中存在的问题, 提出解决这些问题的方法。

1 新型空地信息传输系统的研制及其适应性研究

在“九五”科技攻关“农业气象灾害防御技术研究———人工增雨农业减灾技术研究”课题执行期间, 为实现飞机空中探测信息与地面气象信息的有效结合, 让地面丰富的气象业务信息、指挥信息、地理信息在空中飞行作业中发挥更好的作用, 更好地把握作业时机, 建立集空中探测、定位等信息与气象业务信息于一体的综合显示平台, 中国气象科学研究院人工影响天气研究所联合北京凌亨科技有限责任公司, 在以往飞机人工增雨空地传输系统的基础上, 利用GPS飞机定位测速系统, 研制和发展了我国新型飞机人工增雨信息空地传输显示系统^[4], 并自2003年, 投入到“十五”科技攻关项目“人工增雨技术研究及示范”的河北、河南、北京和青海的综合外场探测试验中。

由于空地信息传输系统在地面仅能进行车载漫游测试, 因此需要在实际增雨作业飞行时(增雨飞机航速在60 ~ 120m s^-1)进行充分的适应性分析研究来提高整套系统的性能。基于北京2004年飞行航迹资料, 结合一段实际飞行数据, 用Matlab对飞行数据品质进行分析。

2 数据分析

执行定检后的夏延-ⅢA云降水探测作业飞机, 在2004年4月29日, 开始了在北京地区的增蓄型飞机人工增雨综合外场探测试验的飞行。利用机载PMS粒子测量系统和AgI末端燃烧器, 结合车载多普勒雷达监测、地面降水自动站观测、卫星遥感监测分析, 对15次降水天气过程进行了21架次的探测作业飞行。于2004年8月15日下午16时完成本年度探测作业飞行试验任务。

探测作业试验飞行期间共收集处理了20架次的飞行探测作业航迹资料, 其中一个架次资料因计算机故障而丢失。通过资料分析整理, 归纳了这套空地信息传输系统GPS资料记录中存在的问题, 并选择其中存在问题较多的2004年8月12日第二架次探测飞行的GPS航迹数据, 进行示例讨论。

2004年8月12日夏延-ⅢA云降水探测作业飞机在北京沙河机场进行了2架次的降水天气系统的探测作业飞行, 其中第二架次飞行探测作业, GPS开始记录时间10:3:36, 飞机离地时间10:9:28, 飞机落地时间12:51:50, 数据记录间隔1s。GPS数据格式为常见的文本格式, 可用Excel打开, 然后在Matlab中利用数据的导入(import)功能得到一个二维的数组。数据每列分别表示时间(北京时), 飞机所处的经(东经)纬(北纬)度, 飞机飞行相对地面海拔高度, 飞机飞行相对地面移动速度, 飞机航向等。表 1给出了这次探测作业起飞10分钟后飞机处于爬升阶段的一个GPS实际数据示例。

利用Matlab分析, GPS航迹数据存在问题可以归纳为7种:1)由于GPS接受天线或系统记录错误, 造成一段时间内一些数据都为0;2)飞行探测作业期间, 在GPS海拔高度分析中, 出现实际飞行中不可能的海拔高度; 3)因笔记本电脑死机造成的航迹记录数据丢失和重复, 部分时间段没有数据记录, 同时存在数据重复出现, 记录数据在一段时间段上不更新; 4)在5月16日以前的飞行中, GPS记录间隔为1 ~ 5s不等, 分析发现与1s一个记录的PM S探测数据记录不一致, 后通过程序改进将其设为1s一个记录, 因此在调试期间的数据分析中, 有些时间段数据记录间隔仍达5s;5)飞机相对地面移动速度单位为km h^-1, 而不是比较直观的m s^-1; 6)除时间记录外, 经纬度、海拔高度、飞机相对地面移动速度及飞机航向均存在小数点位置后有过多尾数的情况; 7)经纬度数据存在60°与百进制的换算问题, 且小数点位置应进2位。

3 处理方法

问题1):通过简单判识后, 直接剔除。

问题2):因为我国目前用于云降水探测飞机(除运-八外)的实际飞行高度均在1万米以下, 这样只要判断飞行高度高于海拔10km以上, 就标识这段数据为不良记录剔出。在对2004年21个架次飞行资料分析中, 只有7月29日11:15:22、8月12日的10:33:59 —10:34:33和12:12: 22 —12:13:14段出现GPS海拔高度记录分别为31524m, 23410m, 和41225m, 这3个数值均与周围高度数据明显不一致, 通过航迹高度分析, 手工剔除即可。由于目前还不清楚造成GPS海拔高度数据异常的原因, 故在此仅只能提醒注意这一问题的出现, 尚不能对这一现象作更多的解析和说明。

对于问题3)和4), 考虑丢失和重复的数据没有统一的格式, 很难直接用程序剔出, 因此利用Matlab把原有时间序列两个相邻记录的差建立一个新的时间序列{a}, 因时间序列都是实数格式的, 当时间序列数据连续变化时, a就是一个不变或小于一定值的值。如GPS每秒一个记录, 则a为1, 如5秒一个记录, 则a为5。这样GPS数据如果重复或超过5s以上间断, 就会造成a的突变, 选取一个合适的临界数值A (1或5), 查找a值大于A的数值位置, 也就能得到记录中时间序列间断的数据段。同理如a等于0, 则对应数据记录重复。采用这种方法, 从这次近2小时飞行数据9870行记录中得到的99个记录断点。

分析这99个记录断点, 对于1s一个记录的航迹数据, 在这次探测作业飞行中, 出现数据断点间隔为5 ~ 66s, 其中超过15s的达51个, 超过40s的达13个。相对于其它探测要素时间分辨率, 可能影响不大, 但对于巡航速度在100m s^-1的夏延-Ⅲ A飞机来讲, 15秒的数据丢失则是在1. 5km的距离上没有探测数据, 这就有可能丢失一个穿越云带的数据, 40s时间飞机则穿过了4km距离, 如果丢失这样尺度的云中航迹数据, 将对云降水物理资料的分析品质造成极大的影响。

对于近10000行数据记录来讲, 这些断点不明显, 因此特别提醒在航迹分析时注意这个问题, 从而检验航迹资料中是否存在这种情况, 可以尽快追踪出导致这种问题出现的原因, 进而从技术上尽早予以解决。

对于问题5), 速度转换只需对这列数值, 除以3.6即可。同样, 对问题6)只需在数据转存时减少其小数点后的尾数位。

对于问题7), 经纬度的换算就需引入转换公式。选取一点(本例中用(E115, N40))为参考点, 以分为单位。a = 11532.194336表示东经115度32.194336分, b = 4052.133789表示北纬40度52.133789分, 利用公式(1)、(2)换算:

(1)

(2)

其中a, b分别表示GPS数据的经度和纬度, Long、Lat分别表示变换后的经纬度, 对于上组经纬度数据, 可以得到(50.960938, 17.895996), 这样就可用转换好的经纬度数据进行航迹分析。

图 1给出了由Matlab对2004年8月12日10:09 —12:56的GPS航迹数据换算处理后, 由Matlab直接生成的飞机立体航迹图。经过处理生成的航迹图, 不再存在断点、异常值和轨迹不连续情况。

4 结语

(1) 我国新型飞机人工增雨空地信息传输系统, 在地面多项测试后, 仍需通过高速移动的飞行平台进行适应性研究分析试验, 来完善和提高其实现数据信息传输和资料存储的功能。系统出现的问题, 主要是所采用的GPS系统处理软件设计还有待完善, 通过航迹数据处理方法的改进, 将可得到质量更好的数据资料。

(2) 通过Matlab对这些航迹数据进行换算和处理, 可直接用图形来显示处理结果, 非常简便。对北京2004年20个探测作业架次的航迹图分析, 所得的航迹图与每个架次的宏观飞行记录对比, 可见分析结果与实际飞行情况完全吻合。

(3) 一些引进设备中的数据格式, 往往需结合具体业务需求进行转换和分析, Matlab是完成这项工作的很好工具。另外通过发展MATLAB用于飞机人工增雨探测作业和大气环境探测的飞行轨迹分析的工具箱, 可获得更为丰富、更为明晰的结果。

(4) 目前对GPS资料的初步处理, 仅仅是针对航迹的分析, 随着工作的进一步深入, 结合其它飞行资料, 可以对飞行探测作业航向、环境风以及飞行速度进行更为系统的研究应用, 利用Matlab的图示功能, 得到4维的航迹图像和更丰富的图形产品, 包括任意一点的飞行速度或者航向。