摘要： 文章针对风廓线雷达实际工程应用中，处理风场数据时的非线性问题，在传统卡尔曼滤波技术的基础上，提出了基于扩展卡尔曼滤波的数据处理方法。选取2011年的风场数据进行计算机仿真分析，结果表明，该方法可以有效去除风场数据中掺杂的噪声干扰，较好的平滑了未滤波曲线，并弥补了经典卡尔曼滤波参数设置方面的不足，滤波效果优于传统的卡尔曼滤波，具有一定的工程应用前景。
Abstract： A data processing method based on Extended Kalman Filter （EKF） is presented in this paper， which focuses on the nonlinear problems of wind data of wind profiler radar. This method is the development of traditional Kalman Filter （KF） in practical engineering applications. To verify the validity of this first—order EKF from different angles， wind data in 2011 is selected for simulation filtering. After processing the measured wind data， the simulated experimental results indicate that this method can effectively remove the noise of wind data. It plays the role of smoothing， and solve the problem of boundary jump， making up for the classical KF parameter settings of the insufficient. EKF is better for data processing of wind profiler radar and have some engineering application prospects.
关键词： 风廓线雷达；数据处理；卡尔曼滤波；一阶扩展卡尔曼滤波
Key words： wind profiler radar；data processing；Kalman Filter；First—order Extended Kalman Filter
中图分类号：[P415.2]文献标识码：A文章编号：1006—4311（2012）28—0014—03
0 引言
风廓线雷达（wind profiler radar，WPR）是一种新型的多普勒测风雷达，以晴空大气作为探测对象，利用大气湍流对电磁波的散射作用对大气风场等物理量进行遥感探测[1]。风廓线雷达在探测的过程中，往往受到各种非气象因素的干扰，从而影响探测数据的精度，导致雷达数据的质量问题，因此对数据进行处理是十分必要的。目前比较成熟的风廓线雷达数据处理方法主要有一致性平均法、最优插值法、WW算法、中值滤波法和滑动平均滤波法等[2]，上述方法对于风廓线雷达的数据处理效果明显，为下一步研究工作提供了可靠的数据支持。但也必须注意到，这些方法存在一定的局限性，如对于风切变的判断，每个方法都有其使用的前提，破坏了这些前提，以上方法就不能够使用了。此外，杨馨蕊[3]尝试将传统的线性卡尔曼滤波（KF）方法应用到风廓线雷达数据处理当中，并有效的去除了一部分高频干扰，达到了一定的滤波效果，为风廓线雷达数据处理开拓了新的思路。
但是必须注意到，在风廓线雷达测风活动中，风向和风速随高度的分布是一个随机过程，在每个高度上的取值都是随机的，因此，虽然各时间点和各高度上的取值之间会有波及，但是很难保证都是线性关系。因此，在实际应用中，需要处理的问题常是非线性的。本文将传统卡尔曼滤波方法进行推广，针对实际工程中风场数据的非线性特性，提出了基于非线性的扩展卡尔曼滤波（EKF）技术对风廓线雷达数据质量进行控制。
1 风廓线雷达EKF基本原理
目前，卡尔曼滤波（KF）技术已应用于气象预报之中，对于风廓线雷达数据处理方面也有了初步的研究，并取得了一定的滤波效果。Kalman滤波是卡尔曼（R. E. Kalman）于1960年提出的，从与被提取信号有关的观测量中通过算法估计出所需信号的一种滤波算法[4]，它不需要保留用过的观测数据，当测得新的观测数据后，可按照一套递推公式（公式参见下文扩展卡尔曼滤波部分）算出新的估计量，不必重新计算。此外它还打破了对平稳过程的限制，可用于对时变随机信号的滤波。图1给出了卡尔曼滤波估计的流程图（图中将k作为当前时刻）。
扩展卡尔曼滤波（EKF）是在卡尔曼滤波的基础上针对非线性系统提出的一种改进方式，比较适合处理风廓线雷达数据处理当中风场的非线性问题，它通过泰勒级数展开的线性化处理方式，对非线性方程作线性化处理，再结合经典的卡尔曼滤波进行滤波估计。在估计过程中利用系统的状态方程和观测方程的统计特性形成滤波算法，更多的利用了系统前后的状态信息，实际上是一种最优估计方法。本文以风廓线雷达风速数据的处理为例，分析扩展卡尔曼滤波技术在雷达数据处理当中的应用，针对风速数据的非线性特性，建立风廓线雷达系统的非线性状态方程为（设当前时刻为k+1）：
X（k+1）=f（k，X（k））+V（k）（1）
式（1）中，X（k）为k时刻的系统状态输入量，V（k）表示过程噪声。在对风廓线雷达风速数据进行处理的过程中，假定没有控制输入，并假定过程噪声V（k）是加性零均值白噪声，其方差为Q（k）。
量测方程为：z（k）=h（k，X（k））+W（k）（2）
其中：量测噪声W（k）也假定是加性零均值白噪声，其方差为R（k）。
假定过程噪声和量测噪声序列是彼此不相关的，并具有初始状态估计■（0|0）和协方差矩阵P（0|0）。和线性情况一样，假定k时刻的估计为：■（k|k）≈E[X（k）|Zk]（3）
它是一个近似的条件均值，其相伴协方差矩阵为P（k|k）。为了达到预测的状态■（k+1|k），对式（1）中的非线性函数在■（k|k）附近进行泰勒级数展开，取其一阶或者二阶项，以产生一阶或者二阶EKF，其中一阶泰勒级数展开式为：X（k+1）=f（k，■（k|k））+fX（k）[X（k）—■（k|k）]+（高阶项）+V（k）（4）
其中：fX（k）=[？荦Xf′（k，X）]′■
=■[f■（X）…f■（X）]■
=■■（5）
是向量f的雅克比矩阵，在状态的最近估计上取值，其中x1、x2、…、xn为n维状态向量X（k）的元素，在风廓线雷达风速数据处理当中，分别代表各观测时间点上、每个高度的风速观测值。
如果泰勒级数展开式中保留到二阶或三阶，则可得到二阶和三阶扩展卡尔曼滤波。Phanenf R.J.对不同阶数的扩展卡尔曼滤波性能进行了仿真分析。仿真结果表明，二阶扩展卡尔曼滤波的性能比一阶的要好，而二阶以上的扩展卡尔曼滤波性能与二阶相比并没有明显的提高，所以超过二阶以上的扩展卡尔曼滤波一般都不被采用。二阶扩展卡尔曼滤波的性能虽然要优于一阶的，但二阶的计算量很大，所以一般情况下只采用一阶扩展卡尔曼滤波。本文采用一阶扩展卡尔曼滤波对数据进行处理。
一阶扩展卡尔曼滤波的公式系包括：
状态的一步预测■（k+1|k）=f（k，■（k|k））（6）
协方差的一步预测
P（k+1|k）=fX（k）P（k|k）f■■（k）+Q（k）（7）
增益K（k+1）=P（k+1|k）h■■（k+1）S—1（k+1）（8）
状态更新方程
■（k+1|k+1）=■（k+1|k）+K（k+1）[z（k+1）—■（k+1|k）]（9）
协方差更新方程
P（k+1|k+1）=[I—K（k+1）hX（k+1）]P（k+1|k）[I+K（k+1）hX（k+1）]′—K（k+1）R[k+1]K′（k+1）（10）
其中：I为与协方差同维的单位矩阵。
利用上述算法，就可以根据当前观测值z（k+1）和前一个估计值■（k|k）来给出当前风速数据真值的估计■（k+1|k+1）。
2 实测数据计算机仿真处理及分析
这里选取日本京都大学所设风廓线雷达在印度尼西亚坤甸地区测取的风速数据，对上述EKF滤波方法进行仿真验证。观测参数如表1所示。
为了获取风廓线雷达上空三维风速信息，至少需要三个不共面的波束。而为了提高探测精度，我们所取探测数据利用的风廓线雷达采用五波束形式，波束方向参见表1所示。其中，1波束为垂直方向，其他四个波束是在方位上均匀分布的倾斜波束，分别指向东、南、西、北四个方向。
在直角坐标系下，将风速分解为u、v、w三个分量。我们将东、西两个方向测得的数据进行合成，得到u分量；将南、北两个方向测得的数据进行合成，得到v分量。垂直波束测得的数据作为w分量。在运用扩展卡尔曼滤波进行数据处理的过程中，选取的距离库范围为：0km—10km，分别将每个距离库的u、v分量作为系统状态输入矩阵X（k）进行运算，可以计算得出滤波后的风速最优化估算值■（k+1|k+1）。
图2为2011年7月17日凌晨2时50分风速u、v分量滤波前后对比图。使用卡尔曼及一阶扩展卡尔曼滤波对风廓线雷达数据进行质量控制时，处理的数据层面为：水平风分量u、v，选取的距离库范围为：0.2km～9.3km，分别对每个距离库的u、v分量进行滤波处理。
从图2中可以看出：在未滤波的情况下，原始测量数据中掺杂较多的噪声干扰。进行传统卡尔曼滤波处理后，明显去除了高频干扰。但我们发现，经过传统卡尔曼滤波处理之后，曲线在2km以下，出现了边界突跳和较大震荡问题，这是由于卡尔曼滤波方法受模型参数的影响很大。而采用一阶扩展卡尔曼滤波后，由于其采用的是非线性化处理方法，滤波模型更适合探测数据表现为非线性的风廓线雷达数据处理问题，因此，既对数据起到了平滑作用，有效去除了包括间歇性干扰和高频噪声干扰在内的高频成分，同时解决了边界突跳问题，得到了比线性卡尔曼滤波更好的滤波效果。
以上是根据风场实测数据，分析在某一时刻不同高度、不同风速分量上的一阶卡尔曼滤波效果。为了验证此方法在风廓线雷达风速数据处理当中的普遍实用性，我们从另一角度进行了滤波实验。取2011年2月1日印度尼西亚坤甸地区测得的风速数据，并将u、v、w三个分量上的数据进行合成。根据图3的几何关系（其中？兹是波束的天顶角，？渍表示方位角，取正北方向？渍=0），可以推出径向速度Vr（？兹，？渍）与风速u、v、w三分量的关系为：
Vr（？兹，？渍）=usin？兹sin？渍+vsin？兹cos？渍+wcos？兹（11）
选取2011年2月1日24小时之内风速数据分量，将其按照公式（11）合成后得到径向速度，并运用扩展卡尔曼滤波对其进行数据处理，输入量为各个时间点同一高度上的径向速度Vr（？兹，？渍），观察滤波后的风速最优化估算值■（k+1|k+1）。图4给出了2011年2月1日印尼坤甸地区上空径向风速滤波对比图。
由图4可以看到，在0～2小时之间，卡尔曼滤波过程出现了突跳问题。在取同一高度不同时间的风速分量，并将各个风速分量合成为径向速度之后，运用一阶卡尔曼技术进行滤波数据处理，得到了和上个实验相似的结果：与传统卡尔曼滤波相比，一阶卡尔曼滤波不仅能够对风廓线雷达全天风速数据的平均值进行有效平稳处理，而且避免了边界突跳问题，滤波效果优于传统卡尔曼滤波方法。
3 小结
本文针对实际风廓线雷达数据处理过程中遇到的风场非线性问题，提出了基于扩展卡尔曼滤波的风廓线雷达数据处理方法。通过实测数据处理，表明该方法在处理雷达数据的过程中能够有效去除噪声干扰，既对数据起到了平滑的作用，同时解决了边界突跳问题，弥补了经典卡尔曼滤波参数设置方面的不足。与线性卡尔曼滤波方法风速处理曲线进行比较，扩展卡尔曼滤波显示出自身的非线性特性，从两个角度分别进行仿真实验，结果表明其在工程应用中的普遍适用性。需要指出的是：扩展卡尔曼滤波结果的好坏与状态噪声和观测噪声的统计特性有关，在EKF的递推滤波过程中，状态噪声和观测噪声的协方差矩阵保持不变，如果这两个噪声协方差矩阵估计不够准确，则容易产生误差累积，导致滤波发散。
参考文献：
[1]何平.相控阵风廓线雷达.北京：气象出版社，2006：13—15.
[2]Weber B.L.， D.B. Wuertz， D.C. Welsh， R.McPeek. Quality controls for profiler measurements of winds and RASS temperatures. J.Atmos. Oceanic Technol.， 1993， 10： 452—464.
[3]杨馨蕊.风廓线雷达数据质量控制方法探讨及其应用.中国气象科学研究院，2009，5：50—58.
[4]D. Simon. Kalman filtering with state constraints： a survey of linear and nonlinear algorithms. Control theory and application， IET， 2010， 4（8）： 1303—1318.
推荐访问： 数据处理 仿真 方法 研究 风廓线
