引言

国内对于海洋表面风场的研究主要利用地基多普勒雷达和风廓线雷达^[1-3]的探测结果，但其监测范围小，无法满足大范围精确观测的需要，因此有必要发展我国星载微波散射计探测海洋表面风场，国内有关风场的探测主要集中于中高层风场探测^[4-5]。目前我国气象部门在国产的星载探测仪器对于海洋表面风场的研究还处于起步阶段，主要是利用美国Seawinds和欧洲ASCAT (Advanced SCATterometer) 的资料开展研究工作。由于海面风场主要由台风气旋风场和非气旋风场构成，非气旋风场 (即去除台风气旋覆盖区域之后的风场，该区域内风向变化不会超过90°) 风向块状模糊去除问题是海洋表面风场反演的一项关键技术，本文利用欧洲ASCAT资料开展了该项关键技术的研究。

ASCAT^[6-7]是一种推帚型扫描方式、扇形波束的散射计，其工作频率为5.5 GHz，对于雨衰不敏感，通过排列与系统两侧45°，90°和135°的各3个辐射波导天线，产生仰角宽度为40°、方位角宽度为0.4°的扇形波束。由于采用推帚式描扫工作方式，地面轨道上的同一风元 (wind vector cell，WVC) 一般会有3个雷达后向散射截面积NRCS (normalized radar cross section) 实测值 (前向、中向、后向)，利用这些不同观测条件下获取的NRCS值和相应的观测参数，即可反演出海面风场。

由于NRCS与风速存在着对数正相关，与相对风向x存在双调和关系，所以在理论上只要有2个以上实测NRCS，就可反演出1个确定的风矢量。但由于大气、云水和雨衰，以及卫星运行状态的不稳定性，得到的NRCS会受到污染，其误差满足平均值为零的高斯分布。这些受污染的NRCS，使得在近似风速条件下，出现风向180°变化的解，通过最大似然法进行风场反演，也往往只能得到一组模糊解。模糊解一般为2~4个，其中只有1个接近真实风矢量，其余的称为伪解。所以，要得到每个风元上的真实风矢量，必须对初步反演的结果进行模糊去除。

模糊去除一般采用圆中数滤波法，但是传统的圆中数滤波法^[8-13]无法去除块状模糊，同时传统的圆中数滤波法还会由于第1风场的缺陷使得滤波迭代毁坏相邻的风矢量值。本文采用一种加强的圆中数滤波法对非气旋风场进行模糊去除，能够有效去除块状模糊问题，该方法定义简单明确，计算量小，且收敛速度较快。在满足圆中数滤波法初始风场真实风矢量比例必须达50%以上的前提条件后，即使不是随机分布也可以得到较好效果 (不集中于某一连续区域)，所以滤波器工作性能并不一定像传统圆中数滤波法一样完全取决于初始风场的质量。

本文提出的加强型圆中数滤波法能够根据非气旋风场风矢量分布特点，可以克服传统圆中数滤波法的苛刻条件 (如各风元模糊风场必须随机分布，不可有块状模糊)，从该风元的几个模糊解中提取出与真实矢量最近的风矢量解，克服块状模糊，取得良好效果。但在业务运行该方法前，应当使用其他星载仪器确定台风云图覆盖区域，在去除台风覆盖区域后，在非气旋风场区域使用本文方法。因此该方法为业务化非气旋洋面风场数据提取给出了新方案。

1.   加强型圆中数滤波算法

加强型圆中数滤波算法首先通过非气旋风场特点进行初始化，能够有效去除块状模糊问题，使得导致圆中数滤波失效并且影响邻域的数据得到修正。然后通过在风场二维空间M行N列中开一定大小的窗口，计算该窗口数据的圆中数，然后在窗口中心风元对应的几个模糊解中找出与圆中数最接近的一个替代当前窗口中心的风矢量，移向下一位置，通过反复迭代，直到风场不再改变或者迭代次数达到预设的最大迭代次数为止。最后对风场中的个别模糊点进行平滑处理。

1.1   非气旋第1风场初始化

在气象卫星刈幅范围内的非气旋风场中，风向难以出现大角度的突变，风向主要集中在以该范围内风向平均值为中心的180°范围内，因此根据非气旋风场风向分布特点，可以对最大似然估计解算^[14-15]的第1风场进行初始化，按照风向将360°分为8个区间：[0°，45°)，[45°，90°)，[90°，135°)，[135°，180°)，[180°，225°)，[225°，270°)，[270°，315°)，[315°，360°)，每个区间的代表角度为θ₁=22.5°，θ₂=67.5°，θ₃=112.5°，θ₄=157.5°，θ₅=202.5°，θ₆=247.5°，θ₇=292.5°，θ₈=337.5°，它们分别代表 8个角度范围的特征值。对第1风场风向进行统计，统计出各个范围中风向的个数，即某个代表角度θ_i区间里有n个角，将这8个区间中具有最多角度个数的区间提取出来，并且提取其左右相邻两边区间的风向数量进行加权求风向平均值。

