阅读说明：本技术 基于多波长激光雷达的云凝结核数浓度垂直廓线反演方法 (Cloud condensation kernel number concentration vertical profile inversion method based on multi-wavelength laser radar ) 是由 吕敏 王志恩 李占清 罗涛 刘�东 于 2020-06-12 设计创作，主要内容包括：本发明公开了一种基于多波长激光雷达的云凝结核数浓度垂直廓线反演方法,该方法包括利用多波长激光雷达探测气溶胶光学参数；利用探空修正相对湿度对气溶胶光学参数的影响；选择双模态对数正态分布的气溶胶模型,建立考虑气溶胶粒子形状的T-matrix散射模型；建立不同类型气溶胶的气溶胶光学参数与谱分布参数之间相关联的查找表；基于激光雷达探测到的气溶胶光学参数通过查找表得到对应的谱分布参数值,通过临界活化半径计算得到云凝结核数浓度垂直廓线的结果。优点：提高不同类型气溶胶和不同相对湿度条件下云凝结核数浓度垂直廓线反演的准确性,实现云凝结核数浓度垂直廓线的自动反演,为改进气候模式中云凝结核数浓度的参数化方案提供更准确的观测资料。(The invention discloses a cloud condensation nucleus number concentration vertical profile inversion method based on a multi-wavelength laser radar, which comprises the steps of detecting aerosol optical parameters by using the multi-wavelength laser radar; correcting the influence of the relative humidity on the optical parameters of the aerosol by utilizing sounding; selecting a bimodal lognormal distribution aerosol model, and establishing a T-matrix scattering model considering the shape of aerosol particles; establishing a lookup table correlating aerosol optical parameters and spectral distribution parameters of different types of aerosols; and obtaining corresponding spectral distribution parameter values through a lookup table based on the aerosol optical parameters detected by the laser radar, and calculating through the critical activation radius to obtain the result of the cloud condensation nucleus number concentration vertical profile. The advantages are that: the accuracy of inversion of the cloud condensation nucleus number concentration vertical profile under different types of aerosol and different relative humidity conditions is improved, automatic inversion of the cloud condensation nucleus number concentration vertical profile is achieved, and more accurate observation data are provided for improving a parameterization scheme of the cloud condensation nucleus number concentration in a climate mode.)

基于多波长激光雷达的云凝结核数浓度垂直廓线反演方法

技术领域：

本发明涉及激光雷达技术领域，尤其涉及基于多波长激光雷达的云凝结核数浓度垂直廓线反演方法。

背景技术：

云凝结核(CCN)是指在云中一定过饱和条件下能够被活化成云滴的部分大气气溶胶粒子，可以改变云的微物理过程，进而影响云的光学性质、云体的生命周期和降水等，产生间接气候效应，是目前气候模式中最重要的不确定性因子之一。随着现代工农业生产的发展和人类活动范围的扩展，气溶胶的时空分布越来越复杂。尽管国内外开展了大量的气溶胶观测试验，但是气溶胶作为CCN所产生的间接气候效应仍具有很高的不确定性。究其原因，主要是CCN的时空分布及其演变的观测数据仍然缺乏，特别是其垂直分布信息很难获得，因此现有气候模式在模拟CCN数浓度的垂直分布上仍存在着较大的误差。如果能够精确得到CCN数浓度的空间分布特征，则能改进气候模式中CCN数浓度的参数化方案，进而降低气候变化预测的不确定性。

在现有观测对CCN数浓度进行直接观测的研究平台包括：地面观测和飞机观测。目前用于直接观测CCN的仪器主要是云凝结核计数器。单独使用该仪器进行CCN的地面观测时，水平空间覆盖范围小，无法反映气溶胶的整体特性和演变过程。此外，最突出的不足之处在于，由于气溶胶垂直分布的不均一性，地面与云底附近的CCN特性存在较大差异。而若将该仪器由飞机携带观测时，尽管可以获得垂直方向以及云底的气溶胶和CCN数浓度，但其观测时飞行成本高且只适用于部分大型外场实验，因此，飞机观测能够获得的资料非常有限。

遥感反演方法是目前间接获取CCN空间分布的主要手段，就卫星遥感观测而言，由于受到垂直分辨率低、云污染以及气溶胶吸湿增长等的影响，很难准确获取云底的CCN数浓度。激光雷达作为一种主动遥感手段，以其精细的时空分辨率、大的垂直探测距离和稳定的连续观测能力在大气探测中得到了愈来愈广泛的应用，为CCN数浓度的垂直分布观测提供了可能性。但目前大部分研究者需要联合地面的其他仪器来间接获得CCN数浓度的垂直分布信息，观测成本高；并且反演过程中使用了仅适用于球形粒子的Mie散射模型，而实际上气溶胶粒子并非严格球形，尤其是沙尘气溶胶；尚未考虑不同相对湿度条件下的气溶胶吸湿性对CCN数浓度反演结果的影响，反演结果存在较大的误差。目前国内尚未有一套成熟的基于多波长激光雷达的光学参数的CCN数浓度垂直廓线的反演方法。

