基于非球形冰晶的星载94/220GHz双频云雷达云微物理参数反演研究

吴琼; 仰美霖; 尹红刚; 商建; WU Qiong; YANG Mei-Lin; YIN Hong-Gang; SHANG Jian

网刊加载中。。。

使用Chrome浏览器效果最佳，继续浏览，你可能不会看到最佳的展示效果，

确定继续浏览么?

复制成功，请在其他浏览器进行阅读

基于非球形冰晶的星载94/220GHz双频云雷达云微物理参数反演研究 PDF

- ORCID：
吴琼 ¹
✉
- ORCID：
仰美霖 ²
- ORCID：
尹红刚 ¹
- ORCID：
商建 ¹
✉

1. 国家卫星气象中心，北京 100081； 2. 北京城市气象研究院，北京 100089

中图分类号： P412.25

最近更新：2021-10-18

DOI：10.11972/j.issn.1001-9014.2021.05.006

摘要

利用94/220GHz双频云雷达的探测优势，结合短柱、三瓣子弹花和星状雪花三种非球形冰晶粒子的散射特性，对双频云雷达云微物理参数反演算法进行了研究。通过建立双频比和中值直径的关系，再通过后向迭代的双频反演算法计算出谱参数的最优解，进而反演云微物理参数并分析反演误差。最后，通过数值试验，对考虑系统噪声和随机噪声后的反演误差进行了比较分析，并提出了仪器信噪比的指标需求。结果表明：双频联合可以较好地反演非球形冰晶的云微物理参数，为了保证0.2 g/m³以下冰水含量的反演误差小于30%，仪器的信噪比不能低于11.39dB，为国内开展星载云雷达指标设计和云微物理参数反演研究提供参考。

关键词

双频反演; 云雷达; 粒子谱参数; 非球形冰晶

引言

星载云雷达是测量全球云的三维分布尤其是垂直结构的重要遥感仪器之一，相比可见、红外以及被动微波仪器，它可以穿云入雨，给数值模式以及气候模式提供云的宏微观参数，提高预报精度。目前，国际上唯一在轨运行的星载云雷达是2006年由美国宇航局NASA（National Aeronautics and Space Administration）和加拿大空间局CSA（Canadian Space Agency）联合发射的Cloudsat卫星上装载的W波段（94GHz）云雷达CPR（Cloud Profiling Radar）。另外，欧空局ESA（European Space Agency）和日本宇宙航空研究开发机构JAXA（Japan Aerospace Exploration Agency）也正计划于2022年发射地球云、气溶胶和辐射测量卫星EarthCare（Earth Cloud，Aerosols and Radiation Explorer），搭载在其上的主载荷之一就是带有多普勒功能的W波段云雷达。目前，关于下一代星载云雷达系统方案的设计还没有最终确定，并且NASA也没有正在资助中的其他云卫星和载荷任务^［

1］。尽管如此，对未来云雷达的一些概念研究还是在一直不断推进，包括考虑增加多频、多极化和扫描等功能，以及使用小卫星或者多星组网的观测方式。关于雷达频点，W波段是测云的首选频点，因为它具有足够的灵敏度可以观测到大部分的云，但是单频反演算法往往依赖于经验的统计关系，反演精度受到一定的制约。研究表明^{［参考文献 2

百度学术}2］，如果仅使用94GHz的雷达资料，对0.01g/m³的冰水含量的反演误差为+85%到-45%。但是，如果有辅助的云滴谱信息，反演误差可以缩小到+50%到-35%。因此，使用双频甚至多频联合来获取云滴谱从而提高云微物理参数的反演精度是未来星载云雷达数据反演与应用的发展方向之一。

国内星载云雷达的研究起步较晚，目前关于云雷达的论证主要集中在W波段和观测卷云微小冰晶粒子的G波段（220 GHz）。在硬件方面，中国电子科技集团公司第三十八研究所、中国航天科技集团公司九院704所以及南京信息工程大学等几家单位都成功研制了94 GHz测云雷达，并已有机载云雷达首次外场校飞试验的计划。另外，中国航天科技集团公司五院和八院于2015年分别完成了220 GHz云雷达原理样机的研制，雷达灵敏度突破了-40 dBZ，但是目前在大功率发射机、天线研制以及外定标等方面还有一定的技术瓶颈。

双频反演算法研究最早可以追溯到20世纪70年代，并随着携带双频降水测量雷达DPR（Dual-frequency Precipitation Radar）的全球降水测量卫星GPM（Global Precipitation Mission）的发射得到了进一步的发展和应用^［

3-7］。因此，国际上对双频联合反演算法的关注主要集中在低频Ku/Ka波段的降水/雪反演上，对高频W/G波段联合反演云微物理参数的研究成果鲜见发表。国内近年才开始关注星载W/G波段的云雷达，如吴举秀等分析了94/220 GHz探测不同冰晶粒子时散射和衰减特性的差异，指出了冰晶粒子的散射和衰减与冰晶粒子尺寸以及形状之间的关系^{［参考文献 8

百度学术}8］。吴举秀还分析了94/220 GHz的双波长比对云滴谱参数变化的敏感性，论证了雷达探测距离和冰水含量以及云层厚度之间的关系^{［参考文献 9

百度学术}9］。此外，吴琼等初步研究了94/220 GHz联合反演云微物理参数的算法，但是仅仅将冰晶粒子假设为球形，没有考虑冰晶粒子的复杂形状对反演的影响^{［参考文献 10

百度学术}10］。因此，本文的主要目的是从我国星载云雷达探测的实际需求出发，基于非球形冰晶粒子，探索星载94/220 GHz云雷达的反演算法，研究算法的可行性、适用范围以及反演误差，并期望从反演的角度为我国未来星载云雷达的仪器指标设计提供一定的参考依据。