对第1风场进行初始化，将每个风元4个模糊解的角度θ_i，j^k，i∈M，j∈N，k∈4按照从最大可能解到最小可能解的顺序依次和该风向平均值θ进行比较，取差值Δθ_{min (i，j)}最小的风向θ_i，j^k作为真解，

并且记录对应的风速值，对各个风元的模糊解进行重新排序，初始化整个M行N列风元后得到第2风场。

1.2   圆中数的迭代计算

传统的圆中数[9，16-18]：圆中数θ满足分布在两半圆中的风向数目相等。在该定义下，会出现多个圆中数，需从中选择与圆平均数最接近的一个作为唯一的圆中数。因此，该定义下的圆中数计算较复杂。对于散射计模糊去除，由于已有第2风场中每个风元都有几个模糊解θ_ij^k可供选择，所以不必要在整个风矢量窗口中进行搜索确定圆中数，然后在模糊解中选择一个与圆中数最接近的模糊解做为窗口中心的风矢量；可以直接从θ_ij^k中选出一个使式 (2) 达到最小的一个模糊解代替窗口中心的风矢量。若选取M×N的二维风矢量窗口进行滤波运算，这里引入一种分区的圆中数定义，将一个矩形区域风场分为9个区，每个区的领域范围不同，导致了迭代取值范围的不同。分区定义见式 (3)。

加强型圆中数计算与传统圆中数计算相比具有如下优点：① 通过初始化第1风场，第2风场的模糊性降低了，克服了传统圆中数滤波方法无法解决的风场块状模糊问题。② 计算过程简单，不需统计圆直方图，也不需计算圆平均数。③ 该定义下的圆中数只有一个，不会出现圆中数多解现象。④ 分区后的圆中数滤波法不会因为边界区域风元取值狭窄问题，使得计算无法进行。

1.3   加强型圆中数滤波的运算过程

① 非气旋风场初始化，对最大似然估计计算的第1风场中的模糊解进行重新排序，得到第2风场。

② 确定滤波器的参数 (窗口大小、权值、最大迭代次数)。

③ 分区圆中数滤波：对于 (i，j) 位置的风元，根据式 (3) 分为9个区域，每个区域有特定的迭代区域。选出其中一个模糊解做为圆中数，代替原窗口中心的风矢量。滤波器移向下一位置，重复上述操作。

④ 迭代：重复步骤③，直到当前风场不再改变或迭代次数达到最大迭代次数为止。

⑤ 对当前风场数据中个别模糊数据进行精细化平滑处理，得到加强型圆中数滤波风场。

2.   ASCAT风场模糊去除

2.1   ASCAT第1风场中模糊解的分布特点

通过排列与系统两侧45°, 90°, 135°的各3个辐射波导天线，产生仰角宽度为40°、方位角宽度为0.4°的扇形波束。由于采用推帚式描扫工作方式，左、右刈幅宽度为550 km，左、右刈幅间隔672 km。其主要观测几何参数及地面轨道情况如图 1所示。

图  1  ASCAT的轨道示意图

Figure  1.  Swath of ASCAT

下载: 全尺寸图片 幻灯片

由于天线固定推帚式扫描模式，每一个星下点风元都会有3次测量值，前向、中向和后向。左、右刈幅对称排布，最大似然估计方法求出的第1风场，每个风元都有多个可能的模糊解，在实际解算中留下最有可能的4个解，作为后面解模糊的风元参数库。

2.2   ASCAT风场初始化方案

传统的圆中数滤波方法不能有效去除风元的连续块状模糊问题，如果不进行初始化而直接对第1风场进行圆中数滤波，不仅不能去除块状模糊，而且这些块状模糊解还会污染相邻区域的真解，使得去模糊彻底失败。

由于第1风场中模糊解正确率超过50%，且非气旋风场的特点是有限区域内风向不会出现大角度突变，因此结合此特点，首先对第1风场进行初始化，得到第2风场。

2.3   ASCAT模糊去除过程

根据上述ASCAT第2风场中模糊解分布的特殊性，对整个区域再次进行圆中数滤波，得到加强型圆中数滤波结果。

2.4   模糊去除方法的适用性讨论

从以上对模糊去除方法的描述，对于在空间上连续、渐变的简单风场来说，该方法显然是适用的。对于复杂风场 (如有锋面)，如果风场突变180°以上时，部分失效。综上所述，本文提出的模糊去除方法可以解决简单风场分布和部分复杂风场分布情况。

3.   试验分析

3.1   数据选取与参数设定

