CN105068063A - 一种基于回波全波形分解的地物特征获取方法 - Google Patents

Abstract

本发明涉及一种基于回波全波形分解的地物特征获取方法，包括以下步骤：1)根据激光高度计的参数获取激光高度计的激光发射信号的数学模型；2)获取目标地物的三维地物模型以及激光高度计发射传输过程的数据误差；3)根据三维地物模型以及激光发射信号模型并考虑激光高度计发射传输过程的数据误差，获取激光发射信号与目标地物相互作用产生回波的数学模型；4)对激光发射信号与目标地物相互作用产生的回波采用非线性最小二乘高斯分解算法进行全波形分解，获得分解后的高斯波；5)根据分解后的高斯波的波形获取目标地物的垂直高度，并估计激光高度计的实际测量精度。与现有技术相比，本发明具有精度高、成本低、适用范围广等优点。

Description

一种基于回波全波形分解的地物特征获取方法

技术领域

本发明涉及测绘、摄影测量与遥感领域，尤其是涉及一种基于回波全波形分解的地物特征获取方法。

背景技术

自上世纪70年代以来，卫星测高技术已经得到迅速发展并取得了巨大的研究成果。精确的测高任务已经改变了人类对地球的认识和观测方式。而高精度的卫星测高使我们有能力并系统的进行与其相关的各种科学研究。2003年1月13日，美国NASA发射了ICESat(Ice，CloudandLandElevationSatellite)，其上搭载了一个激光测高系统GLAS(GeoscienceLaserAltimeterSystem)，到2009年止，已经采集了密集的全球范围高精度高程数据，在陆地上主要应用于测定全球地形，作为地形图和数字高程模型的参考基准。然而城市区域特有的地物类型，导致传统的全波形分解方法无法完整准确的提取信息，需要在城区星载激光测距雷达仪器设计方案基础上，建立激光测距雷达的回波模拟模型，并从波形数据中提取地物特征，进行精度评价。

发明内容

本发明的目的就是为了克服上述现有技术存在的缺陷而提供一种精度高、成本低、适用范围广的基于回波全波形分解的地物特征获取方法。

本发明的目的可以通过以下技术方案来实现：

一种基于回波全波形分解的地物特征获取方法，包括以下步骤：

1)根据激光高度计的参数获取激光高度计的激光发射信号的数学模型；

2)获取目标地物的三维地物模型以及激光高度计发射传输过程的数据误差；

3)根据三维地物模型以及激光发射信号模型并考虑激光高度计发射传输过程的数据误差，获取激光发射信号与目标地物相互作用产生回波的数学模型；

4)对激光发射信号与目标地物相互作用产生的回波采用非线性最小二乘高斯分解算法进行全波形分解，获得分解后的高斯波；

5)根据分解后的高斯波的波形获取目标地物的垂直高度，并估计激光高度计的实际测量精度。

所述的步骤1)中激光发射信号波形的数学模型为高斯函数：

$p_T\left(t\right)=A_T\exp \{-\frac{{(t-t_{C,T})}^2}{2{s_T}^2}\}$

其中，p_T(t)为高斯函数的峰值，t_C,T为高斯函数的中心位置，s_T为高斯函数的标准差。

所述的步骤2)中通过数字表面模型表示地表和地物表面的起伏和分布三维地物模型。

所述的步骤2)中数据误差包括激光高度计发射和接收不平行造成的误差、脉冲信号传播时间的测量误差、震动误差、反光镜的旋转误差、脉冲零点误差、大气中水分子、悬浮物和大气的折射、散射造成的测量误差以及大气温度、湿度和气压变化造成的测量误差。

所述的步骤3)中激光发射信号与目标地物相互作用产生回波的数学模型为：

$p_R\left(t\right)\≈\underset{i=1}{\overset{N}{\Σ}}\{\left({\mathrm{\ϵ}}_i\frac{4p_T\left(t\right)}{{\mathrm{\πR}}_i^2{\mathrm{\β}}_T^2}\right)*{\mathrm{\σ}}_i\left(t\right)\}\·\frac{1}{4{\mathrm{\πR}}_0^2}\·\frac{{\mathrm{\πD}}_R^2}{4}\·{\mathrm{\η}}_{atm}\·{\mathrm{\η}}_{sys}$