本文先基于非球形水晶的散射特性数据库并假设粒子谱服从伽马分布，建立了双频比和粒子谱参数中值直径之间的关系曲线，在此基础上，对若干种不同形状的冰晶粒子，利用后向迭代的反演算法实现了双频94/220 GHz云微物理参数的反演，并进行了误差分析。接着，通过比较不同形状的冰晶粒子在系统噪声和随机噪声影响下冰水含量反演精度的差异，初步分析了粒子形状对抗噪能力的影响。最后，通过调节粒子谱参数，统计了冰水含量和反演误差之间的关系，并最终提出了仪器信噪比的期望指标。

1 非球形冰晶散射特性数据库

在毫米波波段，计算非球形冰晶粒子散射较常用的方法有T矩阵、时域有限差分法FDTD（finite difference time domain）以及离散偶极子近似法DDA（discrete dipole approximation）。T矩阵法虽然计算比较准确并且速度很快，但是目前只能解决类似椭球这种具有轴对称性粒子的散射问题。FDTD虽然适用于各种形状，但是当尺度因子超过20时，计算效率会大大降低^［

11］。DDA同样也适用于任意形状、非均匀和各向异性粒子的散射计算，但是计算速度很慢。为此，国外有学者基于DDA方法开发了非球形粒子的单散射特性数据库。其中，Liu的数据库^{［参考文献 12

百度学术}12］虽然只有10个频点，但是涵盖了目前星载云和降水雷达进行在轨观测和预研的5个主要频点（13.4、35.6、94、220以及340 GHz）。因此，文中选择Liu开发的数据库来进行星载94/220 GHz云雷达的反演算法研究。该数据库主要提供了11种非球形冰晶粒子的等质量球半径、尺寸参数、最大尺度以及在5种温度条件下的单散射消光截面、单散射散射截面、不对称因子和相函数。关于非球形粒子形状的基本信息如表1所示。文中主要基于其中的短柱hexs、三瓣子弹花ros3以及星状雪花sstr这三种形状进行了反演。

表1 数据库中包含的频点、温度以及粒子形状

Table 1 Characteristics of ice particles in Liu’s database

频率（GHz）	13.4，35.6，85.5，94，118，150，166，183，220，340
温度（⁰C）	0，-10，-20，-30，-40
形状	编号	英文全称	英文简称
长柱	0	hexagonal long column	hexl
短柱	1	hexagonal short column	hexs
厚柱	2	hexagonal block column	hexb
厚板	3	hexagonal thickplate	hexf
薄板	4	hexagonal thin plate	hexp
三瓣子弹花	5	3-bullet rosette	ros3
四瓣子弹花	6	4-bullet rosette	ros4
五瓣子弹花	7	5-bullet rosette	ros5
六瓣子弹花	8	6-bullet rosette	ros6
星状雪花	9	sector snowflake	sstr
树枝状雪花	10	dendrite snowflake	sden

2 基本理论

采用Gamma分布来描述冰晶粒子的谱分布^［

13］

N (D) = N_{0} D^{μ} e x p (\frac{- [3.67 + μ] D}{D_{0}})

，

（1）

其中， $N_{0}$ 为数浓度， $D$ 是粒子直径， $D_{0}$ 为中值直径， $μ$ 是形状参数。在粒子谱已知的情况下，等效雷达反射率因子 $Z_{e}$ 、冰水含量 $I W C$ 以及等效粒子半径 $r_{e}$ 可分别表示为^［

14］：

Z_{e} = \frac{λ^{4}}{π^{5} {|K_{w}|}^{2}} \int_{0}^{\infty} N (D) σ_{b} (D) d D

，

（2）

I W C = \frac{π ρ}{6} \int_{0}^{\infty} D^{3} N (D) d D

，

（3）

r_{e} = \frac{\int_{0}^{\infty} r^{3} N (r) d r}{\int_{0}^{\infty} r^{2} N (r) d r}

，

（4）

其中， $λ$ 代表波长， $σ_{b}$ 代表后向散射截面， $K_{w} = (m^{2} - 1) / (m^{2} + 2)$ ， $m$ 是水的复折射指数。0度时，220 GHz的 ${|K_{w}|}^{2}$ 取值为0.46，94GHz的 ${|K_{w}|}^{2}$ 取值为0.68，35GHz的 ${|K_{w}|}^{2}$ 取值为0.88。 $ρ$ 是粒子的密度，文中统一取值为0.916g/cm³， $r$ 是粒子半径。

雷达波在传播过程中会受到水凝物（雨、雪、云冰等）和气体的衰减，其中，水凝物的衰减系数定义如下：

