CN114355367A - 一种基于星载单光子激光雷达数据测量浅海水深的方法 - Google Patents

Abstract

一种基于星载单光子激光雷达数据测量浅海水深的方法，首先分析星载单光子激光雷达在不同环境下浅海区域的光子信号与噪声高程特点，构建光子点云的初始滤波模型；随后基于初滤波结果分析不同环境下光子点云的空间密度分布，构建水面和水下光子点云分离模型，实现水面有效光子信号和海浪波廓线提取；再根据不同环境的水下光子点云在不同水深的空间密度分布特点，构建随水深变化的自适应光子密度滤波模型与水深提取方法，分析星载单光子激光雷达的水深测量能力和精度；最后基于星载单光子激光雷达的水深探测结果与遥感反射率光谱之间的关系，使用星载多源遥感数据开展主‑被动融合的高精度区域化浅海水深反演，实现境内外高精度浅海水深数据提取。

Description

一种基于星载单光子激光雷达数据测量浅海水深的方法

技术领域

本发明涉及航天遥感测绘技术领域，具体是一种基于星载单光子激光雷达数据测量浅海水深的方法。

背景技术

我国有300多万平方公里的辽阔海域，大陆岸线1.8万公里，岛屿岸线1.4万公里，海岛1.1万余个，水深在50米以内的浅海海域面积达50万平方公里。目前我国获取全球高程数据主要依赖卫星遥感技术，由于高程精度制约，难以满足海图要素修编、海基精确武器远程打击、海礁军事设施建设规划等需求。鉴于手段有限、地形复杂等因素，目前我国针对海岛礁、滩涂、海岸线浅海地形测量效率比较低下，尤其一些偏远远洋区域无法进行探测。

单光子模式激光测高仪具有灵敏度高、重复频率高、测量精度高、重量轻、体积小等诸多优势，其代表着新一代天基激光雷达的未来发展趋势。2018年9月，美国国家航空航天局已经发射第一颗对地观测单光子多波束激光雷达卫星ICESat-2，ICESat-2作为ICESat-1的后续卫星其主要搭载了具有6波束的高级地形测量激光高度计系统(AdvancedTopographic Laser Altimeter System，ATLAS)，10kHz的激光脉冲重复频率，沿轨激光脚点间隔仅为0.7m；该载荷实现了ICESat-1没有实现的多波束推扫功能，主要目的是获取测量地球两极的冰盖形貌，首次实现了星载单光子激光雷达的全球对地观测。目前ATLAS单光子激光雷达数据已经向全球发布，并展示出较强的浅海水深探测潜力。

与传统的激光测深相比，单光子激光测高仪解决了能量和重复频率之间的矛盾问题，采用光子计数的单光子激光测高仪在设计和数据处理方法上有了较大改变，获取有效信号时，它不再专注于用高能量发射获取高信噪比的波形，而是采用概率探测的模式同时记录了噪声和信号光子，这也对数据处理方法提出了新的挑战，通过改进传统的数据处理方法，在低信噪比的光子云中也能做到有效信号的提取。

目前国内外关于单光子激光测高仪的点云数据去噪方法主要分为以下三类：基于栅格的点云数据去噪方法、基于局部统计信息的点云数据去噪方法和基于聚类的点云数据去噪方法。根据单光子计数测量实验和研究，发现光子计数激光雷达中噪声主要来至于背景光，光子计数激光雷达设备中死时间的长短，如何有效地消除噪声提取目标信号，是保障和提高光子计数激光测量精度的关键。美国罗切斯特理工大学卡尔森成像科学中心提出利用基于密度的空间聚类模型的滤波算法，能够有效滤除光子计数激光雷达点云数据中的噪声。美国麻省理工大学从探测概率优化等理论问题出发，提出一种基于栅格划分和设定计数阈值的目标信号提取算法。

目前我国没有发射相关的单光子激光雷达卫星，机载单光子激光雷达设备和相关数据也十分稀少，基于单光子激光雷达测量的相关研究主要采用国外开源数据，处理对象主要都集中在陆地和植被等领域，针对海岛礁和海岸线区域的浅海水深的研究更是凤毛麟角。除此之外，我国浅海区域大、分布广，很多浅海区域环境复杂。受湿度影响，裸露的浅滩、礁石、泥沙等表面反射率大致在3％左右，远低于传统的陆地激光雷达测绘时20％左右反射的典型应用条件和环境。基于低信噪比的光子计数激光测量技术作为一种高灵敏度的新型测绘体制，可以更好地适应浅海区域的探测条件和环境。

