水深是海洋信息的重要组成部分之一, 其中浅海水深是海洋地形、海上航行和浅海自然资源开发利用最基础同时也是最重要的海洋空间地理信息之一。人工探测、船载声呐以及激光点云等传统手段对水深探测步骤复杂并且消耗大量人力、物力、财力。而遥感水深反演可以大尺度、迅速、低成本^[1]获取水深数据, 其中光学遥感受大气条件、水体杂质等影响较大, 微波遥感则不受这些因素干扰, 因此微波遥感探测水下地形在海洋科学中具有重要的意义。

合成孔径雷达(Synthetic Aperture Radar, SAR)^[2], 又称主动式微波成像雷达, 依靠电磁散射进行距离探测和目标信息提取。自20世纪50年代SAR技术出现后, 在科研、经济、军事领域发挥了重要作用。据现有理论体系, 可将SAR水下地形成像过程概括为: 海流与水下地形相互调制改变表层流速, 然后风致表面短波与海面表层流相互调制改变表层短波空间分布, 最后雷达波和表层短波相互调制改变归一化雷达截面积(Normalized Radar Cross Section, NRCS)^[3]。

SAR浅海水下地形遥感研究可以追溯到1969年, 第一颗海洋雷达卫星Seasat首次向人们展示了SAR浅海水下地形探测的巨大潜力^[4]。1978年, 海洋卫星发射并传回了大量包含浅海水下特征的雷达影像, 自此, 将SAR用于浅海地形研究的潜质引起了广泛的关注^[5]。美军在1982年开展的雷达海上观测研究中, 对浅海地形的成像理论进行初步的探究^[6-7]。Alpers和Hennings提出的SAR浅海水下地形的成像理论模型为SAR浅海水下地形定量化研究迈出了重要的一步^[8]。杨劲松等^[9]对文氏海浪谱采取完全非线性积分法对SAR图像通过仿真计算, 得到了合成孔径雷达在浅海信息反演中的最佳雷达参数、海况的等信息。夏长水等^[10]利用SAR仿真模型对渤海塘沽海区的水深进行反演。杨俊钢^[11]结合海浪波高频谱解析表达式和成像机理, 设计了多源多时相的浅海地形探测方法, 通过最小迭代函数反演台湾岛浅滩不同地形下的水深。范开国等改进传统仿真模型的缺陷, 基于MS4程序和成像机理建立新的仿真模型, 对风速风向与图像的特征关系提出了新的认识^[12]。黄韦艮等^[13]基于浅海水下地形SAR成像机理提出的水深探测新技术对台湾岛浅滩海域进行了遥感探测实例研究。赵现斌等^[14]提出了一种全极化的合成孔径雷达反演参数设计方法, 对辐射分辨率和信噪比等技术参数的最佳取值范围进行研究。郑全安等^[15]介绍了SAR在探测浅海地形的研究进展、发展前沿和几个具有代表性的雷达探测地形的实例, 为研究更复杂的地形打下了基础。Bian等^[16]提出了一种基于涌浪参数和水深之间的线性色散关系的多时相合成孔径雷达水深探测方法。

传统的SAR水深遥感反演模型中所需的参数过多, 且初始精度不高, 往往需要加入迭代模型进行至少10多次的迭代运算才能完成精度提升^[17], 而机器学习法可以简化反演的过程同时提高结果的精度。其中梯度提升决策树模型在遥感的监测中远优于常见的几种模型, 将其用于水深反演保持高精度的同时可以降低算法的迭代次数提高反演效率。

