激光雷达系统数据处理方法
本发明属于用于气象应用的激光雷达领域，特别是一种用于气象应用的激光雷达系统数据处理方法。
背景技术：
大气探测激光雷达是当今大气和气象观测的一个重要工具。利用特定的激光雷达可以获取距地面数百米至数公里处的风廓线分布。激光雷达系统的数据处理算法会极大程度地影响激光雷达系统在探测精确度、时间分辨率和有效数据获取率等指标上的性能表现。评价激光雷达系统数据处理方法的指标，主要包括处理的所需的时间、处理结果的精确度、低信噪比下结果的正确率、对扫描方向数量的需求。目前已有的数据处理算法，主要是dswf、fswf、mfas算法。dswf算法速度快，适合实时探测的激光雷达数据处理，但是在低信噪比下结果的正确率很低；fswf、mfas算法在低信噪比下结果的正确率都较高，但是需要大量的搜索或迭代运算，不适用于实时探测的激光雷达数据处理。随着激光器、探测器和数字化仪等系统组件性能的提升，激光雷达系统在单位时间内获取的回波数据数量有了大幅增长。对于风廓线的实时观测应用，亟需一种能够兼具计算耗时较少、在低信噪比条件下正确率高——这两项优点的数据处理方法。
在背景技术部分中公开的上述信息仅仅用于增强对本发明背景的理解，因此可能包含不构成在本国中本领域普通技术人员公知的现有技术的信息。
技术实现要素：
针对现有技术中存在的问题，本发明提出一种激光雷达系统数据处理方法，利用激光雷达系统参数和探测获取的回波信号频谱结果，快速计算风廓线，且在低信噪比下仍有较高的正确率。
本发明的目的是通过以下技术方案予以实现：
一种激光雷达系统数据处理方法包括以下步骤：
第一步骤中，激光雷达系统采集时域信号，按照距离门划分时域信号并计算得到回波信号频谱，对回波信号频谱进行处理得到第一峰值频率和信噪比，利用第一峰值频率、信噪比以及激光雷达系统参数计算有效风速估计的误差下限和第一径向风速估计值；
第二步骤中，对回波信号频谱进行处理，并计算出搜索频谱峰值位置时所被限定的频带范围，利用被限定的频带范围和处理后的频谱计算第二峰值频率，利用第二峰值频率以及激光雷达系统参数计算第二径向风速估计值；
第三步骤中，利用所述误差下限和所述第二径向风速估计值，以高斯分布函数计算第一径向风速估计值对应的权重系数；
第四步骤中，利用多个径向上的第一径向风速估计值、所述权重系数和激光雷达系统参数，使用加权正弦波拟合算法计算矢量风速结果。
所述的方法中，回波信号频谱可以为功率谱、幅度谱，或是多次脉冲回波的频谱的累加或平均。
所述的方法中，回波信号频谱包括大气中目标微粒的多普勒回波频谱以及噪声的频谱。
所述的方法中，第一步骤中，所述回波信号频谱的处理是频谱噪声基线的平坦化，利用平坦化后的频谱计算第一峰值频率和信噪比，基于所述回波第一峰值频率、信噪比以及激光雷达系统参数计算有效风速估计的误差下限和第一径向风速估计值。
所述的方法中，第一步骤中，计算有效风速估计的误差下限采用的参数包括激光波长、对回波信号进行模数转换的采样率、回波信号频谱累加(或平均)脉冲数、距离门内采样点数、归一化激光脉冲频谱宽度和信噪比。
所述的方法中，第一步骤中，计算所述第一径向风速估计值的参数包括激光波长、激光器出射激光与本振激光频率差值、发射系统与接收系统的夹角和第一峰值频率。
所述的方法中，第二步骤中，所述回波信号频谱的处理包括依次进行谱拼接、二维低通滤波、一维二阶差分、标准化、二值化和连通域查找，搜索频谱峰值位置时所被限定的频带范围是根据连通域查找的结果划定的频带范围，利用一维二阶差分或标准化的频谱计算在限定的频带范围内的第二峰值频率.
所述的方法中，第二步骤中，计算所述第二径向风速估计值的参数包括激光波长、激光器出射激光与本振激光频率差值、发射系统与接收系统的夹角和第二峰值频率。
所述的方法中，第三步骤中，高斯分布函数以第一径向风速估计值为输入自变量，利用有效风速估计的误差下限计算标准差，以第二径向风速估计值为期望，得到的概率密度值为权重系数的计算结果。
所述的方法中，所述激光雷达系统参数至少包括激光波长、激光器出射激光与本振激光频率差值、发射系统与接收系统的夹角、对回波信号进行模数转换的采样率、回波信号频谱累加脉冲数或平均脉冲数、距离门内采样点数、归一化激光脉冲频谱宽度和激光探测指向的方位角与俯仰角。
和现有技术相比，本发明具有以下优点：
本发明在处理数据时花费时间少，且在低信噪比条件下的结果仍具较高的正确率。
附图说明
通过阅读下文优选的具体实施方式中的详细描述，本发明各种其他的优点和益处对于本领域普通技术人员将变得清楚明了。说明书附图仅用于示出优选实施方式的目的，而并不认为是对本发明的限制。显而易见地，下面描述的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。而且在整个附图中，用相同的附图标记表示相同的部件。
在附图中：
图1是本发明的激光雷达系统数据处理方法的步骤示意图；
图2是本发明一个具体实施例的平坦化和查找到的峰值的示例结果；
图3是本发明的一个实施例的第一步骤的第一径向风速估计值和有效风速估计的误差下限的示例结果；
图4是本发明的一个实施例的第二步骤的谱拼接的示例结果；
图5是本发明的一个实施例的第二步骤的二维低通滤波的示例结果；
图6是本发明的一个实施例的第二步骤的一维二阶差分的示例结果；
图7是本发明的一个实施例的第二步骤的二值化、连通域查找和搜索频谱峰值位置时所被限定的频带范围的示例结果；
图8是本发明的一个实施例的第二步骤的径向风速估计值的示例结果；
图9是本发明的一个实施例的第三步骤的权重系数的示例结果；
图10是本发明的一个实施例的第四步骤的矢量风速的示例结果与其他两种常见算法的对比结果。
以下结合附图和实施例对本发明作进一步的解释。
具体实施方式
下面将参照附图图1至图10更详细地描述本发明的具体实施例。虽然附图中显示了本发明的具体实施例，然而应当理解，可以以各种形式实现本发明而不应被这里阐述的实施例所限制。相反，提供这些实施例是为了能够更透彻地理解本发明，并且能够将本发明的范围完整的传达给本领域的技术人员。
需要说明的是，在说明书及权利要求当中使用了某些词汇来指称特定组件。本领域技术人员应可以理解，技术人员可能会用不同名词来称呼同一个组件。本说明书及权利要求并不以名词的差异来作为区分组件的方式，而是以组件在功能上的差异来作为区分的准则。如在通篇说明书及权利要求当中所提及的“包含”或“包括”为一开放式用语，故应解释成“包含但不限定于”。说明书后续描述为实施本发明的较佳实施方式，然所述描述乃以说明书的一般原则为目的，并非用以限定本发明的范围。本发明的保护范围当视所附权利要求所界定者为准。
为便于对本发明实施例的理解，下面将结合附图以具体实施例为例做进一步的解释说明，且各个附图并不构成对本发明实施例的限定。
一个具体的激光雷达系统数据处理方法实施例中，该方法包括以下步骤：
第一步骤s1中，激光雷达系统采集时域信号，按照距离门划分时域信号并计算得到回波信号频谱，对回波信号的频谱进行处理得到频谱参数，利用频谱参数，即第一峰值频率、信噪比，以及激光雷达系统参数计算有效风速估计的误差下限σv和第一径向风速估计值
第二步骤s2中，对回波信号频谱进行处理，并计算出搜索频谱峰值位置时所被限定的频带范围，利用被限定的频带范围和处理后的频谱计算第二峰值频率，利用频谱参数以及激光雷达系统参数计算第二径向风速估计值
第三步骤s3中，利用所述误差下限和所述第二径向风速估计值，以高斯分布函数计算第一径向风速估计值对应的权重系数w；
第四步骤s4中，利用多个径向上的第一径向风速估计值
本发明通过对回波信号频谱进行处理，计算出搜索频谱峰值位置时所被限定的频带范围。限定后的范围是依据信号自身分布确定的，既没有减小风速探测范围，又提高了低信噪比下查找到有效的风速估计的几率。利用被限定的频带范围和处理后的频谱计算频谱参数，利用频谱参数以及激光雷达系统参数计算第二径向风速估计值；利用所述误差下限和所述第二径向风速估计值，以高斯分布函数计算第一径向风速估计值对应的权重系数，而无需迭代或查找，计算速度快。这一方法相比于其他常用算法兼具了计算耗时较少、在低信噪比条件下正确率高的优点。
在一个具体实施例中，所述激光雷达系统参数至少包括激光波长λ、激光器出射激光与本振激光频率差值faom、发射系统与接收系统的夹角φ、对回波信号进行模数转换的采样率fs、回波信号频谱累加脉冲数n、距离门内采样点数m、归一化激光脉冲频谱宽度wn和激光探测指向的方位角θi和俯仰角
在一个具体实施例中，激光雷达系统为大气探测激光雷达系统。
在一个具体实施例中，回波信号频谱指，利用雷达采集的时域信号，按照距离门划分后，利用周期图法计算得到的频谱结果。具体形式是累加后的功率谱。包括大气中目标微粒的多普勒回波频谱以及噪声的频谱。
在一个具体实施例中，径向风速估计值是沿激光探测指向的风速，
在一个具体实施例中，搜索频谱峰值是搜索在指定范围里的频谱的最大值所在的频谱的峰，频谱的峰是最大值所在的频谱突起，峰值频率是频谱的峰所在的频率，多普勒频移是峰值频率相对于出射激光与本振激光频率差值的移动量。
在一个具体实施例中，有效的风速估计是查找到的频谱的峰是大气中目标微粒的多普勒回波频谱时，对应的径向风速估计结果，无效的风速估计是查找到的频谱的峰是噪声频谱时，对应的径向风速估计结果。在实际中很难分辨在低信噪比区域得到的频谱的峰是噪声谱中的起伏，还是大气中目标微粒的多普勒回波频谱。
有效风速估计的误差下限是指，即使径向风速估计值是有效的，也不代表准确的风速值，在计算过程中引入的有效的风速估计值与真值统计误差在理论上大于某个值，而实际系统的统计误差只会大于此值，这个值就是误差下限。对于利用周期图法计算的频谱，误差下限可通过下式计算
式中snr是当前距离门的信噪比，λ是激光波长，fs是采样率，n是累加脉冲数，m是单个距离门内的采样点数，wn是归一化谱宽。
在一个具体实施例中，第一步骤s1中，所述回波信号频谱的处理是进行频谱噪声基线的平坦化。平坦化时使用回波信号频谱中的预采集信号的频谱作为先验噪声基线，并从后续距离门的功率谱中扣除。利用回波信号频谱计算的参数至少包括信噪比snr和第一峰值频率fpeak1。对一个距离门的结果平坦化前后的对比如图2所示，可以看出，未平坦化的功率谱中包含明显的噪声基线，如果直接查找最大值，频谱的峰将会在边界处出现。图中结果的第一峰值频率fpeak1＝77.8mhz。信噪比是频谱的峰的面积与噪声面积之比。
在一个具体实施例中，第一步骤s1中，有效风速估计的误差下限的计算需要的系统参数包括对回波信号进行模数转换的采样率fs、回波信号频谱累加脉冲数n、距离门内采样点数m、归一化激光脉冲频谱宽度wn。
在一个具体实施例中，第一步骤s1中，有效风速估计的误差下限的计算需要的根据回波信号频谱计算的参数包括信号的信噪比snr。这里使用的信噪比是宽带信噪比。
在一个具体实施例中，第一步骤s1中，计算有效风速估计的误差下限使用如下的公式：
其中，k是距离门序号，i是方位角序号。
在一个具体实施例中，第一步骤s1中，第17个距离门上的有效风速估计的误差下限的计算结果如图3所示。在图中，横坐标代表激光探测指向的方位角，在每一个方位角上都获取了32个距离门的结果，图中只展示了每个径向第17个距离门的结果。在方位角90°上获取的结果的效风速估计的误差下限最大，但是与其余结果的差距不大，只根据误差下限难以将这一结果排除。
在一个具体实施例中，第一步骤s1中，计算所述第一径向风速估计值的激光雷达系统参数包括激光波长λ、激光器出射激光与本振激光频率差值faom、发射系统与接收系统的夹角φ，计算所述第一径向风速估计值的采用的频谱参数包括第一峰值频率fpeak1。
在一个具体实施例中，第一步骤s1中，计算所述第一径向风速估计值使用如下的公式：
在一个具体实施例中，第一步骤s1中，第17个距离门上的不同方位角上的第一径向风速估计值的结果如图3所示。可以看出，虽然第一径向风速估计值是公认的当前距离门的无偏估计，但是在方位角90°上获取的结果计算的径向风速估计与其余方向相差较大。
在一个具体实施例中，第二步骤s2中，对回波信号频谱进行的处理包括谱拼接、二维低通滤波、一维二阶差分、标准化、二值化和连通域查找。第二步骤中的处理是对逐个径向的功率谱依次进行操作的，因此以下的处理步骤中省略了代表径向的索引i。单个距离门的累加功率谱记为psdraw(j)，
与psdraw进行卷积：
psdg(j，k)＝hg(s，t)*psdraw(j，k)；
结果如图5所示。在经过二维低通滤波后，频谱中的噪声有了明显的减弱，在原本被噪声淹没的距离门上，也能观察到连贯的频谱的峰，说明经过这一处理，有机会在低信噪比的区域计算出有效的径向风速估计。一维二阶差分是逐个距离门地依次进行二阶差分，可以记为使用hd＝[-1，2，-1]t与psdg进行卷积：
psdd(j，k)＝hd*psdg(j，k)；
结果如图6所示，可以看出，一维二阶差分的步骤中没有扣除先验噪声基线，但是得到的结果中已经没有了明显的噪声基线，同时在频谱的峰两侧形成了“深谷”，突出了频谱的峰。标准化按照均值为0，标准差为1的标准，逐距离门地进行；二值化使用大津法计算整体的阈值并进行二值化。二值化的结果如图7所示，深灰色部分是二值化结果为0的部分，其余代表的结果为1。
在一个具体实施例中，第二步骤s2中，搜索频谱峰值位置时所被限定的频带范围是根据连通域查找的结果划定的，对回波信号频谱进行处理后，利用一维二阶差分或标准化后的频谱在限定的频带范围内计算第二多普勒频移。在一个具体实施例中，被限定的频带范围是被限制在与第一个距离门的结果有相交的连通域内，同时连通域分叉后的结果从搜索带宽中剔除。被限定的频带范围的结果展示在图7中，以“*”标注的中灰色区域为被限定的频带范围，可以看出，经过二值化和连通域查找，第二峰值频率查找的限定的频带范围相比于原始的搜索范围有了很大的减小，这就显著提高了查找到有效的风速估计的几率，同时，这一区域是根据频谱的峰的分布在原始的搜索范围中确定的，也就没有影响激光雷达的风速探测范围。
在一个具体实施例中，第二步骤s2中，径向风速估计值的计算需要的系统参数包括激光波长λ、激光器出射激光与本振激光频率差值faom、发射系统与接收系统的夹角φ，计算所述第一径向风速估计值的采用的频谱参数包括第二峰值频率fpeak2。在一个具体实施例中，第二峰值频率fpeak2是使用一维二阶差分后的功率谱的信号峰在搜索范围内查找的。
在一个具体实施例中，第二步骤s2中，计算所述第二径向风速估计值使用如下的公式：
在一个具体实施例中，第二步骤s2中，第17个距离门上的径向风速估计值的结果如图8所示。可以看出经过第二步骤s2的处理，处在方位角90°上获取的结果计算的径向风速估计与其余方向的差距较小。
在一个具体实施例中，第三步骤s3中，高斯分布函数以第一径向风速估计值为输入自变量，利用有效风速估计的误差下限计算标准差，以第二径向风速估计值为期望，使用的权重系数的计算公式为：
在一个具体实施例中，第三步骤s3中，第17个距离门上的权重系数结果如图9所示。当方位角为90°时，以第二径向风速估计为期望，与第一径向风速估计值的差大于10m/s，远大于半宽度，因此权重系数为0。
在一个具体实施例中，第四步骤s4中，使用8个径向上的结果进行加权正弦波拟合计算得到矢量风速估计值。需要的系统参数包括激光探测指向的方位角θi和俯仰角
在一个具体实施例中，第四步骤s4中，第17个距离门上的矢量风速估计值被重新分解到各个径向上，以便于理解算法的输出与输入径向风速估计值的关系，并解读各个算法的性能表现。分解后的结果如图10所示。dswf算法的结果明显地受到错误结果的影响，拟合结果出现了向上的偏移，与大部分径向风速估计值之间的差距都较大。本算法的计算结果和fswf算法的计算结果则接近大部分第一径向风速估计值，可以认为合成的矢量风速估计值较为正确。
在一个具体实施例中，算法的用时使用处理2000个径向结果的用时来评价，dswf算法需要花费2.959s，fswf算法需要花费359.191s，而本算法需要使用7.421s。可以看出，本发明在处理数据时花费时间少，且在低信噪比条件下的结果仍具有较高的正确率。
- 还没有人留言评论。精彩留言会获得点赞！