简介:这份PDF文献聚焦激光雷达探测大气气溶胶的数据处理研究,面向大气科学、环境监测及遥感方向的学习者与科研人员,帮助理解米散射激光雷达的系统构成与反演算法原理。资源包内含1个PDF文件,大小约180KB,内容源自期刊论文,涵盖激光发射单元、接收望远镜、Si:APD单光子计数器探测以及回波信号分析等关键环节。文中详细介绍了云和气溶胶消光系数及衰减后向散射系数的反演方法,包括斜率法、Klett法与Fernald法等常用策略,并给出24小时连续观测的处理结果,便于读者对照激光雷达方程理解重叠因子修正与参数计算流程。目前已有208人学习,适合需要夯实激光雷达数据处理基础、撰写相关论文或开展大气探测实验的研究者参考借鉴。
1. 激光雷达气溶胶数据处理:从原始回波到可信廓线的关键一跃
拿到一台激光雷达,最让人头疼的往往不是硬件调试,而是后面那一长串原始数据。尤其是做大气气溶胶测量时,回波信号里混着背景光、几何重叠因子、距离平方衰减,还有各种噪声,直接画出来的曲线根本没法看。激光雷达测量大气气溶胶的数据处理研究,核心就是解决从原始光子计数或模拟回波,到最终气溶胶消光系数、后向散射系数廓线的全过程。这套流程决定了你能不能从一台几万块的激光雷达里,榨出真正有物理意义的数据。适合谁看?做大气环境监测的、搞气象观测的、以及刚接手激光雷达数据的研究生。如果你手里有数据但不知道怎么处理,或者处理完发现负值满天飞、边界值乱跳,那这篇笔记就是给你写的。我会把整个数据处理框架拆开,从信号预处理到反演算法,再到踩过的坑,一步步讲清楚。
2. 原始回波信号预处理:把噪声和背景光先摁住
激光雷达原始信号通常有两种形式:模拟信号和光子计数信号。不管哪种,第一步都是扣除背景噪声。背景噪声的来源很多,太阳光、探测器暗电流、大气分子散射都会贡献。常见做法是取远距离末端的一段信号做平均,作为背景基线。但这里有个细节:如果大气边界层很高,远端信号可能还没完全衰减到背景水平,这时候直接扣背景就会把有效信号也扣掉。我一般会先画一遍距离校正信号,看看远端是否平坦,再决定背景段的位置。
2.1 背景扣除与距离平方校正
背景扣除之后,紧接着就是距离平方校正。激光雷达方程里信号强度随距离平方衰减,所以要把接收到的信号乘以距离的平方,才能还原出大气的后向散射特性。这一步看起来简单,但距离的起点选错,整个廓线形状都会变。距离起点应该是激光出射点到接收望远镜光轴的几何交点,而不是简单的望远镜位置。很多商用激光雷达会在元数据里给出这个值,如果没有,就需要用重叠因子反推。
import numpy as np def preprocess_lidar_signal(raw_signal, distance, bg_start_idx, bg_end_idx): """ 激光雷达原始信号预处理 raw_signal: 原始回波信号数组 distance: 对应的距离数组,单位米 bg_start_idx: 背景段起始索引 bg_end_idx: 背景段结束索引 """ # 计算背景均值 bg_mean = np.mean(raw_signal[bg_start_idx:bg_end_idx]) # 扣除背景 signal_bg_corrected = raw_signal - bg_mean # 避免负值影响后续对数运算 signal_bg_corrected = np.maximum(signal_bg_corrected, 0) # 距离平方校正 range_corrected = signal_bg_corrected * distance ** 2 return range_corrected, bg_mean这段代码里,bg_start_idx和bg_end_idx的选择直接决定背景扣除的效果。通常我会选最远端的 10% 到 20% 的数据点,但前提是确认这段信号已经衰减到接近零。如果远端信号还有明显起伏,说明背景段选得太近,需要往后挪。np.maximum那一步是为了防止扣除背景后出现负值,虽然物理上不应该有负值,但实际数据里噪声会导致负值出现,直接取对数会报错。
2.2 重叠因子校正与信号平滑
几何重叠因子是激光雷达近距离信号失真的主要原因。发射光束和接收视场在近距离没有完全重合,导致近场信号被低估。校正方法有实验法和理论计算法。实验法一般用水平均匀大气假设,通过比较不同距离的信号斜率来反推重叠因子。理论计算则需要知道激光发散角、望远镜视场角、两者间距等参数。我一般先用实验法快速评估,如果重叠因子在几百米内就接近 1,那后续反演可以忽略这段;如果重叠区域延伸到 1 公里以上,就必须做校正。
信号平滑是另一个容易翻车的地方。平滑窗口太宽,会把气溶胶层的精细结构抹掉;窗口太窄,噪声又压不住。常见做法是用滑动平均或者小波变换。滑动平均简单,但会引入相位偏移;小波变换能保留突变特征,但参数不好调。我的经验是,先用 5 到 9 点的滑动平均试一下,如果气溶胶层边界变模糊了,就换小波。平滑之后一定要检查信噪比,如果平滑后信噪比还是低于 3,那这段数据基本不可用。
3. 气溶胶消光系数反演:Klett 法和 Fernald 法怎么选
预处理完的信号,下一步就是反演消光系数。激光雷达方程里有两个未知数:消光系数和后向散射系数,直接求解是欠定的。所以需要假设两者之间的关系,也就是激光雷达比。Klett 法和 Fernald 法是两种最常用的反演方法。Klett 法假设后向散射系数和消光系数成幂律关系,适合气溶胶为主的情况;Fernald 法把分子散射和气溶胶散射分开处理,需要知道分子消光系数,适合边界层以上气溶胶较少的场景。
3.1 Klett 反演法的参数设置与边界值选择
Klett 法的核心公式里,边界值的选择至关重要。边界值通常选在远端,假设那里大气均匀,消光系数已知。如果边界值选得太远,信号噪声会放大;选得太近,又可能把气溶胶层截断。我一般会先画距离校正信号的对数曲线,找一段斜率稳定的区域作为边界。边界值的大小可以参考大气能见度或者太阳光度计数据,如果没有,就用经验值,比如 1e-5 到 1e-4 每米。
def klett_inversion(range_corrected, distance, lidar_ratio, boundary_value, boundary_idx): """ Klett 法反演气溶胶消光系数 range_corrected: 距离平方校正后的信号 distance: 距离数组 lidar_ratio: 激光雷达比,典型值 30-70 sr boundary_value: 边界处的消光系数 boundary_idx: 边界点索引 """ # 取对数 log_signal = np.log(range_corrected) # 计算积分项 integral = np.cumsum(range_corrected) * (distance[1] - distance[0]) # 边界处的积分值 integral_boundary = integral[boundary_idx] # 反演消光系数 extinction = range_corrected / (range_corrected[boundary_idx] / boundary_value - 2 * lidar_ratio * (integral - integral_boundary)) return extinction这里lidar_ratio的取值直接影响反演结果。气溶胶类型不同,激光雷达比差别很大。城市气溶胶通常在 40 到 60 之间,沙尘气溶胶可以到 50 以上,海洋气溶胶偏低。如果不知道具体类型,先用 50 试,然后根据反演出的消光系数廓线是否合理来调整。boundary_value和boundary_idx需要配合使用,边界点一般选在信号信噪比还不错的远端,比如 3 到 5 公里处。
3.2 Fernald 法分离分子与气溶胶散射
Fernald 法的思路是把分子散射和气溶胶散射分开。分子消光系数可以通过标准大气模型计算,比如美国标准大气。气溶胶消光系数则通过迭代求解。Fernald 法对边界值的要求比 Klett 法更敏感,因为分子散射的贡献在远端占比更大。如果边界值选得不对,反演出的气溶胶消光系数会出现负值。
def fernald_inversion(range_corrected, distance, molecular_extinction, molecular_backscatter, lidar_ratio_aerosol, boundary_value, boundary_idx): """ Fernald 法反演气溶胶消光系数 molecular_extinction: 分子消光系数廓线 molecular_backscatter: 分子后向散射系数廓线 lidar_ratio_aerosol: 气溶胶激光雷达比 """ # 分子激光雷达比 lidar_ratio_molecular = 8 * np.pi / 3 # 计算分子后向散射与消光比 ratio_molecular = molecular_backscatter / molecular_extinction # 迭代求解 extinction_aerosol = np.zeros_like(distance) extinction_aerosol[boundary_idx] = boundary_value for i in range(boundary_idx - 1, -1, -1): # 简化迭代公式,实际使用时需要根据具体文献调整 numerator = range_corrected[i] * np.exp(-2 * (lidar_ratio_aerosol - lidar_ratio_molecular) * molecular_extinction[i] * (distance[i+1] - distance[i])) denominator = range_corrected[boundary_idx] / boundary_value + 2 * lidar_ratio_aerosol * numerator extinction_aerosol[i] = numerator / denominator return extinction_aerosolFernald 法的迭代方向是从边界点向近端推进,所以边界点的选择决定了整个廓线的基准。如果边界点选在气溶胶层内部,反演结果会严重失真。我一般会选在气溶胶层以上、分子散射为主的区域,比如 5 到 8 公里。molecular_extinction和molecular_backscatter可以用标准大气模型算,也可以用地基微波辐射计或者探空数据。如果没有实测数据,用标准大气也能凑合,但精度会打折扣。
4. 数据处理中的避坑指南:那些让你白干一整天的细节
做激光雷达数据处理,最怕的不是算法复杂,而是细节没注意,结果全错。下面这几条是我和同行们踩过的坑,每条都按现象、原因、解决来写。
4.1 避坑一:背景扣除后信号出现大面积负值
现象:扣除背景后,远端信号变成负值,距离平方校正后负值更明显。原因:背景段选得太近,把还有效的信号当成了背景。或者探测器饱和导致远端信号被截断,背景均值算出来偏高。解决:先画原始信号,确认远端是否平坦。如果远端有起伏,把背景段往后挪。如果探测器饱和,需要换用低增益通道或者加衰减片重新测量。
4.2 避坑二:反演出的消光系数出现负值
现象:Klett 或 Fernald 反演后,某些距离上的消光系数是负数。原因:边界值选得太小,或者激光雷达比设得太大。也可能是信号预处理时平滑过度,把真实信号抹掉了。解决:先检查边界值,用太阳光度计或者能见度数据校准。如果边界值没问题,把激光雷达比调小 10% 到 20% 再试。平滑窗口不要超过 9 点,否则会引入虚假的负值。
4.3 避坑三:重叠因子校正后近场信号反而更差
现象:做了重叠因子校正,近场信号反而出现异常峰值或者凹陷。原因:重叠因子曲线是用水平大气假设反推的,如果实际大气不均匀,反推的重叠因子就不准。或者校正时距离起点没对齐,导致校正曲线偏移。解决:先用水平均匀天气的数据做重叠因子,比如清晨或者阴天。校正时确保距离起点和重叠因子曲线的起点一致。如果近场信号还是不对,干脆把 500 米以内的数据标记为不可用,不要强行校正。
4.4 避坑四:不同时间的数据拼在一起出现断层
现象:把不同时间测的廓线拼成时间序列,发现相邻时刻消光系数跳变。原因:每次测量的背景噪声不一样,或者激光能量有波动。也可能是反演时边界值每次都在变。解决:每次测量都单独算背景,不要用固定值。激光能量波动可以用监测通道归一化。边界值尽量固定,如果必须变,记录变化原因。拼时间序列前,先做一致性检查,把明显异常的时刻剔除。
4.5 避坑五:信噪比低的数据强行反演
现象:反演出的廓线噪声极大,看不出任何气溶胶层结构。原因:原始信号信噪比太低,可能是天气不好、激光能量下降或者探测器老化。解决:先算信噪比,如果低于 3,这段数据直接放弃。不要试图用强平滑来救,平滑只会把噪声变成虚假的结构。如果必须用,就做时间平均,比如 10 分钟平均,但要注意气溶胶层的变化尺度。
5. 从廓线到应用:气溶胶边界层高度提取与验证技巧
处理完消光系数廓线,下一步通常是提取气溶胶边界层高度。这是激光雷达数据最常用的应用之一。方法有很多,梯度法、小波变换法、曲线拟合法。梯度法最简单,找消光系数梯度最大的点。但梯度法对噪声敏感,容易把气溶胶层内部的波动当成边界层顶。小波变换法抗噪能力强,但尺度参数不好选。我一般先用梯度法快速看一下,再用小波变换法验证。
import pywt def boundary_layer_height(extinction, distance, wavelet='db4', scale=10): """ 用小波变换提取气溶胶边界层高度 extinction: 消光系数廓线 distance: 距离数组 wavelet: 小波基 scale: 尺度参数 """ # 对消光系数做小波变换 coeffs = pywt.cwt(extinction, scales=[scale], wavelet=wavelet) # 取模极大值对应的距离 modulus = np.abs(coeffs[0]) blh_idx = np.argmax(modulus) return distance[blh_idx]这里scale参数决定了小波变换的尺度。尺度太小,会把气溶胶层内部的细节当成边界层顶;尺度太大,边界层顶的位置会偏移。我一般会试 5 到 20 之间的几个尺度,看哪个尺度下提取的高度和探空数据最接近。如果没有探空数据,就看时间序列是否连续,如果边界层高度在一天内变化平滑,说明尺度选得合适。
验证方法也很重要。我习惯用两种方式交叉验证:一是和微波辐射计或者探空数据对比,二是看时间高度图上的结构是否合理。如果提取的边界层高度在时间序列上出现频繁跳变,那多半是算法参数没调好。另外,气溶胶边界层高度和云底高度容易混淆,如果消光系数廓线在边界层以上还有明显的峰值,那可能是云层,需要区分开。
最后说一个我自己的习惯:每次处理完数据,我都会把原始信号、距离校正信号、消光系数廓线、边界层高度画在一张图上,从头到尾看一遍。如果中间哪一步出现异常,图上会很明显。这个习惯帮我省了很多后悔药。希望帮到你。
本文还有配套的精品资源,点击获取