发明内容

本发明所要解决的技术问题是提供一种基于星载单光子激光雷达数据测量浅海水深的方法，可以有效解决上述背景技术中提出的问题。

为解决上述问题，本发明所采取的技术方案是：一种基于星载单光子激光雷达数据测量浅海水深的方法，包括以下步骤：

S1、获取ICESat-2卫星测得的原始数据，包括单光子激光数据、高分遥感数据和潮汐数据；其中在单光子激光数据中获取光子点云的密度结构数据，在高分遥感数据获取高分多光谱影像；

S2、单光子激光点云数据预处理：基于光子点云的密度结构数据进行预处理，通过距离选通技术和初始滤波算法确定目标点或区域的区域范围，减少大量的噪声光子点云数量；

S3、多光谱影像得太阳耀斑消除：采用线性回归模型建立可见光波段与近红外波段的对应关系，来消除可见光波段影像中的耀斑分量；

S4、水面于水下光子分离：将单光子激光雷达数据按高程值从上至下，以Δh为间隔将高程划分为M个切片，据每个切片单元统计光子数量，利用正态分布曲线对高程切片中心高程与切片单元中的光子信号数量之间的关系进行拟合；以切片中心高程为中心，获取每个水面光子的高程；

S5、自适应不同水深光子密度滤波：基于单光子激光雷达探测技术理论以及水体中光子信号的数据特点，确定水面与水下信号点的分界线，并自适应获取初始椭圆搜索区域的长轴、短轴以及密度滤波的初始最小点数阈值；分析随信号光子高程值的变化，高程切片内光子数量的变化特征，确定搜索椭圆的长轴、短轴以及密度滤波的最小点数阈值与水深之间的变化关系；

S6、遥感影像的水陆分离：利用近红外影像和边缘提取算法，开展水陆分离和岸线边缘提取，采用水陆分离模板用于识别陆地和水体，水陆边缘的识别采用估算摄影时刻前后一段时间内的平均海面高度；

S7、水深反演模型构建：通过搜索水深值相近的临近样本点作为解算数据集，从而将整个目标区域划分成多个小区域，每个小区域分别用不同水深样本点进行反演，获得更高的反演精度。

作为本发明的进一步优选方案，所述步骤S2确定目标区域中有效光子的大概高程范围，用高斯分布函数进行局部拟合，根据高斯峰值、幅宽度等参数为基础，结合先验高程阈值，得到高斯分布函数拟合关系式：

关系式中，G(h)表示高斯分布函数，A_maxpho表示波形振幅，σ为半波宽，h为当前点高程值，h_u表示探测获取最大极值所对应的光子高程值；从而确定目标区域中有效光子信号的初始范围。

作为本发明的进一步优选方案，所述步骤S2目标区域的有效光子S_signal由多个数据子区R_g(g＝1,L,M)中有效光子的总和构成，关系式为：

关系式中，R_g(g＝1,L,m)表示在沿轨飞行方向划分为m个分段，g为沿轨分段下标值，n表示在当前沿轨分段内光子点的总合数，i为当前沿轨分段内某点点号，

表示当前分块区域的有效光子数，

为当前分段的有效光子数。

作为本发明的进一步优选方案，所述步骤S4中以Δh为间隔将高程划分为M个切片得到关系式：

根据每个切片单元统计光子数量得到关系式：

关系式中，h_i为每个高程切片中心高程，x_i为切片单元中的光子信号数量。

作为本发明的进一步优选方案，所述步骤S4中正态分布曲线对高程切片中心高程h_i与切片单元中的光子信号数量x_i之间得到关系式：

关系式中，x_maxpho表示包含光子信号最多的高程切片中的光子数量，h_maxpho表示切片中心高程；以切片中心高程h_maxpho为中心，通过阈值确定水面光子高程分界线，并得到水面光子高程下限H_min和上限H_max。

