朱星源，苏洁,，宋梅，杨茜,，梁韵,
(1.中国海洋大学 海洋与大气学院，山东 青岛 266100；
2.中国海洋大学 物理海洋教育部重点实验室，山东 青岛 266100；
3.中国科学院西北生态环境资源研究院 冰冻圈科学国家重点实验室，甘肃 兰州 730000；
4.中国科学院南海海洋研究所热带海洋环境国家重点实验室，广东 广州 510301)
渤海位于37°～41°N 之间，由辽东湾、渤海湾、莱州湾和中央浅海盆地组成。三面被陆地包围，属于半封闭海区，与外界海水的热量交换较少，在冬季西伯利亚高压、太平洋副热带高压等因子的影响下，结冰现象显著[1−4]。入冬之后，随着负积温的累积和冷空气的侵袭，特别是强寒潮的暴发和延续，海冰面积会不断扩大，厚度也会随之增加[5−8]。海冰会对海上的交通运输和生产活动产生影响，甚至会引发自然灾害，如1969 年渤海大冰封，持续的冰封不仅导致海上航行受阻，“海二井”石油平台还直接被海冰推倒[9]；
2009–2010 年1 月中下旬的渤海冰情为近40 年同期最重，冰情等级为4 级，对环渤海地区产生了严重影响，直接经济损失达到了63.18 亿元[10]；
近10 年渤海最重冰情在2012–2013 年冬季，冰情等级为3.5 级，水产养殖受灾面积超过2.292×104hm2，直接经济损失达3.22 亿元[11]。为了了解渤海海冰的变化趋势，研究海冰灾害的成因并进行预测，对渤海海冰的准实时监测是非常必要的。
渤海海区范围较小，微波数据的分辨率较低，不足以进行有效的监测，采用更高分辨率的可见光数据进行海冰反演不失为一种合理有效的方法。渤海海冰密集度可见光数据反演的研究已有不少[12–16]，但作为海上冰情重要指标的海冰厚度，由于反演参数的获取较为困难，进展相对缓慢。Grenfell[17]提出的海冰厚度与反照率呈现指数关系的计算公式是目前较常用的可见光冰厚反演算法；
谢锋等[18]将此算法应用于AVHRR（Advanced Very High Resolution Radiometer）数据，对辽东湾海域冰厚进行反演；
Yuan 等[19]根据不同区域海冰光谱特性将渤海划分为5 个区域，也采用此算法基于AVHRR 数据分区域进行了渤海海冰厚度的计算；
Su 和Wang[20]基于Grenfell[17]和谢锋等[18]的结论，使用MODIS（Moderate-resolution Imaging Spectroradiometer）数据得出渤海冰厚计算值，但并未进行与实测数据的定性分析验证；
Liu 等[21]则将该算法应用于GOCI（Geostationary Ocean Color Imager）卫星数据，反演渤海海冰厚度。在这类反照率与海冰厚度指数关系模型的应用中发现，反照率不仅与海冰厚度有关，也受到海水中泥沙悬浮物的影响[22]，然而渤海海区受黄河泥沙输入影响，再加上渤海海区内水动力环境相当复杂，导致泥沙悬浮物浓度在时间和空间尺度上变化剧烈[23]，这也影响到了反照率与海冰厚度指数关系算法的准确性[18–19,24]。除海冰厚度与反照率指数关系的算法外，吴龙涛等[14]简单地使用MODIS 不同波段的分段线性法计算冰厚，其结果与石油平台获得的冰厚相比偏厚；
Ning 等[25]通过海冰在MODIS 和TM（Thematic Mapper）卫星不同波段的光谱特性对海冰进行分类，进而大致估算海冰厚度；
Yuan 等[26–27]从海冰的光辐射传输过程出发，提出了一种光学遥感海冰厚度半经验模型，该模型考虑到了其中各参数的时空异质性，并利用MODIS 数据多通道的反射率反演海冰厚度，但是存在一定局限性，容易低估高反射率区域的海冰厚度。
冰水分离是反演海冰厚度的重要步骤，近年来人们做了大量的研究。吴奎桥等[28]结合了MODIS 通道1 反射率和通道32 亮温进行冰水识别；
Su 和Wang[20]使用双通道比值阈值判别法识别海冰；
Su 等[29]利用灰度共生矩阵纹理分析方法进行渤海海冰探测；
Zhang 等[30]介绍了一种分类与回归树的方法检测海冰；
Su 等[31]提取了海冰的温度特征和纹理特征，然后使用支持向量机进行海冰范围估计；
由于渤海是我国的内海，黄河口附近海区注入大量泥沙、水体较为浑浊，这导致目前常用的方法大多存在悬浮泥沙误判的问题；
通过海冰纹理识别冰区的方法虽然能排除高浓度泥沙海区，但是对于平整薄冰的检测能力不足。Li和Yang[32]提出了一种基于多种海冰特征的线性分解法用于冰水分离，在提取海冰形状特征时，创造性地使用了Canny 算子对海冰的边缘和裂缝进行识别来区分泥沙与海冰[33]，为提取海冰范围提供了一种新思路，遗憾的是该方法没有充分完善Canny 算子的提取结果，不同水色海区交界区域会被误判。另外，通过线性分解模型提取海冰范围也存在不易提取平整薄冰等问题。
综上所述，这些研究为获得渤海海冰厚度可见光遥感数据奠定了一定的基础，但仍存在一些问题：在冰水分离方面，Canny 算子识别海冰边缘和裂缝作为一种提取海冰范围新思路，使用时仍需解决不同水色海区交界区域被误判和不易提取裂缝较少的平整薄冰区等问题；
在冰厚计算方面，泥沙悬浮物提高了海水反照率，影响了反照率与海冰厚度指数关系算法的准确性。为此，本文将针对Canny 算子提取海冰时存在的问题进行改进和完善，构造出一套基于Canny算子的完整的自动化海冰范围提取方法；
而对于海冰厚度计算，为了降低泥沙影响和提高计算结果准确性，本文将基于渤海海区的物理特征，通过试验确定海冰厚度与反照率指数关系算法中的相关参数，实现对算法的改进，并且将最终的反演结果与平台实测数据进行对比和分析。
2.1 MODIS 遥感数据
海冰厚度反演使用的数据为1 000 m 分辨率的MODIS L1B 数据和MODIS 03 地理信息数据，由美国宇航局提供并且已经过几何校正。中分辨率成像光谱仪MODIS 搭载在TERRA 和AQUA 卫星上，可以覆盖从可见光到近红外（0.405～14.385 μm）的36 个通道，最高空间分辨率可达250 m[34]，可以在每日上午、下午各获得一次渤海海区观测数据，具有免费获取、时空分辨率和光谱分辨率高等特点，使其可以有效地进行渤海海冰的监测。
本文主要使用的MODIS 通道和相对应的光谱范围如表1 所示，其中波段1、3 和4 的辐射率用于获取真彩图和灰度图，波段1、6 的反射率用于云剔除，波段31、32 的亮温用于冰水判别，波段1、2、3、4、5 和7 的反射率用于计算海冰厚度。数据预处理包括辐射定标和太阳高度角订正。为了避免天气原因造成的影响，本文选取的MODIS 数据均在晴空条件下。
表1 MODIS 部分波段光谱范围Table 1 Spectral range of MODIS partial bands
2.2 渤海海上石油平台实测数据
本文采用的实测数据为渤海海上石油平台的观测数据，包括2013–2014 年冬季和2015–2016 年冬季的JZ20-2 平台（40.500°N,121.352°E），2020–2021 年冬季的JZ20-2 和JZ9-3 平台（40.664°N,121.462°E）冰厚观测数据，观测时间为当地时间早上8 点至下午7 点，每小时进行一次海冰最大厚度和平均厚度观测；
2009–2010 年冬季和2012–2013 年冬季数据引自Zeng等[24]和Karvonen 等[35]的文章，来源于渤海海上石油平台JX1-1（39.977 5°N,121.058°E）、JZ25-1S（40.243°N,121.025°E）、JZ20-2 和JZ9-3 平台的海冰厚度现场观测数据，测量方法为目视观测，测量时间为当地时间下午6 点。平台具体位置如图1 所示。本文在使用实测数据时，将全部实测数据随机分为训练数据集和测试数据集，分别占总数据量的63%和37%，训练数据集用于海冰厚度算法的改进，测试数据集用于验证和分析算法改进效果。
图1 获取实测数据的石油平台位置Fig.1 Positions of oil platform observing the measured data
本文首先基于MODIS 数据中云的光谱特征进行云剔除；
继而使用Canny 算子提取海冰裂缝和边缘，再通过进一步处理实现冰水分离；
在计算海冰厚度时对海冰厚度与反照率指数关系模型中参数设置进行了分析和优化，使其更加符合渤海海区的水文特征，其中，将海水反照率参数由固定值变为随海域实际情况而改变的动态数值，利用渤海海上石油平台实测数据获取了更为准确的海冰衰减系数估计值。具体反演流程及改进情况如图2 所示。
图2 渤海海冰厚度反演算法流程Fig.2 Flow chart of Bohai Sea ice thickness retrieval algorithm
3.1 云剔除
本文选取渤海冰区晴空下的MODIS 数据进行反演，而水区上空的云在冰水分离阶段容易被误判为冰，还会影响算法中需要的海水反照率的计算，因此，首先要对渤海海区进行云剔除。云和海冰在可见光波段都有较高的反射率，而海冰在近红外波段反射率显著下降，因此本文采用的云剔除方法是：构建band 1和 b and 6反 射率r的归一化指数R1_6[36]，
对于式（1）结果的频数分布曲线，Bayes 分类准则认为将云峰和水冰峰之间的谷值作为阈值[37]，误判概率最小，剔除低于阈值的像元。
3.2 冰水判别方法的改进
3.2.1 海冰范围提取原理
冰区有不规则裂缝和边缘，而水区表面平滑，这是冰区和水区之间表面特征的重要差异，Canny 边缘检测算法可以提取冰区的边缘以及内部裂缝[32]，其重要优势是双阈值检测，可以有效地减少渤海清、浊水的边缘误判。本文利用Canny 算子的优势，针对引言中分析的算法缺陷，加入一系列处理步骤，包括真彩图灰度化、高斯模糊、空洞填充、灰度图二值化和表面温度（ST）自动化阈值判别法，提出一个可以自动化进行冰水判别并提取渤海海冰范围的新方案。
以2010 年1 月22 日的MODIS L1B 渤海海区数据为例（图3），具体分析冰水分离过程。
图3 2010 年1 月22 日渤海MODIS 真彩图Fig.3 MODIS true color image of the Bohai Sea on January 22,2010
首先将真彩图进行灰度化（图4a），使其可以进行基于单色图像的Canny 算子裂缝和边缘提取，提取效果如图4b 所示。由图可见，绝大多数的海冰被密集的裂缝覆盖，但是仍有部分平整薄冰裂缝较为稀疏，如果只是提取裂缝密集的区域，则会遗漏部分平整薄冰，因此需要分别采用不同方式提取裂缝密集冰和平整薄冰；
同时，部分裂缝较为密集的区域位于清水浊水交界处或者海冰边缘外侧，不属于冰区；
更重要的一点是图4b 由线状的裂缝组成，无法对裂缝密集的区域直接提取，综上所述，进行海冰提取还需要解决3 个问题：（1）通过线变面将密集的裂缝形成面状区域，使裂缝密集区作为区域得以直接提取；
（2）通过缺失冰区填充得以提取裂缝较为稀疏的薄冰区；
（3）去除非冰区的裂缝密集区。
3.2.2 提取裂缝密集区
第一步，将裂缝密集区作为区域进行整体提取。Li 和Yang[32]通过计算以每个像元为中心的圆形区域内裂缝的密度来给每个像元赋值，从而将线图变为面图，但是这种方法最终给出的图像较为模糊，不易区分裂缝密集区的边界。本文采用的方法是对图像进行高斯模糊，高斯模糊的优势在于增大了距离中心像元较近像元的权重，中心像元的权重最高，因此裂缝密集区边界更为清晰。对图4b 进行高斯模糊处理后，得到图4c，可以看到裂缝密集的区域亮度较高，再对图4c 进行二值化处理即可提取裂缝密集区域（图4d）。
3.2.3 提取裂缝稀疏薄冰区
第二步，提取裂缝较为稀疏的薄冰区。由于Canny算子对海冰边缘更为敏感，导致容易忽视海冰内部半封闭或封闭的平整薄冰，如图4d 的红圈区域。为了避免出现这种情况，我们具体做法为进行一次膨胀算法，将半封闭缺失薄冰变为封闭缺失薄冰，二值图空洞填充后进行一次腐蚀算法，即可完成对内部缺失薄冰的再提取，结果如图4e 所示。
3.2.4 误判像元去除
第三步，去除误判为冰的水像元。由图4e 可见提取到的区域存在水像元并且边缘过于光滑，不利于3.3 节改进冰厚计算方法时所需的海冰外缘线的判断，所以要对误提取的部分进行修正。图4f 为基于图4e 的裂缝区灰度图，由图可见水像元的灰度值明显较低，可以使用最大类间法进行低灰度像元去除，而最大类间法在对像元进行分类时，倾向于在类内方差较大的那一类像元（海冰像元）中选择阈值，为了避免造成海冰范围损失，采用约束灰度范围的改进型最大类间方法[38]。结果如图4g 所示，消除了部分水像元而且冰区的边缘更加符合对海冰外缘线的人工视觉判断。
3.2.5 表面温度自动化阈值判别法
解决了以上3 个问题后，冰区（图4g 白色区域）仍包括了一些泥沙浓度较高的浊水或清–浊水交界处的像元，将其去除才能完成海冰范围提取。海冰和海水在温度上的差异可以作为一个重要的判据[32,39]，可以有效地修正冰水分离时泥沙悬浮物造成的海冰误判，使用MODIS 数据31、32 波段亮温计算出图4g 中白色区域的表面温度（ST）如图4h[40]，由于海冰ST 往往低于海水ST，因此可以通过设置合适的阈值，剔除ST 高于阈值的部分。
关于阈值的取值，本文参考强度比算法的思路，提出了基于ST 频数比例的提取方法。强度比算法最早用来处理可见光航拍图像，以每个像元与相邻像元灰度值的差值大于临界值作为频数统计的判断依据[41]，但是此法难以适用于渤海这种情况较为复杂的大范围区域。为了实现对阈值进行较为准确自动化选取，本文的算法如下：
以2010 年1 月22 日的数据为例，（1）对整个渤海海区的ST 值域进行频数统计，记为 ϵ(k)，k为份数，间隔为0.02 K。（2）对图4h 的识别为冰区像元的ST 进行频数统计（图5a），记为 δ(k)。（3）最后计算二者的比值χ(k)=δ(k)/ϵ(k)，将 χ (k)称为粗糙度，画出粗糙度曲线如图5b 所示。由于海冰在之前的一系列提取流程中被保留下来，表现为低温区的粗糙度较高而接近1，而随着ST 升高达到某个值时粗糙度迅速下降，粗糙度迅速下降的区间就是阈值所在的区间，这里，我们取粗糙度为0.4 时对应的ST 作为阈值，所对应的阈值为272.2 K，该阈值随渤海ST 分布的变化而变化。将ST 小于阈值的像元剔除，最终的冰水分离结果如图4i 所示，白色部分即为最终提取到的海冰范围。
图4 2010 年1 月22 日渤海冰水分离过程Fig.4 Process diagram of ice-water separation of the Bohai Sea on January 22,2010
图5 频数分布（a）和粗糙度分布（b）Fig.5 Frequency distribution (a) and roughness index distribution (b)
3.3 海冰厚度算法及改进
Grenfell[17]和 Allison 等[42]在极地调查中发现，当海冰厚度从2 cm 增加到9 cm 时，海冰反照率由0.11 增加到0.24。Grenfell[17]提出反照率与海冰厚度呈指数关系模型为
式中，h为海冰厚度，单位为m；
αmax为无限厚冰反照率，一般取0.7；
系数为海水反照率；
µα为海冰的反照率衰减系数，谢锋等[18]的估计值为1.209。
Su 和Wang[20]使用上述模型基于MODIS 数据在渤海海区进行海冰厚度计算。将式（2）中的海水反照率 αsea设 为0.06。计算反照率 α(h)时使用Liang[43]提出的反照率 α经验公式，在MODIS 的6 个波段之间建立 了函数关系为
计算海冰厚度时一般将海水反照率 αsea设为常数[18,20–21]，但渤海的泥沙悬浮物浓度在时间和空间尺度上变化剧烈[23]，将海水反照率 αsea设为固定数值的做法与实际情况明显不符，也会影响冰厚计算结果的准确性，因此有必要在不同的时间、地点根据实际情况设置不同海水反照率。而冰下的海水由于被海冰覆盖，无法进行直接观测。为了解决这个问题，比较可行的办法是将与海冰相邻的无冰海区的海水反照率外推到有冰区。下面仍以2010 年1 月22 日渤海数据为例具体分析外推过程。
图6a 为当日的渤海海区反照率分布，红线为基于形态学提取的海冰外缘线[44]。为了排除外缘线以外可能存在的冰水混合像元误判，经试验，本文将海冰外缘线进一步外推4 个像元得到图6a 中的蓝线，将蓝线所在的宽度为3 个像元的狭长条带区域定义为最邻近海冰海水区，根据空间相关性理论，距离近的事物关联更紧密[45]，理论上这个区域内的海水反照率比较接近海冰区内的海水反照率。将最邻近海冰海水区像元作为插值节点，冰区像元作为被插值点，使用反距离加权插值法对反照率进行插值，结果如图6b 所示，由图可见插值后的海水反照率分布基本符合渤海泥沙分布特点。在本研究中使用上述算法将海水反照率 αsea由固定值0.06 变为随海域实际情况而改变的范围在0.05～0.13 之间的动态数值。
图6 渤海海区反照率（a）和海水反照率插值（b）Fig.6 The albedo of Bohai Sea (a) and the albedo of sea water interpolation (b)
谢锋等[18]在使用AVHRR 数据反演海冰厚度时，对衰减系数 µα的估算方法为：将式（2）变形为式（4），海冰厚度h使用结（融）冰度日法通过实测气象数据估算[46]，α(h)由卫星观测给出，αsea和 αmax分别取0.1 和0.7，然后反推出 µα的值，因此这个方法并没有采用冰厚实测数据参与计算，并且结（融）冰度日法本身存在误差，同时他们也指出卫星数据来源不同，定标方式不同，µα的取值也不同，再加上本文已经对算法中的αsea做 了优化，不再为常数，所以对 µα进行重新计算是十分有必要的。
为了得到相对准确的衰减系数 µα的估计值，本文随机抽取部分冰厚实测数据h作为训练数据集（具体占比见2.2 节），将对应日期的MODIS 相关波段数据代入式（3）得到整个海区的反照率。再通过最邻近海冰海水区反照率插值算法，获得实测数据所在位置冰下海水反照率 αsea。然后将冰厚实测数据h，海水反照率 αsea以 及海冰反照率 α(h)代入式（4），最终反推出 µα。
由式（4）可知，衰减系数 µα由比值运算推出，当海冰实测厚度h较小时，分子分母数值均较小，即使是微小的偏差都会导致最终结果产生较大误差，最终反推出的衰减系数 µα也更容易发散。如图7 所示，实测厚度较小的薄冰衰减系数 µα较为分散，随冰厚增大 µα分布越来越集中。实际上，海冰的衰减系数主要由海冰类型，海冰内部卤水泡、叶绿素等物质的吸收系数及其体积分数决定[17,22]，但是通过卫星无法获取具体的相关参数，所以在实际应用中一般将衰减系数 µα设为固定值。为了得到衰减系数 µα相对准确的确定值，需要对反推出的 µα进行进一步处理。目前图7 中所有数据点衰减系数平均值为2.19，为了减小薄冰（实测冰厚小于6 cm）带来的较大偏差，仅让厚冰的衰减系数（图7 中黑点）参与计算，算得其平均值为1.85（图7中蓝色实线），标准差为0.58，然后进行质量控制，仅对处于平均值±标准差以内的数值进行平均（图7 中两条红线之间），最终算得 µα的估计值为1.74（图7 中蓝色虚线）。
图7 训练数据集实测冰厚与衰减系数散点图Fig.7 Scatter plot between attenuation coefficient and measured sea ice thickness of training data sets
3.4 冰厚计算敏感性试验
本文通过敏感性试验研究指数关系模型中各参数对反演结果产生的影响。图8 反映了海冰厚度反演结果对海冰反照率 αice、海水反照率 αsea和衰减系数µα的 敏感性。如图所示，当海冰反照率 αice与海水反照率 αsea相等时，海冰厚度反演结果为0，随着海冰反照率 αice的增加，冰厚反演结果呈现指数上升趋势。在海冰反照率 αice不变的情况下，越低的海水反照率αsea对应越大的冰厚反演结果。那么在泥沙浓度较高的水区（如莱州湾，渤海湾），其实际的海水反照率αsea较高（0.09～0.13），如将其设为以往的固定值0.06，反演的海冰厚度往往偏大，甚至无冰海区，其反演结果也可以达到5 cm 左右。而本文使用的根据实际情况设置海水反照率的做法，会显著减小高浓度泥沙区海冰厚度的反演结果，对于反照率为0.15 的海冰，计算冰厚时将海水反照率 αsea由0.06 增大到0.1，会使反演结果由8.7 cm 减小到5 cm（µα设为1.74）。另外，本文利用实测数据和MODIS 数据重新估算了渤海海区海冰的反照率衰减系数，取常数为1.74，与之前使用的1.209 相比[18]，冰厚反演结果整体缩小30.5%。
图8 模型中各参数对反演结果的影响Fig.8 The influence of each parameter in the model on the retrieval result
仍以2010 年1 月22 日为例，图3 为当日渤海海区MODIS 真彩图。使用本文提出的冰水分离算法提取海冰范围，然后使用指数关系模型改进前与改进后算法，以及Yuan 等[26]的半经验模型法进行海冰厚度计算，最终结果如图9 所示。从冰水分离效果来看，本文所提取海冰范围与当日真彩图（图3）的海冰范围一致，外缘线也与真彩图的目测外缘线基本吻合。从海冰厚度反演结果来看，改进后的反演结果（图9b）与之前相比有所减小，在泥沙悬浮物浓度高的区域尤为明显，比如渤海湾和莱州湾大面积几乎透明的薄冰（图3），改进前算法（图9a）对这些薄冰的反演结果多处于0.1～0.15 m 之间，而改进后算法结果大多在0～0.05 m 的区间内。再将算法改进前后结果与Yuan等[26]算法结果（图9c）进行对比分析，空间平均绝对偏差分别为0.045 m 和 0.027 m。
图9 渤海海冰厚度反演结果对比Fig.9 Comparison of retrieval results of Bohai Sea sea ice thickness
由于与Yuan 等[26]半经验模型法的相似性高只代表与同类型反演算法的比较情况，无法作为误差的检验标准。因此本文利用渤海海上石油平台实测数据来检验算法改进前后的结果；
同时为了检验本文对反演算法中海水反照率 αsea与衰减系数 µα的改进效果，分别改变其中一个参数而另一个参数保持不变，将所得结果与实测数据进行比较，4 个试验的结果见表2。需要指出的是在进行误差分析时使用的数据集为测试数据集，与3.3 节反推衰减系数 µα时使用的训练数据集不同。
表2 海冰厚度测试数据集实测数据与反演结果Table 2 Measured data of test data sets and retrieval results of Bohai Sea sea ice thickness
续表 2
如表2 所示，T0（改进前算法厚度）与平均实测厚度之间的平均误差、平均绝对误差和均方根误差分别为6.66 cm、7.05 cm 和8.25 cm；
T1（改进后算法厚度）各项误差分别为0.49 cm、2.74 cm 和3.75 cm，与T0 相比分别降低了93%、61%和67%；
其中仅改变海水反照率 αsea的T2 平均误差、平均绝对误差和均方根误差分别降低45%、33%和30%；
仅改变衰减系数µα时的T3 误差降低的幅度比T2 更大，各项误差分别降低61%、47%和43%。4 个试验结果与平均实测厚度的相关系数分别为0.434、0.485、0.485 和0.434，其中T0 与T3，T1 与T2 相关系数相同，这是由于衰减系数 µα与海冰厚度反演结果之间呈反比例关系，因此仅衰减系数 µα不同不会影响相关性。我们同时也将算法结果与最大实测海冰厚度进行了比较，结果显示T0 高于最大厚度，T1 低于最大厚度，且相关性更强（0.480），值得注意的是T2 各项误差均低于T0 和T1，相关系数与T1 一致。考虑到过境时间，上午星TERRA过境时间为当地时间10:30 前后，下午星AQUA 过境时间为当地时间13:30 前后，此时海冰厚度比较接近一天之中的平均厚度，因此在与实测数据进行比较时，仍以平均实测冰厚为准。为了进行对比，本文还计算了Yuan 等[26]半经验模型法与实测数据的误差，结果显示与平均实测厚度的各项误差低于T0 高于T1，但与最大实测厚度最为接近。
最后，选取2020–2021 年冬季渤海石油平台的观测数据，与反演结果进行比较。由图10 可见，总体来看JZ20-2 平台实测冰厚和反演结果相较于JZ9-3 平台波动起伏较大，这是由于JZ20-2 平台位置更加靠近海冰边缘（图11），极易受到海冰范围增加和缩小的影响；
改进前的算法结果与平均实测数据相比偏大，有时甚至高于最大实测厚度（如2021 年1 月8 日JZ20-2 平台），而改进后算法的整体厚度为3 种方法中最小，除1 月5 日、13 日低于JZ20-2 平台平均实测厚度，1 月8 日高于最大实测厚度以外，多介于平均厚度和最大厚度之间。接下来，通过真彩图和海冰厚度反演结果（图11）对以上提到的3 个特殊时间节点进行具体分析：由图11a 可见，1 月5 日JZ20-2 平台周边的海冰比较松散，海冰的漂移会相对更加明显，平台观测时间与卫星观测时间的不一致可能是导致反演有偏差的原因；
另外，密集度较小的海区海冰可能偏薄，而本反演算法对厚度小于6 cm 海冰的反演结果尚不是十分准确。1 月13 日真彩图显示JZ20-2 平台周边基本没有海冰（图11d），海冰边缘在平台以北，这与冰厚反演结果一致但与平台实测不同（图11h），其原因也与平台观测与卫星观测不同步有关。1 月8 日的JZ20-2 平台3 种算法反演结果一致大于最大实测冰厚，由图11c 可见海冰表面明显偏白，可能有积雪覆盖，同样情况的还有1 月7 日真彩图（图11b），积雪会明显提高海冰反照率进而影响反演结果，导致JZ9-3 平台1 月7 日和8 日冰厚反演结果均显著高于降雪之前的1 月5 日反演结果。但由于缺少石油平台降雪量数据支持，通过真彩图对降雪的判断是一种猜测。
图10 2020–2021 年冬季反演结果与实测数据的时间序列比较（JZ20-2 平台和JZ9-3 平台）Fig.10 Time series comparison between measured data and retrieval results in 2020–2021 winter (JZ9-3 and JZ20-2)
本文利用MODIS 遥感数据实现了渤海海冰范围的自动化提取，并由渤海的实际情况出发，针对海冰厚度指数关系模型中海水反照率 αsea与 衰减系数 µα参数的设置进行了反演算法的改进，并将改进前后算法结果及前人算法结果与实测数据进行检验和比对，主要结论如下：
（1）基于Canny 边缘检测算子，加入真彩图灰度化、高斯模糊、空洞填充、灰度图二值化和ST 自动化阈值判别法等步骤，提出了一种渤海海冰的自动化提取算法。该算法在一定程度上克服了Canny 边缘检测算子对平整薄冰区检测能力不足的弊端，可以有效提取高浓度泥沙海区的薄冰，并且通过应用最大类间法和粗糙度法实现了冰水分离过程中的阈值自动化提取。
（2）基于渤海的水文特征，对反照率与海冰厚度指数关系模型中的参数重新设置，实现对算法的改进。选取了海冰外缘线外侧的狭长条带内海水像元作为插值节点，向海冰像元点所在位置进行插值，获得冰覆盖区海水反照率 αsea，进而将海冰厚度指数关系模型中的海水反照率 αsea由固定值变为随海域实际情况而改变的动态数值。使用海上石油平台的冰厚实测数据和MODIS 数据重新反推衰减系数 µα，为了减小偏差，反推过程中排除了实测厚度6 cm 以下的薄冰，确定其估计值为1.74。
（3）利用海上石油平台实测数据中的测试数据集进行验证，结果显示：改进后的海冰厚度指数关系模型的反演结果与实测数据间误差明显缩小，平均绝对误差由7.05 cm 缩小到2.74 cm，相关系数由0.434 提高到0.485；
但反演结果较最大实测厚度偏薄，多介于平均厚度与最大厚度之间。
需要说明的是，计算冰下海水反照率 αsea时，由于本文使用的插值算法为反距离加权插值法，理论上最近的插值节点与被插值点之间距离太远会增大插值结果与真实值之间的差异，最终冰厚反演结果也会受到影响，尤其是重冰期远离海水的海冰，插值误差造成的影响更加明显，因此还需要进一步考虑如何提高αsea取值的合理性。本文利用海冰平均厚度实测数据对衰减系数 µα进行了重新反推计算，所得衰减系数µα为1.74，本文也尝试采用了最大实测厚度进行反推，所得衰减系数 µα为1.18，代入算法中使得反演结果高于前者47.5%，但由于篇幅所限，未在正文中展开讨论；
计算时排除了6 cm 以下的薄冰，因此薄冰反演结果的准确性有所下降，并且最终结果不是通过理论推导所得，而是经验性的估计值，实际上海冰的衰减系数并不是固定不变的，而是上下波动的，如何通过卫星判断不同的海冰类型，进而获得更加准确的衰减系数 µα值，是未来继续优化反演算法的关键之一。
海冰表面覆盖积雪会显著提高海冰反照率，导致冰厚反演结果偏大，而且积雪厚度越大，对结果造成的影响越大。在未来的研究中，需要先区分冰上积雪和厚冰，再进一步对冰厚反演算法进行修正。同时也建议平台及岸站观测部门在监测海冰实况时，增加对降雪情况的记录。
猜你喜欢 衰减系数海冰渤海 渤海大学作品精选声屏世界(2022年17期)2022-12-18末次盛冰期以来巴伦支海-喀拉海古海洋环境及海冰研究进展海洋通报(2021年3期)2021-08-14渤海竞渡帆高举——记渤海轮渡集团党委书记、总经理于新建水上消防(2020年3期)2020-07-25结合时间因子的校园论坛用户影响力分析方法研究智能计算机与应用(2017年5期)2017-11-08落水洞直径对岩溶泉流量影响的试验研究南水北调与水利科技(2017年4期)2017-07-29基于SIFT-SVM的北冰洋海冰识别研究电子技术与软件工程(2016年24期)2017-02-23辽代“斡鲁朵”内的渤海人移民通化师范学院学报(2016年11期)2017-01-15HT250材料超声探伤中的衰减性探究中国测试(2016年3期)2016-10-17累积海冰密集度及其在认识北极海冰快速变化的作用中国学术期刊文摘(2016年8期)2016-02-13基于TerraSAR-X全极化数据的北极地区海冰信息提取自然资源遥感(2014年3期)2014-02-27