本文发表于《Remote Sensing of Environment》期刊第250卷，文章编号112047，于2020年8月14日在线发表。该研究由山东科技大学测绘与空间信息学院马跃、武汉大学电子信息学院李松、山东科技大学海洋科学与工程学院刘振、武汉大学测绘遥感信息工程国家重点实验室杨必胜、山东科技大学测绘与空间信息学院阳凡林、新南威尔士大学科学学院王效华等学者合作完成。研究得到了国家科技重大专项、国家自然科学基金、国家重点研发计划、中国博士后科学基金等项目的支持。
本文属于类型a，即单项原始研究论文。
该研究聚焦于卫星遥感水深反演领域。水深（bathymetry）数据对于海洋、海岸带生态系统及相关研究至关重要。传统浅水水深测量主要依赖船载单波束或多波束回声测深仪（echo sounder）以及机载水深激光雷达（airborne bathymetric lidar），但这些方法在偏远或敏感区域往往难以实施。卫星多光谱影像已被广泛用于基于经验模型（empirical model）和物理模型（physics-based model）的浅水水深估计，然而经验模型通常需要实测水深数据作为先验信息，这限制了其在偏远礁区等无实测数据区域的应用。ICESat-2（Ice, Cloud, and land Elevation Satellite-2）卫星搭载了全球首台星载光子计数激光雷达ATLAS（Advanced Topographic Laser Altimeter System），其绿色激光可穿透水柱，在非常清洁的水体中最大探测深度约40米，能够提供沿轨分布的精确水深点。本研究旨在利用ICESat-2水深点替代实地辅助水深点来训练经典经验模型，结合Sentinel-2多光谱影像，实现仅凭卫星遥感数据即可生成浅水水深图。
研究的工作流程主要包含三个步骤。第一步是ICESat-2信号光子检测。ICESat-2的ATL03原始数据包含大量噪声光子，尤其在白天太阳背景噪声很强。官方提供的置信度（confidence）参数在海底信号光子检测中表现不佳，会丢失许多海底信号光子并保留大量噪声光子。为此，研究提出了一种改进的基于密度的噪声应用空间聚类（DBSCAN）算法来从ATL03原始数据光子中检测信号光子。该算法将每连续10000个原始光子作为一个沿轨片段进行计算。由于ICESat-2在海洋区域的垂直范围仅为60米，因此不进行垂直分段，而是直接利用全部光子计算给定半径ra内的信号与噪声期望光子数sn1。然后利用垂直范围内最低5米层的光子计算噪声期望光子数sn2。根据白天和夜晚太阳背景噪声的显著差异，白天邻域半径ra设为1.5米，夜晚设为2.5米。自适应阈值minpts由公式minpts = (sn1 - sn2)/ln(sn1/sn2)计算，且最小值为3。利用该DBSCAN方法，从每个沿轨片段中提取海面和海底信号光子。
第二步是海底光子的水深误差改正。ICESat-2的ATL03数据已改正了大气延迟、固体潮和系统指向偏差等误差，但未改正水深相关误差，包括水柱折射效应、水面折射效应和水面波动效应。研究首先从检测出的海面光子中计算局部平均海面lm和均方根（RMS）波高，将高程低于局部平均海面减去3倍RMS波高的光子识别为海底光子。然后针对每个激光脉冲，利用射线追踪方法计算相应的沿轨水深。先计算未改正水柱折射的初步激光距离rp，再通过rc = rp·n1/n2 ≈ 0.75rp进行水柱折射改正，其中n1为大气折射率（约等于1），n2为水的折射率（约等于1.334）。接着考虑水面波动影响，将当前水位偏差δd（等于局部平均海面与当前海面之差）加到改正后的激光距离上，得到水深r = (n1/n2)rp + (lm - lc)。此外，研究还首次改正了由海面坡度引起的折射偏差。根据斯涅尔定律（Snell’s law），折射角θ2 = asin(n1·sinθ1)/n2，水平误差δx = sin(θ1 - θ2)·r，垂直误差δz = tan[(θ1 - θ2)/2]·sin(θ1 - θ2)·r。最后利用OTPS2潮汐模型将水深归算到平均海平面（MSL）。最终水深公式为d = {1 - tan[(θ1 - θ2)/2]·sin(θ1 - θ2)}·[(n1/n2)rp + (lm - lc)] - δht。
第三步是基于ICESat-2水深点的卫星水深反演。研究使用了两种经典经验模型：线性波段模型（linear band model）和波段比值模型（band ratio model）。线性波段模型利用对数变换后的多光谱波段反射率作为预测因子计算水深，波段比值模型则利用高吸收波段与低吸收波段的反射率比值与水深建立对数线性关系。Sentinel-2影像首先经过大气校正，使用C2RCC（Case 2 Regional Coast Colour）处理器，然后使用去耀斑（deglint）算子进行太阳耀斑校正。ICESat-2水深点被用作先验测量来训练模型参数。
研究选取了两个研究区：南海南沙群岛永乐环礁（Yongle Atoll）和巴哈马东南部的阿克林岛与长岛之间的潟湖。在每个研究区分别使用四期Sentinel-2影像和多个日期的ICESat-2数据。在永乐环礁附近，使用了2018年10月22日、2019年2月22日和2019年4月12日三个日期的ICESat-2轨迹数据以及2019年2月24日、2020年3月10日、2020年3月20日和2020年3月25日四期Sentinel-2影像。验证数据为2012年底机载Optech Aquarius激光雷达系统采集的甘泉岛附近实测水深。在阿克林岛与长岛附近的潟湖区域，训练ICESat-2数据来自2019年2月11日、2019年3月12日和2019年9月2日，验证ICESat-2数据来自2018年11月12日和2019年6月3日，四期Sentinel-2影像分别获取于2019年2月6日、2019年10月14日、2019年12月23日和2020年1月27日。
研究结果显示，改进的DBSCAN方法在白天和夜晚均能较好地检测海底信号光子，明显优于ATL03官方结果。在永乐环礁附近，2018年10月22日白天采样的数据中，官方结果未能检测到水深超过5米的大部分水下信号光子，而本研究方法检测到的海底光子与光学影像反映的地形吻合良好。在22米水深的阿克林岛与长岛附近潟湖，检测到的海底光子最大深度约25米。经过误差改正后，各研究区的局部海面RMS波高在0.14至0.17米之间，水面平均坡度约为8°至11°。水平误差约为0.04倍水深，垂直误差约为0.0008倍水深，对于20米水深水平误差约80厘米、垂直误差约1.6厘米，可忽略不计。
在模型训练方面，永乐环礁附近线性波段模型的四期影像r²为0.95至0.97，均值为0.96，训练RMSE为0.96至1.23米，均值为1.08米；波段比值模型的r²为0.93至0.95，均值为0.94，训练RMSE为1.17至1.41米，均值为1.30米。在阿克林岛与长岛附近潟湖，线性波段模型的平均r²为0.93，训练RMSE均值为0.63米；波段比值模型的平均r²为0.93，训练RMSE均值为0.52米。在验证方面，永乐环礁附近甘泉岛区域线性波段模型四期影像的平均r²为0.91，平均RMSE为1.44米；波段比值模型的平均r²为0.85，平均RMSE为1.85米。在阿克林岛与长岛附近潟湖，线性波段模型四期影像的平均r²为0.90，平均RMSE为1.18米；波段比值模型的平均r²为0.89，平均RMSE为1.24米。
在空间一致性分析方面，将四期水深图的均值与单期水深图进行比较。永乐环礁附近线性波段模型四期空间差异分布的r²均值分别为0.97和0.95，RMSE均值分别为0.77米和0.87米。阿克林岛与长岛附近潟湖线性波段模型四期空间差异分布的r²均值为0.98，RMSE均值为0.35米；波段比值模型的r²均值为0.97，RMSE均值为0.48米。
研究结论表明，利用ICESat-2水深点替代实地辅助水深点训练经验模型是可行的，水深反演精度良好，两个模型在两个研究区的RMSE均低于或接近最大水深的10%。多期数据的一致性好，平均r²为0.97。线性波段模型在精度和一致性方面略优于波段比值模型，可能是因为线性模型使用了三个Sentinel-2波段，而波段比值模型只使用了两个波段。阿克林岛与长岛附近潟湖的反演结果略优于永乐环礁，原因可能是永乐环礁的机载验证数据采集于2012年底，与Sentinel-2影像获取时间存在时间差，以及两区域水质条件不同。
该研究的主要创新点在于：第一，提出了一种改进的DBSCAN方法从ICESat-2原始噪声数据中有效提取水深信号光子，性能优于官方置信度参数；第二，首次分析并改正了ICESat-2数据中由水面波动效应和水面坡度折射效应引起的水深误差；第三，利用ICESat-2水深点作为先验测量，使基于卫星多光谱影像的经验模型不再依赖于当地实地先验测量，这对于偏远礁区、岛屿周边和内陆水域等难以获取实测数据的区域具有重要意义。未来，借助免费开放的ICESat-2数据和Sentinel-2影像，该方法可扩展到更大尺度，实现全球范围沿海、岛礁周边和内陆水域浅水水深图的获取。