北京S波段天气雷达夜间晴空回波产生原因

滕玉鹏 陈洪滨 马舒庆 李思腾 吴东丽 周 燕

1)(中国科学院大气物理研究所， 北京 100029)2)(中国科学院大学地球科学学院， 北京 100049)3)(中国气象局气象探测中心， 北京 100081)4)(北京城市气象研究院， 北京 100089)5)(中国农业科学院植物保护研究所， 北京 100193)

引 言

S波段天气雷达在云雨天气观测中起重要作用。日常应用中，尤其在春季、夏季、秋季夜间，经常能观测到明显的晴空回波。对于夜间晴空回波的产生原因与分类目前主要存在两种观点[1]：一种认为夜间晴空回波主要是晴空大气回波，以大气湍流造成的折射指数起伏的布拉格(Bragg)散射为主；另一种则认为是生物体散射回波，因昆虫和鸟类的生物迁飞迁徙所致。

具体而言，气象研究人员通常认为晴空回波的多普勒速度和谱宽是湍流运动所致，并将其应用于大气风场分析及湍能耗散等临近预报中[2]。但有研究表明，鸟类和昆虫尤其大型昆虫会污染雷达的多普勒速度[3-4]。相反，虽然生物研究人员根据夜间晴空回波推测生物移动，尤其是昆虫的跨区域移动[5]，但错误的判识结果会影响农业和林业部门后续的研判与决策。

研究人员对晴空回波机制存在不同观点。魏鸣等[6]、黄琴等[7]通过分析南京地区晴空回波的演变特征及大气垂直结构，认为南京地区上空晴空回波是由大气湍流造成的布拉格散射所产生。而美国国家大气研究中心(NCAR)的Wilson等[8]则认为生物体散射是佛罗里达州和科罗拉多州地区上空晴空回波的主要产生原因。由于我国研究应用中多将晴空回波视为杂波进行滤除，以排除其对降水回波的干扰[9-13]，因此，针对晴空回波产生机制的讨论鲜见报道。

北京作为我国政治中心及北方经济中心，对气象保障，尤其是短时临近预报要求较高。同时，北京地处渤海与燕山山脉之间，是连接东北生物迁徙迁飞的陆路通道，对迁徙迁飞生物的监测关系我国粮食安全和生态保护。故对于北京地区晴空回波的分析具有现实意义。

本文重点着眼于生物季节性迁飞迁徙特征，结合L波段探空数据，讨论造成北京地区夜间晴空回波的主要原因。同时，分析大气层结对晴空回波强度的影响，提出造成晴空回波日变化的原因。

1 晴空回波产生机制

晴空回波解译的分歧在于无云雨时电磁波散射原因。由生物体散射产生的晴空回波，以下简称为生物回波；而由大气湍流的布拉格散射造成的大气晴空回波，称为湍流回波。

1.1 生物回波

有关生物体的散射目前尚未建立确切、通用的模型，这是因为生物体表面外形和介电常数相较于降水粒子更为复杂多变[14]。由于昆虫个体发育状况不同，且在空中飞行姿态存在差异，通常研究人员使用经验公式，或将生物体近似成椭球体模拟生物富含水分的躯干部分[15]。基于椭球体模型，生物体差分反射率平均可达5～10 dB。同时，基于局地各向同性湍流概念，研究者认为湍流回波的差分反射率应为零，进而使差分反射率成为美国研究者利用双偏振天气雷达识别生物回波的主要依据[15-17]。

不仅如此，实际应用中X波段，甚至C和S波段的生物体散射都进入Mie散射区，使生物体散射的定量分析变得更加困难。但鉴于鸟类飞行速度远高于昆虫，速度分析仍是区别鸟类与昆虫的一种经典方法。

1.2 大气晴空回波

2 数据处理

2.1 数 据

北京地区S波段天气雷达位于北京市观象台(CINRAD/SA)，本文所用数据为2018年3月1日—10月18日每6 min 1次的雷达基数据。由于降水回波会严重干扰晴空回波的研究结果，故除特别注明外，存在降水回波的雷达数据全部被人工剔除。

相同时间段的X波段天气雷达基数据用于人工识别和剔除降水天气。X波段天气雷达数据选自北京市X波段双偏振多普勒天气雷达网，组网的雷达设备分别位于密云、昌平、顺义、通州和房山，每部设备每3 min完成1次体扫。

