基于相干多普勒激光雷达的斜程湍流参数反演方法研究
1 引言
大气湍流作为大气中一种广泛存在的运动形式, 具有尺度覆盖范围大、变化速度快的特点。大气湍流对航空安全的影响不可忽视[1], 起飞或降落轨迹上的湍流 (斜程湍流) 会导致飞机颠簸, 严重的会造成飞机失控, 引发重大安全事故。飞机起飞或降落轨迹上大气湍流的准确探测, 不仅能够降低由大气湍流导致的飞机事故率, 而且对于提高机场运行效率和研究不同空间位置的大气湍流特征具有积极意义。
高时空分辨率和高精度的风场数据是低空斜程湍流反演的基础。传统的超声风速计安放高度通常为10 m, 有限的探测高度限制了其在斜程湍流探测中的应用[2];天气雷达是航空气象保障中重要的气象观测设备, 在湿度较大的天气下可以发挥作用, 但天气雷达空间分辨率差、数据刷新率低且在晴空无云的天气情况下探测能力不足[3]; 风廓线雷达只能够获得其所在位置上方的风廓线, 无法获取较大范围内的风场情况, 因此无法获得飞机起降轨迹上的风场数据[4]。近年来, 相干多普勒测风激光雷达凭借时空分辨率高、轻便、扫描模式多样等特点成为晴空条件下航空气象保障的可靠设备。相干多普勒测风激光雷达能够根据不同观测需求设置不同的扫描模式, 进而获取低空风场。
针对不同扫描模式采用对应的速度结构函数法反演大气湍流参数是相干激光雷达观测大气湍流的难点和关键。近年来, 相关研究已陆续开展。2002年, Frehlich和Cornman[5]利用相干多普勒激光雷达模拟数据的径向速度估计获取了模拟湍流速度场的空间统计特性, 并利用径向速度计算得到的结构函数计算出湍能耗散率 ε 和积分尺度。2005年, Smalikho等[6]、Frehlich等[7]利用脉冲相干激光雷达距离-高度-显示 (RHI) 扫描模式的观测结果从多普勒谱宽与高度结构函数两个方向反演大气湍流参数, 并对实验结果采用数值模拟, 证实了这两种方法的可靠性。2008年, Frehlich[8]采用固定俯仰角、改变方位角的平面-位置-显示 (PPI) 扫描模式, 使用纵向结构函数方法和横向结构函数方法反演边界层内的湍流参数, 并将反演结果同超声风速计的探测结果进行了对比, 数据的一致性较好。2011年, Chan和Lee[9]对PPI扫描径向风场数据进行子扇区划分, 使用速度结构函数计算得到每个子扇区内的 ε1/3, 进而获取了 ε1/3 的空间分布情况。2017年, Smalikho和Banakh[10]使用速度-方位-显示 (VAD) 扫描模式下的方位角结构函数反演大气湍流参数, 并将该方法的适用范围由对流边界层扩展到稳定边界层。2017至2019年, 中国海洋大学翟晓春等[11,12]采用VAD、PPI和RHI等多种激光雷达观测模式研究了大气边界层内的湍流垂直结构特征, 分析了不同粗糙度下垫面影响下的风机尾流与大气湍流的相互作用特点, 并且探究了大气湍流对飞机尾涡演化过程的影响。
而针对飞机起降轨迹上的斜程湍流, 香港天文台开展了一系列飞机下滑道区域内的湍流观测研究。2007年, Kwong和Chan[13]证明了相干多普勒测风激光雷达能够在近实时的情况下生成沿飞机滑行路径的 ε 剖面。2010年, Chan[14]利用方位角结构函数和速度结构函数进行了下滑道上的低空湍流反演与预警, 获取了 ε, 两种计算方法的结果一致, 并与飞行数据以及风切变与湍流预警系统 (WTWS) 的结果一致。2014年, Hon和Chan[15]将速度结构函数法应用到针对飞机起降阶段沿滑行路径区域扫描的下滑道扫描模式, 获得了沿滑行路径的 ε1/3 剖面, 与飞机的飞行员报告匹配较好。2010年, Chan[14]利用方位角结构函数和速度结构函数进行了下滑道上的低空湍流反演与预警, 获取了 ε, 两种计算方法的结果一致, 并与飞行数据和风切变与湍流预警系统 (WTWS) 结果一致。然而对于斜程空间范围内的湍流参数反演研究依旧处于起步阶段, 同时缺乏对于斜程空间湍流的连续观测。
本文利用安放于某机场的相干多普勒测风激光雷达, 获取下滑道扫描模式下的高时空分辨率风场数据, 将获取的径向风速数据进行5°方位角和60 m距离的扇区划分, 进一步基于Kolmogorov局地均匀各向同性湍流理论, 并使用横向结构函数法估算斜程湍流的 ε1/3 分布情况。具体分析了某机场南段跑道2018年12月15日凌晨 01:38―01:50 以及 02:22―02:35 两个时间段的观测数据, 获得的 ε1/3 斜程空间分布情况。最后, 在同区域、同时段内, 使用相同数据计算了与 ε1/3相同量纲的风切变强度, 并将二者进行了对比验证。
2 理论方法
大气湍流在产生、发展和消散的过程中会形成各种尺度的湍涡, 某些小尺度湍涡具有各向同性的特性, 可以进行统计模拟, 因此在实际研究中, 小尺度湍涡的研究是目前湍流研究领域的关注点。在Kolmogorov局地均匀各向同性湍流理论的假设中, 当雷诺数足够高时, 存在湍流强度仅由 ε 确定的惯性子区。国际民用航空组织 (ICAO) 规定, 使用 ε1/3 来衡量大气湍流强度: 0.3~0.5 m2/3/s为中度湍流, 大于0.5 m2/3/s为严重湍流。
速度结构函数法可利用相干多普勒激光雷达的径向风速数据估算满足Kolmogorov局地均匀各向同性理论假设下的 ε 值, 是当前相干多普勒激光雷达获取大气湍流 ε 的常用方法。该方法将由实测数据计算的速度结构函数与模型结构函数进行拟合, 获得 ε 值。在拟合过程中, 需要根据Kolmogorov局地均匀各向同性理论对拟合出的湍流特征尺度进行限制, 以确保其处于惯性子区内。
2.3 Kolmogorov局地均匀各向同性理论
根据Kolmogorov定律, 当雷诺数足够大时, 湍流获得充分发展。在湍涡的发展过程中, 大尺度湍涡的能量传递给次级的湍涡, 显然是各向异性的。但在串级传输的过程中, 存在一部分在统计特性上是各向同性的小尺度湍涡, 即在局地均匀各向同性区域内存在一个仅由 ε 确定的惯性子区, 其尺度
2.4 速度结构函数法
速度结构函数法需要对实测数据计算的速度结构函数同模型结构函数进行拟合来获取湍流参数, 在实际湍流研究中, 主要关注的是惯性子区内湍涡的贡献。速度结构函数根据数据的扩展方向分为横向结构函数和纵向结构函数两种。在本研究中激光雷达下滑道扫描模式的方位角分辨率高, 因此使用横向结构函数计算方法可以更加精确地反演大气湍流参数。横向结构函数
式中
近地面层区域的模型结构函数可以简化表示为
式中
当
式中
当
联立
由于激光雷达光束始终存在一定的脉冲宽度, 实际测量中探测体积平均效应无法避免, 这将影响结构函数计算。Frehlich 等[7]基于上述理论推导出探测体积平均效应下的结构函数模型。
当激光雷达的横向分辨率
式中
式中
式中
2.5 斜程湍流反演算法
图 2. 2018年12月15日01:41:17的平均风场分布图 (a) 和风速脉动分布图 (b)
Fig. 2. Mean radial velocity (a) and radial velocity fluctuation (b) of wind field at 01:41:17 on 15 December 2018
为获取倾斜下滑道区域内的不同空间位置的湍流参数并提高运算效率, 采用空间区域划分方法将扫描区域划分为多个子扇区 (每个子扇区尺度为 5°方位角、2个距离库), 计算每个子扇区内的横向结构函数并使用自协方差法进行系统随机误差校正。后将计算的结构函数与模型结构函数进行最小二乘法拟合, 获得大气湍流参数 ε。
图 3. 2018年12月15日01:44:28距激光雷达210 m处的结构函数估计
Fig. 3. Structure function estimation of turbulence at 210 m from the lidar at 01:44:28 on 15 December 2018
3 结果与分析
3.1 实验设置
中国海洋大学于2018年11月至2019年8月在某机场开展湍流观测实验。机场所在区域于秋冬季节盛行西北风, 且风速较大, 跑道端周围多建筑与高大树木, 下垫面情况复杂, 易产生地形诱导风切变。实验期间, 将青岛镭测创芯科技有限公司研制的Wind3D 6000-AP相干测风激光雷达 (详细参数见
表 1. Wind3D 6000-AP相干测风激光雷达技术指标
Table 1. Technical specifications of Wind3D 6000-AP coherent wind lidar
|
实验中激光雷达的观测模式为下滑道扫描模式。本次实验中激光雷达安放位置距离跑道边缘线约为70 m, 以保证激光光束与下滑道夹角不大于30°, 减小侧风分量的影响, 使激光雷达径向风速可以近似代替飞机真实经历的风速(真实风速在跑道方向上的风速分量)。下滑道扫描模式下的激光雷达扫描参数为: 俯仰角范围3°~4°, 方位角范围148°~174°, 方位角的扫描速度约为0.2 °/s~0.3 °/s, 距离分辨率为30 m, 风速测量精度约0.1 m/s, 风向测量精度约5°。
3.2 数据分析和讨论
3.2.1 稳定大气条件下的观测数据
安放在机场跑道南端的激光雷达于2018年12月15日01:38―01:50获取的径向风速分布如
图 5. 2018年12月15日01:38―01:50径向风速分布图
Fig. 5. Distribution of radial wind velocity during 01:38―01:50 on 15 December 2018
将每次扫描反演得到的不同空间位置的大气湍流参数按照空间位置作图, 如
图 6. 2018年12月15日01:38―01:50 ε1/3 的空间位置分布图
Fig. 6. Spatiotemporal distributions of ε1/3 during 01:38―01:50 on 15 December 2018
3.2.2 波动大气状态案例分析
2018年12月15日凌晨02:22―02:35, 跑道南端的激光雷达共进行了9次扫描, 反演获取了约13min的径向风速, 如
图 7. 2018年12月15日02:22―02:35径向风速分布图
Fig. 7. Distribution of the radial wind velocity during 02:22―02:35 on 15 December 2018
图 8. 2018年12月15日02:22―02:35 ε1/3 的空间位置分布图
Fig. 8. Spatiotemporal distributions of ε1/3 during 02:22―02:35 on 15 December 2018
3.2.3 湍流特征参数与风切变强度对比印证
为验证湍流特征参数反演算法的准确性, 采用经验证的下滑道风切变识别法[18]计算的风切变强度值进行对比印证。下滑道风切变识别法通过逆风廓线的增大或减小来判断飞机在某一特定方向上是否遭遇风切变。该方法可以进行风切变强度值
式中
对
调整后的风切变强度值
图 9. 2018年12月15日01:42:53逆风廓线 (a) 与ε1/3空间分布 (b) 以及02:30:51逆风廓线 (c) 与ε1/3 空间分布 (d)
Fig. 9. Headwind profile (a) and spatiotemporal distributions of ε1/3 (b) at 01:42:53, headwind profile (c) and spatiotemporal distributions of ε1/3 (d) at 02:30:51 on 15 December 2018
图 10. 2018年12月15日 01:38―01:50 (a) 和02:18―02:38 (b) 时段的ε1/3 与风切变强度对比
Fig. 10. Comparison of ε1/3 and I1/3 during 01:38―01:50 (a) and 02:18―02:38 (b) on 15 December 2018
4 结论
应用高时空分辨率相干多普勒激光雷达观测数据, 采用横向结构函数法估算出凌晨两个时间段 (01:38―01:50、02:22―02:35)下滑道扫描模式的斜程湍流参数, 分析了观测时段内大气湍流参数的斜程空间分布情况, 捕捉到了较为稳定大气状态下的轻微大气扰动。将估算的大气湍流强度与使用相同数据、不同计算方法计算的风切变强度值进行比对, 二者的变化趋势具有较好的一致性。然而风切变强度只代表下滑道方向上的风切变强度水平, 而空间平均后的ε1/3代表一次下滑道扫描范围内的整体湍流强度水平, 因此两者在时间上的变化趋势无法完全对应。
斜程湍流参数的观测反演研究具有广阔的应用前景, 如在气象观测中获取复杂地形条件下的湍流信息,校正天文望远镜在斜程观测中因大气湍流造成的观测误差等。在后期的实验中计划在激光雷达的扫描区域
内设置多个安装超声风速计的测风杆, 利用超声风速计的数据对激光雷达反演的湍流参数进行验证对比, 修正大气湍流反演算法, 进一步提高激光雷达的大气湍流探测反演精度。
[1] International Civil Aviation Organization. Safety report [R]. 2020.
[2] 刘艳华, 李富余, 张宏升, 等. 超声风速仪与三轴风速仪测风的比较研究[J]. 气象水文海洋仪器, 2003, 20(3): 7-16.
Liu Y H, Li F Y, Zhang H S, et al. The comparison between sonic-anemometer and three-component propeller anemometer[J]. Meteorological Hydrological and Marine Instrument, 2003, 20(3): 7-16.
[3] 李柏, 古庆同, 李瑞义, 等. 新一代天气雷达灾害性天气监测能力分析及未来发展[J]. 气象, 2013, 39(3): 265-280.
Li B, Gu Q T, Li R Y, et al. Analyses on disastrous weather monitoring capability of CINRAD and future development[J]. Meteorological Monthly, 2013, 39(3): 265-280.
[4] Armijo L. A theory for the determination of wind and precipitation velocities with Doppler radars[J]. Journal of the Atmospheric Sciences, 1969, 26(3): 570-573.
[5] Frehlich R, Cornman L. Estimating spatial velocity statistics with coherent Doppler lidar[J]. Journal of Atmospheric and Oceanic Technology, 2002, 19(3): 355-366.
[6] Smalikho I, Köpp F, Rahm S. Measurement of atmospheric turbulence by 2-μm Doppler lidar[J]. Journal of Atmospheric and Oceanic Technology, 2005, 22(11): 1733-1747.
[7] Frehlich R, Meillier Y, Jensen M L, et al. Measurements of boundary layer profiles in an urban environment[J]. Journal of applied meteorology and climatology, 2006, 45(6): 821-837.
[8] Frehlich R. Doppler lidar measurements of winds and turbulence in the boundary layer[J]. IOP Conference Series: Earth and Environmental Science, 2008, 1: 012017.
[9] Chan P W, Lee Y F. Application of short-range lidar in wind shear alerting[J]. Journal of Atmospheric and Oceanic Technology, 2012, 29(2): 207-220.
[10] Smalikho I N, Banakh V A. Measurements of wind turbulence parameters by a conically scanning coherent Doppler lidar in the atmospheric boundary layer[J]. Atmospheric Measurement Techniques, 2017, 10(11): 4191-4208.
[11] Zhai X C, Wu S H, Liu B Y. Doppler lidar investigation of wind turbine wake characteristics and atmospheric turbulence under different surface roughness[J]. Optics Express, 2017, 25(12): A515-A529.
[12] 翟晓春. 相干多普勒激光雷达的湍流探测与ALADIN机载样机A2D风场反演方法研究 [D]. 青岛: 中国海洋大学, 2019.
ZhaiX C. Atmospheric Turbulence Detection Using Coherent Doppler Lidar and Wind Field Retrieval of ALADIN Airborne Demonstrator(A2D) [D]. Qingdao: Ocean University of China, 2019.
[13] KwongK M, ChanP W. LIDAR-based turbulence intensity calculation along glide paths [C]. 14th Coherent Laser Radar Conference, 2007.
[14] Chan P W. LIDAR-based turbulence intensity calculation using glide-path scans of the Doppler Light Detection And Ranging (LIDAR) systems at the Hong Kong International Airport and comparison with flight data and a turbulence alerting system[J]. Meteorologische Zeitschrift, 2010, 19(6): 549-563.
[15] Hon K K, Chan P W. Application of LIDAR-derived eddy dissipation rate profiles in low-level wind shear and turbulence alerts at Hong Kong International Airport[J]. Meteorological Applications, 2014, 21(1): 74-85.
[16] Kolmogorov A N. Dissipation of energy in the locally isotropic turbulence[J]. Proceedings of the Royal Society of London Series A: Mathematical and Physical Sciences, 1991, 434(1890): 15-17.
[17] MoninA S, YaglomA M, LumleyJ L. Statistical Fluid Mechanics: Mechanics of Turbulence [M]. Array Cambridge, MA: MIT Press, 1971.
[18] Zhang H W, Wu S H, Wang Q C, et al. Airport low-level wind shear lidar observation at Beijing Capital International Airport[J]. Infrared Physics & Technology, 2019, 96: 113-122.
[19] Zhang H W, Liu X Y, Wang Q C, et al. Low-level wind shear identification along the glide path at BCIA by the pulsed coherent Doppler lidar[J]. Atmosphere, 2021, 12(1): 50.
[20] 刘晓英, 吴松华, 张洪玮, 等. 基于相干多普勒测风激光雷达的不同成因类型的低空风切变观测[J]. 红外与毫米波学报, 2020, 39(4): 491-504.
Article Outline
陈晓敏, 张洪玮, 孙康闻, 吴松华. 基于相干多普勒激光雷达的斜程湍流参数反演方法研究[J]. 大气与环境光学学报, 2023, 18(1): 1. Xiaomin CHEN, Hongwei ZHANG, Kangwen SUN, Songhua WU. Inversion methods of slant turbulence parameters based on coherent Doppler lidar[J]. Journal of Atmospheric and Environmental Optics, 2023, 18(1): 1.