k_{h y d r o m e t e r} = 0.434 \int_{0}^{\infty} σ_{t} (D) N (D) d D

，

（5）

$σ_{t}$ 代表衰减截面。考虑衰减以后，雷达实际测量的反射率因子用下式表示^［

15］：

Z_{m} (r) = Z_{e} (r) e x p [- 0.2 I n (10)] \int_{0}^{r} k (s) d s = Z_{e} (r) A (r)

，

（6）

其中， $A (r)$ 是衰减因子， $k$ 是每一层衰减系数的总和，主要由以下几部分组成：

k = k_{o^{2}} + k_{H_{2} O} + k_{h y d r o m e t e r}

，

（7）

$k_{H_{2} O}$ 代表水汽衰减、 $k_{o^{2}}$ 表氧气衰减。如果将整个路径上的衰减系数进行积分，可以得到整个路径的积分衰减 $P I A$ ，其中 $r_{s}$ 代表地表。

P I A = 10 l o g_{10} A (r_{s}) = 2 \int_{0}^{r_{s}} k (s) d s

（8）

3 后向迭代的双频反演算法

3.1　算法描述

后向的双频反演算法一共有三种，分别是差分衰减法、表面参考技术法和基于非表面参考技术的后向迭代法^［

4］。因为冰晶粒子的散射不满足差分衰减法要求的瑞利散射条件，也无法提前获知表面参考技术法要求的路径积分衰减信息，所以文中基于非表面参考技术的后向迭代法来进行反演。

把公式（1）代入公式（2），再把 $N_{0}$ 提取出来，可得：

Z_{e} (r) = N_{0} (r) I_{b} (D_{0} (r))

，

（9）

其中，

I_{b} (D_{0} (r)) = \frac{λ^{4}}{π^{5} {|K_{w}|}^{2}} \int_{0}^{\infty} D^{μ} e x p (\frac{- [3.67 + μ] D}{D_{0} (r)}) σ_{b} (D) d D

（10）

由公式（10）可知， $I_{b}$ 只和粒子谱参数 $D_{0}$ 相关，和 $N_{0}$ 无关，因此双频联合反演的基本思路是通过两个频率雷达反射率因子的比值消去粒子谱参数中的一个未知参数 $N_{0}$ ，使得该比值只与另一个粒子谱参数 $D_{0}$ 相关^［

16-17］。根据定义和公式（6），两个频点雷达反射率因子的比值表示为：

\begin{array}{l} D F R (r) (d B) = 10 l o g 10 [\frac{Z_{e 1} (r)}{Z_{e 2} (r)}] \\ = 10 l o g 10 [\frac{Z_{m 1} (r) A_{2} (r)}{Z_{m 2} (r) A_{1} (r)}] \end{array}

，

（11）

其中，下标中的1、2分别代表双频组合中的低频和高频，该函数关系或者查找表需要在反演前就预先建立好。实际反演中，如果已知双频比，就能反演出 $D_{0}$ ， $D_{0}$ 可以表示为：

D_{0} (r) = F^{- 1} (D F R (r))

（12）

在此基础上，联合方程（6）和 $D_{0}$ ，可以进一步求解 $N_{0}$ ， $N_{0}$ 的表达式为：

N_{0} (r) = \frac{Z_{m} (r)}{A (r) I_{b} (D_{0} (r))}

（13）

雷达波在传输过程中有一定程度的衰减，因此必须先对测量的雷达反射率因子进行衰减订正后才能使用 $D F R$ 来进行反演。双频反演的衰减订正和单频反演^［

15］不同，单频的衰减订正独立于反演，一般是完成整条廓线的衰减订正后再做反演。而双频的衰减订正和反演基本上是同步进行的。反演前，预先给整层路径的衰减因子

A (r_{s})

设定一个初值，从回波底开始，联合公式（11）~（13）先进行回波底的衰减订正并完成

D_{0}

和

N_{0}

的反演。此时，回波底以上其他距离库的衰减因子可以统一表述为：

A (r_{j}) = A (r_{s}) e x p

[

0.2 I n (10) r \sum_{p = j}^{s} N_{0} (r_{p}) I_{t} (D_{0} (r_{p}))]

，

（14）

其中，

I_{t} (D_{0} (r)) = 0.434 \int_{0}^{\infty} D^{μ} e x p (\frac{- [3.67 + μ] D}{D_{0} (r)}) σ_{t} (D) d D

，

（15）

$j$ 代表距离库号。在此基础上，采用逐库订正的方式往上进行衰减订正直到回波顶，同时反演出 $D_{0}$ 和 $N_{0}$ 的廓线并计算 $A (r_{s}) 。$ 如果反演得到的衰减因子 $A (r_{s})$ 满足公式（16）的约束条件，则退出计算，如果不满足，则使用新计算的 $A (r_{s})$ 代替 $A (r_{s}) 的初值$ ，进行迭代计算。

|\frac{A_{i, n + 1} (r_{s}) - A_{i, n} (r_{s})}{A_{i, n + 1} (r_{s})}| < ε_{i}

，

（16）