作为本发明的进一步优选方案所述步骤S4中获取每个水面光子的高程，从而形成水面光子数据集SP和水下光子数据集WP，得到等式：

等式中，

表示任意点i的高程，i＝1,2,L n，n为光子总个数，H_min为水面光子高程下限，H_max为水面光子高程上限。

作为本发明的进一步优选方案，所述步骤S7基于局部水深约束的多波段水深反演模型，所述模型有效地将水深变化信息作为约束，融入到多波段水深反演模型，得到公式：

式中，Z表示对应位置水深值，N表示波段总数，L(λ_i)是传感器收到的i波段辐射能量(浅水)，L_∞(λ_i)为对应的深水区能量；a₀…a_N为经验回归参数，这些参数根据水深样本点与遥感辐射量进行回归分析解算得出；由于采用小范围的精确模型参数代替全局粗略模型参数，从而获得更高的反演精度。

与现有技术相比，本发明提供了一种基于星载单光子激光雷达数据测量浅海水深的方法，具备以下有益效果：

该发明能够更好地适应浅海区域的探测条件和环境；实现水面有效光子信号和海浪波廓线提取；分析星载单光子激光雷达的水深测量能力、精度和可靠性；使用星载多源遥感数据开展主-被动融合的高精度区域化浅海水深反演，实现境内外高精度浅海水深数据提取。

附图说明

图1为本发明方法流程图；

图2为光子点云预处理原理和数据切片示意图：

图3为不同数据切片中的高斯拟合与光子滤波示意图；

图4为基于线性回归的耀斑消除原理示意图；

图5为光子高程方向分块后光子数量统计分块中心光子高度示意图；

图6为随着光子高度变化的搜索椭圆示意图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。

参照图1，本发明提供一种基于星载单光子激光雷达数据测量浅海水深的方法，基本思路为：首先分析星载单光子激光雷达在不同环境下浅海区域的光子信号与噪声高程特点，构建光子点云的初始滤波模型；随后基于初滤波结果分析不同环境下光子点云的空间密度分布，构建水面和水下光子点云分离模型，实现水面有效光子信号和海浪波廓线提取；再根据不同环境的水下光子点云在不同水深的空间密度分布特点，构建随水深变化的自适应光子密度滤波模型与水深提取方法，分析星载单光子激光雷达的水深测量能力、精度和可靠性；最后基于星载单光子激光雷达的水深探测结果与遥感反射率光谱之间的关系，使用星载多源遥感数据开展主-被动融合的高精度区域化浅海水深反演，实现境内外高精度浅海水深数据提取。

作为本发明的一个具体实施例:

步骤一：获取ICESat-2卫星测得的原始数据，包括单光子激光数据、高分遥感数据和潮汐数据；其中在单光子激光数据中获取光子点云的密度结构数据，在高分遥感数据获取高分多光谱影像；

步骤二、单光子激光点云数据预处理：ICESat-2卫星搭载的单光子激光测高仪的距离窗为6km，且噪声点随机分布在整个距离窗，噪声点分布区域十分广泛，达到几千米。而一定时间内的目标点的分布区域一般只有几百米，甚至几十米。因此，需要进行光子点云的预处理，及通过距离选通技术和初始滤波算法确定目标点或区域的大概区域范围，减少大量的噪声光子点云数量，参照图2所示；

单光子激光测高仪的点云数据为剖面式分布，单光子激光测高仪的点云数据为剖面式分布，其沿飞行轨迹为线状分布，因而可以将三维点云数据转变为二维点云数据表示，即表示为离散点的水平距离值(沿轨距离值)或时间，以及垂直距离值(光子高程)；

根据目标区域中的卫星轨迹，在沿轨方向将光子数据用N个纵向切片进行切分，从而形成N+1个数据子区域，如图2所示。根据区域内的点云数据建立高程统计分布图，即在每个数据子区中将点云数据按照等高程区间划分为M个高程片，统计每个高程片单光子事件的个数。统计每个分区中的光子数量，对高程统计分布直方图，构建分布曲线；为了减少探测数据随机性对滤波结果的影响，可采用均值低通滤波器对直方图曲线进行平滑滤波；