ASCAT散射计每条轨道包括1629行、42列个风元数据，每个风元大小为50 km×50 km。本文选取了ASCAT 2009年9月15日轨道号为15092的部分数据 (地理范围是介于9°~29°N，142.8°~153°E，跨度范围：1100 km长，550 km宽，该区域无雨) 和2011年2月27日部分数据 (地理范围是介于3°S~11°N，171.2°~178.2°W，跨度范围：800 km长，550 km宽，该区域无雨区居多)，两个时段数据对模糊去除方法进行验证。由于ASCAT工作在C波段，该波段对云雨导致回波信号衰减并不灵敏，因此没有区分该区域内晴空、云水和雨区的单独影响，而是进行统一的反演考虑。因此滤波窗口大小选7×7较为适宜^[5]。故在本试验中，圆中数滤波器参数设定如下：窗口大小为7×7，最大迭代次数为100次。这里没有考虑权重，因为在窗口较小的情况下权重影响不大。

3.2   2009年9月15日数据分析

为了进行对比，先对所选数据全部用第1风场中模糊解初始化风场进行圆中数滤波，然后采用加强型圆中数滤波方案进行滤波。试验结果如图 2~图 4所示。各图中的箭头方向表示风向，箭头长度表示风速，长度越长表示风速越大，反之风速越小。

图  2  根据2009年9月15日轨道数据得到的第1风场

Figure  2.  The most likely wind vector field after wind retrieval from orbit data on 15 Sep 2009

下载: 全尺寸图片 幻灯片

图  3  根据2009年9月15日轨道数据经过传统圆中数滤波直接处理过后的结果

Figure  3.  Wind vector field after traditional circle median filtering wind retrieval from orbit data on 15 Sep 2009

下载: 全尺寸图片 幻灯片

图  4  根据2009年9月15日轨道数据经过加强型圆中数滤波方案滤波结果

Figure  4.  Wind vector field after enhanced circle median filtering wind retrieval from orbit data on 15 Sep 2009

下载: 全尺寸图片 幻灯片

图 2为经最大似然法算出的第1风场。从图 2可以看出，中间区域出现块状风向反向和不连续现象，右上角区域出现少量风向不连续现象，整个区域有半数以上的风向与真实风向一致，满足传统圆中数滤波条件，即该区域中的伪解可以由圆中数滤波器去除。但是块状模糊问题由传统圆中数滤波算法无法解决，只能解决右上角少量不连续现象。

图 3为图 2的第1风场直接经过传统圆中数滤波过后的结果。与图 2相比，图 3右上角区域的少量风向模糊得到有效去除。但是，中间区域大量连续模糊没有得到改善，出现了所谓的块状模糊，该区域内的风向与周围风向相差160°~180°。故用传统圆中数滤波方法在该区域失效。

图 4为采用本文所提出的加强型圆中数滤波方案进行滤波后的结果。滤波迭代次数为9次，每次风矢量改变数目依次为662，221，160，151，129，96，65，26，0。与图 3该区域进行比较可以发现，图 4中间区域的块状模糊得到有效去除，右上角少量不连续也得到纠正，整个风场风向分布较为连续，达到了加强型圆中数滤波去除风矢量模糊的目的。具体风矢量数据分析见数据偏差分析。

为了对滤波结果精度进行定量分析，将滤波结果风场与欧洲EUMETSAT对同一区域反演的风场逐分辨单元进行偏差、绝对偏差统计，统计数据 (样本量为1638) 如表 1所示。

表  1  本研究结果与EUMETSAT L2风场数据偏差和绝对偏差统计

Table  1.  Statistics on biases between filtered wind vector field and L2 wind vector field from EUMETSAT

统计量	偏差		绝对偏差
	风速/(m·s-1)	风向/(°)	风速/(m·s-1)	风向/(°)
最小值	-10.09	-48.4	0	0
最大值	3.31	57.1	10.09	57.1
平均值	-1.5858	1.0424	1.859	8.426
方差	3.3437	155.17	2.4035	57.1

下载: 导出CSV 

| 显示表格

表 1中的各项指标基于偏差、绝对偏差进行统计。对于偏差，风速偏差满足正态分布，风速误差平均值为-1.5858 m·s^-1，测量风速比真实风速普遍偏小，主要是由于对主动探测的回波功率没有做水汽、云水和氧气衰减修正，使得实际收到回波功率比真实值偏小，由于回波功率与风速存在对数正相关，因此测得的风速偏小，本文重点是解决风向块状模糊问题，没有做大气、云水和雨衰订正，因此风速偏小。风向偏差满足正态分布，风向误差平均值为1.0424°，风向偏差较小，因为块状模糊得到较好解决。对于绝对偏差，平均偏差远小于风场反演精度要求的最大误差值 (风速误差小于2 m·s^-1，风向误差小于20°)。考虑到L2反演精度比较稳定，L2风场与海面真实风场的偏差以及本文滤波结果与L2的偏差均符合正态分布，可以推断绝大部分分辨单元的误差应该在反演精度要求的范围之内。试验表明，模糊解之间风速差异很小，而风向相差较大。加强圆中数滤波去除模糊解的目的是从几个模糊解序列中选出一个最可能接近真实解的风矢量，并不改变模糊解的值，因此滤波精度主要体现在风向偏差上，而表 1中平均风向偏差确实很小 (只有8.426°)。以上分析从定量的角度证明了本文采用的滤波算法是比较可靠的。