梯度提升决策树(Gradient Boosting Decision Tree, GBDT), 属于决策树模型的其中一个分支, 其核心思想是把残差学习作为基础, 优化梯度方向形成大量回归树, 所有回归树结果累加得最终模型, 通过逐次迭代来减少残差同时获取较高精度的结果, 适用于处理各种非线性的数据。由Friedman^[18]提出的基于Boosting的梯度提升算法, 已经被广泛运用于各种领域。Jeany等^[19]将GBDT用于跟踪算法中, 获得了较高的精度。Xu等^[20]提出一种基于多层梯度提升树的高光谱遥感图像的分类方法。陈训来等^[21]基于GBDT模型用新的气象雷达的三维拼图和观测雨量资料建立降水预测模型。同年, 拱健婷等^[22]提出一种新型快速无损药材鉴定方法, 将GBDT模型结合近红外光谱对药材的伪品进行判别。金秀等^[23]发现GBDT算法可以集成多种模型的优势。Wei等^[24]将高光谱影像用于河道水污染检测, 发现随机森林法虽然和梯度提升树有相似的精度, 但是计算更加复杂。李怡静等^[25]把高分一号影像和实测水质数据结合GBDT算法搭建了遥感水质参数反演模型, 证明该方法可以实现高精度的遥感信息监测。孟然等^[26]将GBDT模型用于水深反演, 发现其精度远高于已有的水深反演模型。同年, Yao等^[27]把Landsat图像和地面实测数据与GBDT模型结合用于北部湾沿岸的叶绿素a的遥感监测中, 探究叶绿素a的变化规律和时空分布情况。

综上所述, 将GBDT算法与遥感探测相结合具有明显的优势。本文提出一种SAR图像水深反演方法, 在GBDT中加入NRCS、基础水深数据、风场数据、流场数据建立水深反演模型, 来计算未知区域的水深。全文结构如下: 第一部分介绍研究区域的特点与具体数据, 第二部分介绍水深反演过程使用的方法和各项模型。第三部分详细分析反演的结果并对其进行精度评价。第四部分对所得结果进行讨论, 并分析影响结果的因素。

1 研究区域与数据 1.1 研究区

本文选取的研究区域为杭州湾和长江口南缘相连的浅海区域(图 1), 杭州湾整体形态呈现喇叭状, 其底部地形从湾口到乍浦地势较为平坦, 乍浦以后的地势按0.1‰~2‰坡度变化向东下降, 水深变化大, 由于独特的地形特征使得这里经常出现暴涨潮和涌潮。其北岸为长江三角洲的南部, 沿岸滩地分布较广, 海岸线曲折绵长, 有利于海上贸易和交通运输, 现已凭借其独特的先天优势成为了优良的港口。以上原因使得西面水域较为浑浊, 呈泥浆色, 东面水域较为清亮, 呈浅蓝色, 同时当前水域时常受大风影响, 使得光学遥感水深反演模型在该区域难以发挥作用。

1.2 数据 1.2.1 Sentinel-1遥感卫星数据

