-
早期的激光测量技术一般只记录若干次回波的位置信息,而不记录整个激光脉冲回波,这样可以节省存储空间,提高数据处理效率。实际上, 全波形的激光数据更能适应不同的任务需求,有利于分析相邻激光脉冲之间的关系,通过对回波的分析和处理,可以有效地提高激光测距的精度[1-2]。全波形激光测距系统最早应用于激光测深系统中[3],激光脉冲由激光发射器出射至到达海底的过程中,与水-气界面、水体、海底相互作用,生成复杂的回波信号[4]。直接从采集的数据中实时获取有效、高精度的水深数据是非常困难的[5],通用的方法是在完成数据采集后,对全波形数据进行后处理获取高精度水深数据。
随着全波形激光测距系统的发展,波形数据处理算法也逐渐成为激光雷达数据处理领域研究的热点[6-8]。由于激光脉冲在水体中传输时信号衰减过程复杂,一般的波形数据分解方法较少考虑激光测深回波在水体中传输的部分或将其与海面反射视为整体来处理[8],这些方法对清澈大洋水体具有一定的适用性,若直接应用于浑浊水体,会导致较大的误差。Abdallah提出了Wa-LiD(Water LiDAR)[9]激光测深波形模拟工具,利用这一工具可生成具有不同水体参数、激光测深系统参数的模拟数据集。Abdallah[10]、Abady[11]分别对水体后向散射回波进行模拟,得到了有意义的结果。
本文将进一步探索函数拟合方法在激光测深回波信号处理中的应用,提出使用指数函数模拟激光测深水体后向散射回波成分的方法,海面回波和海底回波则采用高斯函数。测深处理结果分别与Abdalla[10]提出的三角形拟合方法、Abady[11]提出的四边形拟合方法进行比较,并分析了不同函数模拟水体后向散射对测深结果的影响。
-
本文使用的激光测深波形数据由Wa-LiD波形模拟工具得到。对于水面反射波形的建模,使用了Cook和Torrance [12]提出的表面反射模型对海面几何进行模拟,同时考虑了激光在海水中传输时的辐射传输规律。激光测深回波信号与地面激光测距回波信号的主要区别在于水体后向散射部分,水深z处激光信号后向散射功率[9]为:
$$ \begin{matrix} {{P}_{\text{bw}}}\left( z \right)= \\ \frac{{{P}_{e}}T_{\rm{atm}}^{^{2}}{{A}_{R}}{{\eta }_{e}}{{\eta }_{R}}F\text{ }\left( 1-{{L}_{s}} \right){{~}^{2}}{{\beta }_{\text{ }\!\!\pi\!\!\text{ }}}\exp \left( \frac{-2kz}{\cos {{\theta }_{w}}} \right)}{{{\left( \frac{~{{n}_{w}}\cdot H+z}{\cos \theta } \right)}^{2}}} \\ \end{matrix} $$ (1) 式中各参数定义可参考文献[9]。
在以往的研究中[9-11, 13],水体衰减系数k通常设置为较大的范围,而βπ设置为固定的数值,很少讨论βπ对波形分解结果的影响。由Wa-LiD模拟工具和激光测深信号形成的物理过程可知,水体后向散射回波是一系列位于不同水深处高斯波形的叠加,而每个高斯波形都是发射激光脉冲与水深z处的水体后向散射功率的卷积。因此,βπ对水体后向散射波形的振幅具有很大的影响。Petzold[14]利用海洋光学仪器对不同水域的体积散射函数进行了实地测量,得到了514 nm处不同角度的体积散射函数,Mobley [15]将Petzold的实验结果整理后制成了表格,给出了4种不同水体的βπ值,分别为2.936×10-4m-1sr-1,5.019×10-4m-1sr-1,1.030×10-3m-1sr-1,5.686×10-3m-1sr-1。激光测深采用波长一般为532 nm,与514 nm具有相似的散射特性,因此本文将分析文献[15]中4种不同的βπ取值对激光测深波形数据分解的影响。
为了使模拟数据具有更广泛的统计学意义,对应4种不同的βπ取值的每组模拟数据涵盖了不同水质条件下的20 000个模拟波形,水深范围为1~10 m,即每个水深处随机选择水质参数生成2 000个波形。用于生成模拟波形的相关水质参数范围和环境参数范围如下:海面粗糙度r为0.1~0.5;海底反射率Rb为0.05~0.2;海表镜面反射率ks为0.6~0.9;漫衰减系数k(m-1)为0.01~1;海水折射率nw为1.33;双向大气损失因子Tatm2为0.9;水体反射太阳辐射量Is(w/m2/sr/nm)为0.025。激光测深系统的相关参数使用HawkEye激光测深系统的系统参数[11]。
-
激光测深回波信号受多种因素影响,使用含有少量参数的数学表达式模拟测深波形,能够有效地减少未知参数数量。虽然这种简化的数学模拟可能会导致一定程度的误差,但是赋予模拟函数的参数以合理的物理意义后再进行验证,可以证明测深数据的处理结果是否可靠可用。
对于水体后向散射部分的函数模拟方法,Abdallah、Abady分别提出了三角形函数(图 1(a))和四边形函数(图 1(b))模拟方法,三角形函数的数学表达式为:
图 1 使用三角形、四边形、指数函数拟合水体后向散射示意图
Figure 1. Simulate Water Backscatter Echo Using Triangle, Quadrilateral, Exponential Function Respectively
$$ ~{{f}_{\text{tri}}}\left( t \right)={{A}_{\text{tri}}}\times \left\{ \begin{align} &0, t\le a \\ &\frac{t-a}{b-a}, a\le t\le b \\ &\frac{c-t}{c-b}, b\le t\le c \\ &0, c\le t \\ \end{align} \right. $$ (2) 式中,Atri为三角形函数的振幅;(a, b, c)为三角形3个顶点的横坐标位置; t为时间。四边形函数的数学表达式为:
$$ {{f}_{\text{squ}}}\left( t \right)=\left\{ \begin{align} &0, t\le {{d}_{\text{squ}}} \\ &{{A}_{\text{squ1}}}\cdot \text{ }\frac{t-{{d}_{\text{squ}}}}{{{e}_{\text{squ}}}-{{d}_{\text{squ}}}}~, \ \ \ {{d}_{\text{squ}}}\le t\le {{e}_{\text{squ}}} \\ &\frac{{{A}_{\text{squ1}}}\cdot {{g}_{\text{squ}}}-{{e}_{\text{squ}}}\cdot {{A}_{\text{squ2}}}+t\left( {{A}_{\text{squ2}}}-{{A}_{\text{squ1}}} \right)}{{{g}_{\text{squ}}}-{{e}_{\text{squ}}}}, \ \ {{e}_{\text{squ}}}\le t\le {{g}_{\text{squ}}} \\ &{{A}_{\text{squ2}}}\cdot \text{ }\frac{{{h}_{\text{squ}}}-t}{{{h}_{\text{squ}}}-{{g}_{\text{squ}}}}, \ \ \ {{g}_{\text{squ}}}\le t\le {{h}_{\text{squ}}} \\ &0, {{h}_{\text{squ}}}\le t~ \\ \end{align} \right. $$ (3) 式中,Asqu1、Asqu2为四边形第二、第三个顶点的纵坐标;(dsqu, esqu, gsqu, hsqu)为四边形4个顶点的横坐标位置。
三角形函数和四边形函数分别使用三角形和四边形模拟水体后向散射回波,并未考虑水体后向散射的物理特征。由光线在水体中辐射传输可知,激光回波功率在海水中传播时的幅值随水深的增加而指数衰减,体散射回波成分随水深(时间)的变化曲线在理想情况下可以用指数函数表示。因此,本文提出使用指数函数(图 1(c))拟合水体后向散射成分,由此得到水体后向散射的数学表达式为:
$$ {{f}_{\text{ex}}}\left( t \right)=\left\{ \begin{align} &0, t\le {{d}_{\text{ex}}} \\ &{{A}_{\text{ex1}}}\cdot \text{ }\frac{t-{{d}_{\text{ex}}}}{{{e}_{\text{ex}}}-{{d}_{\text{ex}}}}~, \ \ \ {{d}_{\text{ex}}}\le t\le {{e}_{\text{ex}}} \\ &\exp ~\left( \frac{\left( {{g}_{\text{ex}}}-t \right)\cdot \ln {{A}_{\text{ex1}}}-\left( {{e}_{\text{ex}}}-t \right)\text{ }\cdot \ln {{A}_{\text{ex1}}}}{~{{g}_{\text{ex}}}-{{e}_{\text{ex}}}} \right), \ \ {{e}_{\text{ex}}}\le t\le {{g}_{\text{ex}}} \\ &{{A}_{\text{ex2}}}\cdot \text{ }\frac{{{h}_{\text{ex}}}-{{t}_{\text{ex}}}}{{{h}_{\text{ex}}}-{{g}_{\text{ex}}}}~, \ \ \ {{g}_{\text{ex}}}\le t\le {{h}_{\text{ex}}} \\ &0, \ \ \ {{h}_{\text{ex}}}\le t \\ \end{align} \right. $$ (4) 式中,Aex1、Aex2为指数函数左右端点的纵坐标;(eex, gex)为指数函数左右端点的横坐标;(dex, hex)为用于模拟水体后向散射波形起始、结束的横坐标位置。
本文将分别利用这3种拟合方法对模拟波形数据采用非线性递归最小二乘[16]方法进行分解,优化海面回波和海底回波位置,并分析不同的拟合方法对分解结果的影响。
-
在处理激光测深回波信号之前,应对波形信号进行平滑处理,平滑函数可使用与激光发射脉冲同宽的高斯脉冲[2]。为了获取水深,激光雷达回波信号应具有两个可识别的峰值,分别对应于海面回波信号和海底回波信号,如图 2所示。图 2中,As和ts是海面回波峰值振幅和位置,Ab和tb是海底回波峰值振幅和位置。对于这两个信号位置的探测,可采用文献[10]中的方法。
图 2 由Wa-LiD工具模拟激光测深回波波形及平滑波形示例
Figure 2. Example of Laser Bathymetry Echo Waveform and Smoothed Waveform Simulated by Wa-LiD
使用该方法对激光测深回波信号进行探测,当探测到峰值结果为两个时,保留该波形,用于下一步的波形分解实验。因此本实验具体包括如下步骤:①使用Wa-LiD工具生成用于实验的多组激光测深回波数据集;②对测深回波信号进行平滑处理,抑制噪声信号;③使用峰值探测方法探测海面回波和海底回波,保留具有两个波峰的回波信号;④使用不同的波形模拟函数对平滑后的波形进行分解,获取激光测深结果,并进一步分析。
-
使用§2.2提出的峰值探测方法对4个模拟数据集提取可探测波形(具有两个峰值的波形),并计算每个水深处可探测波形占模拟波形比例。
图 3给出了每个模拟数据集在不同水深处可探测波形的百分比。由图 3可以看出,4个不同的数据集中波形探测概率随水深具有相似的变化规律,即均在水深2 m时探测概率达到最大,之后随水深的增加而逐渐减少。在水深1 m处,4个数据集探测概率均较小,这主要是由于水深较浅时海面回波和海底回波产生了叠加,使得海底回波难以探测。在βπ取值较小的前3个数据集中,探测概率随水深呈指数衰减。这与海底回波振幅随水深增加而指数衰减的趋势一致。这种探测概率的衰减趋势也说明了前3个数据集中未被成功探测的波形主要是由于海底回波较弱而湮没在噪声中造成的。而模拟数据集4表现出了与前3个数据集不同的探测概率变化趋势,且在浅水区域探测概率较前3个数据集有大幅下降。一个较为合理的解释是在该数据集中,采用的βπ数值较大,使水体后向散射幅值显著增加,海底回波与水体后向散射回波叠加之后,使得海底回波难以被识别。
图 3 不同水深条件下可探测波形百分比
Figure 3. The Percentage of Waveforms that can be Detected at Different Water Depth
图 4给出了海底回波未被识别的一个示例波形。生成该波形的水质参数为βπ=5.686×10-3m-1sr-1,ks=0.849,k=0.362 m-1,r=0.417,Rb=0.080,水深(垂直深度)3 m。生成该波形的参数取自图 3中模拟数据集4,为了便于显示,模拟波形中未加入噪声。从图 4中可以看出,水体后向散射与海底波形叠加之后,海底回波波峰消失,使用峰值探测方法无法从模拟波形中提取出海底回波波峰位置。因此,模拟数据集4中未被成功探测的波形有两类,第一类为海底回波较弱湮没在噪声中无法探测;第二类为海底回波与水体后向散射叠加而使海底回波“湮没”在水体后向散射中难以被识别。第二类未被成功探测波形主要出现在水深较浅的区域,随着水深增加,水体后向散射呈指数衰减,其“尾端”对海底回波的影响逐渐减小,这也是图 3中随着水深增加,模拟数据集4可探测波形百分比与前3个数据集之间差距缩小的原因。同样的,从图 2还可以看出,同一水深处βπ越大,探测概率越低,这也是由于第二类未被成功探测的波形随βπ增大而增多导致的。
图 4 海底回波与水体后向散射回波叠加示例波形
Figure 4. An Example Waveform Generated by Superimposing Seabed Echo Waveform and Water Backscatter Waveform
对于每一个探测到具有两个波峰的波形,使用§2.2提出的3种波形拟合方法, 对平滑后的激光测深回波信号采用非线性递归最小二乘方法进行分解,优化海面回波和海底回波位置。根据分解结果计算模拟数据真实水深与不同拟合方法获取深度值之间的误差,利用数据集不同深度误差统计结果计算偏差和标准差。图 5以直方图的形式给出了4个数据集不同水深处误差和标准差的大小。图 5中正值为标准差,负值为偏差。从图 5中可以看出:
图 5 4个数据集不同水深处误差统计量(偏差、标准差)直方图
Figure 5. Error Statistics(Bias, Standard Deviation) Histogram at Different Water Depths for Four Datasets
1) 随着体积散射函数β(θ)取180°时数值βπ的增加,3种拟合方法在不同水深处水深测深结果误差随βπ的增加而增加。
2) 所有数据集就偏差而言呈现相同的变化规律:在水深小于3 m的区域,3种拟合方法分解结果之间差别不大;在水深大于3 m的区域,使用指数函数拟合获取的测深结果偏差最小,其次是四边形拟合方法。其中在βπ较小的前两个数据集中,指数函数拟合结果优势较为明显,偏差相对于四边形拟合方法有大幅提高。在第三个数据集部分水深也有同样的现象。
3) 就方差而言,4个数据集仍保持一致, 指数函数拟合结果的方差均为最小,其中前3个数据集中,三角形拟合方法的方差小于四边形拟合方法,在第四个数据集中,水深超过5 m的区域,三角形拟合方法的方差大于四边形拟合方法。通过对波形分解结果的检查发现,在某些参数情况下,四边形拟合结果误差要远大于三角形和指数函数,这是造成前3个数据集中四边形拟合方差均大于三角形拟合方差的原因。图 6中示例波形取自数据集2,即βπ=5.019×10-4m-1sr-1,其他参数为ks=0.831,k=0.282 m-1,r=0.281,Rb=0.177,水深(垂直深度)7 m。图 6(a)显示了四边形拟合结果,可以看出,拟合后的海底回波波峰和振幅均与实际波形有较大偏差;图 6(b)为对图 6(a)中海底回波的放大,图中标记出了不同拟合方法得到的实际海底回波波峰的位置。拟合误差为三角形法-4.26 cm,四边形法-18.04 cm,指数函数法-1.29 cm。而当βπ较大时,三角形拟合又出现较“异常”的误差,造成第四个数据集中部分水深三角形拟合方法的方差大于四边形拟合方法。这种情况如图 7所示,此波形取自数据集4,即βπ=5.686×10-3m-1sr-1,其他参数为ks=0.801,k=0.240 m-1,r=0.143,Rb=0.198,水深(垂直深度)7 m。图 7(a)显示了三角形拟合结果,可以看出,拟合后海底高斯脉冲完全偏离了实际位置,脉冲宽度也被放大;图 7(b)为对图 7(a)中海底回波的放大,图中同样标记出了不同拟合方法得到的实际海底回波波峰的位置。拟合结果为:三角形拟合误差-241.37 cm,四边形拟合误差-57.71 cm,指数函数拟合误差-15.27 cm。
4) 三角形拟合方法、四边形拟合方法拟合结果的偏差和标准差随水深变化均具有较大的变化,而指数函数拟合方法结果几乎不受水深的影响,是3种方法中拟合结果最稳定的,这非常有利于对测深结果误差进一步校正。
通过对模拟数据集波形拟合结果的分析可以看出,不论是三角形拟合方法还是四边形拟合方法,在对水体后向散射回波信号建模时都会产生误差较大的情况,这主要是因为这两种函数使用直线拟合水体后向散射随水深的变化趋势,没有考虑水体后向散射回波功率随水深呈指数衰减这一趋势。指数函数拟合方法从辐射传输规律角度出发,遵循激光在海水中传输的物理过程,因此更能够适用于不同的水质条件和水深。
-
本文在原有的激光测深回波信号拟合方法的基础上,结合信号形成的物理过程,提出使用指数函数模拟激光测深水体后向散射回波成分的方法。在对模拟激光测深回波预处理之后,本文使用峰值探测方法检测海面和海底回波,发现模拟数据集4表现出了与前3个数据集不同的探测概率变化趋势,且在浅水区域探测概率较前3个数据集有大幅下降。这是由于前3个数据集未被成功探测的波形因海底回波较弱而湮没在噪声中,而模拟数据集4中除了海底回波较弱造成部分海底信号未被成功探测之外,还有大量的未被成功探测的回波信号是由于海底回波与水体后向散射叠加而使海底回波“湮没”在水体后向散射回波中造成的。
对于每一个探测到的具有两个波峰的波形,本文分别使用三角形、四边形、指数函数拟合方法对海面回波和海底回波位置进行优化,发现测深结果误差随βπ的增加而增加。就偏差而言,指数函数拟合获取的测深结果最优,其次是四边形拟合方法。测深结果误差的标准差在不同的βπ取值下呈现不同的结果,进一步对波形分解案例的检查发现,不论是三角形拟合方法还是四边形拟合方法,在对水体后向散射回波信号建模时都会产生误差较大的情况。这也正是造成测深误差统计结果不稳定的原因。而指数函数拟合方法从辐射传输规律角度出发,遵循水体后向散射回波功率随水深呈指数衰减这一趋势,因此更能够适用于不同的水深和水质条件。
The Impact of Different Fitting Functions for Water Backscatter Waveforms on the Accuracy of Laser Sounding
-
摘要: 为了提高激光测深数据处理结果的精度,在原有的激光测深回波信号拟合方法的基础上,结合信号形成的物理过程,提出使用指数函数模拟激光测深水体后向散射回波成分的方法。利用Wa-LiD(Water LiDAR)工具生成了多组模拟数据集以对比和检验提出的新方法。对于每个可探测回波,海面回波和海底回波均使用高斯函数进行拟合,水体后向散射波形分别使用三角形、四边形、指数函数进行拟合。实验结果显示,指数函数测深误差统计量明显优于三角形和四边形拟合方法,且更适用于拟合不同水质和水深条件下水体后向散射回波信号。同时检验了体积散射函数取180°时数值对测深回波信号探测结果的影响。Abstract: Considering the physical process of laser sounding echo signal formation, a fitting method using exponential function to simulate water column contribution in Lidar waveforms is proposed on the basis of the existing Lidar waveforms fitting method. In order to compare and test the new method, pluralities of sets of simulated data sets are generated for further analyses. For each detectable echoes, the sea surface echoes and the seabed echoes are fitted by Gauss function, and the backscattering echoes are fitted by triangle, quadrilateral and exponential function. The results show significantly more accurate laser sounding estimates for the use of exponential function method compared with the use of a triangular function or a quadrilateral function to fit the column contribution. Moreover, exponential function is more suitable for fitting water column contribution in Lidar waveforms of different water quality and water depth. Meanwhile, the paper examines the impact of the volume scattering function at 180° on laser sounding estimates.
-
-
[1] Jutzi B, Stilla U. Laser Pulse Analysis for Reconstruction and Classification of Urban Objects[J].International Archives of Photogrammetry, Remote Sensing and Spatial Information Sciences, 2003, 34(3/W8):151-156 http://www.researchgate.net/publication/26919955_Laser_pulse_analysis_for_reconstruction_and_classification_of_urban_objects [2] Hofton M A, Minster J B, Blair J B. Decomposition of Laser Altimeter Waveforms[J].IEEE Transactions on Geoscience and Remote Sensing, 2000, 38(4):1989-1996 doi: 10.1109/36.851780 [3] Roncat A, Wagner W, Melzer T, et al. Echo Detection and Localization in Full-Waveform Airborne Laser Scanner Data Using the Averaged Square Difference Function Estimator[J]. Photogrammetric Journal of Finland, 2008, 21(1):62-75 https://www.researchgate.net/publication/248359607_Echo_Detection_and_Localization_in_Full-Waveform_Airborne_Laser_Scanner_Data_using_the_Averaged_Square_Difference_Function_Estimator [4] Tulldahl H M, Steinvall K O. Analytical Waveform Generation from Small Objects in LiDAR Bathymetry[J].Applied Optics, 1999, 38(6):1021-1039 doi: 10.1364/AO.38.001021 [5] Guenther G C, Cunningham A G, Larocque P E, et al. Meeting the Accuracy Challenge in Airborne LiDAR Bathymetry[J]. Proc. Earsel Symp. Workshop on LiDAR Remote Sensing of Land & Sea, 2000(1):1-27 http://agris.fao.org/agris-search/search.do?recordID=AV20120133452 [6] Guenther G C, Mesick H C. Analysis of Airborne Laser Hydrography Waveforms[C]. International Society for Optics and Photonics Orlando Technical Symposium, Orlando, 1988 http://proceedings.spiedigitallibrary.org/proceeding.aspx?articleid=1252318 [7] Wong H, Antoniou A. Characterization and Decomposition of Waveforms for LARSEN 500 Airborne System[J]. IEEE Transactions on Geoscience and Remote Sensing, 1991, 29(6):912-921 doi: 10.1109/36.101370 [8] Allouis T, Bailly J S, Pastol Y, et al. Comparison of LiDAR Waveform Processing Methods for very Shallow Water Bathymetry Using Raman, Near-Infrared and Green Signals[J]. Earth Surface Processes and Landforms, 2010, 35(6):640-650 https://www.researchgate.net/publication/229804649_Comparison_of_LiDAR_Waveform_Processing_Methods_for_Very_Shallow_Water_Bathymetry_Using_Raman_Near-Infrared_and_Green_Signals [9] Abdallah H, Baghdadi N, Bailly J S, et al. Wa-LiD:A New LiDAR Simulator for Waters[J]. IEEE Geoscience and Remote Sensing Letters, 2012, 9(4):744-748 doi: 10.1109/LGRS.2011.2180506 [10] Abdallah H, Bailly J S, Baghdadi N N, et al. Potential of Space-Borne LiDAR Sensors for Global Bathymetry in Coastal and Inland Waters[J]. Selected Topics in IEEE Journal of Applied Earth Observations and Remote Sensing, 2013, 6(1):202-216 doi: 10.1109/JSTARS.2012.2209864 [11] Abady L, Bailly J S, Baghdadi N, et al. Assessment of Quadrilateral Fitting of the Water Column Contribution in LiDAR Waveforms on Bathymetry Estimates[J]. IEEE Geoscience and Remote Sensing Letters, 2014, 11(4):813-817 doi: 10.1109/LGRS.2013.2279271 [12] Cook R L, Torrance K E. A Reflectance Model for Computer Graphics[J]. ACM Transactions on Graphics (TOG), 1982(1):7-24 [13] Wang C, Li Q, Liu Y, et al. A Comparison of Waveform Processing Algorithms for Single-Wavelength LiDAR Bathymetry[J].ISPRS Journal of Photogrammetry and Remote Sensing, 2015(101):22-35 https://www.sciencedirect.com/science/article/pii/S0924271614002718 [14] Mobley C D. Light and Water:Radiative Transfer in Natural Waters[M]. New York:Academic Press, 1994 [15] Petzold T J. Volume Scattering Functions for Selected Ocean Waters[R]. San Diego: Scripps Institution of Oceanography, 1972 [16] Marquardt D W. An Algorithm for Least-Squares Estimation of Nonlinear Parameters[J]. Journal of the Society for Industrial & Applied Mathematics, 1963, 11(2):431-441 doi: 10.1137/0111030 -