确定目标区域中有效光子的大概高程范围，即探测分布曲线中的峰值点，用高斯分布函数进行局部拟合，根据高斯峰值、幅宽度等参数为基础，结合先验高程阈值，从而确定目标区域中有效光子信号的初始范围，如图3所示；高斯分布函数拟合用公式(1)表示，

式中，G(h)表示高斯分布函数，A_maxpho表示波形振幅，描述了矩形内的光子数量，σ为半波宽，h为当前点高程值，h_u表示探测获取最大极值所对应的光子高程值；

整个目标区的有效光子S_signal由多个数据子区R_g(g＝1,L,m)中有效光子的总和构成，如公式(2)所示；

式中，R_g(g＝1,L,m)表示在沿轨飞行方向划分为m个分段，g为沿轨分段下标值，n表示在当前沿轨分段内光子点的总合数，i为当前沿轨分段内某点点号，

表示当前纵向分段区域的有效光子数，

为当前横向(沿轨)分段的有效光子数。

步骤三，多光谱影像的太阳耀斑消除：遥感影像的辐射处理主要是对海洋遥感影像中存在的太阳耀斑进行提取和消除，以及影像的增强，从而获取清晰的水下纹理，提高海洋影像的匹配精度。太阳耀斑是影响海洋遥感图像辐射质量的主要因素，参考图4，由于水体对近红外波段(NIR)的强烈吸收，可认为NIR波段的辐射强度仅仅由大气散射和太阳耀斑组成，大气校正之后NIR波段残留的辐射强度即为太阳耀斑部分。本申请，拟采用线性回归模型建立可见光波段(VIS)与近红外波段(NIR)的对应关系，参考图4，来消除可见光波段影像中的耀斑分量。

步骤四，水面与水下光子分离：将单光子激光雷达数据按高程值从上至下或从下至上，以Δh为间隔将高程划分为M个切片如公式(2)，对每个切片单元统计光子数量；以每个高程切片中心高程h_i为横轴，切片单元中的光子信号数量x_i为纵轴，构建二者之间的关系，如式(3)所示；

利用正态分布曲线高程切片中心高程h_i与切片单元中的光子信号数量x_i之间的关系进行拟合，拟合关系记作f(x)。如式(4)所示，其中x_maxpho表示包含光子信号最多的高程切片中的光子数量，h_maxpho表示该切片的中心高程，参考图5所示；

以切片中心高程h_maxpho为中心，通过阈值确定水面光子高程分界线，并得到水面光子高程下限H_min和上限H_max；将水面高程上限下限值作为边界，获取每个水面光子的高程，形成水面光子数据集SP和水下光子数据集WP，如式(5)所示。

式中，

表示任意点i的高程，i＝1,2,L n，n为光子总个数。

步骤五，自适应不同水深光子密度滤波：基于单光子激光雷达探测技术理论以及水体中光子信号的数据特点，为克服光子计数数据处理存在的问题，本发明一种单光子激光雷达自适应阈值水深提取算法；该方法以椭圆为滤波搜索窗口，随着水深的增加，椭圆滤波窗口会自适应的进行椭圆窗口放大，从而对不同深度处的光子进行滤波处理，获取最后的有效光子信号数据，参考图6所示；

确定水面与水下信号点的分界线，并自适应获取初始椭圆搜索区域的长轴、短轴以及密度滤波的初始最小点数阈值；通过分析随着信号光子高程值的变化，高程切片内光子数量的变化特征，确定搜索椭圆的长轴、短轴以及密度滤波的最小点数阈值与水深之间的变化关系。基于椭圆的空间密度聚类滤波算法与采用圆形的空间滤波窗口，能更好对空间密度分布不均匀光子点云数据进行空间滤波，同时对密度变换的方向性具有较好的敏感性；

本发明提出了一种单光子激光雷达自适应阈值水深提取算法，使用的椭圆滤波窗口邻域的形状由两点p和q的距离函数的选择决定，用dist(p,q)表示；这里两个参数是点的eps邻域，由dist(p,q)≤eps定义，以及该eps邻域中的最小点数Min_pts；