发明内容

:

针对现有技术存在的问题，本发明的目的在于提供一种基于多波长激光雷达的云凝结核数浓度垂直廓线反演方法，采用适用于球形和非球形粒子的散射模型，考虑气溶胶吸湿增强因子对反演结果的影响，提高不同类型气溶胶和不同相对湿度条件下CCN数浓度反演的准确性，弥补现有云凝结核数浓度垂直廓线观测资料的匮乏。

为解决上述技术问题，本发明通过以下方案进行实施；

一种基于多波长激光雷达的云凝结核数浓度垂直廓线反演方法，包括以下步骤：

S1：利用多波长激光雷达观测反演气溶胶的光学参数，包括：355和532nm波长上的消光系数(α₃₅₅，α₅₃₂)、355、532和1064nm波长上的后向散射系数(β₃₅₅，β₅₃₂，β₁₀₆₄)、532nm的激光雷达比和退偏振比；

S2：利用步骤S1获得的激光雷达比和退偏振比判断气溶胶的类型，主要类型包括三种:城市工业气溶胶(Type 1)，生物质燃烧气溶胶(Type 2)以及沙尘气溶胶(Type 3)；

S3：利用探空数据获得大气相对湿度廓线结合步骤S1中的气溶胶光学参数观测值(β₃₅₅，β₅₃₂，β₁₀₆₄，α₃₅₅，α₅₃₂)，计算气溶胶光学参数的吸湿增强因子和考虑相对湿度后的光学参数修正值；

S4：选择双模态对数正态分布作为三种类型气溶胶的谱分布模型，建立考虑气溶胶粒子形状的T-matrix散射模型，结合地面Aeronet长期观测数据建立不同类型气溶胶的谱分布参数与气溶胶光学参数之间相关联的查找表；

S5：根据步骤S2、步骤S3和步骤S4计算获得的气溶胶类型、气溶胶光学参数修正值和查找表，反演获取气溶胶谱分布参数的垂直分布，在此基础上进一步计算得到不同过饱和比下云凝结核数浓度的垂直廓线。

进一步的，所述步骤S1包括：

S11：对多波长激光雷达探测到的回波信号数据进行预处理，包括背景噪声订正和重叠因子订正；

S12：采用Raman-Mie方法反演355、532nm的气溶胶消光系数和后向散射系数(α₃₅₅，α₅₃₂，β₃₅₅，β₅₃₂)，其中消光系数和后向散射系数的反演公式如下：

其中为距离z处的N₂的分子数密度；和表示距离z处大气分子在波长λ₀和λ_R上的消光系数；和为距离z处气溶胶在波长λ₀和λ_R上的消光系数；和为距离z和z_c处大气分子在波长λ₀上的后向散射系数，z_c是标定高度；和为激光雷达系统接收的来自距离z和z_c处的N₂的弹性后向散射回波信号；和为激光雷达系统接收的来自距离z和z_c处的N₂的拉曼后向散射回波信号。

采用Fernald方法反演1064nm的后向散射系数反演公式为：

其中S₁和S₂分别是气溶胶和大气分子的激光雷达比，P₁₀₆₄(z)激光雷达系统接收的在距离z的1064nm波长上的回波信号，为大气分子在1064nm波长上的后向散射系数。

S13：根据步骤S12中得到的光学参数，计算激光雷达比R以及退偏振比δ的公式分别如下：

其中和为气溶胶在波长λ上、高度z处的消光系数和后向散射系数，其中波长λ可为355或532nm；下标p和s分别表示与发射激光偏振方向平行和垂直的两个方向；P_s(z,λ)和P_p(z,λ)分别为激光雷达在高度z处接收到的大气后向散射回波功率的平行偏振分量和垂直偏振分量；k(λ)＝k_p(λ)/k_s(λ)，为两个通道的增益常数比，即标定因子，主要受激光雷达系统的退偏振效应和两个通道的探测效率等因素的影响。

进一步的，所述步骤S2的气溶胶类型判断依据为：城市工业气溶胶(Type 1)，生物质燃烧气溶胶(Type 2)以及沙尘气溶胶(Type 3)的532nm的激光雷达比和退偏振比范围分别为45-70sr和0-6％，55-70sr和2-10％，30-50sr和15-30％。