${\mathrm{\ϵ}}_i=\frac{2}{\mathrm{\π}ab}\·\exp \{-\frac{1}{2}\[{\left(\frac{2x_i}{a}\right)}^2+{\left(\frac{2y_i}{b}\right)}^2\]\}\·\frac{{\mathrm{\πR}}_i^2{\mathrm{\β}}_T^2}{4}$

${\mathrm{\σ}}_i\left(t\right)=\frac{4\mathrm{\π}}{{\mathrm{\Ω}}_i\left(t\right)}{\mathrm{\ρ}}_i\left(t\right)A_{S,i}\left(t\right)$

其中，p_R(t)为地物回波波形，p_T(t)为发射脉冲功率，R_i为地物与激光高度计之间的距离，R₀为所有地物与激光高度计的近似距离，β_T为发射脉冲的波束宽度，ρ_i(t)为地物反射率，σ_i为后向散射截面，A_S,i(t)为地物受到脉冲能量照射的面积，Ω_i(t)地物散射立体角，D_R为接收天线孔径直径，η_atm表示大气对信号传输的影响参数，η_sys表示硬件系统对信号传输的影响参数，a、b为激光测距雷达激光光束的照射大脚印内椭圆的长半轴和短半轴，x_i、y_i分别为第i个地物的坐标值，ε_i为激光测距雷达发射的脉冲能量呈椭圆高斯分布引入的参数。

所述的步骤4)中分解后的高斯波的波形数学模型为：

$w\left(t\right)=\underset{m=1}{\overset{N_p}{\Σ}}{A}_m\exp \[-\frac{{(t-t_m)}^2}{2{\mathrm{\σ}}_m^2}\]+\mathrm{\ϵ}$

其中，N_p是回波中峰值个数，A_m是第m个峰值的大小，ε是噪声影响而产生的偏差，t_m是峰值的位置，σ_m是第m个峰值的标准差。

所述的步骤5)中估计激光测距雷达的实际测量精度估计方法为：

采用分解后的高斯波的波形数学模型对脉冲信号进行波形拟合，并确定回波波形的中心位置，估计精度RMSE的表达式为：

$RMSE=\sqrt{\frac{\mathrm{\Σ}{(h-h_{set})}^2}{n}}$

h＝c*(t₁-t₂)/2

其中，h为通过回波波形拟合计算得到的地物的垂直高度，h_set为建模时设定的地物垂直高度，即真值，n为仿真次数，c为光速，t₁和t₂为相邻脉冲的时刻。

与现有技术相比，本发明具有以下优点：

一、精度高：通过仿真模拟城市不同地表目标回波波形，对城市激光高度计的系统设计参数提供重要的理论依据。

二、成本低：激光系统和技术的开发费用高昂，需要一个准确的系统仿真模型来指导工程实践、降低费用，并为数据处理、地物信息提取算法开发提供必要的模拟数据源。

三、建立具有城市特征场景的模型。

四、适用范围广：建立模型数据库，存储各种目标模型，以及实验数据库，不同的目标模型对应的波形特征以及实验结果，方便用户使用。

五、回波数学模型考虑到了脚印内里的脉冲能量从中间向外沿逐渐减弱这一特性，即地物在脚印内的位置不同，接收到的能量是不同的。

六、利用回波数学模型模拟时，考虑到了大气、硬件系统对信号传输的影响，使模拟的回波信号更接近真实环境中的回波信号。

七、精确描述、控制城市地物场景，模拟出回波波形，有利于分析地物面积位置等微小变化是如何影响到激光回波信号的。

附图说明

图1为本发明的方法流程图。

图2为实施例1中地面反射率为0.65时的波形图。

图3为实施例1中地面反射率为0.06时的波形图。

图4为实施例2中地面反射率为0.65时的波形图。

图5为实施例2中地面反射率为0.06时的波形图。

图6为实施例3中地面反射率为0.65时的波形图。

图7为实施例3中地面反射率为0.06时的波形图。

图8为实施例4中地面反射率为0.65时的波形图。

图9为实施例4中地面反射率为0.06时的波形图。

具体实施方式

下面结合附图和具体实施例对本发明进行详细说明。

如图1所示，本方法主要有以下三步组成：

激光发射信号模拟

激光测距雷达的发射波形可近似为半峰宽度(FWHM，FullWidthHalfMagnitude)为8ns的高斯函数。

$p_T\left(t\right)=A_T\exp \{-\frac{(t-t_{C,T})}{2{s_T}^2}\}$