3.3   2011年2月27日数据分析

为了进行对比，先对所选数据全部用第1风场中模糊解初始化风场进行圆中数滤波，然后采用加强型圆中数滤波方案进行滤波。试验结果如图 5、图 6所示。各图中的箭头方向表示风向，箭头长度表示风速，长度越长表示风速越大，反之风速越小。

图  5  根据2011年2月27日轨道数据经过传统圆中数滤波直接处理过后的结果

Figure  5.  Wind vector field after traditional circle median filtering wind retrieval from orbit data on 27 Feb 2011

下载: 全尺寸图片 幻灯片

图  6  根据2011年2月27日轨道数据经过加强型圆中数滤波方案滤波结果

Figure  6.  Wind vector field after enhanced circle median filtering wind retrieval from orbit data on 27 Feb 2011

下载: 全尺寸图片 幻灯片

图 5为第1风场直接经过传统圆中数滤波过后的结果。该图上部大部分区域风场较好，但是左下角区域出现了块状风场模糊问题。右边下半部分区域出现少量甚至小块状连续模糊问题，该区域内的风向变化与周围风向相差较大。故用传统圆中数滤波方法在下半部分区域失效。

图 6为采用本文所提出的加强型圆中数滤波方案进行滤波后的结果。滤波迭代次数为11次，每次风矢量改变数目依次为853，521，363，251，117，72，57，37，21，7，0。与图 6该区域进行比较，可以发现，图 5左下角区域的大面积块状模糊得到有效去除，右下角不连续也得到纠正，整个风场风向分布较为连续，达到了加强型圆中数滤波去除风矢量模糊的目的。具体风矢量数据分析见数据偏差分析。为了显示差异，图 5和图 6的图像箭头归一化参数不同。

为了对滤波结果精度进行定量分析，将滤波结果风场与欧洲EUMETSAT对同一区域反演的风场逐分辨单元进行偏差、绝对偏差统计，统计数据 (样本量为1155) 如表 2所示。

表  2  本研究结果与EUMETSAT的L2风场数据偏差和绝对偏差统计

Table  2.  Statistics on biases between filtered wind vector field and L2 wind vector field from EUMETSAT

统计量	偏差		绝对偏差
	风速/(m·s-1)	风向/(°)	风速/(m·s-1)	风向/(°)
最小值	-8.09	-30.4	0	0
最大值	4.17	51.7	8.03	71.1
平均值	-1.932	2.1457	1.651	10.317
方差	5.1473	175.14	3.2156	78.42

下载: 导出CSV 

| 显示表格

表 2中的各项指标基于偏差、绝对偏差进行统计。对于偏差，风速偏差满足正态分布，风速误差平均值为-1.932 m·s^-1，测量风速比真实风速普遍偏小，主要原因同表 1。风向偏差满足正态分布，风向误差平均值为2.1457°，风向偏差较小，因为块状模糊得到良好解决。对于绝对偏差，平均偏差远小于风场反演精度要求的最大误差值 (风速误差小于2 m·s^-1，风向误差小于20°)。考虑到L2反演精度比较稳定，L2风场与海面真实风场的偏差以及本文滤波结果与L2的偏差均符合正态分布，可以推断绝大部分分辨单元的误差应该在反演精度要求的范围之内。试验表明，模糊解之间风速差异很小，而风向相差较大。加强圆中数滤波去除模糊解的目的是从几个模糊解序列中选出一个最可能接近真实解的风矢量，并不改变模糊解的值，因此滤波精度主要体现在风向偏差上，而表 2中平均风向偏差确实很小 (只有10.317°)。以上分析从定量的角度证明了本文采用的滤波算法比较可靠。

4.   结语

本文采用加强型圆中数洋面风场模糊去除方法是在传统圆中数定义基础上，结合散射计第1风场中模糊解在地面轨道的空间分布特性优化得到的。该方法与传统圆中数滤波法相比，主要特点是能够有效去除非气旋风场块状模糊问题，对圆中数的定义简单，计算量小，易收敛；并且经理论分析和试验验证，对于去除台风覆盖区域的非气旋风场分布和部分复杂风场分布情况，该洋面风场反演方法比传统的圆中数滤波方法更加具有通用性和普适性。