其中， $i$ 代表雷达频点， $n$ 代表迭代计算的次数， $r_{s}$ 代表距离库底， $A (r_{s})$ 代表从距离库顶到距离库底衰减因子的积分， $ε$ 代表约束量。反演出粒子谱参数 $N_{0}$ 和 $D_{0}$ 以后，通过公式（3）和（4）就能得到云微物理参数冰水含量和有效粒子半径。

3.2　双频比 $D F R$ 和中值直径的关系

双频比 $D F R$ 和中值直径的关系是双频反演的核心。一个单调变化的 $D F R$ 和 $D_{0}$ 关系曲线对反演比较有利，否则一值多解的问题会给反演带来很大的不确定性^［

5］。因此，了解不同形状非球形粒子的双频比

D F R

和中值直径

D_{0}

的关系是开展反演的前提。

图1给出了hexs、hexb、hexp、ros3、sstr和sden共6种粒子的 $D F R$ 和 $D_{0}$ 的关系曲线，每幅图中的三条曲线分别代表DDA计算的94/220 GHz的 $D F R$ 随 $D_{0}$ 的变化曲线、DDA计算的35/94 GHz的 $D F R$ 随 $D_{0}$ 的变化曲线以及将非球形粒子等效成等质量的球形粒子后用Mie散射计算的94/220 GHz的 $D F R$ 随 $D_{0}$ 的变化曲线。对hexb而言，非球形的94/220 GHz的 $D F R$ 曲线随着 $D_{0}$ 先快速增大，到0.06 cm以后再缓慢减小，Mie散射计算的曲线与DDA计算的曲线分布趋势基本一致，但是Mie散射计算的 $D F R$ 明显偏小。因此，这种情况下如果用计算球形粒子的Mie散射来计算非球形粒子的话，会低估 $D F R$ 的值。此外，35/94 GHz的 $D F R$ 曲线先小幅增大再大幅单调递增，这和94/220GHz的 $D F R$ 曲线具有很大的差异。从反演的角度来讲，在 $D_{0}$ 小于0.06 cm的时候，因为94/220 GHz的 $D F R$ 明显大于35/94 GHz的 $D F R$ ，并且94/220 GHz的 $D F R$ 曲线以一个较大的斜率单调递增，对 $D_{0}$ 的变化更加敏感，此时用94/220 GHz来探测hexb具有较大的优势。但是当 $D_{0}$ 大于0.06cm以后，94/220 GHz的 $D F R$ 开始逐渐减小，一个 $D F R$ 对应两个 $D_{0}$ ，出现不适定问题。而35/94GHz的 $D F R$ 尽管数值上没有94/220 GHz的 $D F R$ 大，但是却始终保持单调递增，这对获得确定的数值解是比较有利的。和hexb相比，hexs最大的区别在于94/220 GHz的 $D F R$ 曲线增加到一定程度以后基本趋于饱和，因此对hexs的探测，94/220 GHz的组合有非常明显的优势。同理，单纯从反演的角度，对于ros3和sstr，用94/220 GHz的组合来探测比较合适，而hexp和sden则使用较低的35/94 GHz的组合更加合适，因为低频的 $D F R$ 曲线不仅单调递增而且幅度比94/220 GHz大，对 $D_{0}$ 的变化更加敏感。此外，hexl以及hexf的曲线变化趋势与hexb很类似，ros4、ros5、ros6则和ros3类似，不再赘述。

图1 6种非球形粒子对应的 $D F R$ 和 $D_{0}$ 曲线

Fig. 1 $D F R$ and $D_{0}$ curves of six non spherical particles

总的来讲， $D F R$ 和 $D_{0}$ 的关系受粒子形状以及频率的影响很大。因为考虑到94/220 GHz对反演hexs、ros3以及sstr这三种形状的冰晶粒子具有优势，因此文中仅基于这三种非球形粒子来开展云微物理参数的反演研究以及敏感性试验。

3.3　模拟数据和反演方案

为了模拟94/220 GHz对云的观测，假设一个高度为4 km、垂直分辨率为250 m的云层，每层的 $D_{0}$ 固定为0.05 cm， $N_{0}$ 固定为3 $\times$ 10⁵cm^-1m^-3。基于该粒子谱和散射数据库，通过联合公式（2）、（5）-（8）可以模拟测量的雷达反射率因子、等效雷达反射率因子和冰晶粒子衰减。其中，测量的雷达反射率因子是反演的输入。下图是在上述粒子谱假设下，模拟的94和220 GHz雷达分别探测sstr和hexs这两种冰晶粒子时测量的雷达反射率因子和等效雷达反射率因子分布图。

从图2（a）可以看出，因为每个距离库假设的粒子谱是固定的，所以对于同一个频点而言，每个距离库的等效雷达反射率因子是相同的，从上到下表现为一条直线。但是由于后向散射、复折射指数以及波长的差异，导致220 GHz的等效雷达反射率因子比94 GHz的等效雷达反射率因子小约4.5 dBZ。此外，消光截面随着雷达频点的增大而增大，使得94 GHz时冰晶粒子的衰减远小于220 GHz时的衰减。当粒子形状由sstr换成hexs后（如图b），模拟的等效雷达反射率因子和衰减都随之变大。也就是说，在粒子谱相同或者冰水含量基本相同的情况下，如果粒子形状不同，雷达接收到的反射率因子是有明显差异的，因此在设计雷达灵敏度的指标时，必须综合考虑冰水含量和粒子形状的影响。比如在图2的个例中，探测0.0293 g/m³的sstr时，雷达的灵敏度只要达到-4.2 dBZ，而探测0.029 3 g/m³的hexs时，雷达的灵敏度则需要提升到-5 dBZ。