进一步的，所述步骤S3需要计算气溶胶光学参数的吸湿增强因子，公式如下：

其中RH是通过探空数据获取的相对湿度，r是粒子的干半径，ξ(RH,r)表示干气溶胶半径为r、相对湿度为RH时的气溶胶光学参数，RH_ref被选取为参考相对湿度，用来表示较干的大气环境。

进一步的，所述步骤S4包括:

S41：选择双模态对数正态分布的气溶胶谱分布模型可以表示为：

其中i＝f,c分别表示气溶胶细模态和粗模态的粒子，N_ti是第i个模态的总的粒子数浓度，r_i ⁿ表示在第i个模态的气溶胶数浓度谱分布的中值半径，lnσ_i表示第i个模态的谱宽。

S42：在利用T-matrix散射模型计算非球形气溶胶粒子的光学参数时，主要把非球形抽象成椭球粒子、圆柱粒子和切比雪夫粒子这三种典型形状。其中椭球粒子的形状用轴比率a/b进行描述，指短轴长a和长轴长b之比；圆柱粒子的形状由轴比率D/L进行描述，是底面直径D与母线长L之比；切比雪夫粒子用T_n(X)来表示，n和X分别指变形程度和变形参数。其中球形粒子的轴比率为1，因此在使用T-matrix模型计算球形粒子的光学参数时，能够将理论公式简化为Mie理论的表达式。

S43：结合每种类型气溶胶的Aeronet多年观测站点数据建立合理的查找，为减少查找表的大小和维度并权衡了算法的计算时间和准确度，查找表由两部分组成，可表示为：

式中，B_f和B_c是分别根据气溶胶的细模态谱分布参数(σ_f,r_f,r)和粗模态谱分布参数(σ_c,r_c,r)的范围预先计算得到的数据矩阵，其中σ_f、r_f和r_c的计算间隔分别为0.01、0.002和0.01μm，σ_c为已知量。计算中气溶胶粒子半径r的范围被限制在0.01-10μm之间，计算时r的对数等距间隔为0.002(无量纲)。气溶胶粒子数浓度N_ti范围为100-10000个/cm³，间隔为10个/cm³。

进一步的，所述步骤S5包括：

S51：通过寻找目标约束方程最小值的方法反演得到气溶胶谱分布参数，其目标约束方程如下：

式中，g_i表示已知的光学参数观测值(β₃₅₅,β₅₃₂,β₁₀₆₄,α₃₅₅,α₅₃₂)，g'_i表示待求的谱分布参数通过T-matrix散射模型计算得到的光学参数模拟值(β’₃₅₅,β’₅₃₂,β’₁₀₆₄,α’₃₅₅,α’₅₃₂)。

S52：根据步骤S2得到的气溶胶类型对应的吸湿性参数κ，计算不同过饱和比S下被活化的气溶胶临界半径的公式如下：

r＝D/2；

式中，S为水的饱和比，可以计算得到过饱和比S_c的值(S_c＝S-1)；D为气溶胶粒子的湿直径，r为气溶胶粒子的湿半径，D_d为气溶胶粒子的干粒径；σ_s/a为0.072，M_w为水分子的摩尔质量，R为统一气体常数，T为温度，单位为K，ρ_w为水的密度。

S53：根据步骤S51和S52的计算结果，CCN数浓度可以通过以下表达式简化计算得到：

其中N_ccn为CCN的数浓度。

有益效果：本发明与现有技术相比的具有如下有益效果：

(1)利用T-matrix散射模型计算光学参数值，考虑了球形与非球形气溶胶粒子的散射特征差异，克服了仅使用Mie散射模型带来的误差，提高了对不同类型气溶胶的CCN数浓度反演结果的准确性。

(2)利用计算得到的吸湿增强因子将实测的气溶胶光学参数校正到干的气溶胶光学参数，修正较高的大气相对湿度情况下气溶胶吸湿性对CCN数浓度反演结果的影响。

(3)本发明可以应用到CCN数浓度的连续观测中，反演方法便于程序化，能够适用于无人值守的自动反演测量，获得CCN数浓度垂直廓线的连续观测数据集，为改进气候模式中云凝结核数浓度的参数化提供更多和更准确的观测资料。

附图说明

图1为基于多波长激光雷达的云凝结核数浓度垂直廓线反演方法的总体流程图；

图2为三种类型气溶胶的归一化谱分布；

图3三种类型气溶胶粒子的临界半径与过饱和比之间的关系图；

