一种基于解析模型的多视场海洋激光雷达数据正则化反演方法

1.本发明属于海洋激光雷达技术领域，尤其是涉及一种基于解析模型的多视场海洋激光雷达数据正则化反演方法。


背景技术：

2.海洋是全球最大的碳储存库，海洋中的浮游植物是初级生产力的重要来源，通过叶绿体的光合作用和线粒体的呼吸作用参与全球碳循环过程。
3.目前的海洋探测手段包括原位仪器、被动遥感和主动遥感，其中主动遥感领域的激光雷达由于其能够实现昼夜探测、不受纬度限制、可获取海洋垂直剖面信息等优势获得广泛关注，成为了一种重要的海洋探测手段。海洋激光雷达利用激光脉冲的飞行时间和能量强度获取海洋的水体参数信息，脉冲的飞行时间反映了海水深度，脉冲的能量强度反映了水体的固有光学特性。在仅考虑单次散射的情况下，激光雷达的回波信号强度将以光束衰减系数c呈指数衰减，然而海水的致密性使得激光雷达能够在水中发生多次散射，造成时间和空间的激光脉冲展宽，导致激光雷达的衰减偏离理想情况，给水体的固有光学特性反演带来了问题。
4.多次散射效应对激光雷达回波信号的影响来自于两个方面，一方面是激光雷达的系统参数设置，另一方面是海水本身的固有光学特性。就激光雷达的系统参数设置而言，接收视场角是一项重要的影响因素。在较小的视场角下，大量多次散射信号将溢出接收视场，使得回波信号的衰减接近于光束衰减系数c，而在较大的视场角下，多次散射信号能够被系统接收，使得回波信号衰减接近于a+bb。通过设置视场角大小，可以获得受到不同程度多次散射效应影响的激光雷达信号，从而实现海水固有光学特性、颗粒物粒径分布以及水中浮游植物叶绿素浓度的反演。
5.现有系列研究，对海水的多次散射效应进行了讨论，提出了海洋垂直剖面信息的探测及反演方法，如公开号为cn107976686a的中国专利文献基于参数灵敏度评价方法，讨论了工作在蓝绿激光波段的海洋激光雷达视场角选择问题，然而该方法仅对多视场海洋激光雷达设置进行了优化，并未提出水体固有光学特性的反演方法。公开号为cn107831485a的中国专利文献基于walker#mclean模型提出了一种船载多视场激光雷达的数据反演方法，然而这种方法采取的相函数固定，无法针对实际后向散射情况进行相函数更换，此外其引入了多次散射系数估计、散射系数近似和最大深度近似等假设，可能影响探测精度，且无法提供颗粒物粒径、叶绿素浓度等参数。
6.因此，亟需提出一种能够应用于多视场海洋激光雷达的多参数高精度反演方法，为激光雷达探测海洋垂直剖面提供更加全面的理论支撑。本发明提出了一种基于解析模型的多视场海洋激光雷达数据反演算法，利用解析模型计算设定水体固有光学特性参数的激光雷达理论回波信号，通过灵敏度模型推算设定参数与实际参数之间的差值，从而实现多视场海洋激光雷达数据到水体固有光学特性参数以及水中浮游植物叶绿素浓度的反演。


技术实现要素：