L波段探空仪数据同样选自北京市观象台内高空气象观测站，高空探测系统采用GFE(L)-1型二次测风雷达和GTS1型电子探空仪[21]。观测站每日施放两次探空气球，分别为07：15(北京时，下同)和19：15。探空数据采用2018年3月1日—10月18日探空原始秒级数据，其中2018年8月16日—9月4日和10月5日的数据因故未能获得。

2.2 预处理

由于雷达基数据采用极坐标形式进行存储，因此需要对其进行处理。极坐标的雷达数据转换为直角坐标时，采用Delaunay三角剖分算法进行填充[22]，同时考虑地球曲率对高度影响，对雷达回波垂直廓线的处理参考平均法[23]。虽然雷达文件内回波的最小数值为-33 dBZ，但为保证数据质量，参考晴空回波实际情况，反射率因子小于-10 dBZ 即视为空值。为避免地物及异常回波，大于25 dBZ的回波也被剔除。降水数据的剔除参考北京X波段双偏振多普勒天气雷达组网拼图。

由探空数据获得大气层结状况，需要计算梯度理查逊数与水平风切变。对气象探空数据各物理量的垂直廓线线性插值后，计算表征湍流的梯度理查逊数Ri[24]。在实际计算中，垂直温度梯度与风速梯度可由对数内插公式确定[25]。对于水平风切变ΔU的计算采用式(1)，

(1)

式(1)中，u1和u2分别为z1和z2高度上水平风速的东西分量，v1和v2分别为z1和z2高度上水平风速的南北分量。

在后续分析中，采用熵对数据降维以表征其特征。图像纹理特征变化对于回波类型分析十分重要[26]，由于图像一维熵不能反映图像的纹理[27]特征，使用灰度共生矩阵[28]，可在一定程度上弥补这一不足。

2.3 数据分类

生物迁徙迁飞活动存在明显季节性和方向性，而对于大气晴空回波，风向变化导致的湍流变化没有生物活动明显。因此，雷达数据可按照生物向北迁徙迁飞或向南迁徙迁飞时间段，将雷达数据分为北迁时段与南迁时段，并根据当日夜间探空数据风向，将当日雷达数据细分为南风天气条件与北风天气条件，进而比较在不同季节、不同风向下，雷达数据差异。

对于北京所处的中国北方地区，通常4月起生物开始向北迁移，8—10月则返回南方[29-36]。以草地螟为例，一年共4次迁飞，春季1次，夏季2次，秋季1次。幼虫于内蒙古南部、陕西北部、河北西部越冬，成虫在春季随西南风迁飞，秋季回迁[30]。为了提高迁徙迁飞时的存活几率，生物往往选择顺风移动以保存体力。对二点委夜蛾(可在北京地区越冬)的观测分析表明，6—8月的高峰期捕获数量与空中925 hPa的偏南气流密切相关[37]。虽然不同种类昆虫的飞行高度各异，但一般成层高度集中于200～1000 m[38-41]。

一般认为，当100～1000 m高度风向符合生物迁飞需求时，生物回波的组合反射率因子将产生差异，且在回波强度、回波分布及各类数据纹理值存在差异。大气晴空回波则不会出现类似情况，虽然不同风向影响气团性质，进而影响晴空回波，但这种影响非常有限。同时，没有证据表明生物活动会影响大气状态，即当存在生物活动时，其雷达探测到的晴空回波必然为生物回波与大气晴空回波两者的叠加，且生物回波随风向变化比大气晴空回波更明显。如果北京地区S波段晴空回波存在随风向的明显变化，则生物回波可能性较大；如果晴空回波随风向没有明显变化，则以湍流回波为主。

综上所述，本研究选择100 m，750 m和1.5 km高度作为表征生物活动的3个高度，即地面、925 hPa 高度和边界层顶，以3—7月为生物向北迁徙迁飞季节(北迁)，以8—10月为生物向南迁徙迁飞季节(南迁)，以19：15探空数据代表生物起飞时风向，对雷达数据进行分类。

3 晴空回波随时间与风向变化

3.1 晴空回波日变化

由于晴空回波多见于3—11月，因此随机选择期间偶数月份的第1日，即2018年4月1日、6月1日、8月1日、10月1日20:00雷达组合反射率因子作为示例展示(图1)。由图1可见，晴空回波以雷达为中心，向远距离延伸且减弱，回波强度一般小于20 dBZ。北部回波较强区域受山脉影响，疑似地物回波但其多普勒速度不为零。

图1 北京S波段天气雷达2018年20:00晴空回波反射率因子 (a)4月1日，(b)6月1日，(c)8月1日，(d)10月1日Fig.1 Beijing S-band weather radar clear-air echo reflectivity factor at 2000 BT in 2018(a)1 Apr,(b)1 Jun,(c)1 Aug,(d)1 Oct