其中，A_T是高斯函数的峰值，t_CT是高斯函数的中心位置，s_T是其标准差。

信号传输过程会使测高数据产生误差。主要包括高度计发射和接受不平行，脉冲信号传播时间的测量误差、震动误差、反光镜的旋转及脉冲零点误差，仪器姿态误差，大气中的水分子、悬浮物以及大气的折射、散射等，将会使测量产生偏差，另外激光脉冲的传播还受到大气温度、湿度和气压的影响。

激光波束与典型特性地物目标相互作用的数学模型

不同地物目标的模拟回波波形对测高仪系统参数设计具有重要意义。激光测距雷达可以精确记录发射脉冲信号从卫星发出经过光斑内地表目标的反射到接受器接收到回波信号的时间，连续的记录激光雷达回波构成了回波波形。激光雷达回波方程如下式所示：

$p_R\left(t\right)\≈\underset{i=1}{\overset{N}{\Σ}}\{\left({\mathrm{\ϵ}}_i\frac{4p_T\left(t\right)}{{\mathrm{\πR}}_i^2{\mathrm{\β}}_T^2}\right)*{\mathrm{\σ}}_i\left(t\right)\}\·\frac{1}{4{\mathrm{\πR}}_0^2}\·\frac{{\mathrm{\πD}}_R^2}{4}\·{\mathrm{\η}}_{atm}.{\mathrm{\η}}_{sys}$

D_R

${\mathrm{\ϵ}}_i=\frac{2}{\mathrm{\π}ab}\·\exp \{-\frac{1}{2}\[{\left(\frac{2x_i}{a}\right)}^2+{\left(\frac{2y_i}{b}\right)}^2\]\}\·\frac{{\mathrm{\πR}}_i^2{\mathrm{\β}}_T^2}{4}$

${\mathrm{\σ}}_i\left(t\right)=\frac{4\mathrm{\π}}{{\mathrm{\Ω}}_i\left(t\right)}{\mathrm{\ρ}}_i\left(t\right)A_{S,i}\left(t\right)$

其中，p_R(t)是地物回波波形，p_T(t)是发射脉冲功率，R_i是地物与发射机之间的距离，R₀为所有地与激光高度计的近似距离，β_T是发射脉冲的波束宽度，ρ_i(t)是地物反射率，σ_i后向散射截面，A是地物受到脉冲能量照射的面积，Ω_i(t)地物散射立体角，D_R是接受天线孔径直径，η_atm表示大气对信号传输的影响，η_sys表示硬件系统对信号传输的影响，a、b是脚印内椭圆的长半轴和短半轴，x_i、y_i分别是地物i离脚印长轴、短轴的垂直距离，即在以脚印中心为原点，以椭圆短轴为x轴，以椭圆长轴为y轴的坐标系下，地物i的坐标值。

激光作为一种光源，其光束截面内的光强分布是不均匀的，即光束波面上各点的振幅是不相等的。参数ε_i是考虑到脉冲能量从脚印中心向周围逐渐衰减，呈椭圆高斯分布而引入的。以上模拟回波波形的过程，是假设建筑物只有屋顶存在能量入射，侧面没有能量入射。

激光测距仪激光回波全波形分析算法

卫星激光测距雷达通过测定激光发射到地面的往返时间来确定其到地面的距离，进而精确测定地物的高度。为了精确测量距离，需要对脉冲信号进行波形拟合，精确确定回波波形的中心位置。模型拟合可以估计平均表面高，脉冲振幅，脉冲宽度及信号噪声。

波形拟合算法有非线性最小二乘最优算法、最大似然估计算法、蒙特卡洛马尔科夫链方法等。发射激光脉冲是高斯波，可以把回波看成是带有偏差的高斯波型。高斯分解实质上是一个曲线拟合的过程，关键是数学模型和参数初始值的选择，因为这两个过程决定了模型对回波的分解程度，以及分解出来的高斯波的位置和形状。

波形模型可写为：

$w\left(t\right)=\underset{m=1}{\overset{N_p}{\Σ}}{A}_m\exp \[-\frac{{(t-t_m)}^2}{2{\mathrm{\σ}}_m^2}\]+\mathrm{\ϵ}$

其中，N_p是回波中峰值个数，A_m是第m个峰值的大小，ε是噪声影响而产生的偏差，t_m是峰值的位置，σ_m是第m个峰值的标准差。可由非线性最小二乘法计算模型参数。提取地物高度信息首先要确定地面回波波形的高斯波个数，而噪声对此影响很大，需要合理有效的滤波算法进行去噪。

