CN116068672B - 一种动力与统计降尺度相结合的山区气温精细化预报方法 - Google Patents
一种动力与统计降尺度相结合的山区气温精细化预报方法 Download PDFInfo
- Publication number
- CN116068672B CN116068672B CN202310201913.9A CN202310201913A CN116068672B CN 116068672 B CN116068672 B CN 116068672B CN 202310201913 A CN202310201913 A CN 202310201913A CN 116068672 B CN116068672 B CN 116068672B
- Authority
- CN
- China
- Prior art keywords
- air temperature
- forecast
- mountain
- wrf
- site
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000013277 forecasting method Methods 0.000 title claims abstract description 15
- 238000000034 method Methods 0.000 claims abstract description 70
- 238000007637 random forest analysis Methods 0.000 claims abstract description 37
- 238000004458 analytical method Methods 0.000 claims abstract description 26
- 230000008569 process Effects 0.000 claims abstract description 17
- 238000001914 filtration Methods 0.000 claims abstract description 11
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 9
- 230000004907 flux Effects 0.000 claims description 50
- 238000009826 distribution Methods 0.000 claims description 39
- 230000008859 change Effects 0.000 claims description 30
- 238000004088 simulation Methods 0.000 claims description 23
- 230000005855 radiation Effects 0.000 claims description 19
- 238000004364 calculation method Methods 0.000 claims description 17
- 238000012937 correction Methods 0.000 claims description 12
- 238000005338 heat storage Methods 0.000 claims description 11
- 238000012876 topography Methods 0.000 claims description 10
- 238000012544 monitoring process Methods 0.000 claims description 8
- 238000000179 transient infrared spectroscopy Methods 0.000 claims description 8
- 230000009467 reduction Effects 0.000 claims description 7
- 238000012360 testing method Methods 0.000 claims description 6
- 230000001419 dependent effect Effects 0.000 claims description 5
- 238000012952 Resampling Methods 0.000 claims description 4
- 238000002790 cross-validation Methods 0.000 claims description 4
- 238000005457 optimization Methods 0.000 claims description 4
- 238000012549 training Methods 0.000 claims description 4
- 230000003044 adaptive effect Effects 0.000 claims description 3
- 230000002776 aggregation Effects 0.000 claims description 3
- 238000004220 aggregation Methods 0.000 claims description 3
- 230000008878 coupling Effects 0.000 claims description 3
- 238000010168 coupling process Methods 0.000 claims description 3
- 238000005859 coupling reaction Methods 0.000 claims description 3
- 238000003066 decision tree Methods 0.000 claims description 3
- 230000003247 decreasing effect Effects 0.000 claims description 3
- 230000000694 effects Effects 0.000 claims description 3
- 238000002474 experimental method Methods 0.000 claims description 3
- 230000010354 integration Effects 0.000 claims description 3
- 238000004445 quantitative analysis Methods 0.000 claims description 3
- 238000012216 screening Methods 0.000 claims description 2
- 238000011160 research Methods 0.000 abstract description 6
- 230000000750 progressive effect Effects 0.000 description 5
- 238000010586 diagram Methods 0.000 description 4
- 238000006243 chemical reaction Methods 0.000 description 2
- 241000132092 Aster Species 0.000 description 1
- 238000009825 accumulation Methods 0.000 description 1
- 230000006978 adaptation Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 238000010801 machine learning Methods 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000002310 reflectometry Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Chemical compound O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01W—METEOROLOGY
- G01W1/00—Meteorology
- G01W1/10—Devices for predicting weather conditions
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation
Landscapes
- Environmental & Geological Engineering (AREA)
- Engineering & Computer Science (AREA)
- Life Sciences & Earth Sciences (AREA)
- Atmospheric Sciences (AREA)
- Biodiversity & Conservation Biology (AREA)
- Ecology (AREA)
- Environmental Sciences (AREA)
- Air Conditioning Control Device (AREA)
Abstract
本发明公开了一种动力与统计降尺度相结合的山区气温精细化预报方法，包括：S1，采用WRF模式动力降尺度；S2，山区地面气象站点观测气温的空间代表性分析；S3，基于随机森林算法的WRF预报气温降尺度；S4，选择一阶自适应卡尔曼滤波方法，对随机森林降尺度得到的参考站点未来24小时预报气温进行订正；S5，基于地表能量平衡理论对山区气温进行优化预报。本发明能够实现复杂地形与下垫面条件下的山区气温30m空间分辨率的精细化预报目标，为山区气象业务应用与地理学相关过程模型研究提供技术支撑。
Description
技术领域
本发明属于精细化天气预报的技术领域，尤其是涉及一种动力与统计降尺度相结合的山区气温精细化预报方法。
背景技术
山区地貌类型复杂、地形条件多变、高程差异巨大，对局地天气特征影响剧烈。在现实生活中，农业、交通、体育、电力等行业对山区精细化气温预报服务的需求非常强烈。为准确刻画地形对局地微气候的影响，开展山地气象及其相关领域研究，需要空间分辨率高达几十米甚至米级的逐时气温预报数据。然而，全球数值预报模式、区域数值天气预报系统的空间分辨率分别为10km、1km左右，难以准确预报与表达山区复杂地形条件下的局地尺度气温时空变异特征。山区气温时空分布及其变化不仅与大尺度大气运动有关，而且与海拔高度、坡度、坡向、遮蔽度等地形要素密切相关，同时受地表下垫面类型与特性等参数的影响巨大。局地尺度山区气温的空间变异甚至超过气候变化的幅度，忽略这种局地差异将给相关研究与应用带来较大误差。因此，高时空分辨率的山区气温准确预报已成为天气预报研究与业务的难点与重点。数值预报模式主要从大气科学角度出发，重点关注大气中的热量、动量、水汽平衡等物理过程、通过求解描述天气演变过程的流体力学和热力学方程组，预报未来一定时段的大气运动状态，但由于模式空间分辨率较低，难以准确模拟预报局地尺度、非均匀复杂地表的近地面气温。统计降尺度方法物理机制不够清晰，不足以定量表达局地尺度上多个地表要素对气温变化的非线性影响，模型稳定性与移植性需要进一步改善；而且，许多降尺度研究针对过去一段时间的再分析资料或模式预报数据进行，其中采用了大量气象站点观测数据参与过程模拟与精度评估，这些方法通常无法用于未来时刻山区气温的精细化预报。
发明内容
解决的技术问题：本发明在数值模式动力降尺度与统计降尺度基础上，结合地表能量平衡原理，综合考虑地表多个要素对气温变化的影响，提出了一种动力与统计降尺度相结合的山区气温精细化预报方法，以实现山区气温30m空间分辨率、逐小时的精细化预报目标。
技术方案：
一种动力与统计降尺度相结合的山区气温精细化预报方法，所述预报方法包括以下步骤：
S1，采用WRF模式动力降尺度；具体包括：
选择WRF模式，导入气象背景场驱动数据和30m空间分辨率的下垫面土地类型、地形、植被覆盖度遥感数据，设定多参数方案的不同组合形式，以站点观测气温数据作为参考，对比评估不同参数方案用于WRF模拟的精度，确定适合于试验区的模式最优参数化方案；开展未来24小时气温预报的WRF模式模拟，输出330m空间分辨率的预报气温空间分布；
S2，山区地面气象站点观测气温的空间代表性分析；具体包括：
通过点-面时空尺度匹配，对山区地面站点观测气温进行空间代表性分析，选择330m分辨率气温分布与对应时刻WRF模式模拟得到的预报气温之间差异最小的一个或者多个站点作为备选参考站点；
S3，基于随机森林算法的WRF预报气温降尺度；具体包括：
在330m空间分辨率下，将WRF模式模拟得到的预报气温作为因变量，地表遥感参数、坡度、坡向、海拔高度参数作为自变量，利用随机森林算法训练降尺度模型，采用K折交叉验证方法，确定随机森林模型中的树节点预选变量个数和决策树数目这两个关键参数，得到最优随机森林模型，进而计算模型拟合残差并采用双线性内插方法重采样为30m分辨率；将30m空间分辨率的各类地表参数输入随机森林模型，在计算结果基础上加入模型拟合残差，得到降尺度后的30m空间分辨率的预报气温；
S4，选择一阶自适应卡尔曼滤波方法，对随机森林降尺度得到的参考站点未来24小时预报气温进行订正；
S5，基于地表能量平衡理论对山区气温进行优化预报；具体包括：
在30m空间分辨率的尺度上，通过分析地形和下垫面特性参数对山区不同地点各通量项以及气温的影响差异，以代表性参考站点降尺度后的气温预报数据作为背景气温，模拟得到预报山区其它地点的气温。
进一步地，步骤S1中，采用WRF模式动力降尺度的过程包括以下步骤：
选择WRF模式，采用五重嵌套网格，垂直方向分为35层，顶层气压为50hPa，气象背景场驱动数据为全球FNL再分析资料，输入30m空间分辨率的包括下垫面土地利用类型、地形、植被覆盖度在内的遥感数据，采用逐3天的分段积分方式，并利用分析场Nudging方法对外层粗网格进行格点逼近，将计算域3的WRF预报结果作为LES的初始场，模式输出时间间隔为60分钟；
围绕辐射、微物理、积云对流、边界层、近地层以及陆面过程参数，以及WRF-LES耦合中的次网格闭合方案、LES计算域的侧边界与顶部边界条件确定，设定多参数方案的不同组合形式，利用WRF模式模拟预报330m分辨率的24小时气温空间分布，并以站点观测气温数据作为参考，对比评估不同参数方案用于WRF模拟的精度，以精度最优为筛选条件，确定适合于试验区的模式最优参数化方案；
在上述参数化方案优化基础上，开展未来24小时气温预报的WRF模式模拟，输出330m空间分辨率的预报气温空间分布。
进一步地，步骤S2中，山区地面气象站点观测气温的空间代表性分析的过程包括以下步骤：
在Landsat-8 TIRS热红外传感器观测时刻，建立气象站点观测气温与TIRS传感器反演地表温度、地形参数、地表遥感指数之间的多元回归模型，估算30m空间分辨率的气温空间分布；
利用30m分辨率的气温数据进行平均聚合，得到330m分辨率气温分布，将其与WRF模式输出格点进行空间配准，并计算其与对应时刻WRF预报气温之间的差值；
针对每个地面气象观测站点位置，比较两种数据之间的差异，选择差异较小的站点作为备选参考站点。
进一步地，在不同季节分别开展山区地面气象站点观测气温的空间代表性分析，并比较不同站点位置的WRF预报、遥感估算差值的平均值、标准差，选择平均值、标准差最小的一个或者多个地面气象站点作为区域空间代表性程度高的站点。
进一步地，步骤S4中，选择一阶自适应卡尔曼滤波方法，对随机森林降尺度得到的参考站点未来24小时预报气温进行订正的过程包括以下步骤：
利用公式（1），在有观测值的t时刻，计算随机森林降尺度得到的预报气温
与站点观测气温/>
之差，得到参考站点位置的模式预报误差/>
：
采用离散逼近法迭代计算模式误差
：
t+1时刻，模式预报订正结果
为：
对于无观测值的预报时段，利用随机森林降尺度预报的订正结果
替代公式（1）中的/>
；
利用有站点观测数据的时段，通过试验设定权重系数w的取值范围为0.01~0.99，变化间隔0.01，分析比较权重系数w的不同取值情况下24小时预报的总误差，确定总误差最小时对应的取值作为权重系数w的最终值，将其用于WRF模式预报气温订正。
进一步地，步骤S5中，基于地表能量平衡理论对山区气温进行优化预报；具体包括：
在30m空间分辨率的尺度上，通过分析地形、下垫面特性参数对山区不同地点各通量项以及气温的影响差异，以代表性参考站点降尺度后的气温预报数据作为背景气温，建立公式（4）用于模拟预报山区其它地点的气温；利用参考气象站点计算其与其它地点之间的通量项差异，并在公式（4）中考虑温度平流差异，在一定程度上可以降低山区地表能量平衡闭合效果差的影响，而且利用参考站点预报气温能够约束求解公式（4），提高气温预报精度；
其中，t为预报时刻；
为山区待模拟计算地点的预报气温；/>
为代表性参考站点降尺度预报气温的卡尔曼滤波订正结果；/>
为各通量变化引起的待模拟计算地点气温变化值，/>
分别为各通量变化引起的参考站点气温变化值；/>
为待模拟计算地点与参考站点因高度不同而产生的气温差异值；/>
为冷暖空气水平运动引起的温度平流差异值。
进一步地，在参考站点位置上，采用公式（5）计算各通量变化引起的气温变化值
：
式中，
是净辐射通量，在WRF模式预报和遥感监测方法基础上，根据山区复杂地形与下垫面对入射短波辐射、发射长波辐射的影响以及周围环境散射反射与辐射的贡献计算得到；G是地表热储存通量，根据净辐射通量/>
与地表热储存通量之间的经验关系模型估算得到；L是潜热通量，m是湍流交换系数，/>
是地表显热通量阻抗系数；
计算出待模拟计算地点的净辐射通量、地表热储存通量、潜热通量、湍流交换系数和地表显热通量阻抗系数后，将其代入公式（5），得到待模拟计算地点气温变化值
。
进一步地，通过不同季节与月份站点观测气温随高度变化规律的定量分析，以及t时刻WRF模式气温预报结果，计算山区不同地点的气温垂直递减率
，并结合模拟计算地点与参考站点之间的高度差异值/>
，采用公式（6）计算得到待模拟计算地点与参考站点因高度不同而产生的气温差异值/>
：
根据参考站点与模拟地点的温度水平分布梯度，以及风矢量数据，采用公式（7）计算得到冷暖空气水平运动引起的温度平流差异值
：
式中，u、v分别是东西、南北方向的风速；
、/>
分别表示东西、南北方向的温度差异。
进一步地，步骤S5中，对于拟模拟山区地点，利用WRF模式预报、随机森林降尺度结果、遥感监测数据，估算或以
表达公式（5）至公式（7）中的各项参数，并以参考站点预报气温/>
为驱动，通过公式（4）迭代求解。
有益效果：
本发明提出了一种动力与统计降尺度相结合的山区气温精细化预报方法，该方法在数值模式预报基础上，采用实验观测、遥感监测、理论模拟相结合的手段，构建由模式动力降尺度提供气温预报初始值、随机森林统计降尺度提升空间分辨率、地表能量平衡理论提高预报精度的层级递进方法。本发明能够实现复杂地形与下垫面条件下的山区气温30m空间分辨率的精细化预报目标，为山区气象业务应用与地理学相关过程模型研究提供技术支撑。本发明的预报方法在河北省张家口市崇礼区进行了应用，结果表明本发明提出的预报方法不仅能够有效提升预报气温的空间分辨率，而且能够显著提高气温预报精度。
附图说明
图1为本发明实施例的动力与统计降尺度相结合的山区气温精细化预报方法流程图；
图2a为WRF模式预报的2020年1月19日2时气温空间分布示意图；
图2b为WRF模式预报的2020年1月19日14时气温空间分布示意图；
图3a为本发明实施例预报的2020年1月19日2时气温空间分布示意图；
图3b为本发明实施例预报的2020年1月19日14时气温空间分布示意图；
图4a为2020年1月19日2时WRF预报与观测气温散点图；
图4b为2020年1月19日2时本发明实施例预报与观测气温散点图；
图4c为2020年1月19日14时WRF预报与观测气温散点图；
图4d为2020年1月19日14时本发明实施例预报与观测气温散点图。
具体实施方式
下面的实施例可使本专业技术人员更全面地理解本发明，但不以任何方式限制本发明。
本发明实施例的目的在于提出一种山区气温精细化预报方法。本山区气温精细化预报方法在大气状况与地表要素对山区气温时空分布的定量影响分析基础上，综合考虑陆-气相互作用过程，构建模式动力降尺度、机器学习统计降尺度、地表能量平衡理论相结合的层级递进方法。
1）WRF模式动力降尺度
选择WRF模式，采用五重嵌套网格，垂直方向分为35层，顶层气压为50hPa，气象背景场驱动数据为美国NCEP/NCAR发布的全球FNL（Final Operational Global Analysis）再分析资料，输入30m空间分辨率的下垫面土地利用类型、地形、植被覆盖度等遥感数据，采用逐3天的分段积分方式，并利用分析场Nudging方法对外层粗网格进行格点逼近，将计算域3的WRF预报结果作为LES的初始场，模式输出时间间隔为60分钟。
围绕辐射、微物理、积云对流、边界层、近地层以及陆面过程等参数，以及WRF-LES耦合中的次网格闭合方案、LES计算域的侧边界与顶部边界条件确定，设定多参数方案的不同组合形式，利用WRF模式模拟预报330m分辨率的24小时气温空间分布，并以站点观测气温数据作为参考，对比评估不同参数方案用于WRF模拟的精度，确定适合于试验区的模式最优参数化方案。在上述参数化方案优化基础上，开展未来24小时气温预报的WRF模式模拟，输出330m空间分辨率的预报气温空间分布。
2）山区地面气象站点观测气温的空间代表性分析
在地表异质性特点显著的复杂山区，为了客观评价动力降尺度、统计降尺度在不同地理位置上预报气温的准确性，以便于后续利用地表能量平衡理论进一步提高气温预报精度，需要通过点-面时空尺度匹配，开展山区地面站点观测气温的空间代表性分析。可通过地面观测、遥感监测、WRF模式模拟之间的尺度转换关系，开展符合研究区实际状况的地面气象站点空间代表性分析。
首先，在Landsat-8 TIRS热红外传感器观测时刻，建立气象站点观测气温与TIRS传感器反演地表温度、地形参数、地表遥感指数之间的多元回归模型，估算30m空间分辨率的气温空间分布；然后，利用30m分辨率的气温数据进行平均聚合，得到330m分辨率气温分布，将其与WRF模式输出格点进行空间配准，并计算其与对应时刻WRF预报气温之间的差值；最后，针对每个地面气象观测站点位置，比较两种数据之间的差异，选择差异较小的站点作为备选参考站点。为了提高参考站点选择的准确性，可以在不同季节分别开展类似计算，并比较不同站点位置的WRF预报、遥感估算差值的平均值、标准差，选择平均值、标准差都较小的地面气象站点作为区域空间代表性程度较高的站点。
3）基于随机森林算法的WRF预报气温降尺度
在330m空间分辨率下，将WRF模式预报气温作为因变量，地表遥感参数、坡度、坡向、海拔高度等参数作为自变量，利用随机森林算法训练降尺度模型，采用K折交叉验证方法，确定随机森林模型中的树节点预选变量个数和决策树数目这两个关键参数，得到最优随机森林模型，进而计算模型拟合残差并采用双线性内插方法重采样为30m分辨率；然后，将30m空间分辨率的各类地表参数输入随机森林模型，在计算结果基础上加入模型拟合残差，得到降尺度后的30m空间分辨率的预报气温。
4）代表性参考站点的预报气温订正
选择一阶自适应卡尔曼滤波方法，对随机森林降尺度得到的参考站点未来24小时预报气温进行订正。利用公式（1），在有观测值的t时刻，计算随机森林降尺度得到的预报气温
与站点观测气温/>
之差，得到参考站点位置的模式预报误差/>
：
采用离散逼近法迭代计算模式误差
：
t+1时刻，模式预报订正结果
为：
对于无观测值的预报时段，公式（1）中的
可利用随机森林降尺度预报的订正结果F(t)替代。一阶自适应卡尔曼滤波方法的关键是确定权重系数w的取值。利用有站点观测数据的时段，通过试验设定w取值范围为0.01~0.99，变化间隔0.01，分析比较不同w取值情况下24小时预报的总误差，确定总误差最小时对应的w值，将其用于WRF模式预报气温订正。
5）基于地表能量平衡理论的山区气温优化预报方法
在30m空间分辨率的尺度上，通过分析地形、下垫面特性等参数对山区不同地点各通量项以及气温的影响差异，以代表性参考站点降尺度后的气温预报数据作为背景气温，建立公式（4）用于模拟预报山区其它地点的气温。利用参考气象站点计算其与其它地点之间的通量项差异，并在公式（4）中考虑温度平流差异，在一定程度上可以降低山区地表能量平衡闭合效果差的影响，而且利用参考站点预报气温能够约束求解公式（4），提高气温预报精度：
其中，t为预报时刻；
为山区待模拟计算地点的预报气温；/>
为代表性参考站点降尺度预报气温的卡尔曼滤波订正结果；/>
为各通量变化引起的待模拟计算地点气温变化值，/>
分别为各通量变化引起的参考站点气温变化值；/>
为待模拟计算地点与参考站点因高度不同而产生的气温差异值；/>
为冷暖空气水平运动引起的温度平流差异值。
进一步地，在参考站点位置上，根据公式（5），利用湍流交换系数m、地表显热通量阻抗系数
计算气温变化值/>
；在待模拟计算地点位置上，计算地表热储存通量、潜热通量、湍流交换系数和地表显热通量阻抗系数后，将其代入公式（5），可得到待模拟计算地点气温变化值/>
。
在WRF模式预报和遥感监测方法基础上，根据山区复杂地形与下垫面对入射短波辐射、发射长波辐射的影响，以及周围环境散射反射与辐射的贡献，计算净辐射通量
；根据净辐射通量/>
与地表热储存通量之间的经验关系模型估算地表热储存通量G，利用湍流模式、整体空气动力学等方法估算空气动力学阻抗、下垫面表面阻抗、潜热通量L、湍流交换系数m和地表显热通量阻抗系数/>
等参数。
通过不同季节与月份站点观测气温随高度变化规律的定量分析，以及t时刻WRF模式气温预报结果，计算山区不同地点的气温垂直递减率
，并结合模拟计算地点与参考站点之间的高度差异/>
，计算公式6中的/>
，根据参考站点与模拟地点的温度水平分布梯度，以及风矢量数据，计算公式中的/>
：
式中，u、v分别是东西、南北方向的风速；
、/>
分别表示东西、南北方向的温度差异。
在近地表风速估算过程中，考虑到相邻像元位置海拔高度差异较小，且相邻位置上方一定高度上（如200m）的风速近似相等，可利用参考站点为中心向外逐一计算近地表风速。以参考站点位置WRF预报风速的一阶卡尔曼滤波时序订正结果为基础，通过已知位置风速、摩擦风速、200m高度风速之间的换算关系推算下一个邻近像元的近地表风速。
在上述关键参数估算基础上，对于拟模拟山区地点，利用WRF模式预报、随机森林降尺度结果、遥感监测等数据，直接估算或以
表达公式（5）-（7）中各项参数，并以参考站点预报气温/>
为驱动，通过方程（4）迭代求解/>
。
实例
本发明实施例在WRF模式动力降尺度、随机森林统计降尺度基础上，结合地表能量平衡理论，开展山区气温层级递进精细化预报，以下为示例的具体实施步骤。技术流程见图1。
1）试验区选择：选择河北省张家口市崇礼区作为试验区，崇礼境内山地面积约占80%，海拔高度从812m到2169m变化，平均海拔高度为1485m，形成了“连绵不断的山连山、沟套沟”的地貌特征，是开展山区气温精细化预报方法检验的理想区域。
2）WRF模式动力降尺度：选择WRF-ARW模式，采用五重嵌套网格，格距分别为27km，9km，3km，1km，0.33km，计算域5完全覆盖河北省张家口市崇礼区，气象背景场驱动数据为美国NCEP/NCAR发布的1°×1°网格全球FNL（Final Operational Global Analysis）再分析资料，设定Thompson微物理参数、YSU行星边界层参数、RRTM长波辐射参数、Dudhia短波辐射参数和Noah陆面过程模式，模拟时段为2020年1月18日12:00至1月20日0:00（spin-up6h），共计36小时，输出2020年1月19日、330m空间分辨率的逐时气温预报空间分布。
3）地面气象站点观测气温的空间代表性分析：选择2020年1月22日11时成像的Landsat8 OLI/TIRS传感器数据，对OLI可见光与近红外波段图像进行大气校正，计算NDVI、NDSI、反照率，选择TIRS第10波段利用单窗算法反演地表温度；利用ASTER DEM数据提取高程，并计算坡度、坡向；根据崇礼区74个自动气象站的经纬度数据，提取对应位置的地表温度、NDVI、NDSI、海拔高度、坡度、坡向等数据作为自变量，以自动气象站11时观测气温作为因变量，建立气温遥感估算模型，并将其重采样至330m分辨率；计算WRF模式预报的11时气温与遥感估算气温之间的差值，选择差值较小的站点作为参考站点，本发明实验中选择“云顶5号”气象站作为参考站点。
4）基于随机森林算法的WRF预报气温降尺度：在330m空间分辨率下，将WRF模式预报气温作为因变量，Landsat8 OLI遥感图像计算出的NDVI、NDSI遥感指数、多波段反射率、反照率、坡度、坡向、海拔高度等参数作为自变量，采用K折交叉验证方法，训练得到最优的随机森林模型，输入30m空间分辨率的各类地表、地形参数，得到降尺度后的30m空间分辨率的预报气温。
5）参考站点预报气温订正与基于地表能量平衡理论的山区气温优化预报：利用公式（1）、（2）、（3）描述的一阶自适应卡尔曼滤波方法，对参考气象站点“云顶5号”的2020年1月19日预报气温进行订正；利用公式（4）、（5）、（6）、（7），根据地表能量平衡理论对随机森林降尺度后的每个像元位置的预报气温进行优化，顺序估算气温垂直递减率、湍流交换系数、地表显热通量阻抗系数、冷暖空气水平运动引起的温度平流差异等中间参数，以及各通量变化引起的待模拟计算地点、参考站点的气温变化，最后计算得出30m空间分辨率的层级递进精细化预报气温空间分布。
6）实验结果分析：将WRF动力降尺度预报气温空间分布与本发明实施例得到的预报气温进行对比，以2020年1月19日2时、14时为例，图2a、图2b、图3a和图3b分别给出了WRF预报气温和本发明实施例预报气温的空间分布。根据图2a、图2b、图3a和图3b，在两个不同时刻，两种预报方法得出的气温空间分布存在较大差异，总体表现为本发明实施例预报气温在空间分布细节上更为精细，最高气温与最低气温之间的差异相对较小。进一步地，利用74个站点气温观测数据，选择地理加权回归方法，获得观测气温的空间分布，并随机选择部分位置作为参考气温，图4a至图4d给出了两个时刻的WRF预报气温、本发明实施例预报气温与站点观测气温的差异。根据图4a和图4b，在2时，与WRF预报气温相比，本发明实施例预报气温与观测气温之间的决定系数R2提高了0.39，均方根误差下降约1.2℃；根据图4c和图4d，在14时，本发明实施例预报气温与观测气温之间的决定系数R2提高了0.21，均方根误差下降约0.8℃；在两个不同时刻上，本发明实施例预报气温与观测气温散点分布更接近于直线y=x。
7）实验结论：本发明实施例结果说明模式动力降尺度提供气温预报初始值、随机森林统计降尺度提升空间分辨率、地表能量平衡理论提高预报精度的层级递进模型，可以较高的精度预报未来24小时、30m分辨率的逐时气温空间分布，能够为复杂下垫面与地形条件下的山区气温精细化预报提供技术支撑。
以上仅是本发明的优选实施方式，本发明的保护范围并不仅局限于上述实施例，凡属于本发明思路下的技术方案均属于本发明的保护范围。应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理前提下的若干改进和润饰，应视为本发明的保护范围。
Claims (9)
1.一种动力与统计降尺度相结合的山区气温精细化预报方法，其特征在于，所述预报方法包括以下步骤：
S1，采用WRF模式动力降尺度；具体包括：
选择WRF模式，导入气象背景场驱动数据和30m空间分辨率的下垫面土地类型、地形、植被覆盖度遥感数据，设定多参数方案的不同组合形式，以站点观测气温数据作为参考，对比评估不同参数方案用于WRF模拟的精度，确定适合于试验区的模式最优参数化方案；开展未来24小时气温预报的WRF模式模拟，输出330m空间分辨率的预报气温空间分布；
S2，山区地面气象站点观测气温的空间代表性分析；具体包括：
通过点-面时空尺度匹配，对山区地面站点观测气温进行空间代表性分析，选择330m分辨率气温分布与对应时刻WRF模式模拟得到的预报气温之间差异最小的一个或者多个站点作为备选参考站点；
S3，基于随机森林算法的WRF预报气温降尺度；具体包括：
在330m空间分辨率下，将WRF模式模拟得到的预报气温作为因变量，地表遥感参数、坡度、坡向、海拔高度参数作为自变量，利用随机森林算法训练降尺度模型，采用K折交叉验证方法，确定随机森林模型中的树节点预选变量个数和决策树数目这两个关键参数，得到最优随机森林模型，进而计算模型拟合残差并采用双线性内插方法重采样为30m分辨率；将30m空间分辨率的各类地表参数输入随机森林模型，在计算结果基础上加入模型拟合残差，得到降尺度后的30m空间分辨率的预报气温；
S4，选择一阶自适应卡尔曼滤波方法，对随机森林降尺度得到的参考站点未来24小时预报气温进行订正；
S5，基于地表能量平衡理论对山区气温进行优化预报；具体包括：
在30m空间分辨率的尺度上，通过分析地形和下垫面特性参数对山区不同地点各通量项以及气温的影响差异，以代表性参考站点降尺度后的气温预报数据作为背景气温，模拟得到预报山区其它地点的气温。
2.根据权利要求1所述的动力与统计降尺度相结合的山区气温精细化预报方法，其特征在于，步骤S1中，采用WRF模式动力降尺度的过程包括以下步骤：
选择WRF模式，采用五重嵌套网格，垂直方向分为35层，顶层气压为50hPa，气象背景场驱动数据为全球FNL再分析资料，输入30m空间分辨率的包括下垫面土地利用类型、地形、植被覆盖度在内的遥感数据，采用逐3天的分段积分方式，并利用分析场Nudging方法对外层粗网格进行格点逼近，将计算域3的WRF预报结果作为LES的初始场，模式输出时间间隔为60分钟；
围绕辐射、微物理、积云对流、边界层、近地层以及陆面过程参数，以及WRF-LES耦合中的次网格闭合方案、LES计算域的侧边界与顶部边界条件确定，设定多参数方案的不同组合形式，利用WRF模式模拟预报330m分辨率的24小时气温空间分布，并以站点观测气温数据作为参考，对比评估不同参数方案用于WRF模拟的精度，以精度最优为筛选条件，确定适合于试验区的模式最优参数化方案；
在上述参数化方案优化基础上，开展未来24小时气温预报的WRF模式模拟，输出330m空间分辨率的预报气温空间分布。
3.根据权利要求1所述的动力与统计降尺度相结合的山区气温精细化预报方法，其特征在于，步骤S2中，山区地面气象站点观测气温的空间代表性分析的过程包括以下步骤：
在Landsat-8 TIRS热红外传感器观测时刻，建立气象站点观测气温与TIRS传感器反演地表温度、地形参数、地表遥感指数之间的多元回归模型，估算30m空间分辨率的气温空间分布；
利用30m分辨率的气温数据进行平均聚合，得到330m分辨率气温分布，将其与WRF模式输出格点进行空间配准，并计算其与对应时刻WRF预报气温之间的差值；
针对每个地面气象观测站点位置，比较两种数据之间的差异，选择差异较小的站点作为备选参考站点。
4.根据权利要求1所述的动力与统计降尺度相结合的山区气温精细化预报方法，其特征在于，在不同季节分别开展山区地面气象站点观测气温的空间代表性分析，并比较不同站点位置的WRF预报、遥感估算差值的平均值、标准差，选择平均值、标准差最小的一个或者多个地面气象站点作为区域空间代表性程度高的站点。
5.根据权利要求1所述的动力与统计降尺度相结合的山区气温精细化预报方法，其特征在于，步骤S4中，选择一阶自适应卡尔曼滤波方法，对随机森林降尺度得到的参考站点未来24小时预报气温进行订正的过程包括以下步骤：
利用公式（1），在有观测值的t时刻，计算随机森林降尺度得到的预报气温
与站点观测气温/>
之差，得到参考站点位置的模式预报误差/>
：
采用离散逼近法迭代计算模式误差
：
t+1时刻，模式预报订正结果
为：
对于无观测值的预报时段，利用随机森林降尺度预报的订正结果
替代公式（1）中的
；
利用有站点观测数据的时段，通过试验设定权重系数w的取值范围为0.01~0.99，变化间隔0.01，分析比较权重系数w的不同取值情况下24小时预报的总误差，确定总误差最小时对应的取值作为权重系数w的最终值，将其用于WRF模式预报气温订正。
6.根据权利要求1所述的动力与统计降尺度相结合的山区气温精细化预报方法，其特征在于，步骤S5中，基于地表能量平衡理论对山区气温进行优化预报；具体包括：
在30m空间分辨率的尺度上，通过分析地形、下垫面特性参数对山区不同地点各通量项以及气温的影响差异，以代表性参考站点降尺度后的气温预报数据作为背景气温，建立公式（4）用于模拟预报山区其它地点的气温；利用参考气象站点计算其与其它地点之间的通量项差异，并在公式（4）中考虑温度平流差异，在一定程度上可以降低山区地表能量平衡闭合效果差的影响，而且利用参考站点预报气温能够约束求解公式（4），提高气温预报精度；
其中，t为预报时刻；
为山区待模拟计算地点的预报气温；/>
为代表性参考站点降尺度预报气温的卡尔曼滤波订正结果；/>
为各通量变化引起的待模拟计算地点气温变化值，/>
分别为各通量变化引起的参考站点气温变化值；/>
为待模拟计算地点与参考站点因高度不同而产生的气温差异值；/>
为冷暖空气水平运动引起的温度平流差异值。
7.根据权利要求6所述的动力与统计降尺度相结合的山区气温精细化预报方法，其特征在于，在参考站点位置上，采用公式（5）计算各通量变化引起的气温变化值
：
式中，
是净辐射通量，在WRF模式预报和遥感监测方法基础上，根据山区复杂地形与下垫面对入射短波辐射、发射长波辐射的影响以及周围环境散射反射与辐射的贡献计算得到；G是地表热储存通量，根据净辐射通量/>
与地表热储存通量之间的经验关系模型估算得到；L是潜热通量，m是湍流交换系数，/>
是地表显热通量阻抗系数；
计算出待模拟计算地点的净辐射通量、地表热储存通量、潜热通量、湍流交换系数和地表显热通量阻抗系数后，将其代入公式（5），得到待模拟计算地点气温变化值
。
8.根据权利要求7所述的动力与统计降尺度相结合的山区气温精细化预报方法，其特征在于，通过不同季节与月份站点观测气温随高度变化规律的定量分析，以及t时刻WRF模式气温预报结果，计算山区不同地点的气温垂直递减率
，并结合模拟计算地点与参考站点之间的高度差异值/>
，采用公式（6）计算得到待模拟计算地点与参考站点因高度不同而产生的气温差异值/>
：
根据参考站点与模拟地点的温度水平分布梯度，以及风矢量数据，采用公式（7）计算得到冷暖空气水平运动引起的温度平流差异值
：
式中，u、v分别是东西、南北方向的风速；
、/>
分别表示东西、南北方向的温度差异。
9.根据权利要求8所述的动力与统计降尺度相结合的山区气温精细化预报方法，其特征在于，步骤S5中，对于拟模拟山区地点，利用WRF模式预报、随机森林降尺度结果、遥感监测数据，估算或以
表达公式（5）至公式（7）中的各项参数，并以参考站点预报气温
为驱动，通过公式（4）迭代求解/>
。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310201913.9A CN116068672B (zh) | 2023-03-06 | 2023-03-06 | 一种动力与统计降尺度相结合的山区气温精细化预报方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310201913.9A CN116068672B (zh) | 2023-03-06 | 2023-03-06 | 一种动力与统计降尺度相结合的山区气温精细化预报方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN116068672A CN116068672A (zh) | 2023-05-05 |
| CN116068672B true CN116068672B (zh) | 2023-06-16 |
Family
ID=86182135
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202310201913.9A Active CN116068672B (zh) | 2023-03-06 | 2023-03-06 | 一种动力与统计降尺度相结合的山区气温精细化预报方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN116068672B (zh) |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107390298A (zh) * | 2017-07-19 | 2017-11-24 | 云南电网有限责任公司电力科学研究院 | 一种复杂山地下垫面大风的模拟方法及装置 |
| CN110941790A (zh) * | 2019-09-27 | 2020-03-31 | 成都信息工程大学 | 基于高分辨率数值的无人机低空飞行气象信息处理方法 |
| CN112100921A (zh) * | 2020-09-16 | 2020-12-18 | 平衡机器科技(深圳)有限公司 | 一种基于wrf和随机森林获取风资源风速的方法 |
| CN112180472A (zh) * | 2020-09-28 | 2021-01-05 | 南京北极光智能科技有限公司 | 一种基于深度学习的大气能见度集成预报方法 |
| CN115032720A (zh) * | 2022-07-15 | 2022-09-09 | 国网上海市电力公司 | 基于随机森林的多模式集成预报在地面气温预报中的应用 |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| KR102591379B1 (ko) * | 2021-07-30 | 2023-10-19 | 대한민국 | 랜덤포레스트 기계학습과 기상 관측 및 모델 예측 자료를 이용한 실시간 시정 산출 방법 |
-
2023
- 2023-03-06 CN CN202310201913.9A patent/CN116068672B/zh active Active
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107390298A (zh) * | 2017-07-19 | 2017-11-24 | 云南电网有限责任公司电力科学研究院 | 一种复杂山地下垫面大风的模拟方法及装置 |
| CN110941790A (zh) * | 2019-09-27 | 2020-03-31 | 成都信息工程大学 | 基于高分辨率数值的无人机低空飞行气象信息处理方法 |
| CN112100921A (zh) * | 2020-09-16 | 2020-12-18 | 平衡机器科技(深圳)有限公司 | 一种基于wrf和随机森林获取风资源风速的方法 |
| CN112180472A (zh) * | 2020-09-28 | 2021-01-05 | 南京北极光智能科技有限公司 | 一种基于深度学习的大气能见度集成预报方法 |
| CN115032720A (zh) * | 2022-07-15 | 2022-09-09 | 国网上海市电力公司 | 基于随机森林的多模式集成预报在地面气温预报中的应用 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN116068672A (zh) | 2023-05-05 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN105403664B (zh) | 一种基于wrf‑chem的大型点污染源大气环境影响评价方法 | |
| Linares-Rodríguez et al. | Generation of synthetic daily global solar radiation data based on ERA-Interim reanalysis and artificial neural networks | |
| Kedia et al. | Estimates of change in surface meteorology and urban heat island over northwest India: Impact of urbanization | |
| Evans et al. | Middle East climate simulation and dominant precipitation processes | |
| Hart et al. | Daily reference evapotranspiration for California using satellite imagery and weather station measurement interpolation | |
| Xu et al. | Spatial and temporal dynamics of urban heat island and their relationship with land cover changes in urbanization process: a case study in Suzhou, China | |
| CN112785024A (zh) | 一种基于流域水文模型的径流计算和预测方法 | |
| CN109657988B (zh) | 基于hasm和欧氏距离算法的烟叶品质分区方法 | |
| Tiwari et al. | Sensitivity of the Himalayan orography representation in simulation of winter precipitation using Regional Climate Model (RegCM) nested in a GCM | |
| Sahu et al. | Assimilating observations to simulate marine layer stratocumulus for solar forecasting | |
| CN112784395B (zh) | 一种预测和模拟河流水体总磷浓度的方法 | |
| Khesali et al. | A method in near-surface estimation of air temperature (NEAT) in times following the satellite passing time using MODIS images | |
| Dong et al. | Simulation of United States mesoscale convective systems using GFDL’s new high-resolution general circulation model | |
| Wu et al. | The impact of multi-sensor land data assimilation on river discharge estimation | |
| CN117493476A (zh) | 一种融合物理机制和人工智能的径流回溯模拟方法及系统 | |
| Nitis et al. | The effects of basic flow and topography on the development of the sea breeze over a complex coastal environment | |
| CN115544706A (zh) | 一种小波和XGBoost模型集成的大气细颗粒物浓度估算方法 | |
| Pattnayak et al. | How good is regional climate model version 4 in simulating the monsoon onset over Kerala? | |
| Gerelchuluun et al. | Air temperature distribution over Mongolia using dynamical downscaling and statistical correction. | |
| CN116068672B (zh) | 一种动力与统计降尺度相结合的山区气温精细化预报方法 | |
| Nykanen | Linkages between orographic forcing and the scaling properties of convective rainfall in mountainous regions | |
| Finsen et al. | Using radar altimetry to update a large-scale hydrological model of the Brahmaputra river basin | |
| CN116258065A (zh) | 一种云天高时空分辨率近地面气温估算方法 | |
| CN107688712A (zh) | 一种基于dem和ndvi的气温降尺度方法 | |
| Nguyen et al. | Observation and simulation of wind speed and wind power density over bac Lieu region |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |