基于运动估计与补偿的雷达拼图预测分析*

秦运龙 张冰松 杨代才 王迎迎

（湖北省气象信息与技术保障中心 武汉 430074）

1 引言

随着复杂强对流气象天气的增多，雷达拼图对于强天气预警预报尤为重要，其及时性和完整性需求也越来越高。在业务传输中，由于设备、系统和网络等原因可能导致雷达拼图漏缺或延迟。如何及时较好地补缺雷达拼图，并提前预测预判雷达拼图的演变在短时强对流预报中具有较强的需求。本文根据雷达拼图特点，提出了一种基于运动估计和运动补偿的方法预测雷达拼图［1～2］，保障雷达拼图时序上的完整性，并可提前预判短时间内雷达拼图演变情况。

2 运动估计

每种雷达拼图产品都是连续动态变化的，每6min 一张图，因此连续雷达拼图的纹理在短时间内是不会发生较大变化，通过纹理特征可以进行运动估计判断雷达图变化与走向。

设雷达拼图实况在时序上的拼图为I，时序拼图序列集为{Ij，j∈(1，2，3…)} ，其中Ij表示第 j 张雷达拼图。将每张拼图划分为多个8*8 像素的宏块，宏块序列集表示为 {MBj(m，n)，j∈(1，2，3…)，m，n∈(1，2，3…)}其中MBj(m，n)表示Ij拼图的第m 行第n 列个宏块。宏块作为拼图变化的最小运动单元，每个宏块前后移动的横向与纵向位移对应一个运动矢量［3］，宏块MBj(m，n)的运动矢量记为MVj(m，n)。

2.1 纹理特征提取

图像的纹理特征是图像的低频信息，表征了图像的具体的内容，通过对每个宏块进行二维离散余弦变换（Discrete Cosine Transform，DCT）可获得宏块的二维频谱信息［4～5］，从而可取的低频的纹理信息数据。设宏块MBj(m，n)对应的8*8 空间像素值为f(x，y)，则其频谱变化为

式中，M=8，N=8，

F(u，v) 为宏块对应8*8 大小的二维频谱信息。通过提取F(u，v)的低频纹理信息来压缩数据量获得宏块特征。其中F(0，0)为直流信息，表征图像的主要纹理信息，其余为交流信息，且频率越低，纹理信息越丰富。二维频谱数组中越靠近左上角频率越低，信息量越大，其信息权重也越高，通过设置W1，W2，W3，三个不同权重系数提取压缩后频谱。

通过计算频谱值平方和得到宏块的特征值T，其中W1＞＞W2＞＞W3。即有宏块MBj(m，n)对应特征 值 为Tj(m，n) ，Tj(m，n) 反应了整个宏块MBj(m，n)的纹理特征。由于每张拼图由多个宏块组成，则计算第j 张拼图Ij得到宏块集MBj特征值集Tj，同理可得到时序拼图序列集为I对应宏块集MB特征值集T［6］。

2.2 运动矢量计算

天气雷达拼图反应的是大气云层及水汽含量等特征变化情况，是一个缓慢连续，有一定继承性的过程，不会出现较大的跳变运动，因此相邻近雷达拼图间宏块移动距离有限。可根据预测的时间序列跨度大小，选择适当的宏块特征值匹配搜索范围，降低计算量，减少其它区域特征值引起的误判匹配［7～9］。

图1 运动估计

如图1，计算从Ij到Ij-1雷达拼图的运动矢量集Tj，从而预测Ij+1各宏块对应的位移。设图像Ij有 M 行 N 列个宏块，搜索范围为 L，以 L=5 为例，即Ij某宏块MBj(m，n)在Ij-1中对应宏块MBj-1(m，n)搜索半径为5 个宏块距离进行特征值匹配，最接近的作为相似度最高的宏块，即计算其中l∈(-5，-4，-3，…，3，4，5) 。则Ij-1与Ij相似宏块为MBj-1(m-m1，n-n1) ，MBj(m，n) ，其中 (m1，n1) 为相似宏块的横向和纵向位移，即为运动矢量MVj(m，n)，MVj(m，n)作为MBj(m，n)）预测下一张雷达图Ij+1对应相似宏块MBj+1(m+m1，n+n1) 的运动矢量。计算所有Ij宏块运动矢量得到运动矢量集MVj。同理从MBj-2(m-m2，n-n2)到MBj(m，n)宏块的位移(m2，n2)即为MBj(m，n)预测未来第二张雷达图对应相似宏块MBj+2(m+m2，n+n2)的运动矢量，依次类推，即可获得宏块的运动变化矢量集［10～11］。随着预测雷达拼图的时间跨度增大可适当增大宏块特征值匹配搜索范围L的值。

3 运动补偿

雷达拼图在连续缓慢移动过程中回波强度也在发生一定的变化，运动补偿则反映了雷达拼图前后回波强度变化情况，因此在预测雷达拼图时，不能简单地将原图宏块的像素值复制到预测图像相应匹配宏块位置，而应该加上一定的运动补偿信息［12～13］。

在获得运动矢量的基础上，可以很迅速地找到匹配宏块对应关系，例如Ij与Ij-1的宏块对应关系MBj-1(m-m1，n-n1)对应MBj(m，n)。则该宏块的运动残差cj(m，n)=MBj(m，n)-MBj-1(m-m1，n-n1)。Ij预测第Ij+1时的运动补偿值rcj取前两个残差值平均值

4 图像预测

设图像Ij有 M 行 N 列个宏块，则Ij-1与Ij相似宏块MBj-1(m-m1，n-n1)、MBj(m，n)，其中 (m1，n1)为相似宏块的横向和纵向位移。则其预测的下一张图对应相似宏块为MBj+1(m+m1，n+n1)，通过宏块MBj(m，n)像素值和对应运动补偿值rc（jm，n）得到预测的相似宏块的像素值［14～15］为MBj+1(m+m1，n+n1)=MBj(m，n)+rcj(m，n)。计算所有相似宏块像素值得到第j+1张拼图的预测拼图Ij+1。

5 实验结果与分析

本文采用长江流域雷达拼图组合反射率产品作为实验数据集，图片大小为768×1024 像素，图片每 6min 生成一个，通过调试，取 W1=20，W2=6，W3=2，L=5。在实际预测实验过程中，考虑到图片边界静止不变，且宏块搜索时需进行边界处理，将图片周边8 个宏块均不做预测，被预测雷达拼图像素值初始化为用于作为预测参考的雷达拼图像素值，在预测计算过程中通过宏块匹配和补偿更新被预测雷达拼图像素值，其预测结果如图2，其中第一行为实况，第二行为预测结果。

图2 雷达拼图预测结果

通过计算预测拼图和原实况拼图灰度值相似度，判断预测拼图的质量，拼图相似度以相关系数表征，计算相关系数r的值。

多次预测得相关系数r 平均值变化情况如图3所示。实验证明拼图相关系数大于0.9时其内容与实况产品基本一致，可以满足预报需求，但在色度上仍有一定偏差。相关系数低于0.9 时，预测误差逐步扩大，不宜提供作预报服务。

图3 预测结果与实况相关系数变化曲线

6 结语

本文提供一种基于运动估计与补偿的雷达拼图预测分析方法，该方法利用雷达拼图在纹理变化上连续性的特点，能够较好补缺和提前预测雷达拼图，实验证明，该方法针对短时间范围内的雷达拼图预测结果与实况数据吻合度较高，但在色度上还存在一定的偏差，但针对时间跨度略大的预测，其偏差较大，还需进一步研究。