图4为通过气溶胶的临界活化半径和谱分布计算CCN数浓度的结构示意图。

具体实施方式

结合实施例和附图来详细说明本发明，但本发明并不仅限于此。

图1是本发明提出的基于多波长激光雷达的云凝结核数浓度垂直廓线反演方法的总体流程图，如图1所示，该实施例中基于多波长激光雷达的云凝结核数浓度垂直廓线反演方法包括以下步骤：

步骤S1：利用多波长激光雷达观测反演气溶胶的光学参数，包括：355和532nm波长上的消光系数(α₃₅₅，α₅₃₂)、355、532和1064nm波长上的后向散射系数(β₃₅₅，β₅₃₂，β₁₀₆₄)、532nm的激光雷达比和退偏振比；

步骤S2：利用步骤S1获得的激光雷达比和退偏振比判断气溶胶的类型，主要类型包括三种:城市工业气溶胶(Type 1)，生物质燃烧气溶胶(Type 2)以及沙尘气溶胶(Type3)；

步骤S3：利用探空数据获得大气相对湿度廓线结合步骤S1中的气溶胶光学参数观测值(β₃₅₅，β₅₃₂，β₁₀₆₄，α₃₅₅，α₅₃₂)，计算气溶胶光学参数的吸湿增强因子和考虑相对湿度后的光学参数修正值；

步骤S4：选择双模态对数正态分布作为三种类型气溶胶的谱分布模型，建立考虑气溶胶粒子形状的T-matrix散射模型，结合地面Aeronet长期观测数据建立不同类型气溶胶的谱分布参数与气溶胶光学参数之间相关联的查找表；

步骤S5：根据步骤S2、步骤S3和步骤S4计算获得的气溶胶类型、气溶胶光学参数修正值和查找表，反演获取气溶胶谱分布参数的垂直分布，在此基础上进一步计算得到不同过饱和比下云凝结核数浓度的垂直廓线。

具体实施步骤S1时，对多波长激光雷达探测到的回波信号数据进行预处理，包括背景噪声订正和重叠因子订正，其中背景噪声订正是将背景信号P_bg从激光雷达测量的原始信号P₀(z)中减去，得到P_s(z)＝P₀(z)-P_bg；重叠因子订正采用在给定垂直方向消光后向散射比的情况下直接对拉曼和弹性散射回波信号进行计算，避免过多假设带来的误差。

所述步骤S1，采用Raman-Mie方法反演355、532nm的气溶胶消光系数和后向散射系数(α₃₅₅，α₅₃₂，β₃₅₅，β₅₃₂)，其中消光系数和后向散射系数的反演公式为：

其中为距离z处的N₂的分子数密度；和表示距离z处大气分子在波长λ₀和λ_R上的消光系数；和为距离z处气溶胶在波长λ₀和λ_R上的消光系数；和为距离z和z_c处大气分子在波长λ₀上的后向散射系数，z_c是标定高度；和为激光雷达系统接收的来自距离z和z_c处的N₂的弹性后向散射回波信号；和为激光雷达系统接收的来自距离z和z_c处的N₂的拉曼后向散射回波信号。

所述步骤S1，采用Fernald方法反演1064nm的后向散射系数反演公式为：

其中S₁和S₂分别是气溶胶和大气分子的激光雷达比(S₁取决于气溶胶粒子的谱分布和折射指数，其值至少跨越一个数量级，一般为10-100sr，S₂为8π/3)，P₁₀₆₄(z)激光雷达系统接收的在距离z的1064nm波长上的回波信号，为大气分子在1064nm波长上的后向散射系数。

具体实施步骤S1时，根据上述步骤中得到的光学参数，计算激光雷达比R以及退偏振比δ的公式分别如下：

其中和为气溶胶在波长λ上、高度z处的消光系数和后向散射系数，其中波长λ可为355或532nm；下标p和s分别表示与发射激光偏振方向平行和垂直的两个方向；P_s(z,λ)和P_p(z,λ)分别为激光雷达在高度z处接收到的大气后向散射回波功率的平行偏振分量和垂直偏振分量；k(λ)＝k_p(λ)/k_s(λ)，为两个通道的增益常数比，即标定因子，主要受激光雷达系统的退偏振效应和两个通道的探测效率等因素的影响。

具体实施步骤S2时，气溶胶类型判断依据为：城市工业气溶胶(Type 1)，生物质燃烧气溶胶(Type 2)以及沙尘气溶胶(Type 3)的532nm的激光雷达比和退偏振比范围分别为45-70sr和0-6％，55-70sr和2-10％，30-50sr和15-30％。