以水面和水下分离后的光子集合为基础，在水面光子集合中获取高程方向上的最大差值ΔSH，在ICESAT-2沿轨方向上以光子计数激光雷达沿轨方向分辨率的长度Δd为分界，计算沿轨方向分辨率Δd与ΔSH的比值并获取均值Ratio_ab，该值即为搜索椭圆的长轴与短轴的比值，从而获取长轴长度a；将水面光子高程下限H_min和上限H_max差值的绝对值作为椭圆短轴b，用公式(6)描述，式中表示整个目标区域按沿轨方向用Δd划分的子数据集合数量；

在水面光子数据集SP中，遍历所有光子数据点，以a和b分别为搜索椭圆的长轴与短轴构建搜索区域，统计搜索区域内的点数p_v(v＝1,…,N_sp)，并计算点数均值，获取初始最小点数Min_pts；

在水下光子数据集WP中，椭圆搜索长轴a_i、短轴b_i、最小点数阈值min_pts与水面光子长轴、短轴以及最小点数阈值初始值以及光子高程相关，该过程描述于公式(8)；

所有数据进行基于椭圆滤波核的密度滤波，其中水面光子椭圆滤波核长轴、短轴以及最小点数阈值取值分别为a，b，min_pts，水下光子椭圆滤波核长轴、短轴以及最小点数阈值随着光子高度变化，其值分别为a_i、b_i、min_pts。通过该自适应滤波方法可以有效地对光子信号噪声较大、水面与水下光子信号分布不均匀的情况，自适应选择滤波阈值，实现自动、快速、高效的有效光子点云提取，自动进行水面与水下有效光子信号的分离并获取连续水深数据。

步骤六，遥感影像的水陆分离：利用近红外影像和边缘提取算法，开展水陆分离和岸线边缘提取，采用水陆分离模板用于识别陆地和水体，确保在后续的双介质浅海测量中仅对水体部分进行光线折射改正，而水陆边缘则主要用于估算摄影时刻前后一段时间内的平均海面高度。

步骤七、水深反演模型构建：利用光学影像反演浅海水深是遥感在海洋领域的重要应用之一，它根据辐射传输方程或经验公式建立水体采样单元电磁波辐射能量与水深的映射关系，利用光学影像的灰度或反射率计算各像元对应位置的未知水深值，获得高密度的水深反演结果；由于理论模型中需要大量的水体参数获取比较困难，因此利用光学影像反演浅海水深时常常采用经验模型；本发明根据前期研究和实验结论，提出一种基于局部水深约束的多波段水深反演模型；该模型有效地将水深变化信息作为约束，融入到多波段水深反演模型，如公式(9)所示；由于采用小范围的精确模型参数代替全局粗略模型参数，从而可以获得更高的反演精度；

式中，N表示波段总数，Z表示对应位置水深值，L(λ_i)是传感器收到的i波段辐射能量(浅水)，L_∞(λ_i)为对应的深水区能量，a₀…a_N为经验回归参数，这些参数由水深样本点与遥感辐射量进行回归分析解算得到；经研究表明，多波段模型用于处理水体底质均一、波段反射比一致的浅海区域水深反演的效果较好；

遥感反演模型中为了更好的将水深变化信息作为约束，提高水深反演精度，因此在反演模型使用的样本中，通过搜索水深值相近的临近样本点作为解算数据集，从而将整个目标区域划分成多个小区域，每个小区域分别用不同水深样本点进行反演。

以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明实施例技术方案的精神和范围。

Claims (7)

1.一种基于星载单光子激光雷达数据测量浅海水深的方法，包括以下步骤：

S1、获取ICESat-2卫星测得的原始数据，包括单光子激光数据、高分遥感数据和潮汐数据；其中在单光子激光数据中获取光子点云的密度结构数据，在高分遥感数据获取高分多光谱影像；

S2、单光子激光点云数据预处理：基于光子点云的密度结构数据进行预处理，通过距离选通技术和初始滤波算法确定目标点或区域的区域范围，减少大量的噪声光子点云数量；

S3、多光谱影像得太阳耀斑消除：采用线性回归模型建立可见光波段与近红外波段的对应关系，来消除可见光波段影像中的耀斑分量；