图2 sstr（a）和hexs（b）94和220GHz雷达反射率因子分布图

Fig. 2 Radar Reflectivity factor profiles of sstr（a）and hexs（b）for 94 and 220GHz

为了进行合理的反演，假设两个频率 $A (r_{s})$ 的初值为1（ $P I A$ 为0），迭代计算的最大次数为100， $ε$ 为0.001。使用三次样条插值函数和 $D F R$ - $D_{0}$ 曲线来反演 $D_{0}$ ，并且设定边界条件以避免外插引入不合理误差。图3是反演的基本流程图。一般情况下，收敛的迭代次数都小于100次。但是在衰减较大或者加上噪声以后，可能会出现不收敛的情况，为了使反演收敛以便分析噪声对反演精度的影响，直接使用已知的 $A (r_{s})$ 来进行双频反演，如图中的虚线框所示。

图3 反演流程图

Fig. 3 Flowchart of retrieval algorithm

3.4　反演结果

基于上述的反演流程以及图2中模拟的测量的雷达反射因子廓线，使用后向迭代的双频反演算法对粒子谱参数和云微物理参数进行了反演，其中，sstr的反演结果如图4所示。

图4 反演值和真值的比较

Fig. 4 Comparison between parameters of the model and those retrieved from iteration method

从图4中可以看出，反演值和真值（模拟值）的廓线基本重合，表明后向迭代的双频反演算法可以用来进行94/220 GHz双频云雷达的反演。反演误差统一用相对误差来表示，表2统计了三种不同冰晶粒子反演廓线的平均相对误差。从数值上看，尽管所有物理量的反演误差都很小，但是hexs相比其他两种粒子的反演误差还是略大一些，这主要是由他们衰减的差异引起的，表明基于非表面参考技术的双频反演算法适用于衰减比较小的情况。随着衰减的增大，衰减订正精度变差，从而使得微物理参数的反演误差也随之变大。另外，和单频反演相比，如果仅仅使用94 GHz的雷达资料，sstr冰水含量的反演误差是3.09%，如果仅仅使用220 GHz的雷达资料，sstr冰水含量的反演误差是6.54%，表明使用双频联合对反演精度有了明显的提升。

表2 云微物理参数的反演误差（%）

Table 2 Retrieval error of cloud microphysical parameters（%）

Type	D₀	Log（N₀）	IWC	r_e
Sstr	1.9828e-04	7.5420e-05	1.9007e-04	1.9421e-04
Hexs	0.007 5	0.002 7	0.004 4	0.007 5
Ros3	9.2164e-04	4.0649e-04	0.0014	9.1078e-04

4 噪声对反演精度的影响分析

由于定标偏差以及接收机噪声等的影响，使得雷达实际接收的回波不像图2中那么理想。这些噪声的存在给反演带来很大的不确定性，直接导致反演精度变差。为模拟更真实的雷达回波从而分析噪声对反演精度的影响，并进一步约束仪器的信噪比，在图2的基础上增加了两类噪声，一是系统噪声，二是随机噪声。

4.1　系统噪声

为考虑系统噪声的影响，在测量的雷达反射率因子廓线上给每个距离库增加或减少一个固定强度的噪声，这里分别取0.5 dB、1 dB和1.5 dB。表3统计了三种冰晶粒子在不同强度系统噪声影响下反演得到的冰水含量的平均误差。结果表明，反演误差随着噪声的增强而增大，并且，增加系统噪声比减少相同强度的系统噪声时反演误差大，这主要是因为增加噪声前后雷达反射率因子（单位为mm⁶m³）的偏差比减少噪声前后雷达反射率因子的偏差大。另外，在相同系统噪声的情况下，hexs的反演误差最大，其次是ros3，sstr的反演误差最小，这和表2中的结果是一致的。需要指出的是，有噪声时的反演误差虽然也受反演算法本身精度的影响，但是最主要的误差来源还是噪声。比如，在+0.5dB噪声时，hexs的反演误差是20.68%，如果不考虑噪声，算法本身引起的误差只有4.82%。

表3 增加（+）或者减少（-）系统噪声时反演误差的比较

Table 3 Comparison between retrieval errors of increasing noise（+） or decreasing noise（-）

Type	+0.5	+1.0	+1.5	-0.5	-1.0	-1.5
Sstr	15.834 3	35.144 4	59.130 8	13.177 6	24.275 0	33.712 7
hexs	40.109 0	119.231 6	326.661 5	20.687 9	35.422 9	46.935 7
Ros3	19.994 1	46.940 0	86.285 4	15.532 0	27.953 8	38.094 1

EarthCare从应用需求的角度提出了最小可探测的冰水含量为0.001 g/m³，并且反演误差要小于30%^［

18］。参照这个标准来约束系统噪声，hexs能够承受的系统噪声小于0.5 dB，sstr和ros3能够承受的系统噪声则可以大于0.5 dB，但必须小于1 dB。

图5（a）是sstr在0.5 dB噪声影响下，冰水含量的反演误差随着高度的散点变化图。反演误差呈现出从距离库顶到距离库底逐渐增大的趋势。事实上，因为衰减的影响，信号在逐渐减弱，因此信噪比在逐渐降低，如图5（b）所示。反演误差最大的20%时对应的信噪比约为1.66 dB。

（a）

（b）

图5 冰水含量反演误差（a）和信噪比分布图（b）

Fig. 5 Retrieval error for ice water content（a）and distribution of SNR（b）

在此情况下，如果降低仪器的信噪比，会增大反演误差。数值试验的结果表明，在整条廓线的反演误差都小于30%的指标要求下，sstr可以承受的最大噪声强度可以达到约0.73 dB，此时的信噪比大约为1.43 dB。同理，hexs对信噪比的需求则要达到8.9 dB，最大噪声强度只有0.25 dB，ros3对信噪比的最低要求为3.19 dB，最大噪声强度为0.51 dB。

4.2　随机噪声

相比系统噪声，随机噪声完全是随机数，没有特定的规律。因此对反演结果的影响比系统噪声来的复杂。受随机噪声的影响，测量的雷达反射率因子廓线会发生明显的抖动，如图6（a）所示，该图在图2的基础上增加了一组均值为0、标准偏差为0.5 dB的高斯随机噪声。图6（b）则为冰水含量的反演误差图，和图6（a）类似，存在明显的抖动。

（a）

（b）

图6 0.5 dB随机噪声影响下的雷达反射率因子廓线（a）和冰水含量反演误差廓线（b）

Fig. 6 The profile of radar reflectivity factor with noise of 0.5 dB（a） and profile of retrieval error for ice water content（b）

多组重复试验的结果表明：就整条廓线的平均误差而言，存在随机噪声时的误差有时比存在系统噪声时的误差大，有时比存在系统噪声时的误差小。但是，每条廓线上反演误差的最大值，一般都是存在随机噪声时的误差大于存在系统噪声时的情况。如图5（a）和图6（b）所示，系统噪声时反演误差的最大值为19%，而随机噪声时反演误差的最大值则为26%。因此，为了满足整条廓线上的反演误差都小于30%的指标，随机噪声对信噪比的要求更高，上述个例中信噪比的最大值为3.78 dB。此外，为了满足受随机噪声影响时反演精度的要求，sstr要求的最低信噪比约为3.5 dB，hexs要求的最低信噪比约为10 dB，ros3要求的最低信噪比约为4.5 dB。

5 敏感性分析

事实上，信噪比是由噪声强度和测量的雷达反射率因子共同决定的，而测量的雷达反射率因子又和粒子谱相关。因此，对仪器设计而言，需要充分考虑冰晶粒子谱分布的合理变化，才能对噪声和信噪比有更好的约束。为此，设计了一组敏感性试验，让 $D_{0}$ 和 $N_{0}$ 在一定的范围内变化，从而统计不同粒子谱情况下反演误差的差异，因为随机噪声的情况比较复杂，这里暂且只分析系统噪声。表4是sstr在0.5 dB系统噪声情况下，冰水含量的平均反演误差。

表4 不同粒子谱参数情况下冰水含量的反演误差

Table 4 Retrieval error of ice water content with different particle size parameters

	D0=0.03	D0=0.04	D0=0.05	D0=0.06
N0=1 $\times$ 10⁵	12.353 4	12.725 0	13.380 8	14.270 7
N0=3 $\times$ 10⁵	12.658 9	13.794 3	15.834 3	18.616 1
N0=5 $\times$ 10⁵	12.967 4	14.895 0	18.423 1	23.206 6
N0=8 $\times$ 10⁵	13.436 2	16.611 0	22.578 7	30.419 0
N0=1 $\times$ 10⁶	13.752 8	17.800 8	25.550 8	35.337 2

当 $D_{0}$ 为0.06 cm， $N_{0}$ 为3 $\times$ 10⁶cm^-1m^-3时，sstr、ros3和hexs测量的雷达反射率因子的最小值分别为-1.37 dBZ，-15.52 dBZ和-42.99 dBZ，其中-42.99 dBZ超出了目前国产220 GHz云雷达的探测能力，意味着对hexs的探测而言，可能只能探测到云顶而无法穿透到云底；而当 $D_{0}$ 为0.02 cm， $N_{0}$ 为3 $\times$ 10⁵cm^-1m^-3时，冰水含量小于0.001 g/m³，低于冰水含量的探测需求。基于上述考虑，将 $D_{0}$ 和 $N_{0}$ 的值约束在表4的范围内。从表中可以看出，不管是在固定 $N_{0}$ 还是固定 $D_{0}$ 的情况下，反演误差都会随着粒子谱参数的增大而增大。这主要是因为随着粒子谱参数的增大，衰减也在相应增加，而前面的分析已经表明反演算法的精度会随着衰减的增大而逐渐变差。如果把反演误差和冰水含量建立关系，可以发现（见图7），反演误差随着冰水含量的变大而变大，只有低于0.16 g/m³的冰水含量在0.5dB系统噪声影响下的平均反演误差小于30%。反过来，如果要使得粒子谱参数最大时（ $D_{0}$ 为0.06 cm， $N_{0}$ 为3 $\times$ 10⁶cm^-1m^-3，此时的冰水含量约为0.2 g/m³）的反演误差也满足精度要求，sstr能够承受的最大系统噪声只有0.35 dB，信噪比则要优于9.2 dB，同理，ros3能够承受的最大系统噪声为0.035 dB，信噪比为11.39 dB。hexs在当前粒子谱下反演不收敛，因此按图3中虚线框中所示，用已知的 $A (r_{s})$ 作为初值进行反演，不再进行迭代计算。此时，hexs对应的信噪比最低约为10.16 dB。实际情况下还需要考虑反演算法引入的误差，因此，对雷达信噪比的要求必然要超过10.16 dB。