具体实施步骤S3时，需要计算气溶胶光学参数的吸湿增强因子，公式如下：

其中RH是通过探空数据获取的相对湿度，r是粒子的干半径，ξ(RH,r)表示干气溶胶半径为r、相对湿度为RH时的气溶胶光学参数，RH_ref被选取为参考相对湿度，用来表示较干的大气环境，一般选取40％。

具体实施步骤S4时，选择双模态对数正态分布的气溶胶谱分布模型可以表示为：

其中i＝f,c分别表示气溶胶细模态和粗模态的粒子，N_ti是第i个模态的总的粒子数浓度，r_i ⁿ表示在第i个模态的气溶胶数浓度谱分布的中值半径，lnσ_i表示第i个模态的谱宽。三种类型气溶胶的归一化谱分布如图2所示。

具体实施步骤S4时，利用T-matrix散射模型计算非球形气溶胶粒子的光学参数模拟值，主要把非球形抽象成椭球粒子、圆柱粒子和切比雪夫粒子这三种典型形状。其中椭球粒子的形状用轴比率a/b进行描述，指短轴长a和长轴长b之比；圆柱粒子的形状由轴比率D/L进行描述，是底面直径D与母线长L之比；切比雪夫粒子用T_n(X)来表示，n和X分别指变形程度和变形参数。其中球形粒子的轴比率为1，因此在使用T-matrix模型计算球形粒子的光学参数时，能够将理论公式简化为Mie理论的表达式。利用T-matrix散射模型计算三种典型形状的气溶胶的光学参数模拟值所用到的公式如下：

其中P(180°)是指通过T-matrix散射模型计算的相函数在散射角为180°的值，i＝1,2代表细模态和粗模态，N_i代表细模态和粗模态的数浓度。所述步骤S1中由激光雷达反演得到的相应波长的后向散射系数和消光系数为光学参数的观测值。

具体实施步骤S4时，结合每种类型气溶胶的Aeronet多年观测站点数据建立合理的查找，为减少查找表的大小和维度并权衡了算法的计算时间和准确度，查找表由两部分组成，可表示为：

式中，B_f和B_c是分别根据气溶胶的细模态谱分布参数(σ_f,r_f,r)和粗模态谱分布参数(σ_c,r_c,r)的范围预先计算得到的数据矩阵，其中σ_f、r_f和r_c的计算间隔分别为0.01、0.002和0.01μm，σ_c为已知量。计算中气溶胶粒子半径r的范围被限制在0.01-10μm之间，计算时r的对数等距间隔为0.002(无量纲)。气溶胶粒子数浓度N_ti范围为100-10000个/cm³，间隔为10个/cm³。反演中考虑的三种类型气溶胶的谱分布参数的典型范围如下表所示：

其中，V_tf/V_tc为细模态气溶胶与粗模态气溶胶的体积浓度之比，m_R和m_I表示复折射指数的实部与虚部的平均值。

所述步骤S5包括：通过寻找目标约束方程最小值的方法反演得到气溶胶谱分布参数，其目标约束方程如下：

式中，g_i表示已知的光学参数观测值(β₃₅₅,β₅₃₂,β₁₀₆₄,α₃₅₅,α₅₃₂)，g'_i表示待求的谱分布参数通过所述步骤S4中的T-matrix散射模型计算得到的光学参数模拟值(β’₃₅₅,β’₅₃₂,β’₁₀₆₄,α’₃₅₅,α’₅₃₂)。

根据步骤S2得到的气溶胶类型对应的吸湿性参数κ，计算不同过饱和比S下被活化的气溶胶临界半径的公式如下：

r＝D/2；

式中，S为水的饱和比，可以计算得到过饱和比S_c的值(S_c＝S-1)；D为气溶胶粒子的湿直径，r为气溶胶粒子的湿半径，D_d为气溶胶粒子的干粒径；σ_s/a为0.072，M_w为水分子的摩尔质量，R为统一气体常数，T为温度，单位为K，ρ_w为水的密度。三种类型气溶胶粒子的临界半径r_c与过饱和比S_c之间的关系计算结果如图3所示。

根据上述步骤的计算结果，CCN数浓度可以通过以下表达式简化计算得到：

其中N_ccn为CCN的数浓度。图4中大于临界活化半径r_c由阴影覆盖的区域即为CCN数浓度。

需要注意的是，以上列举的仅是本发明的若干个具体实施例。显然，本发明不限于以上实施例，还可以有许多变形。本领域的普通技术人员能从本发明公开的内容直接导出或联想到的所有变形，均应认为是本发明的保护范围。