S4、水面于水下光子分离：将单光子激光雷达数据按高程值从上至下，以Δh为间隔将高程划分为M个切片，据每个切片单元统计光子数量，利用正态分布曲线对高程切片中心高程与切片单元中的光子信号数量之间的关系进行拟合；以切片中心高程为中心，获取每个水面光子的高程；

S5、自适应不同水深光子密度滤波：基于单光子激光雷达探测技术理论以及水体中光子信号的数据特点，确定水面与水下信号点的分界线，并自适应获取初始椭圆搜索区域的长轴、短轴以及密度滤波的初始最小点数阈值；分析随信号光子高程值的变化，高程切片内光子数量的变化特征，确定搜索椭圆的长轴、短轴以及密度滤波的最小点数阈值与水深之间的变化关系；

S6、遥感影像的水陆分离：利用近红外影像和边缘提取算法，开展水陆分离和岸线边缘提取，采用水陆分离模板用于识别陆地和水体，水陆边缘的识别采用估算摄影时刻前后一段时间内的平均海面高度；

S7、水深反演模型构建：通过搜索水深值相近的临近样本点作为解算数据集，从而将整个目标区域划分成多个小区域，每个小区域分别用不同水深样本点进行反演，获得更高的反演精度。

2.根据权利要求1所述的一种基于星载单光子激光雷达数据测量浅海水深的方法，其特征在于，所述步骤S2确定目标区域中有效光子的大概高程范围，用高斯分布函数进行局部拟合，根据高斯峰值、幅宽度等参数为基础，结合先验高程阈值，得到高斯分布函数拟合关系式：

关系式中，G(h)表示高斯分布函数，A_maxpho表示波形振幅，σ为半波宽，h为当前点高程值，h_u表示探测获取最大极值所对应的光子高程值；从而确定目标区域中有效光子信号的初始范围。

3.根据权利要求1所述的一种基于星载单光子激光雷达数据测量浅海水深的方法，其特征在于，所述步骤S2目标区域的有效光子S_signal由多个数据子区R_g(g＝1,L,M)中有效光子的总和构成，关系式为：

关系式中，R_g(g＝1,L,m)表示在沿轨飞行方向划分为m个分段，g为沿轨分段下标值，n表示在当前沿轨分段内光子点的总合数，i为当前沿轨分段内某点点号，

表示当前分块区域的有效光子数，

为当前分段的有效光子数。

4.根据权利要求1所述的一种基于星载单光子激光雷达数据测量浅海水深的方法，其特征在于，所述步骤S4中以Δh为间隔将高程划分为M个切片得到关系式：

根据每个切片单元统计光子数量得到关系式：

关系式中，h_i为每个高程切片中心高程，x_i为切片单元中的光子信号数量。

5.根据权利要求4所述的一种基于星载单光子激光雷达数据测量浅海水深的方法，其特征在于，所述步骤S4中正态分布曲线对高程切片中心高程h_i与切片单元中的光子信号数量x_i之间得到关系式：

关系式中，x_maxpho表示包含光子信号最多的高程切片中的光子数量，h_maxpho表示切片中心高程；以切片中心高程h_maxpho为中心，通过阈值确定水面光子高程分界线，并得到水面光子高程下限H_min和上限H_max。

6.根据权利要求1所述的一种基于星载单光子激光雷达数据测量浅海水深的方法，其特征在于，所述步骤S4中获取每个水面光子的高程，从而形成水面光子数据集SP和水下光子数据集WP，得到等式：

等式中，

表示任意点i的高程，i＝1,2,L n，n为光子总个数，H_min为水面光子高程下限，H_max为水面光子高程上限。

7.根据权利要求1所述的一种基于星载单光子激光雷达数据测量浅海水深的方法，其特征在于，所述步骤S7基于局部水深约束的多波段水深反演模型，所述模型有效地将水深变化信息作为约束，融入到多波段水深反演模型，得到公式：

式中，Z表示对应位置水深值，N表示波段总数，L(λ_i)是传感器收到的i波段辐射能量，L_∞(λ_i)为对应的深水区能量；a₀…a_N为经验回归参数，这些参数根据水深样本点与遥感辐射量进行回归分析解算得出；由于采用小范围的精确模型参数代替全局粗略模型参数，从而获得更高的反演精度。

Images