图7 sstr在0.5 dB系统噪声情况下反演误差随着冰水含量的变化

Fig. 7 The variation of retrieval error of sstr with ice water content under 0.5 dB system noise

敏感性分析的结果表明，仪器的信噪比与冰水含量的大小以及粒子形状有关，冰水含量越大，对信噪比的要求越高，并且仅仅由粒子形状引起的信噪比差值就能达到2 dB。因此，对低于0.2 g/m³的冰水含量而言，信噪比至少应该优于11.39 dB才能使得三种不同粒子形状的反演误差同时满足指标要求。如果冰水含量的探测范围大于0.2 g/m³，信噪比也要相应的提升。因此，在信噪比的设计中，需要综合考虑探测需求、粒子形状、灵敏度以及仪器的工业制造水平。

6 结论与讨论

利用94/220 GHz双频云雷达探测冰云的优势，结合后向迭代的反演算法和DDA的非球形散射特性数据库对sstr、hexs和ros3这三种典型非球形冰晶粒子的粒子谱参数以及云微物理参数进行了定量反演。通过数值模拟，评估了系统噪声和随机噪声对反演精度的影响，并探讨了仪器信噪比的指标需求。此外，通过调节粒子谱参数，分析了系统噪声影响下冰水含量和反演误差之间的关系，在此基础上，最终从反演的角度提出了信噪比的参考指标。

通过正演分析，选择sstr、hexs和ros3这三种冰晶粒子来进行94/220 GHz的联合反演研究。结果表明：基于非表面参考技术的后向迭代反演算法在衰减较小时的反演精度较高，随着衰减的增大，会发生不收敛的情况。实际应用中，非常有必要联合表面参考技术来提高衰减订正的精度，从而保证算法收敛。此外，反演的精度和粒子的形状有关，其中，sstr的精度最高，其次是ros3，最差的是hexs，在根本上这是由粒子本身消光特性的差异引起的。而且在上述的模拟试验中发现，如果将sstr误判成hexs，引入的反演误差将会超过50%。因此，反演的前提条件是要准确区分粒子的形状，这仅仅依靠双频是不够的，还需要在未来考虑增加多极化的探测功能。

通过设计噪声的数值试验发现，不管是系统噪声还是随机噪声，反演误差随着噪声的增强而增大。在相同反演误差的约束下，sstr的抗噪能力最强，对信噪比的要求最低，而hexs的抗噪能力最差，对信噪比的要求最高。相比于系统噪声，存在随机噪声时往往对信噪比有更高的要求。另外，粒子的形状会影响冰水含量的探测范围，如果雷达的灵敏度为-40 dBZ，那么hexs能够探测的最大冰水含量约为0.2 g/m³，而ros3和sstr能够探测的最大冰水含量则明显要大于0.2 g/m³，并且反演误差呈现出随着冰水含量的增大而增大的特征。为了满足这三种形状下0.2 g/m³的冰水含量反演误差小于30%的指标要求，信噪比最好要优于11.39 dB。

事实上，文中的结论主要基于比较理想的模拟数据而来，还未考虑真实大气场景的复杂性，包括粒子谱分布廓线的垂直变化以及大气衰减等；另外只分析了三种形状的冰晶粒子，还未考虑更复杂的混合冰云，将在未来的工作中进行进一步的探讨。希望本文的结果可以为国内开展星载云雷达的指标论证、仪器设计和云微物理参数定量反演研究提供参考。

Reference

BATTAGLIA A， KOLLIAS P， DHILLON R， et al . Space-borne cloud and recipation radars：status，challenges and ways forward［J］. Review of Geophysics. 2020，58，e2019RG000686. 10.5194/amt-11-5965-2018 [百度学术]

BROWN P R A， ILLINGWORTH A J， HEYMSFIELD A J， et al.The role of spaceborne millimeter-wave radar in the gloable monitoring of ice cloud ［J］. Journal of applied meteorology，1995，34（11）： 2346-2366. 10.1175/1520-0450(1995)034<2346:trosmw>2.0.co;2 [百度学术]

ROSE C R，CHANDRASEKAR V. A systems approach to GPM Dual-frequency retrieval ［J］. IEEE Transaction on Geoscience and Remote Sensing，2005，43（8）： 1816-1826. 10.1109/tgrs.2005.851165 [百度学术]