续图1

同时，分析2018年4月1日、6月1日及10月1日14：00至次日11：00雷达组合反射率因子变化(图2)， 因8月2日清晨起雷达探测到阵风锋回波，弃用该段数据。由图2可知，这3段数据皆显示典型的晴空回波，且无降水等其他类型回波干扰。

S波段天气雷达晴空回波在日落后迅速增加并在20：00前后保持稳定，日出后消散，同时伴随着回波的增强、稳定与减弱。这与赵海军[2]统计的2008年和2013年5月—9月山东临沂雷达晴空回波结果较为一致。这种夜间回波强度的稳定性说明，每日午夜前后的单一时刻数据可表征当日夜间晴空回波特征。

图2 2018年北京S波段天气雷达监测回波变化(a)4月1日14：00—2日11：00组合反射率因子，(b)6月1日14：00—2日11：00组合反射率因子，(c)10月1日14：00—2日11：00组合反射率因子Fig.2 Composite reflectivity of Beijing S-band weather radar(a)from 1400 BT 1 Apr to 1100 BT 2 Apr in 2018,(b)from 1400 BT 1 Jun to 1100 BT 2 Jun in 2018,(c)from 1400 BT 1 Oct to 1100 BT 2 Oct in 2018

续图2

3.2 晴空回波逐月变化

由于晴空回波强度在每日夜间较稳定，选取3月1日—10月18日每日23：00雷达数据，分析晴空回波垂直强度分布(图3)。在0～300 m高度上，晴空回波强度大体为10～15 dBZ，5月中旬、8月中旬—10月中旬的回波略强于其他时间。同时4月下旬—5月末、8月中下旬—10月中旬，晴空回波强度在1.5 km高度附近较强(约为15 dBZ)。

图4为3月1日—10月18日每日21:00至次日02:00每小时组合反射率因子平均值。由图4可见，夜间晴空回波有2个峰值期(5月中旬和9月初)和1个相对谷值期(6月底—7月初)。

为了更好地显示平均回波强度连续变化，图3中保留一些弱降水时段数据，仅剔除大于25 dBZ 回波。在图4中，考虑到弱降水回波的边缘区域可能包含一定水滴，会使组合反射率因子偏大，因此在图4中剔除超过20 dBZ的回波。有研究认为，降水前昆虫会被气流带到空中形成晴空回波[42]，本研究分别对晴空回波组合反射率因子出现大值时刻进行分析。根据图4，分别选择5月10日23:00，6月7日02:00，7月7日21:00和9月15日21:00这4个呈现组合反射率因子峰值时刻的雷达数据，其中6月7日与7月7日存在降水过程。对6月7日与7月7日雷达数据分析表明：这两日存在与降水回波谱宽、纹理相似[43]，且与晴空回波存在明显差异的降水外围回波，判断当日组合反射率因子偏大是空中水滴所致。而对5月10日和9月15日分析结果表明：其所处时段风向与生物迁飞方向相反，晴空回波不可能由生物迁徙迁飞所产生。

图3 北京S波段天气雷达2018年3月1日—10月18日每日23:00反射率因子垂直分布Fig.3 The vertical distribution of Beijing S-band weather radar reflectivity at 2300 BT from 1 Mar to 18 Oct in 2018

图4 北京S波段天气雷达2018年3月1日—10月18日每日21:00—次日02：00逐小时组合反射率因子Fig.4 Hourly composite reflectivity of Beijing S-band weather radar from 2100 BT to next 0200 BT during 1 Mar-18 Oct in 2018

3.3 回波强度与风向

通过2.3节分析可知，风向是判断晴空回波中生物回波的最主要因素，因此，晴空回波反射率因子与风向的相对变化，可以作为确定晴空回波是否由生物活动产生的依据。对晴空回波在3—7月(北迁)与8—10月(南迁)的反射率因子的平均值进行统计，结果如表1所示。

由表1可知，在100 m高度，风向不影响晴空回波强度。在750 m高度，虽然不同时段风向会导致晴空回波强度出现差异，但差异不明显。在1.5 km 高度，北迁时段北风天气条件下晴空回波要弱于南风天气条件，南迁时段北风天气条件下晴空回波要强于南风天气条件。但无论南迁还是北迁时段，南风天气条件下晴空回波强度相当。考虑平均值的局限性，分别针对不同高度统计反射率因子在不同时段、不同风向下出现日数(图5)。