7.针对现有激光雷达在海洋水体固有光学特性参数探测的局限性，本发明提供了一种基于解析模型的多视场海洋激光雷达数据正则化反演方法，能够利用多视场海洋激光雷达数据反演海洋水体微物理特性参数。
8.在已知多视场海洋激光雷达的视场大小情况下，假定一组水体固有光学特性参数基于解析模型进行仿真，获取多个视场角下理论信号与实际信号之间的差值，利用灵敏度模型可将信号差值转换为水体固有光学特性参数差值，经过多次迭代计算后可获得一组参数能够保证所有视场下的实际信号与仿真信号的误差小于设定阈值，此时的设定参数即为水体的有效光学特性反演值。
9.一种基于解析模型的多视场海洋激光雷达数据正则化反演方法，包括以下步骤：
10.(1)假定激光雷达接收视场角分别为fov1，fov2，fov3，根据准小角度单次散射近似理论，计算深度为z处的水体信号如下：
[0011][0012]
其中，fov
x
代表视场角，x＝1，2，3，c0为系统常数，n为折射率、h为激光雷达系统距水面的高度，τ(z)为单次散射条件下的光学厚度，等于光束衰减系数c在深度上的积分，为探测器有效区域的傅里叶变换，b为后向散射分数，v为空间频率，b代表散射系数，ξ代表水深的积分变量，p(z,p)为前向散射相函数的hankel变换；设置水体相函数为ff相函数，计算对应的hankel变换为
[0013][0014]
其中，代表归一化的ff相函数；
[0015]
(2)利用水体分层模型将水体设置为分层水体，单层水体的深度为δz，第j层的深度zj＝j
×
δz，探测位置第m层处的深度z＝zm＝m
×
δz，考虑频谱间隔为δv，第i处的频率vi＝i
×
δv，最大频率v
max
＝n
×
δv；假定初始条件下，每一层的吸收系数、散射系数、后向散射系数均为a0、b0、b
b0
，后向散射分数b＝b
b0
/b，在初始设置下仿真激光雷达的回波信号，计算公式如下：
[0016][0017]
(3)对于三视场海洋激光雷达而言，fov1，fov2，fov3下仿真得到的激光雷达信号分别为f
fov1
，f
fov2
，f
fov3
，计算仿真信号与实际激光雷达信号的差值为δf
fov1
，δf
fov2
，δf
fov3
；
[0018]
(4)计算水深z处的水体参数灵敏度模型，定义z处受z'处水体参数影响的灵敏度为x＝a,b,bb，分别代表吸收系数a、散射系数b、后向散射系数bb的灵敏度，算得
吸收系数灵敏度散射系数灵敏度和后向散射系数灵敏度
[0019]
(5)利用仿真信号与实际信号的差值δf
fov1
，δf
fov2
，δf
fov3
和三个视场下的水体参数灵敏度模型计算设置的水体参数a0、b0、b
b0
的矫正值δa、δb、δbb，计算公式如下：
[0020][0021]
其中下标m代表了水体的总层数；
[0022]
(6)通过矫正值δa、δb、δbb更新设定的水体固有光学特性参量a、b、bb；
[0023]
(7)将更新后的水体固有光学特性参量重新进行仿真，获得fov1，fov2，fov3三个视场对应的仿真信号f
fov1
，f
fov2
，f
fov3
，计算仿真信号与实测信号之间的相对误差是否小于设定的误差阈值δ
threshold
；
[0024]
若小于设定误差阈值，则认为矫正后的a、b、bb符合精度要求，a、b、bb即为多视场海洋激光雷达的反演结果，若大于设定误差，则重复步骤(2)-步骤(7)，直至相对误差满足阈值要求；
[0025]
(8)依据后向散射分数b和相函数模型，对颗粒物粒径分布斜率μ进行反演，颗粒物粒径分布斜率计算公式如下：
[0026]
μ＝3-2v
[0027]
n＝1.01+0.1542(μ-3)
[0028][0029][0030]
其中n代表水体的相对折射率，v与δ为中间变量；
[0031]
(9)假定水体为一类水体，利用反演得到的a进一步推演叶绿素浓度chl。
[0032]
进一步地，步骤(2)中，所述水体分层模型的水深范围依据所有视场下的有效探测范围确定。
[0033]
步骤(3)中，仿真信号与实际激光雷达信号的差值计算公式如下：
[0034]
δf
fovx
＝f
fovx-f
fovx,exp
[0035]
其中，x＝1，2，3，分别对应fov1，fov2，fov3三种视场角。
[0036]
步骤(4)中，所述吸收系数灵敏度的公式为：
[0037][0038]
散射系数灵敏度的公式为
[0039][0040]
后向散射系数灵敏度的公式为
[0041][0042]
步骤(5)中，水体参数灵敏度模型计算设置的水体参数a0、b0、b
b0
的矫正值时，采用基于奇异值分解的正则化算法对矩阵求逆。
[0043]
步骤(6)中，通过矫正值δa、δb、δbb更新设定的水体固有光学特性参量a、b、bb公式如下：
[0044]
a＝a0+δa
[0045]
b＝b0+δb
[0046]bb
＝b
b0
+δbb[0047]
步骤(7)中，仿真信号与实测信号之间的相对误差计算公式如下：
[0048][0049][0050][0051]
式中，f
fov1,exp
、f
fov2,exp
、f
fov3,exp
分别表示fov1，fov2，fov3三个视场对应的实测信号；
[0052]
步骤(9)中，叶绿素浓度chl的反演公式如下：
[0053][0054]
其中，aw代表纯水的吸收系数，a
c*
代表无量纲的叶绿素吸收系数，λ代表激光器的激光波长。
[0055]
与现有技术相比，本发明具有以下有益效果：
[0056]
1、本发明能够利用多视场海洋激光雷达的三通道数据，基于多次散射效应在不同视场下的强弱差距实现海水固有光学特性a、b、bb的反演，并可以通过一类水体生物光学模型进一步推算海水中的叶绿素浓度。
[0057]
2、本发明采用准单次小角度散射近似理论，能够仿真考虑多次散射效应的激光雷达回波信号，解决了传统单视场激光雷达无法评估多次散射效应强弱的问题，从而获取更精准的模拟信号，保证了迭代反演的精度。
[0058]
3、本发明能够在不依赖原位仪器的情况下，高效准确的获取包括浮游植物叶绿素浓度、颗粒物粒径分布斜率在内的多种海水信息，有助于实现快速全面的海洋水体监测。
[0059]
4、本发明适用于多种平台的多视场海洋激光雷达系统，能够在多种激光雷达工作高度的情况下实现高精度反演。
附图说明
[0060]
图1为本发明提出的基于解析模型的多视场海洋激光雷达数据正则化反演方法的原理图；
[0061]
图2为本发明提出的基于解析模型的多视场海洋激光雷达数据正则化反演方法的计算流程图；
[0062]
图3为在水面处归一化的三视场海洋激光雷达的实际探测信号；
[0063]
图4为分层的海水模型及坐标系示意图；
[0064]
图5为本发明实施例的激光雷达信号反演结果；
[0065]
图6为本发明实施例根据激光雷达探测到的水体参数反演得到的颗粒物粒径分布斜率；
[0066]
图7为本发明实施例根据激光雷达探测到的水体参数反演得到的叶绿素浓度。
具体实施方式
[0067]
下面结合附图和实施例对本发明做进一步详细描述，需要指出的是，以下所述实施例旨在便于对本发明的理解，而对其不起任何限定作用。
[0068]
如图1所示，本发明的中心思想是通过数次迭代保证在一定的水体参数设置条件下，基于解析模型的激光雷达仿真信号能够在所有视场角下均与实际测得的激光雷达信号之间的误差小于设定阈值，可认为此时的水体参数设置即为水体的光学参数反演值。
[0069]
本发明具体的计算流程如图2所示，基于准单次小角度散射近似理论，构建考虑海水多次散射效应的激光雷达回波参数模型进行不同视场的海水回波信号仿真，计算初始光学特性条件下各视场的回波信号，通过仿真信号与实际测量信号之间的误差进行初始光学特性参数矫正，利用矫正后的参数再次进行信号仿真，并对相对误差进行阈值判断，当误差小于阈值时即认为实现了多视场海洋激光雷达的信号反演。
[0070]
图3展示了一组在0m处归一化的三视场海洋激光雷达实测信号，激光雷达的实际接收视场角fov1，fov2，fov3分别为40mrad，80mrad，200mrad。由于不同视场下，激光雷达的有效探测深度不同，进行反演时需要考虑激光雷达的最低有效探测深度。
[0071]
通过以下步骤实现基于解析模型的多视场海洋激光雷达数据正则化反演：
[0072]
第一步，计算三视场海洋激光雷达三个视场角对应通道的信号强度，根据准单次小角度散射近似理论，可计算深度为z处的水体信号如式(1)所示
[0073][0074]
其中，fov
x
代表视场角，x＝1，2，3，c0为系统常数，n为折射率、h为激光雷达系统距水面的高度，τ(z)为单次散射条件下的光学厚度，等于光束衰减系数c在深度上的积分，为探测器有效区域的傅里叶变换，b为后向散射分数，可由bb/b计算而来，v为空间频率，b代表散射系数，ξ代表水深的积分变量，p(z,p)为前向散射相函数的hankel变换。设置水体相函数为ff相函数，可计算对应的hankel变换为如式(2)所示
[0075][0076]
其中，代表归一化的ff相函数。
[0077]
第二步，根据图3信号的有效探测范围设置水体范围，考虑到浅水处的探测器非线性效应，信号反演范围为0～18m。设置水体为分层水体，海水的分层模型如图4所示，单层水体深度δz为0.92m，第j层的深度zj＝j
×
δz，探测位置第m层处的深度z＝zm＝m
×
δz，m取为20，考虑频谱间隔δv为0.1，第i处的频率vi＝i
×
δv，最大频率v
max
＝n
×
δv，n取为2500。设置初始条件下，每一层的吸收系数a0、散射系数b0、后向散射系数b
b0
分别为0.1、0.1、0.03，后向散射分数b＝b
b0
/b，在初始设置下仿真激光雷达的回波信号，其解析表达式如式(3)所示
[0078][0079]
第三步，基于准单次小角度近似的解析理论，计算不同视场角下的激光雷达信号，计算方式如式(3)所示。对于三视场海洋激光雷达而言，假定其三个视场40mrad，80mrad，200mrad下仿真得到的激光雷达信号分别为f
fov1
，f
fov2
，f
fov3
，由于设定的水体参数初始值与实际水体参数之间具有差异，仿真信号与实际信号之间也有一定差距，计算仿真信号与实际激光雷达信号的差值为δf
fov1
，δf
fov2
，δf
fov3
，其计算方式如式(4)所示
[0080]
δf
fovx
＝f
fovx-f
fovx,exp
ꢀꢀꢀ
(4)
[0081]
其中，x＝1，2，3，分别对应40mrad，80mrad，200mrad三种视场角。信号差值体现了设置的水体参数与实际水体参数之间的差距，可基于灵敏度模型对设置水体参数进行矫正，使其贴近于实际水体参数。
[0082]
第四步，计算水深z处的水体参数灵敏度模型，定义z处受z'处水体参数影响的灵敏度为x＝a,b,bb，分别代表吸收系数a、散射系数b、后向散射系数bb的灵敏度，可算得吸收系数灵敏度如式(5)所示
[0083][0084]
散射系数灵敏度如式(6)所示
[0085][0086]
后向散射系数灵敏度如式(7)所示
[0087][0088]
灵敏度模型反映了不同层水体参数特性对于每层信号的影响程度，同时受到水体参数数值以及水深影响，需要考虑二者影响实时更新对应于设置水体参数的全水深灵敏度模型。
[0089]
第五步，利用仿真信号与实际信号的差值δf
fov1
，δf
fov2
，δf
fov3
和三个视场下的水体参数灵敏度模型计算设置的水体参数a0、b0、b
b0
的矫正值δa、δb、δbb，计算公式如式(8)所示
[0090][0091]
由于存在矩阵求逆的过程，在矩阵求解过程中可能会出现异常值，影响反演精度，因此需要采取基于奇异值分解的正则化方法避免异常值的出现，保证δa、δb、δbb的值处于合理范围。
[0092]
第六步，通过误差δa、δb、δbb更新设定的水体固有光学特性参量a、b、bb，计算公式如式(9)所示
[0093][0094]
误差δa、δb、δbb反映了设置水体参数与实际水体参数之间的差距，矫正水体参数差异能够减小仿真信号与实际信号差异，然而由于灵敏度模型随水体参数变化而变化，单次矫正无法保证设置水体参数与实际参数之间的差距满足设定条件，因此要进行多次迭代计算保证仿真信号与实际信号之间的差距小于设定阈值。
[0095]
第七步，将更新后的水体固有光学特性参量重新进行仿真，获得fov1，fov2，fov3三
个视场对应的仿真信号f
fov1
，f
fov2
，f
fov3
，计算仿真信号与实测信号之间的相对误差是否小于设定阈值，相对误差的计算公式如式(10)所示
[0096][0097]
判断相对误差是否小于设定的误差阈值δ
threshold
，此处的δ
threshold
设定为0.001，当δf
fov1
、δf
fov2
、δf
fov3
均小于设定误差时，认为矫正后的a、b、bb符合精度要求，a、b、bb即为多视场海洋激光雷达的反演结果。图5展示了图4所示的信号对应的反演结果，其中红色实线代表了原位测量的水体光学参数a
t
、b
t
、b
bt
，蓝色三角线代表了经过多次迭代后反演得到的水体光学参数a、b、bb，可见两者吻合效果良好，证明了反演算法的合理性。
[0098]
第八步，颗粒物粒径是表征海水内部浮游植物的重要参数，浮游植物的众多生理特性都与颗粒物粒径分布紧密相关。依据后向散射分数b和相函数模型，可对颗粒物粒径分布斜率μ进行反演，其计算公式如式(12)所示。
[0099][0100]
其中n代表水体的相对折射率，v与δ为中间变量。依据实施例的反演结果计算的颗粒物粒径分布斜率μ如图6所示。
[0101]
第九步，在海洋中，浮游植物是重要的初级生产力来源，对于全球生态循环具有重要意义，而叶绿素浓度是体现浮游植物生物量的重要参数，已有研究指出吸收系数a与叶绿素浓度之间具有定量关系，依据吸收系数可以推演叶绿素浓度。由于探测水体非常清澈，可以认为水体为一类水体，采用一类水体生物光学模型利用反演得到的a进一步计算叶绿素浓度chl，叶绿素浓度计算公式如式(11)所示。
[0102][0103]
实施例中的海洋激光雷达发射波长λ为532nm，因此取纯水的吸收系数aw为0.051m-1
，a
c*
代表无量纲的叶绿素吸收系数，在532nm波长下等于0.474，依据实施例的反演结果计算的叶绿素浓度如图7所示。
[0104]
在本实施例中，通过本发明的多视场海洋激光雷达数据正则化反演方法，实现了
水体固有光学特性a、b、bb、颗粒物粒径分布斜率μ、叶绿素浓度垂直剖面特性的反演，与原位仪器测量结果吻合良好，证明了本发明的有效性。
[0105]
以上所述的实施例对本发明的技术方案和有益效果进行了详细说明，应理解的是以上所述仅为本发明的具体实施例，并不用于限制本发明，凡在本发明的原则范围内所做的任何修改、补充和等同替换，均应包含在本发明的保护范围之内。