基于梯度算子的高分三号SAR风场反演方法
陈克海1 张金兰1 解学通2 邓 丹3
（1.广东工贸职业技术学院 测绘遥感信息学院，广东 广州 510510；2.广州大学 地理科学与遥感学院，广东 广州 510006；3.广东国图勘测地理信息有限公司，广东 广州 510180）
摘 要：鉴于高分三号（GF-3）合成孔径雷达（SAR）图像具有黑白相间条纹特征，提出一种基于梯度算子的风场反演方法。对SAR 图像重新采样后，利用Sobel 梯度算子求解海面风向值，接着使用扩展后的CMOD5散射模型推算出风速值。为了测试该方法，对GF-3 SAR 图像进行风场反演实验，并使用欧洲中尺度天气预报中心（ECMWF）再分析风场数据和美国国家浮标中心（NDBC）浮标风矢量数据对反演得到的SAR 风场质量进行评估。实验结果表明：GF-3 SAR 风场相比于ECMWF 再分析风场，风速、风向均方根误差分别为1.79 m/s和11.95°；相比于NDBC 浮标数据，风速、风向绝对偏差分别为0.24 m/s和9.03°。因此，本文方法能有效地从GF-3 SAR图像中反演出较高精度的风场，研究结果可为GF-3 SAR风场反演技术提供参考。
0 引言
目前海面风场遥感监测传感器主要有散射计、高度计和合成孔径雷达（synthetic aperture radar，SAR）。其中，高度计垂直向下观测，能观测星下线附近区域的风速，却无法观测风向值。散射计一般具有两个旋转的波束，这两个波束的极化方式和入射角不同，在较大范围内对地重复观测，在同一地面单元上可产生多个观测值，通过反演算法反演出风速、风向值。但是，散射计一般在25 km×25 km 的范围内只有一个风速、风向反演值，其风场数据仅适用于全球环境研究与应用［1-2］。SAR 空间分辨率高，最高可达到米级，能观测到海面上中小尺度风场结构，为区域环境研究和应用提供数据，是目前海面上中小尺度风场观测的热点方向之一［3］。
SAR 能捕捉到海面风场引起的黑白相间的风条纹，由于风条纹走向基本平行于海面风向，可以利用风条纹信息从SAR 图像中提取海面风向［4-5］。目前SAR 图像风场反演一般分两步走，第一步根据SAR 图像风条纹走向求解出风向值，第二步在已知的后向散射模型函数条件下，根据图像后向散射系数和已求解得到的风向值推算出风速值。可见，SAR 图像风场反演关键在于风向值求解。目前，一般采用频谱分析方法求解SAR 图像风向值。首先利用傅里叶变换等方法将SAR 图像从空间域转换到频率域，然后计算功率谱，并基于功率谱上波峰走向垂直于风条纹走向的假设，从波峰走向推算出海面风向［6-8］。
频谱分析方法具有较高的数学理论依据，但是存在一些缺点，比如频谱分析法要求图像中海洋区域必须完整，区域内不能出现陆地、海岛或者船只等其他物体，否则无法开展傅里叶变换及后续处理。另外，频谱分析法要求图像区域要足够大，否则会降低精度。借鉴在图像边缘检测技术中被普遍采用的灰度梯度算子在提取边缘方向时对处理区域的形状和大小没有特别要求，在SAR 图像风场反演算法中引入灰度梯度算子，用于计算图像各像素灰度值梯度方向。在统计意义上，风条纹走向应垂直于分布概率最大的梯度方向。据此，可以利用梯度算子推算出海面风向。为了测试本文方法，对高分三号（GF-3）SAR图像进行风场反演实验，反演得到的SAR 风场与欧洲中尺度天气预报中心（European Centre for Medium-Range Weather Forecasts，ECMWF）再分析风场数据和美国国家浮标中心（National Data Buoy Center，NDBC）浮标风矢量进行对比，分析SAR风场质量 。
1 数据来源
1.1 GF-3 SAR图像
我国GF-3 卫星于2016 年8 月10 日成功发射［9］，其上搭载的SAR 工作在C 波段。它具有多个成像模式，本实验采用其中的标准条带成像模式，数据等级为L1A，文件保存为GeoTiff 格式。每幅图像约覆盖10 000 km2，空间分辨率约为5 m×7 m。在L1A 文件中，像素行号沿着卫星前进方向递增，而列号由近端向远端递增。每个L1A 文件元数据文件包含定位、观测参数等基础数据，利用这些基础数据可以进行辐射校正、空间定位和观测几何参数计算等预处理。
1.2 辅助数据
本文使用ECMWF 再分析风场数据和NDBC浮标风矢量等辅助数据来评估风场反演精度。
（1）ECMWF 再分析风场数据是欧洲空间局事后整合全球卫星数据和各种实测数据优化得到的高精度数值模拟产品，特点是数据量多、精度较高，被广泛用于海洋卫星产品质量评估［10-11］。选用的ECMWF 再分析风场的空间、时间分辨率分别为0.125°×0.125°和6 h，通过在空间维度上采用二维线性插值和在时间维度上采用三次样条函数插值，可得到任意时刻、位置的风速和风向值。
（2）NDBC 浮标风矢量数据由美国国家浮标中心提供［12］。被选用的浮标位于离岸50 km以上的海域，具有连续的观测能力，每10 min 记录一次风速、风向测量值，通过插值得到任意时刻的风速、风向值。
2 反演方法
2.1 风向求解
在开始风场反演之前，以100 m×100 m 的分辨率对SAR 图像重新进行采样，并将图像划分为25 km×25 km 的子图［13］。风场反演技术是在子图上进行。
2.1.1 计算灰度梯度值
子图各个像素在水平、垂直方向上的灰度梯度值Gx和Gy为
式中，σoij 为窗口内第i 行、第j 列的后向散射系数；Sxij和Syij为分别为在水平、垂直方向上的Sobel算子值，具体数值见图1。在图像空间上移动Sobel 算子，计算子图各个像素在水平、垂直方向上的灰度梯度值Gx和Gy。
图1 Sobel算子
2.1.2 计算灰度梯度方向值
根据梯度值Gx和Gy，计算子图每个像素灰度梯度方向值D为
2.1.3 计算风向可能解
根据各像素灰度梯度方向值，绘制灰度梯度方向概率密度曲线。由于噪声的存在，该曲线存在较大起伏，使用均值滤波器进行平滑处理。对于风条纹而言，风条纹走向的灰度梯度方向必然对应于概率密度曲线峰值，由此可以推算出海面风向ϕ为
式中，φ 为平滑处理后概率密度曲线峰值对应的灰度梯度方向值。由式（4）可知，由灰度梯度值计算得到的风向不可避免地存在180°模糊，即存在两个可能风向，分别为ϕ和ϕ+180°。
以上海面风向是在图像坐标系上计算得到的，需转化到地理坐标系上。对于GF-3 SAR，观测方向垂直于卫星方向，列号沿着观测方向递增，可见图像坐标系x轴正方向为观测方向；行号沿着卫星前进方向递增，可见y 轴正方向为卫星前进方向。根据观测几何关系，SAR 风向ϕ 可从图像坐标系转换到地理坐标系，公式为
式中，Φ 为雷达观测方位角，正负号由雷达观测方向而定，当雷达观测方向为左视时选择正号，为右视时选择负号。
2.1.4 去除风向模糊
根据ECMWF 再分析风场数据，从两个模糊解中挑选出与ECMWF 风向最接近的一个作为最终风向解。
2.2 风速求解
GF-3 SAR 工作在C 波段。CMOD5 为目前比较成熟的C 波段散射模型，但是只有垂直（VV）极化模型，缺少水平（HH）极化模型［14］。对于重置（HH）极化模型值，可使用极化比系数由VV 极化模型值转化得到［15］。
式中，σH0 H和σV0 V分别为HH、VV 的后向散射系数模型值；u 和ϕ分别为风速、风向值；θ为雷达入射角；α为极化比系数，这里取值1。
图2为扩展后的CMOD5散射模型曲线图（入射角为40°），其中图2（a）、2（b）分别对应垂直（VV）和HH 极化，图中各曲线从下到上分别表示风速为5、10、15、20、25、30 m/s 时不同相对风向的后向散射系数模型值。
图2 扩展后的CMOD5散射模型曲线图
从子图中计算出风向之后，利用扩展后的CMOD5散射模型，搜索出与该区域后向散射系数测量值均值最接近的风速值，作为该子图的风速值。
2.3 反演质量评估
使用ECMWF 再风场数据和NDBC 浮标风矢量数据，对反演得到的SAR 风场质量进行评估。对于ECMWF 再风场数据，以风速、风向均方根误差为评估参数；对于NDBC浮标风矢量数据，以风速、风向绝对偏差为评估参数。
设在第i 个匹配数据中，SAR 风速、风向分别为si和di，ECMWF或NDBC风速、风向分别为si0和di0，则以上评估参数可表示为：
（1）风速绝对偏差Esi
（2）风向绝对偏差Edi
（3）风速均方根误差Rs
（4）风向均方根误差Rd
以上风速、风向绝对偏差反映SAR 风速、风向值与参考值的绝对偏差值；风向、风速均方根误差综合反映SAR 风速、风向值与参考值之间的偏差，是衡量反演精度的主要参考指标之一。
3 实验结果
以产品号为2131355 的GF-3 SAR 图像为实验数据，先将图像划分为25 km×25 km 的子图，然后使用本文方法进行风场反演实验，反演得到的SAR 风场与ECMWF 再分析风场和NDBC 浮标矢量数据进行对比分析。
现在以该SAR 图像中第2 行第3 列的子图反演情况为例。图3（a）为该子图反演得到的风向和相应ECMWF 风向。图中灰度图为SAR 图像，粗、细箭头分别表示SAR 风向和ECMWF 风向。显然，该子图具有明显的明暗相间的风条纹，风条纹走向与ECMWF 风向基本一致，这说明该区域的ECMWF 风向非常接近于实际风向，可用于评估反演后的SAR 风场质量。从子图反演结果来看，SAR 风向非常接近于ECMWF 方向，可见该子图风场反演效果很好。
图3 GF-3 SAR子图反演结果示例
为了说明本文算法，这里给出该子图风向计算的中间结果。图3（b）为该子图对应的灰度梯度方向概率密度曲线，图中细曲线为原始的灰度梯度方向概率密度曲线，抖动较大，经平滑处理后为粗曲线。平滑后概率密度曲线峰值对应的灰度梯度方向垂直于海面风向。
图4（a）为整个SAR 图像反演后风场，图中粗箭头为各子图反演得到的SAR 风向值，细箭头为插值后风场，右上角最粗箭头为浮标风向，背景灰度值表示不同的风速值。图4（b）为对应的ECMWF风场图，细箭头为ECMWF风场。为了方便对比，将图4（a）SAR 风向和浮标风向也显示在图4（b）上。从图4（b）可以看出，各子图SAR风向非常接近于ECMWF 风向和浮标风向。经统计，与ECMWF 风场相比，SAR 风速、风向均方根误差分别为1.79 m/s 和11.95°；与浮标最近的子图SAR 风向跟浮标相比，SAR 风速、风向绝对偏差分别为0.24 m/s和9.03°。
图4 GF-3 SAR图像风场反演结果
4 结束语
本文针对SAR 图像提出了一种基于灰度梯度算子的风场反演方法，并使用GF-3 SAR 图像进行实验，反演得到的SAR 风场与ECMWF 再分析风场对比，风速、风向均方根误差分别为1.79 m/s 和11.95°；与NDBC 浮标风矢量对比，风速、风向绝对偏差分别为0.24 m/s 和9.03°，可见该方法可以有效从SAR 图像中反演出较高精度的风场。由于缺乏足够的GF-3 SAR 图像作为实验数据，本文仅作为SAR 海面风场反演的一种探索性研究。另外，由于风条纹的产生与海洋大气边界层的稳定性等因素有关，并非所有风场都能产生风条纹，这意味着并非所有SAR 图像都存在风条纹。对于没有风条纹的SAR 图像，需采用其他方法进行风场反演。
参考文献
［1］ 邹巨洪，林明森，张毅，等.海洋二号卫星微波散射计面元匹配［J］.遥感学报，2017，21（6）：825-834.
［2］ 王志雄，邹巨洪，林明森，等.HY-2卫星微波散射计海面风场地球 物理模式函数研究［J］.遥感学报，2023，27（4）：871-880.
［3］ 王超，潘广东.航天飞机成像雷达海面风矢量观测研究——以1994 年4 月南中国海试验区为例［J］.遥感学报，2000，4（1）：51-54.
［4］ ETLING D，BROWN R A.Roll vortices in the planetary boundary layer：A review［J］.Boundary-Layer Meteorology，1993，65（3）：215-248.
［5］ ALPERS W，BRÜMMER B.Atmospheric boundary layer rolls observed by the synthetic aperture radar aboard the ERS-1 satellite［J］.Journal of Geophysical Research Oceans，1994，99（C6）：12613-12621.
［6］ HORSTMANN J，KOCH W.Measurement of ocean surface winds using synthetic aperture radars［J］.IEEE Journal of Oceanic Engineering，2005，30（3）：508-515.
［7］ 杨劲松，黄韦艮，周长宝，等.合成孔径雷达图像的近岸海面风场反演［J］.遥感学报，2001，5（1）：13-16.
［8］ 张雷，石汉青，杜华栋，等.用小波分析反演星载SAR 图像海面风向［J］.遥感学报，2014，18（1）：215-230.
［9］ 张庆君.高分三号卫星总体设计与关键技术［J］.测绘学报，2017，46（3）：269-277.
［10］ DEE D P，UPPALA S M，SIMMONS A J，et al.The ERAInterim reanalysis：Configuration and performance of the data assimilation system［J］.Quarterly Journal of the Royal Meteorological Society，2011，137（656）：553-597.
［11］ 赵碧，丁静，刘建强，等.海洋白帽的高空间分辨率光学遥感估算分析［J］.遥感学报，2023，27（1）：92-103.
［12］ MEINDL E A，HAMILTON G D.Programs of the national data buoy center［J］.Bulletin of the American Meteorological Society，1992，73（7）：985-993.
［13］ CHEN K，XIE X，LIN M.An adaptive GaoFen-3 SAR wind field retrieval algorithm based on information entropy［J］.IEEE Access，2020，8：204494-204508.
［14］ 艾未华，陈冠宇，戈书睿，等.双频共面合成孔径雷达海面风场反演［J］.遥感学报，2019，23（3）：408-417.
［15］ LIU G H，YANG X F，LI X F，et al.A systematic comparison of the effect of polarization ratio models on sea surface wind retrieval from C-band synthetic aperture radar［J］.IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing，2013，6（3）：1100-1108.
引文格式：陈克海，张金兰，解学通，等.基于梯度算子的高分三号SAR 风场反演方法［J］.北京测绘，2023，37（9）：1218-1222.
［基金项目］ 国家自然科学基金（41876204）；广东工贸职业技术学院科技项目（2023-ZKT-03）
［作者简介］陈克海（1979—），男，广东汕头人，博士，讲师，研究方向为水色遥感、SAR和散射计数据处理。
E-mail: 15404247@qq.com
［通信作者］ 解学通，E-mail：xtxie2013@163.com
转自：“测绘学术资讯”微信公众号
如有侵权，请联系本站删除！