本文使用的哨兵1号(Sentinel-1)卫星影像数据来自地球数据网站(https://search.asf.alaska.edu), 选取2018年1月至2018年11月共计23景数据, 极化方式为VV、VH, 空间分辨率为10 m×10 m, 时间分辨率为12 d。大量的研究表明, 由于海表的特殊性, 同极化的回波信号强于交叉极化, 对海洋的研究通常采用VV极化。

1.2.2 全球水深数据

水深数据来自美国国家海洋大气局(https://www.ngdc.noaa.gov)网站中的ETOPO1数据, 其空间分辨率为1′。该数据结合了海洋水深和陆地地形, 其元数据建立在部分地区和全球数据集的基础上, 也可用于冰面和基岩的研究中, 将该数据加入水深遥感反演模型前需要根据遥感数据的生成时间进行潮汐校正, 便于反演结果更接近真实数据, 提高结果的可信度。

1.2.3 风场和流场数据

海表风场数据采用遥感系统网站(https://remss.com)提供的2018年1月17日至2018年12月19日的CCMP风矢量产品。CCMP是由微波传感器的海表风反演和再分析的组合, 使用该方法创造的海表风场数据可以保证研究区域内的风场平滑过渡, 其时间分辨率为6小时, 空间分辨率为0.25°×0.25°。

海表流场数据由HYCOM网站(https://www.hycom.org)提供的2018年1月17日到2018年11月13日存储在netCDF中的数据, 共记录了40个不同深度层次的海流数据, 空间分辨率为0.08°×0.08°。研究表明, 流速对于合成孔径雷达观测水下地形是至关重要的, 流速越大越易观测地形, 当流向与坡度方向平行时观测效果最佳。

2 方法与反演模型 2.1 图像去噪

SAR图像中固有噪声往往会最大程度接近或超过要提取的信息, 按分辨率的损失情况把噪声的处理方法分为两种: 一种是保留原始分辨率的情况下, 利用滤波算法抑制噪声, 另一种是以损失分辨率为前提对数据进行多视操作抑制斑点噪声^[28]。为保留SAR数据高分辨率的特点, 本文选用中值滤波平滑图像去除斑点噪声。根据SAR的光谱图, 利用IDL (Interactive Data Language)设置不同窗口的中值滤波对图像去噪(将窗口为3×3的去噪图像简称为滤波3), 根据表 1中图像与水深数据的相关系数, 判定最佳去噪方式。

在探索最佳滤波窗口时, 对未去噪和已去噪的数据随机提取浅海的小块区域样本, 建立机器学习的水深反演模型。不同数据反演相关系数从高到低依次为: 滤波9、滤波3、滤波5、原始、滤波7、滤波11。最高为82.25%, 最低为76.00%。图 2可以看出, 从原始数据到滤波7的精度变化不大, 经过滤波3、滤波5, 精度相较原始图像有所提升, 滤波7精度相较原始图像有所下降, 但是整体上精度变化波动较小, 未见明显的提升, 到滤波9精度提升较大(如图 3所示), 从原始的76.81%提升到82.25%, 滤波11精度又明显下降, 即反演精度并不会随着滤波窗口的增加而增加, 前期会在一个小范围内提升, 当到达最佳滤波窗口时, 精度提升较为明显, 此后继续扩大滤波窗口, 精度便会下降。

2.2 SDWI水体指数模型

雷达水体指数提取法^[29](Sentinel-1 Dual-Polarized Water Index, SDWI), 是遥感数据水陆分割常用的阈值法, 即对图像采取一定预处理算法之后, 假设K为图像像元值, T为阈值。当K > T, 该像元为水体; 当K < T, 为非水体, 对图像进行二值化, 将非水体部分的信息设为0, 水体部分设为1, 从而精确提取水体信息。SDWI算法的基本原理为找出同极化和交叉极化数据之间的初始关系, 将两种极化影像相乘之后再扩大10倍来增加水陆像元之间的差异, 在水体信息部分不断拟合、筛选、对比, 最终获取研究区最优的水体信息分割方法, 公式(1)如下:

$ {K_\text{SDWI}} = \ln\left( {10 \cdot {K_\text{VV}} \cdot {K_\text{VH}}} \right) - 8 , $

(1)

式中, K_SDWI、K_VV和K_VH分别为雷达水体指数计算结果、同极化和交叉极化, 这三张图像的像元值。

2.3 梯度提升决策树水深反演模型

科学家对SAR浅海水深反演进行了多年的研究, 提出的反演模型大致可以分为两类, 一是基于AH模型改进后的算法, 二是基于Navier-Stokes方程的算法。前者需要加入大量实测海洋物理数据, 但结果还需要在迭代模型中进行循环运算, 才能满足一般的精度要求。后者需要解算大量复杂的方程。本文通过对大量SAR图像处理与分析, 利用梯度提升决策树算法对包含初步水深信息的SAR图像、风场数据、流场数据进行机器学习训练, 从而计算未知区域的水深(如图 4所示)。

梯度提升决策树是一种基于决策树模型的集成算法, 可以在极少的调参时间下获得较高的精度, 不管是离散型还是连续型数据都可以灵活处理。在数据预测和分析中使用效果较好, 包含二叉分类树和二叉回归树。通过模型迭代运算, 将分散的弱水深学习器提升组合成一个强水深学习器, 提高水深反演精度, 每次迭代关系密切, 需要上一次模型计算的结果输入到下一次的模型中, 该算法中需要选取最优的迭代次数和最大特征数, 避免在水深精度评价中出现过拟合的状况。

首先建立弱水深学习器, 如公式(2)所示:

$ {f_0}(x) = \arg \min \sum\limits_{i = 1}^N {L({y_i},c)} . $

(2)

对已有水深样本$ i{\text{ = 1,2,}}...{\text{,}}N{\text{,}} $计算残差, 如公式(3)所示:

$ {r_{im}} = - {\left[ {\frac{{\partial L({y_i},f({x_i}))}}{{\partial f({x_i})}}} \right]_{f(x)}} = {f_{m - 1}}(x) . $

(3)

将上一个步骤计算得到的残差(负梯度)作为新的真实值, 把数据$ \left( {{x_i},{r_{im}}} \right),i = 1,2,...,N $当作下一决策树的训练数据从而得到新的回归树$ {f_m}(x) $, 对应的叶子节点区域为$ {R_{jm}},j = 1,2,...,J $。其中回归树的子叶节点个数为$ J $.

子叶区域的拟合值为, 如公式(4)所示:

$ {\gamma _{jm}} = \underbrace {\arg \min }_\gamma \sum\limits_{{x_i} \in {R_{jm}}} {L({y_i},{f_{m - 1}}({x_i}) + \gamma )} . $

(4)

对公式更新得到强水深学习器, 如公式(5)所示:

$ {f_m}(x) = {f_{m - 1}}(x) + \sum\limits_{j = 1}^J {{\gamma _{jm}}I(x \in {R_{jm}})} . $

(5)

其最终生成的水深学习器, 如公式(6)所示:

$ f(x) = {f_M}(x) = {f_0}(x) + \sum\limits_{m = 1}^M {\sum\limits_{j = 1}^J {{\gamma _{jm}}I(x \in {R_{jm}})} } . $

(6)

2.4 评价指标

为更加准确描述反演结果的精确性, 本文采用相关系数(Correlation Coefficient, CORR, 记为C_corr)、平均绝对误差(Mean Absolute Deviation, MAD, 记为D_ma)以及均方根误差(Root Mean Square Error, RMSE, 记为E_rms)来检验水深反演的精度, 表达式分别为公式(7)~公式(9), 如下:

$ {C_{corr}} = \frac{{\sum\limits_{i = 1}^m {(x - {x_i})} }}{{\sum\limits_{i = 1}^m {(x - \overline {{x_i}} )} }} , $

(7)

$ {D_{ma}} = \frac{{\sum\limits_{i = 1}^m {|({x_i} - x)|} }}{m} , $

(8)

$ {E_{rms}} = \sqrt {\frac{{\sum\limits_{i = 1}^m {|{{(x - {x_i})}^2}|} }}{m}} . $

(9)

3 结果分析与精度评价 3.1 最佳风速

研究表明, SAR图像上条纹的明暗情况与风速具有较强的相关性^[30], 在低风速下SAR图像偏暗, 高风速下则偏亮^[12]。因此风场是雷达成像过程中至关重要的一个因素, 风速过低使得后向散射信号太弱无法被雷达系统探测, 风速过高则会出现白沫效应和波浪破碎, 使得NRCS出现饱和现象, 直接掩盖水下地形对信号的影响^[31]。本文在选择水深反演的最佳风速时, 参考盛琳等^[32]在利用非线性模型对礁盘的水深反演中使用的逐步回归法, 该方法的优势是可以迅速、运算次数最少同时最大程度保持变量的特性。分析反演初始精度与风速的关系, 可以快速选取最佳风速。

水深反演的相关系数与风速的关系如图 5所示。两者之间并非简单的线性关系, 反演的精度随风速的变化一开始是先增加后减小; 随后会出现一个精度变化滞缓区间, 在该区间内, 精度并不会随着风速的增加有明显的变化, 一旦风速超过该区间, 精度就继续下降。当风速小于3.78 m/s, 时, 精度随着风速的增加而增加, 当风速在3.78~6 m/s时, 精度随着风速的增加而减小, 当风速在6~7.5 m/s时, 精度几乎不变, 当风速大于7.5 m/s时, 精度急速下降, 最佳风速约为3.78 m/s。

3.2 最佳迭代次数

为研究GBDT算法在水深反演中的精度, 预先反复试验, 确定算法中各项核心参数。对雷达、风场、流场数据随机采样, 选定内符合精度进行精度验证, 将采样点一部分作为训练集加入到GBDT算法中建立水深反演模型, 另一部分作为样本集用于精度检验。样本集所得的水深数据与实际水深进行对比, 利用相关系数、平均绝对误差、均方根误差作为精度评定标准, 并将结果绘制统计图表, 如表 2所示。

GBDT模型中最重要的参数之一就是迭代次数(迭代1次简称迭代1, 以此类推)。选取最优迭代次数不仅可以减小算法的运行时间, 还可以提高算法的精度。将样本集逐次增加迭代次数, 记录3种评价指标并进行分析。由图 6可知, 经过第1次迭代后, 样本集的大部分反演结果都非常接近甚至等于原始水深, 但仍然存在个别数据点的反演值与实际值差距较大, 经过第2次、第3次迭代后, 偏差较大的水深点仍然存在, 第4次迭代后, 偏差值较大的水深点全部消失, 第5次迭代的散点图和第4次的散点图未见明显变化, 且各项精度统计数据也几乎相等。

由图 7可知, 第1次迭代后, 相关系数和平均绝对误差都取得了较好的精度, 但是均方根误差达到了3 m; 第2次迭代后, 相关系数有所提升, 但变化幅度小, 平均绝对误差明显减小, 但均方根误差却增大, 远超3 m; 第3次迭代后, 相关系数继续提升, 平均绝对误差和均方根误差均下降, 第4次迭代和第5次迭代后, 3项指标的数据几乎一致。综上所述, 综合考虑算法运行的时间和计算机的硬件情况, 最佳迭代次数为4。

3.3 分段水深反演精度分析

为探究反演结果中误差分布情况, 将水深划分为0~10 m、10~20 m、20~30 m、30~40 m、40~50 m, 表 3记录了3种评价指标, 并将结果分水深制作散点图(图 8)和折线统计图(图 9)。

相关系数、均方根误差和平均绝对误差在0~10 m依次为0.83、0.95 m、1.53 m; 在10~20 m相关系数提升了0.07, 均方根误差和平均绝对误差分别增加0.35 m、0.19 m; 在20~30 m, 除相关系数轻微下降, 另外两种指标的精度明显增加, 均方根误差和平均绝对误差分别下降0.55 m、0.76 m, 即当前水深区间的反演精度明显优于0~10 m和10~20 m的水深; 在30~40 m, 相关系数下降0.03, 均方根误差和平均绝对误差分别增加0.5 m、0.7 m, 与20~30 m相比, 精度下降较为明显; 在40~50 m, 相关系数减少0.11, 均方根误差和平均绝对误差大幅增长。

综上所述, 在分段水深中, 相关系数随水深段的增加先增大后减小, 从高到低依次为10~20 m、20~30 m、30~40 m、0~10 m、40~50 m, 最高为0.9, 最低为0.73。在40~50 m的深水区域相关系数最低, 是因为流场与水下地形相互作用反映的水深信息在雷达影像上的表现力有限, 随着水深的增加, 雷达影像获取的水深信息逐渐减少, 因此在较深水域相关系数下降趋势较为明显。平均绝对误差随水深的增加先增加后减小再增加曲折变化, 从低到高依次为20~30 m、0~10 m、30~40 m、10~20 m、40~50 m, 其中误差最低为0.75 m, 最高为1.89 m。均方根误差随水深的变化情况和平均绝对误差一致, 从低到高依次为20~30 m、0~10 m、30~40 m、10~20 m、40~50 m, 其中误差最低为0.96 m, 最高为2.26 m。从3种精度指标来看, 各水深段的反演精度都比较高, 不管是水深探测深度还是精度都远大于光学遥感水深反演模型。

图 9中相关系数在0~10 m水深区间出现精度反常低值, 猜测可能是由波浪破碎导致。作为近海表面的一种常见现象, 波浪破碎是一种强烈的波与波相互作用现象多出现在岸边, 海浪由于海水底部的阻尼作用使得波峰前倾卷倒, 破碎时有大量空气被卷入水中, 致使波峰出现白浪, 海表产生泡沫、空中溅起大量水滴, 次表层海水形成气泡云。波浪破碎在一定程度上会影响着雷达信号的强度大小从而对获取的NRCS产生干扰, 最终影响水深反演的精度。

3.4 总体水深反演精度分析

为了探究不同深度区间内的反演精度, 将水深值按照0~10 m、0~20 m、0~30 m、0~40 m、0~50 m进行划分, 表 4中记录了3种评价指标, 并将结果分水深制作散点图(图 10)和折线统计图(图 11)。

当水深区间第一次扩大后, 相关系数会增加, 虽然均方根误差和平均绝对误差的数值也随之增加, 但该区间的精度仍然较好。在0~30 m, 相关系数较前一个分段增加了0.1, 均方根误差和平均绝对误差分别下降了0.01 m、0.02 m, 即该区域的精度轻微上升; 在0~40 m, 反演精度和0~30 m水深区间一致; 在0~50 m, 除了相关系数增加, 另外两项指标的精度都明显下降, 并且此时的均方根误差已接近1.6 m, 即反演精度较差。

综上所述, 随着水深增加, 反演的相关系数也会提升, 但是深水区域的实际值和反演值的平均绝对误差和均方根误差都会增加。到达50 m左右水域的均方根误差已经接近1.6 m, 反演值与实际数值的偏差较大, 即雷达数据提取深水域地形信息的能力较弱, 不可以仅根据相关系数的增加而盲目扩大反演水深区间, 此次反演的最佳水深区间为0~40 m。

3.5 反演水深空间分布特征

由于研究区内水深像元点的数据量达到千万级别, 在将反演水深可视化时对计算机的配置具有较高的要求, 可利用金字塔降采样算法降低反演的数据量。而风场和流场数据都存在一定的缺失情况, 为提高反演的精度, 在提取训练集和样本集时都排除了缺失数据所在的区域, 使得反演的图像会存在空洞, 利用插值算法填补缺损像元得到反演水深图像图 12(b)。将反演水深与实际水深图 12(a)进行空间分布特征对比分析。

该区域的水深从杭州湾海岸线开始由浅及深阶梯增加, 浅水区占据大量海域面积。从120°45′E到122°15′E都是深度为10 m左右的浅水区域, 几乎占据整片海域的1/3, 同时近岸浅水域内还存在个别零星分布的深水域。122°15′E到122°40′E内的区域水深是10~20 m, 裹挟了几个零星的浅水域, 几乎包围了整个舟山群岛。该区域向东是水深为20~30 m的区域, 从122°50′E向东继续前进, 经度每增加几分, 水深就会增加10 m。

两幅图像的水深由近岸向深水的延伸趋势几乎完全一致, 都是梯度增加。图 12(a)和图 12(b)在0~10 m和10~20 m水深的分界线明显, 边缘过渡清晰; 10~20 m与20~30 m水深边界略模糊, 但仍可看见较清晰的分界线; 20~30 m与40~50 m的水深界线逐渐模糊, 需要仔细分辨, 再往东水深界线逐渐融合在一起, 肉眼难以分辨。反演水深图能较为完整地展现出研究区内的地形特点。图 12(a)在0~10 m水深海域内的30°45′N附近有一小块明显的深水区域, 东西方向狭长分布, 而反演水深图中该特征也被提取出来, 其分布的范围和形状相似, 只是中间夹杂着一些浅水域像元点。在0~10 m水深的最南面靠近海岸线的区域, 也存在一小块深水区域, 在图 12(b)中也体现了该特征。由此可见, 反演结果在0~40 m水深区间内, 能够较好的反应研究区内的实际水深分布情况, 比较符合当前区域的水下地形特征。

4 讨论 4.1 研究结果与前人结论的对比

本文以杭州湾和长江口南缘相连的浅海区域为例, 探究最佳风速、最佳迭代次数和不同水深段的精度。

通过参考前人所用非线性模型中的逐步回归法得到研究区内反演的最佳风速约为3.78 m/s。此次最佳风速研究结果和前人结论对比见表 5。有三组研究人员采用M4S系统仿真研究SAR浅海地形成像过程中的最佳风速。研究人员和最佳风速分别是: 黄韦艮等^[34], 3~9 m/s; 范开国等^[12], 4~6 m/s; 范开国等^[35], 4 m/s。即使是采取相同的方法, 研究区不同, 所得的最佳风速区间也不完全一样, 但是都在中等风速区间内, 并且较为接近。而Calkoen等^[33]使用BAS水深估算系统获取的最佳风速区间为3~5 m/s。和M4S仿真系统相比, 研究所得的最佳风速仍大致相同。本文的最佳风速为3.78 m/s, 风速结果符合黄韦艮等^[34]、Calkoen等^[33]所得的结论, 并且与范开国等^[35]人所得的最佳风速区间只相差0.22 m/s。

综上所述, 在采用相同的模型估算最佳风速时, 根据研究区域的不同, 所得的结论仍会存在小范围的差异, 但是整体结果都非常接近, 使用不同模型估算最佳风速时, 结果虽然不一样, 但也是十分接近。由此可知, 本文所得的结论具有较高的准确性和较好的适用性。

通过表 6可以看出利用GBDT模型训练不同数据到达最佳精度时的迭代次数差距较大, 少的仅为几次, 多的能达到上千次。计算机在大多数情况下进行迭代运算时消耗的时间和硬件配置的要求都是成倍增长的, 在保证数据达到提取要求的前提下, 迭代次数越少越好。本文使用的NRCS、风场、流场数据和水深建立的水深反演模型仅仅需要进行4次迭代就能达到较好的反演效果。

现有的遥感水深反演方法主要分为两类: 光学、微波。光学遥感的优势是可以直接穿透水层到达水底并且返回能量, 通过计算穿透水层消耗能量与水深之间的关系, 从而得出研究区的水深, 并且现有的光学遥感水深反演理论模型较多^[49], 算法较为成熟, 在众多学者的研究改进下, 其精度也在不断提升, 但是仍存在着无法忽视的重大缺陷: 1)受天气和大气条件影响较大。阴雨天气下, 研究区内数据含云量较高, 目前的云处理技术一般是利用高通滤波弱化一部分干扰或是直接掩膜去除云区域。无论是哪种方法, 都会对水深反演造成严重的影响。2)对研究区内水质情况有严格的要求。如表 7所示, 光学数据的研究区多选取在水质条件好、人类活动少的地方。3)水深探测能力较弱。从探测深度一栏可知, 反演的大部分水深区间处于10~25 m。

然而利用合成孔径雷达反演水深恰好可以弥补这些不足。并且SAR数据在水深反演的潜力仍在探索中。本文最佳探测深度为40 m, 远超光学遥感的探测能力。由精度一列可知: 王小珍^[45]、田庆久等^[41]的反演结果的相关性均低于80%, 而本文的相关性高达97%。表 7中其余研究人员在20 m水深处的均方根误差最低为0.95 m, 最高为1.8 m。而本文反演所得的40 m水深的均方根误差仅为1.56 m。综上所述, 利用GBDT模型进行水深反演能保持较强的水深探测能力的同时也具有较高的反演精度。期待在将来的研究中可以改进模型, 进一步提高水深探测能力。

4.2 潮汐数据对结果的影响

虽然SAR数据具有分辨率高的特点, 但正因如此利用该数据参与反演的数据量远大于低分辨率的遥感数据, 当研究区选择过大就会使得运算时间过长, 因此本文的研究区仅选择了杭州湾附近的浅滩。

基础水深数据是由重力卫星反演后的数据, 在不考虑潮汐影响下获得的水深, 若利用该数据与实际水深进行对比, 需要查找离研究区最近的潮汐检测站数据。中国浅海沿岸设有许多潮汐观测站, 因此同一个研究区附近可能有多个检测站。目前对海洋数据的潮汐校正方法按数据量分类, 可以分为单一站点校正和多站点校正。单一站点校正常用于小研究区, 选取离研究区最近或者位于研究区正中心的站点进行校正, 多站点校正常用于大尺度海洋遥感中, 即选取多个潮汐点, 不同的算法进行矫正。如崔丹丹等^[50]利用潮汐分带技术校正遥感图像。

本文数据来自大鱼潮汐网(https://www.chaoxibiao.net), 该网站数据仅记录每个站点的整点数据, 而遥感图像的生成时间不能与此完全对应。观测数据如图 13所示, 整点潮汐之间的数学关系接近于一阶线性函数, 可据此得出遥感图像对应时间的潮汐数值。不同站点潮汐数据对应的潮高基准面和潮高数据不同, 需要计算两者的差值得出与海平面的距离用于潮汐校正。本文的研究区仅为杭州湾附近的浅海区域, 且选取的最深水域仅为50 m, 因此选用小区域的单一站点进行校正便可。

4.3 流场数据对水深结果的影响

虽然研究区内水深反演图和实际水深图相比已达到较高的相似性, 能较好展示该区域内的水下地形细节, 但是仍有个别区域出现了较为明显的误差, 如图 14所示。同一列图表示相同的水域, 第一行是实际水深图中的3处细节, 对应第二行反演水深图的3处细节。

流场数据对反演结果影响最大, 主要是以下两点: 1)流场数据只存在于水体区域, 而流场数据的分辨率远低于SAR图像, 因此两张图像在海陆边缘处的分界线存在一些差异, 使得0~10 m水域的反演结果较低。2)流场数据存在缺失像元, 而GBDT水深反演模型中该参数是必不可少的, 故在缺少流场数据的区域, 无法使用模型进行水深反演, 只能利用插值算法获得近似水深数据。第一点可以利用水体指数建立掩膜区间, 从而降低边缘差异对反演结果的影响。从图 14可以看出, 该方法处理后的边缘与实际水深边缘无明显差异。而第二点, 使反演结果图与实际水深图中出现明显的差异, 例如: (a1)和(a2)的右上角实际为50~70 m水深, 但是(a2)中表现为40~50 m水深; (b1)和(b2)的右下角在实际水深图中表现为边界圆润的椭圆状深水域, 而反演水深图中并未展现出这一细节; (c1)和(c2)图的左上居中区域的反演水深图中该区域有一小块三角状的浅水区, 但实际水深图中并无这一细节, 且在右下角实际水深图中明显存在60~ 70 m的水深, 而反演水深图中也未出现这一细节。由此可知, 流场数据缺失对反演结果的影响更大。