LE M，CHANDRASEKAR V.An algorithm for drop-size distribution retrieval from GPM dual-frequency precipitation radar［J］. IEEE Transaction on Geoscience and Remote Sensing， 2014，52（11）， 7170-7185. 10.1109/tgrs.2014.2308475 [百度学术]

SETO S，IGUCHI T. Applicability of the iterative backward retrieval method for the GPM dual-frequency precipitation radar［J］. IEEE Transaction on Geoscience and Remote Sensing，2011， 49（6）： 1827-1838. 10.1109/tgrs.2010.2102766 [百度学术]

SETO S，IGUCHI T， OKI T.The basic preformance of a precipitation retrieval algorithm for the global precipitation measurement mission’s single/dual-frequency radar measurements［J］.IEEE Transaction on Geoscience and Remote Sensing，2013， 54（12）： 5239-5251. 10.1109/TGRS.2012.2231686 [百度学术]

SETO S，IGHCHI T.Intercomparison of attenuation correction methods for GPM dual-frequency precipitation radar［J］ Journal of atmospheric and oceanic technology， 2014，32（5）：915-926. 10.1175/jtech-d-14-00065.1 [百度学术]

WU Ju-Xiu，DOU Fang-Li，AN Da-Wei， et al. Sensitivity of dual wavelength reflectivity ratio of 94/220GHz space-borne radar to cloud parameters with non-spherical ice crystals［J］. Acta Meteorologica.Sinica，吴举秀，窦芳丽，安大伟，等.94/220 GHz星载雷达双波长比对非球形冰晶云参数敏感性分析.气象学报）2019， 77（3）： 529-540. [百度学术]

WU Ju-Xiu，DOU Fang-Li，AN Da-Wei， et al.Simulation of scatterting characteristics of non-spherical ice crystal with 94/220GHz millimeter-wavelrngth［J］.J.Infrared Millim.Waves，吴举秀，窦芳丽，安大伟，等.非球形冰晶在94/220 GHz毫米波的散射特性模拟计算.红外与毫米波学报）2016， 35（3）： 377-384. [百度学术]

WU Qiong，YANG Mei-Ling，DOU Fang-Li， et al. A study of cloud parameters retrieval algorithm for spaceborne millimeter wavelength cloud radar［J］.Acta Meteorologica.Sinica，吴琼，仰美霖，窦芳丽，等.星载双频云雷达的云微物理参数反演算法研究.气象学报 2018. 76（1）：160-168. 10.1109/igarss.2018.8517398 [百度学术]

WANG JH， GE JX，WEI M， et al. Research progresses on thepretical computation andexperimental measurement of scattering properties of ice particles［J］.Computing Technology and Automation（王金虎，葛俊祥，魏鸣，等.卷云冰晶粒子散射特性的理论计算与实验测量研究进展.计算技术与自动化2013，32（3）： 128-131. 10.2528/pierm14033106 [百度学术]

LIU G S.A database of microwave single-scattering properties for nonspherical ice particles［J］. Bull. Amer.Meteor. Soc. 2008，89（10）： 1563–1570. 10.1175/2008bams2486.1 [百度学术]

KOSAREV A L，MAZIN I P.An empirical model of the physical structure of upper-layer clouds［J］. Atmos. Res. 1991，26（3）： 213-228. 10.1016/0169-8095(91)90055-2 [百度学术]

Meneghini R， Kozu T. Spaceborne Weather Radar［M］ .Boston.London： Artech House，1990，pp：127-135. 10.1109/igarss.1990.688513 [百度学术]

IGUCHI T，KOZU T，MENEGHINI R， et al.Rain-profiling algoruthm for the TRMM pricipitation radar［J］. Journal Of Applied Meteorology， 2000， 39：2038-2052. 10.1175/1520-0450(2001)040<2038:rpaftt>2.0.co;2 [百度学术]

MARDIANA R， IGUCHI T，TAKAHASHI N.A dual-frequency rainprofiling method without the ues of a surface reference technique. IEEE Transaction on Geoscience and Remote Sensing. 2004，42（10）：2214-2225. 10.1109/tgrs.2004.834647 [百度学术]

WU Qiong， SHANG Jian，YANG Mei-Ling， et al. Airborne Dual-frequency Rian Retrieval Analysis without the use of Surface Reference Technology［J］. Remote Sensing Technology and Application，（吴琼，商建，仰美霖，等.基于非表面参考技术的机载双频降水反演结果分析.遥感技术与应用） 2014， 29（3）： 378-385. [百度学术]

European Space Agency Mission Experts Division，http：//www.esa.int/esapub/sp/sp1279/sp1279_1_EARTHCARE.pdf. 10.5270/esa-a346352 [百度学术]

您是第位访问者

主管单位：中国科学院

主办单位：中国科学院上海技术物理研究所，中国光学学会

地址：上海市玉田路500号电话：021-25051553

51La

首页

学报简介

编委会

征稿简则

版权声明

开放获取

出版道德声明

相关下载

联系我们

常见问题

English Version

基于非球形冰晶的星载94/220GHz双频云雷达云微物理参数反演研究 PDF

摘要

关键词

引言

1 非球形冰晶散射特性数据库

2 基本理论