实施例1：

平地不同反射率：

平地返回的是单个回波，对地面反射率从0.65～0.06的共8种回波分别进行了波形拟合。模拟600次求平均的结果表明，测距精度随着地面反射率降低而减小，最高精度为0.068±0.05cm，最低精度为0.73±0.56cm。拟合的参数值见表1，拟合波形见图2和图3。

表5.1平地不同反射率回波拟合参数值

实施例2：

坡度：

对不同坡度的回波信号进行波形拟合，测距精度随着坡度增加而减小，随反射率的减小而减小。模拟600次求平均的结果表明，坡度5度，测距精度为0.11±0.08cm(ρ＝0.65)，和1.46±0.82cm(ρ＝0.06)；坡度15度，测距精度为0.43±0.34cm(ρ＝0.65)和2.94±2.18cm(ρ＝0.06)，拟合参数结果见表2，波形拟合结果见图4和图5。

表2不同坡度回波拟合参数值

实施例3：

一个建筑物距足印中心不同远近地物：

地物到足印中心的距离越远(0米，5米，10米，12米)，即在足印光斑内所占面积越小，波形拟合精度越差，测距误差越大。模拟600次求平均的结果表明，当距离足印中心0m时，测距误差约0.15±0.11cm(ρ＝0.65)和0.91±0.68cm(ρ＝0.06)；当距离足印中心10m时，测距误差约1.24±0.96cm(ρ＝0.65)和8.04±5.93cm(ρ＝0.06)；当距离足印中心12m时，测距误差约3.3±2.3cm(ρ＝0.65)，当地面反射率为0.06时，地物的回波能量很小为1.69mv，小于噪声阈值(4mv)，和噪声混淆在一起，很难将其探测拟合出来，因此如果回波能量小于噪声阈值，则不进行拟合。600次模拟中，只有6次探测出了地物微弱的回波，测距误差较大，约为17.07cm。拟合参数值见表3，Np1是地物相关参数，Np2是地表的相关参数，波形拟合结果见图6和图7。

表3地物到足印中心的距离不同的回波拟合参数值

实施例4：

一个建筑物距足印中心不同远近地物：

足印内有多个地物，波形拟合精度与地物高度、面积以及在足印中所处的位置有关。模拟600次求平均结果表明，精度最好的为3个地物中的Np1，测距误差为0.5394±0.3802cm(ρ＝0.65)和3.434±2.312cm(ρ＝0.06)，因为对应高10m面积为10.1*9.1m2的建筑物，返回的能量较强；如果Np能量接近噪声阈值(4mv)，则测距精度较差，例如地面反射率为0.06时，3个地物的Np2测距误差为9.494±6.415cm，Np3为8.954±6.756cm，2个地物的Np1测距误差为9.454±7.236cm，Np2为8.554±6.31cm，1个地物的Np1，测距误差为8.8267±6.6961cm。拟合精度较差是因为其回波信号较弱，拟合时峰值的确定以及开始结束位置的确定不精确，导致拟合结果较差。拟合参数值见表4。波形拟合结果见图8和图9。

表4足印内不同地物个数的回波拟合参数值

本发明使用非线性最小二乘高斯分解算法，求得回波信号的每个回波的中心位置、峰值、标准差和噪声偏差。多个地物的拟合波形是多个回波的叠加。综合所有实施例，回波信号越强，波形拟合的越精确，测距误差越小，波形拟合使用合适的初始值，得到的结果越精确，测距误差越小。

Claims (7)

1.一种基于回波全波形分解的地物特征获取方法，其特征在于，包括以下步骤：

1)根据激光高度计的参数获取激光高度计的激光发射信号的数学模型；

2)获取目标地物的三维地物模型以及激光高度计发射传输过程的数据误差；

3)根据三维地物模型以及激光发射信号模型并考虑激光高度计发射传输过程的数据误差，获取激光发射信号与目标地物相互作用产生回波的数学模型；

4)对激光发射信号与目标地物相互作用产生的回波采用非线性最小二乘高斯分解算法进行全波形分解，获得分解后的高斯波；

5)根据分解后的高斯波的波形获取目标地物的垂直高度，并估计激光高度计的实际测量精度。

2.根据权利要求1所述的一种基于回波全波形分解的地物特征获取方法，其特征在于，所述的步骤1)中激光发射信号波形的数学模型为高斯函数：