4.4 ETOPO1数据的精确性

本文采用ETOPO1水深数据构建模型, 该数据海洋部分是用卫星测高技术获得, 深海区域具有较高的准确性, 而部分浅海的水深则通过参考实测数据进行修改。据文献记载, 该数据在经度105~120°E水深小于200 m的沿海地区容易出现地形虚假起伏变化且广西沿海个别区域的数据也会出现异常平滑的状况^[51]。因此将数据应用于浅海区域进行水深反演时, 最好是能借助实测水深点对研究区的数据进行修正, 可以大大提高反演水深的准确性。本次研究并未找到该区域的实测数据, 希望之后的研究中可以增强数据的可信度和可用性。

5 结论

本次SAR浅海水深反演中, 根据后向散射模型对图像进行处理, 从分段和总体两个方面用相关系数、均方根误差和平均绝对误差对反演精度进行分析, 选取最佳风速、迭代次数和数据处理方式, 利用机器学习法进行水深反演。结果如下:

(1) 最佳风速约为3.78 m/s。当风速小于3.78 m/s, 时, 精度随着风速的增加而增加, 当风速在3.78~6 m/s时, 精度随着风速的增加而减小, 当风速在6~7.5 m/s时, 精度几乎不变, 超过该区间后, 精度直线下降。

(2) 最佳迭代次数为4。前三次迭代中, 三项精度指标均有所变化, 第四次迭代和第五次迭代后, 三项指标的数据几乎一致, 精度不再提升。

(3) 在分段水深中, 受近岸波浪破碎影响, 使得最浅水深段精度反常下降, 但各项数据均远大于光学遥感水深反演精度。

(4) 最佳反演的水深区间为0~40 m。相关系数会随着水深的增加而增加, 但是均方根误差和平均绝对误差在0~50 m水深范围有较高的误差, 不可根据相关系数的增加而扩大反演水深区间。

(5) 该区域水深从杭州湾海岸线开始由浅及深阶梯增加。在0~40 m水深的反演结果能够较好地反映研究区内的实际水深分布情况, 比较符合当前区域的水下地形特征。

本文还存在着一些不足之处, 未考虑波浪破碎对水深反演的干扰, 此后的研究中应该考虑建立NRCS与波浪破碎参数之间的关系, 更加精确地获取水深信息。