表1 晴空回波北迁时段(3—7月)与南迁时段(8—10月)反射率因子平均值(单位: dBZ)Table 1 The mean value of reflectivity factor of clear air echo in northward migration period(Mar-Jul) and southward migration period(Aug-Oct)(unit：dBZ)

图5 2018年3月1日—10月18日北京S波段天气雷达反射率因子在不同时段、不同风向下出现日数(a)100 m高度，(b)750 m高度，(c)1.5 km高度Fig.5 The number of days of Beijing S-band radar reflectivity factor appearing in different wind directions and different periods from 1 Mar to 18 Oct in 2018(a)at altitude of 100 m,(b)at altitude of 750 m,(c)at altitude of 1.5 km

图5中南迁时段南风天气时反射率因子均大于7 dBZ，整体稍强于其他时段风向。在100 m高度(图5a)，相同时段、不同风向条件下反射率因子大小较为接近，南迁时段反射率因子略大于北迁时段，最大反射率因子出现在南迁时北风天气。在750 m高度(图5b)与1.5 km高度(图5c)，最大反射率因子均出现在南迁时段南风天气。图5中未体现某一时段某一风向出现明显不同于其他风向的特征，这与生物活动在不同时段选择特定风向进行迁徙迁飞的特点不相符。

3.4 回波图像特征与风向

图像纹理特征刻画了图像的局部模式与排列规则，表现图像本身特征。通过不同时段、不同风向下对S波段天气雷达反射率因子、速度和谱宽3个产品的图像纹理值的对比分析，可进一步确认晴空回波产生原因。

有关图像纹理特征的研究方法中，灰度共生矩阵是一种通过研究灰度的空间相关特性描述纹理的常用方法。图6使用图像熵方法对灰度共生矩阵进行降维处理，使二维矩阵降至一维，方便比较。由图6可知，在100 m高度，灰度共生矩阵自5月起，熵值总体呈稳定趋势，多普勒速度的熵相较于反射率因子和谱宽的熵显得更加离散。在750 m与1.5 km 高度，反射率因子与谱宽的灰度共生矩阵的熵值在3—4月与6—7月存在两次峰值，多普勒速度的共生矩阵的熵值则在5月中旬后变得离散。不同时段、不同风向的灰度共生矩阵的熵值无明显特征(图略)。总之，S波段天气雷达数据在不同生物迁飞时段、不同风向下，雷达产品图像无明显差异。

图6 2018年3月1日—10月18日每日23:00北京S波段天气雷达雷达数据灰度共生矩阵在100 m，750 m和1.5 km高度的熵Fig.6 The entropy of gray level co-occurrence matrix of Beijing S-band radar data at altitudes of 100 m,750 m and 1.5 km at 2300 BT from 1 Mar to 18 Oct in 2018

续图6

3.5 晴空回波与大气状态

由3.3节和3.4节可知，晴空回波在生物南迁时段与北迁时段反射率因子大小与雷达产品图像的灰度共生矩阵熵值，在不同风向天气条件下，未展现与生物迁飞相互对应的特征，可以确定北京S波段天气雷达夜间晴空回波总体上属于由大气湍流产生大气晴空回波。结合北京市观象台的夜间探空数据，可以初步分析晴空回波与大气层结的关系。在分析过程中，参考湍流动能收支方程，选择温度、相对湿度和水平风速3种常规探测数据，与理查森数、温度垂直递减率和水平风切变3种表征大气层结稳定度的数据进行研究(图7)。

图7中，理查森数(图7a)、温度(图7b)、相对湿度(图7c)和水平风速(图7d)并未呈现与晴空回波的反射率因子相似的变化趋势。而在图7e和图7f中，温度垂直递减率和水平风切变在6月初—8月中旬的1 km以上高度出现谷值，与晴空回波反射率因子在7月初出现谷值的趋势相似。为了方便比较，对1～2 km高度上温度垂直递减率与水平风切变的平均值变化进行统计(图8)。

由图8a可见，温度垂直递减率的变化趋势呈“M”形，在5月下旬升至高点后，开始进入下行区间直至8月初，其后再次抬升并于10月初再次达到高点，同时，4月初与8月初温度垂直递减率相近。这些特征与晴空回波反射率因子的强度变化趋势一致。由图8b可见，水平风切变呈“V”形变化趋势，8月初达到最小值，5月中旬与8月下旬水平风切变相近，且8月下旬至10月中旬水平风切变保持稳定。这种变化特征也与晴空回波的相对变化存在相似。因此，初步认为晴空回波受温度垂直递减率与水平风切变的影响。