$p_T\left(t\right)=A_T\exp \{-\frac{{(t-t_{C,T})}^2}{2{s_T}^2}\}$

其中，p_T(t)为高斯函数的峰值，t_C,T为高斯函数的中心位置，s_T为高斯函数的标准差。

3.根据权利要求1所述的一种基于回波全波形分解的地物特征获取方法，其特征在于，所述的步骤2)中通过数字表面模型表示地表和地物表面的起伏和分布三维地物模型。

4.根据权利要求1所述的一种基于回波全波形分解的地物特征获取方法，其特征在于，所述的步骤2)中数据误差包括激光高度计发射和接收不平行造成的误差、脉冲信号传播时间的测量误差、震动误差、反光镜的旋转误差、脉冲零点误差、大气中水分子、悬浮物和大气的折射、散射造成的测量误差以及大气温度、湿度和气压变化造成的测量误差。

5.根据权利要求1所述的一种基于回波全波形分解的地物特征获取方法，其特征在于，所述的步骤3)中激光发射信号与目标地物相互作用产生回波的数学模型为：

$p_R\left(t\right)\≈\underset{i=1}{\overset{N}{\Σ}}\{\left({\mathrm{\ϵ}}_i\frac{4p_T\left(t\right)}{{\mathrm{\πR}}_i^2{\mathrm{\β}}_T^2}\right)*{\mathrm{\σ}}_i\left(t\right)\}\·\frac{1}{4{\mathrm{\πR}}_0^2}\·\frac{{\mathrm{\πD}}_R^2}{4}\·{\mathrm{\η}}_{atm}\·{\mathrm{\η}}_{sys}$

${\mathrm{\ϵ}}_i=\frac{2}{\mathrm{\π}ab}\·\exp \{-\frac{1}{2}\[{\left(\frac{2x_i}{a}\right)}^2+{\left(\frac{2y_i}{b}\right)}^2\]\}\·\frac{{\mathrm{\πR}}_i^2{\mathrm{\β}}_T^2}{4}$

${\mathrm{\σ}}_i\left(t\right)=\frac{4\mathrm{\π}}{{\mathrm{\Ω}}_i\left(t\right)}{\mathrm{\ρ}}_i\left(t\right)A_{S,i}\left(t\right)$

其中，p_R(t)为地物回波波形，p_T(t)为发射脉冲功率，R_i为地物与激光高度计之间的距离，R₀为所有地物与激光高度计的近似距离，β_T为发射脉冲的波束宽度，ρ_i(t)为地物反射率，σ_i为后向散射截面，A_S,i(t)为地物受到脉冲能量照射的面积，Ω_i(t)地物散射立体角，D_R为接收天线孔径直径，η_atm表示大气对信号传输的影响参数，η_sys表示硬件系统对信号传输的影响参数，a、b为激光测距雷达激光光束的照射大脚印内椭圆的长半轴和短半轴，x_i、y_i分别为第i个地物的坐标值，ε_i为激光测距雷达发射的脉冲能量呈椭圆高斯分布引入的参数。

6.根据权利要求1所述的一种基于回波全波形分解的地物特征获取方法，其特征在于，所述的步骤4)中分解后的高斯波的波形数学模型为：

$w\left(t\right)=\underset{m=1}{\overset{N_p}{\Σ}}{A}_m\exp \[-\frac{{(t-t_m)}^2}{2{\mathrm{\σ}}_m^2}\]+\mathrm{\ϵ}$

其中，N_p是回波中峰值个数，A_m是第m个峰值的大小，ε是噪声影响而产生的偏差，t_m是峰值的位置，σ_m是第m个峰值的标准差。

7.根据权利要求1所述的一种基于回波全波形分解的地物特征获取方法，其特征在于，所述的步骤5)中估计激光测距雷达的实际测量精度估计方法为：

采用分解后的高斯波的波形数学模型对脉冲信号进行波形拟合，并确定回波波形的中心位置，估计精度RMSE的表达式为：

$RMSE=\sqrt{\frac{\mathrm{\Σ}{(h-h_{set})}^2}{n}}$

h＝c*(t₁-t₂)/2

其中，h为通过回波波形拟合计算得到的地物的垂直高度，h_set为建模时设定的地物垂直高度，即真值，n为仿真次数，c为光速，t₁和t₂为相邻脉冲的时刻。