温度垂直递减率与水平风切变分别从热力因素与动力因素表征湍流动能的趋势，这代表晴空回波与湍流动能相关。但需要注意，湍流动能与晴空回波的关系还应考虑其他因素。如在冬季，北京S波段天气雷达夜间观测不到晴空回波。而这或许因为北京冬季的低气温与低湿度天气，会使大气折射率变得很小，进而无法产生足够强的后向散射。

图8 北京地区2018年3月1日—10月18日每日20:00 1～2 km高度温度垂直递减率(a)和水平风切变(b)平均值Fig.8 The mean value of 2000 BT temperature vertical decline rate(a) and horizontal wind shear(b) of height from 1 km to 2 km in Beijing from 1 Mar to 18 Oct in 2018

根据相关理论，夜间陆地风切变往往是产生湍流的唯一来源，且逆温层的存在使热力湍流很难产生[44]，即机械湍流占主导地位。本研究认为，热力因素在一定程度上会影响晴空回波的强度。热力因素影响晴空回波强度的机制有两种可能：一是热力因素通过改变大气折射率起伏，进而影响晴空回波强度；二是局地热力不均匀引发小尺度的热力湍流，这些热力湍流导致雷达电磁波产生后向散射，进而改变晴空回波强度。

3.6 晴空回波的推测

基于实际业务应用中S波段天气雷达晴空回波表现出的特征，得到如下推测：日出后，随着地表不断升温，湍流混合加剧，气象要素垂直分布在混合作用下变得趋向均匀，大气折射率起伏减弱，晴空回波消失。随着日落，地表逐渐降温，湍流能量呈准平衡状态，折射率起伏进一步减弱。当太阳辐射完全消失，风切变成为湍流能量唯一来源，混合作用减弱，折射率起伏增强，晴空回波反射率因子迅速增大并达到稳定，此时折射率起伏主要受水汽影响[17]。

目前对于晴空回波及边界层湍流的认识，还存在一些不明确之处。虽然Ottersten[45]早在1969年便提出大气湍流散射公式，该公式基于湍流的均匀各向同性性质，而当代研究表明湍流存在间歇性，在一定程度上否定了均匀各向同性湍流，即大气湍流散射公式适用范围需重新考虑。同样需要重新考虑的还有Zrnic等[16]1998年提出的基于湍流局地均匀各向同性假设的大气晴空回波与生物回波的识别方法。

对于晴空回波与湍流的深入研究，需要充分利用不同频段天气雷达谱数据[46]。因为雷达谱数据中既包含湍谱信息，也包含生物活动特征，可进行更加细致的分析。同时，随着我国天气雷达双偏振升级改造工作的逐步推进，结合模糊逻辑、机器学习等方法的雷达偏振数据的使用，也将为研究工作的开展提供便利[47-49]。

4 结 论

本文利用生物在特定时段迁徙迁飞存在定向性的特点，通过对2018年3月1日—10月18日北京地区S波段雷达数据100 m，750 m和1.5 km高度的反射率因子与雷达产品图像的灰度共生矩阵的熵值，结合北京市观象台夜间探空数据，得到以下主要结论：

1) 北京地区晴空回波反射率因子在日落后增大，夜间保持稳定，日出前减弱；回波的反射率因子自3月起逐渐增大， 5月中旬达到峰值，而后减弱，7月初达到谷值，8月上旬再次增强，并在9月起维持较大值。在垂直方向上，反射率因子在500 m高度以上，逐月变化较500 m高度以下明显，晴空回波高度可超过2 km。

2) 北京地区晴空回波在生物南迁时段与北迁时段，其反射率因子大小与雷达产品图像的灰度共生矩阵的熵值，在不同风向下，未展现与生物迁飞相互对应的特征。

3) 1～2 km高度温度垂直递减率平均值在5月下旬升至高点而后减小，8月初再次抬升并于10月再次达到高点，4月初与8月初温度垂直递减率相近；1～2 km高度水平风切变平均值在8月初达到最小值，5月中旬与8月下旬水平风切变相近，且8月下旬—10月中旬保持稳定。温度垂直递减率与水平风切变的变化趋势与晴空回波一致。

本研究认为北京2018年夜间晴空回波总体上属于大气晴空回波，且受温度垂直递减率与水平风切变共同影响。需要注意的是，对于晴空回波产生原因的争论存在已久，观点不同，本文结论仅针对北京2018年S波段天气雷达的晴空回波而言。不同地区、不同时间结论会存在差异，不能一概而论，同时，也不能否认天气雷达在生物监测方面的潜力。