CN112507558B - 一种基于主控因子修正的参照作物腾发量预测方法 - Google Patents
一种基于主控因子修正的参照作物腾发量预测方法 Download PDFInfo
- Publication number
- CN112507558B CN112507558B CN202011473915.6A CN202011473915A CN112507558B CN 112507558 B CN112507558 B CN 112507558B CN 202011473915 A CN202011473915 A CN 202011473915A CN 112507558 B CN112507558 B CN 112507558B
- Authority
- CN
- China
- Prior art keywords
- value
- follows
- forecast
- control factor
- calculation formula
- 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
- 238000012937 correction Methods 0.000 title claims abstract description 45
- 238000000034 method Methods 0.000 title claims abstract description 33
- 230000000630 rising effect Effects 0.000 title abstract description 13
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 claims abstract description 16
- 238000001914 filtration Methods 0.000 claims abstract description 6
- 238000004364 calculation method Methods 0.000 claims description 79
- 230000002354 daily effect Effects 0.000 claims description 33
- 230000005855 radiation Effects 0.000 claims description 19
- 230000035945 sensitivity Effects 0.000 claims description 17
- 238000005259 measurement Methods 0.000 claims description 12
- 238000012216 screening Methods 0.000 claims description 7
- 238000010219 correlation analysis Methods 0.000 claims description 5
- 238000010206 sensitivity analysis Methods 0.000 claims description 5
- 238000012821 model calculation Methods 0.000 claims description 4
- 230000003203 everyday effect Effects 0.000 claims description 3
- 230000004907 flux Effects 0.000 claims description 3
- 239000002689 soil Substances 0.000 claims description 3
- 230000000875 corresponding effect Effects 0.000 description 8
- 230000009286 beneficial effect Effects 0.000 description 7
- 238000004458 analytical method Methods 0.000 description 5
- 238000010586 diagram Methods 0.000 description 5
- 230000000694 effects Effects 0.000 description 4
- 230000008859 change Effects 0.000 description 3
- 230000007423 decrease Effects 0.000 description 3
- 230000000052 comparative effect Effects 0.000 description 2
- 239000012141 concentrate Substances 0.000 description 2
- 230000001276 controlling effect Effects 0.000 description 2
- 238000011161 development Methods 0.000 description 2
- 238000011156 evaluation Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000012544 monitoring process Methods 0.000 description 2
- 238000004451 qualitative analysis Methods 0.000 description 2
- 230000009467 reduction Effects 0.000 description 2
- 230000002596 correlated effect Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 230000002262 irrigation Effects 0.000 description 1
- 238000003973 irrigation Methods 0.000 description 1
- 239000003621 irrigation water Substances 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 238000012876 topography Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/08—Thermal analysis or thermal optimisation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Air Conditioning Control Device (AREA)
Abstract
本发明公开了一种基于主控因子修正的参照作物腾发量预测方法。该方法包括收集预测地历史实测气象数据；并利用PMT模型计算出历年参照作物腾发量；采用基于多参数法确定PMT模型中的主控因子；采用基于贡献率的卡尔曼滤波法对收集的天气预报值中的PMT模型主控因子预报值进行修正；将修正完成的PMT模型主控因子预报值和其它预测地天气预报值代入PMT模型中进行计算，得到基于主控因子修正的参照作物腾发量预测结果。本发明为了解决现有技术中中长期参照作物腾发量预测偏差较大的问题，采用了基于贡献率的卡尔曼滤波法对气象预报数据中的主控因子进行按月修正，提高了参照作物腾发量预测的准确性，为中长期农业需水量准确预报提供了技术支持，实用性较强。
Description
技术领域
本发明涉及一种参照作物腾发量预测方法，尤其涉及一种基于主控因子修正的参照作物腾发量预测方法。
背景技术
我国属于农业大国，农业灌溉用水量占比大，当下水资源短缺已成为制约农业发展的瓶颈，要实现灌溉农业的可持续发展，需要对农业需水量进行精准计算以及预报。现有两类作物需水量的计算方法，一是直接计算出作物需水量，二是通过计算参照作物腾发量(ET0)来间接计算作物需水量，直接计算法多依赖大量实测资料，受限于监测费用与监测人数等因素，难以大规模推广应用，故现阶段用参照作物腾发量来计算作物需水量较为普遍。
现有技术目前普遍采用Penman-Monteith(P-M)公式计算的ET0，但利用实测气象资料计算获得ET0难以预测中长期ET0变化情况，利用天气预报信息去估算ET0已经被证实是可行的，但现有技术多集中于利用短期天气预报解析数据直接代入ET0预报模型并对其预报精度进行评价，而采用未经修正的气象预报解析数据进行中长期ET0的预测偏差较大，精确度较低。
发明内容
针对现有技术中的上述不足，本发明提供了一种基于主控因子修正的参照作物腾发量预测方法。
为了达到上述发明目的，本发明采用的技术方案为：
一种基于主控因子修正的参照作物腾发量预测方法，包括以下步骤：
S1、收集预测地PMT模型计算所需历史实测气象数据；
S2、将步骤S1收集的预测地历史实测气象数据代入PMT模型中计算出历年参照作物腾发量；
S3、根据步骤S2得到的历年参照作物腾发量采用基于多参数法确定PMT模型中的主控因子；
S4、收集预测地PMT模型计算所需天气预报值；
S5、采用基于贡献率的卡尔曼滤波法对步骤S4中收集的天气预报值中的PMT模型主控因子预报值进行修正；
S6、将步骤S5中修正完成的PMT模型主控因子预报值和步骤S4收集的其它预测地天气预报值代入PMT模型中进行计算，得到基于主控因子修正的参照作物腾发量预测结果。
本发明具有以下有益效果：利用各主控因子月贡献量差异，采用按照月份划分修正气象预报值中的PMT模型主控因子的方法，提高了参照作物腾发量的预测精度，具有实用性。
优选地，步骤S1具体包括收集过去30年以上的长系列预测地日最高气温、日最低气温和风速的实测值。
该优选方案具有以下有益效果：收集过去30年以上的长系列气象信息保证了本发明后续确定主控因子和对主控因子进行修正的过程有足够的基数进行运算，确保了确定和修正主控因子结果的准确性。
优选地，步骤S2包括以下分步骤：
S21、根据日最高气温和日最低气温计算出平均温度T，计算公式为：
T＝(Tmax+Tmin)/2
其中，Tmax为日最高气温，Tmin为日最低气温；
S22、根据日最高气温和日最低气温计算出饱和水汽压es，计算公式为：
S23、根据平均温度T计算出饱和水汽压曲线斜率Δ，计算公式为：
S24、根据日最低气温计算出实际水汽压ea，计算公式为：
S25、根据待测地纬度值和当天序数J计算当天天文辐射Ra，计算公式为：
其中，Gsc为太阳常数，δ为太阳倾角，dr为日地相对距离，ωs为日落时角；
根据日最高气温，日最低气温和天文辐射Ra计算出太阳短波辐射Rs，计算公式为：
Rs＝Kr(Tmax-Tmin)0.5Ra
其中，Kr为地区调整系数；
根据日最高气温，日最低气温、天文辐射Ra和太阳短波辐射Rs计算出作物表面上净辐射量Rn，计算公式为：
其中，Z为预测地高程；
S26、根据不同高程的风速值计算出2m高处风速值u2，计算公式为：
其中，uz为z米高处风速值。
S27、根据步骤S21到S26的计算结果通过PMT模型计算历年参照作物腾发量，计算公式为：
其中，ET0为参照作物腾发量，G为土壤热通量，γ为湿度计常数。
该优选方案具有以下有益效果：Penman-Monteith(P-M)模型计算的ET0被应用在任何气候区都是最精确的，PMT模型与PM模型的区别在于通过日最高温度和日最低温度来计算太阳短波辐射而非通过天气类型，采用PMT模型作为预测模型既保证了准确率也保证了数据获取简单方便。
优选地，步骤S3包括以下分步骤：
S31、将日最高温度，日最低温度和风速及对应ET0进行相关性分析，计算出每天各值和ET0的相关系数r，计算公式为：
其中，xi为控制因子，即日最高温度，日最低温度或风速的值，yi为ET0值，为xi的平均值，/>为yi的平均值；
S32、将日最高温度，日最低温度和风速及对应ET0进行敏感性分析，计算出每天各值和ET0的敏感系数Sx，计算公式为：
其中，ΔET0为ET0的变化量，x为各控制因子的值，Δx为x的变化量；
S33、根据步骤S32计算的各控制因子的敏感系数，计算出各控制因子对ET0的贡献率RCx和贡献量Gx，计算公式为：
其中，Trendx为某控制因子x的气候倾向率；|avx|为某控制因子x的多年均值；为ET0的多年日均值。
S34、根据步骤S31、步骤S32和步骤S33的结果，按照筛选规则确定PMT模型中的主控因子，筛选规则为：
Main＝Max[abs(r)or abs(Sx)or abs(RCx)or abs(Gx)]
其中，Main为主控因子，Max[]为取最大值对应的气象因子，abs()为取绝对值。
该优选方案具有以下有益效果：采用筛选规则，通过相关性分析、敏感性分析、贡献量分析和贡献率分析来共同确定PMT模型的主控因子，确保了主控因子确定结果的准确性。
优选地，步骤S4具体包括收集预测地中长期天气预报风速等级，并解析获得PMT模型计算所需预报风速值，收集预测地中长期天气预报日最高温度预报值，日最低温度预报值。
该优选方案具有以下有益效果：获取的值为天气预报中易于获取的值，预测所需数据获取方法简单方便。
优选地，步骤S5包括以下分步骤：
S51、对步骤S34确定的主控因子进行1至30d预报结果精度计算，确定预报精度明显下降时段；
S52、对步骤S51确定的预报精度明显下降时段的主控因子预报值进行修正，依据主控因子在该时段各个月份的贡献率取值范围分成4组，表示为：
其中ωG为权重系数；
以步长0.1在各分组0到1范围内进行权重系数试算，并计算相应的订正结果，通过各分组比较不同权重系数下订正后的预报值精度来确定该分组的最佳权重系数，分组规则为对该时段每个月份在0到1范围内以步长0.1来调整权重系数进行计算，确定各月最佳权重系数，计算公式为：
其中，t为当前预报时间，P(t)为订正误差值，b(t-1)为前一日主控因子预报值和实测值误差，P(t-1)为前一日订正误差值；
S53、根据步骤S52确定的最佳权重系数对PMT计算模型中的主控因子进行修正，修正公式为：
F(t)＝f(t)-P(t)
其中，F(t)为修正预报值；f(t)为预报值。
该优选方案具有以下有益效果：通过逐月对权重系数进行调整的方式提高了主控因子预报值的精度，进而提高了通过PMT模型预测参照作物腾发量的预测精度。
优选地，步骤S51包括以下分步骤：
S511、计算回归系数b，计算公式为：
其中，oi为第i个实测数据，pi为第i个预测数据，i＝1,2，…，n；n为统计样本数；
S512、计算决定系数R2，计算公式为：
其中，为实测数据组的平均值，/>为预测数据组的平均值；
S513、计算均方根误差RMSE，计算公式为：
S514、计算相对误差RE，计算公式为：
S515、计算平均绝对误差MAE，计算公式为：
S516、计算一致性指数d，计算公式为：
S517、根据步骤S511到步骤S516的计算结果，确定预报精度明显下降时段。
该优选方案具有以下有益效果：通过回归系数b、决定系数R2、均方根误差RMSE、相对误差RE、平均绝对误差MAE和一致性指数d来评价主控因子预报值和主控因子实测值的误差，保证了精度评价结果和精度下降阶段选取的准确性。
附图说明
图1是本发明一种基于主控因子修正的参照作物腾发量预测方法流程示意图；
图2是本发明实施例中控制因子与参照作物腾发量的相关系数示意图；
图3是本发明实施例中控制因子与参照作物腾发量的敏感系数示意图；
图4是本发明实施例中控制因子与参照作物腾发量的贡献率示意图；
图5是本发明实施例中控制因子与参照作物腾发量的贡献量示意图；
图6是本发明实施例中修正前PMT模型预报精度计算结果示意图；
图7是本发明实施例中主控因子与参照作物腾发量的贡献率箱型示意图；
图8是本发明实施例中修正后的PMT模型及修正前的PMT模型预测ET0的均方根误差RMSE对比示意图；
图9是本发明实施例中修正后的PMT模型及修正前的PMT模型预测ET0的平均绝对误差MAE对比示意图；
图10是本发明实施例中修正后的PMT模型及修正前的PMT模型预测ET0的相对误差RE对比示意图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本发明，并不用于限定本发明。
请参照图1，本发明提供了一种基于主控因子修正的参照作物腾发量预测方法，包括以下步骤：
S1、收集预测地PMT模型计算所需历史实测气象数据；
本发明实施例中，步骤S1具体包括收集过去30年以上的长系列预测地日最高气温、日最低气温和风速的实测值。
实测最高气温、最低气温和风速资料来源于国家气象信息中心地面气候资料日值数据集。
S2、将步骤S1收集的预测地历史实测气象数据代入PMT模型中计算出历年参照作物腾发量；
本发明实施例中，步骤S2包括以下分步骤：
S21、根据日最高气温和日最低气温计算出平均温度T，计算公式为：
T＝(Tmax+Tmin)/2
其中，Tmax为日最高气温，Tmin为日最低气温；
S22、根据日最高气温和日最低气温计算出饱和水汽压es，计算公式为：
S23、根据平均温度T计算出饱和水汽压曲线斜率Δ，计算公式为：
S24、根据日最低气温计算出实际水汽压ea，计算公式为：
S25、根据待测地纬度值(采用弧度单位表示)和当天序数J(当天为一年中的第几天)，计算当天天文辐射Ra，计算公式为：
其中，Gsc为太阳常数，值为0.082MJ/(m2·min)，δ为太阳倾角，其值为dr为日地相对距离，其值为/>ωs为日落时角，其值为/>
然后根据日最高气温，日最低气温和天文辐射Ra计算出太阳短波辐射Rs，计算公式为：
Rs＝Kr(Tmax-Tmin)0.5Ra
其中，Kr为地区调整系数，内陆地区和沿海地区分别取0.16和0.19；
最后根据日最高气温，日最低气温、天文辐射Ra和太阳短波辐射Rs计算出作物表面上净辐射量Rn，计算公式为：
S26、根据不同高程的风速值计算出2m高处风速值u2，计算公式为：
其中，uz为z米高处风速值。
根据天气预报中风力等级的预报信息，参照GB/T 35227-2017《地面气象观测规范风向和风速》，可将风力等级转化为一个具体的范围值，将其中不同高程处所测风速转化为2m处数值用以代入PMT模型。
S27、根据所述步骤S21到S26的计算结果通过PMT模型计算历年参照作物腾发量，计算公式为：
其中，ET0为参照作物腾发量，G为土壤热通量，这里取0，γ为湿度计常数，值为0.66hPa/℃。
S3、根据步骤S2得到的历年参照作物腾发量采用基于多参数法确定PMT模型中的主控因子；
本发明实施例中，步骤S3包括以下分步骤：
S31、将日最高温度，日最低温度和风速及对应ET0进行相关性分析，计算出每天各值和ET0的相关系数r，计算公式为：
其中，xi为控制因子，即日最高温度，日最低温度或风速的值，yi为ET0值，为xi的平均值，/>为yi的平均值；
皮尔逊相关系数r是最常用的相关系数，本文用于衡量各气象要素与ET0间的线性关系，定性分析PMT模型中ET0的主控因子，r的取值在[-1，1]之间，绝对值越接近于1，表明相关性越强，绝对值越接近于0则表明相关度越弱。若r>0，变量之间为正相关；若r<0，变量之间为负相关。
S32、将日最高温度，日最低温度和风速及对应ET0进行敏感性分析，计算出每天各值和ET0的敏感系数Sx，计算公式为：
其中，ΔET0为ET0的变化量，x为各控制因子的值，Δx为x的变化量，敏感系数Sx绝对值的大小反映了各控制因子敏感程度的大小；
采用敏感系数描述各气象要素对ET0的影响大小，定性分析PMT模型中ET0的主控因子，敏感系数定义为潜在蒸散发变化与气象要素变化率之比，敏感系数正或负分别表示ET0随着气象因子的增加而增加或减少，其绝对值的大小反应敏感程度的大小。
S33、根据步骤S32计算的各控制因子的敏感系数，计算出各控制因子对ET0的贡献率RCx和贡献量Gx，计算公式为：
其中，Trendx为某控制因子x的气候倾向率；|avx|为某控制因子x的多年均值；为ET0的多年日均值。
S34、根据步骤S31、步骤S32和步骤S33的结果，按照筛选规则确定PMT模型中的主控因子，筛选规则为：
Main＝Max[abs(r)or abs(Sx)or abs(RCx)or abs(Gx)]
其中，Main为主控因子，Max[]为取最大值对应的气象因子，abs()为取绝对值。
S4、收集预测地PMT模型计算所需天气预报值；
本发明实施例中，步骤S4具体包括收集预测地中长期天气预报风速等级，并解析获得PMT模型计算所需预报风速值，收集预测地中长期天气预报日最高温度预报值，日最低温度预报值。
预报周期为1～30d的日最高气温、日最低气温、风力等级等公共天气预报数据来源于中国天气网。
S5、采用基于贡献率的卡尔曼滤波法对步骤S4中收集的天气预报值中的PMT模型主控因子预报值进行修正；
本发明实施例中，步骤S5包括以下分步骤：
步骤S5包括以下分步骤：
S51、对步骤S34确定的主控因子进行1至30d预报结果精度计算，确定预报精度明显下降时段；
本发明实施例中，步骤S51包括以下分步骤：
S511、计算回归系数b，计算公式为：
其中，oi为第i个实测数据，pi为第i个预测数据，i＝1,2，…，n；n为统计样本数；
S512、计算决定系数R2，计算公式为：
其中，为实测数据组的平均值，/>为预测数据组的平均值；
S513、计算均方根误差RMSE，计算公式为：
S514、计算相对误差RE，计算公式为：
S515、计算平均绝对误差MAE，计算公式为：
S516、计算一致性指数d，计算公式为：
S517、根据步骤S511到步骤S516的计算结果，确定预报精度明显下降时段。
S52、对步骤S51确定的预报精度明显下降时段的主控因子预报值进行修正，依据主控因子在该时段各个月份的贡献率取值范围分成4组，表示为：
其中ωG为权重系数；
以步长0.1在各分组0到1范围内进行权重系数试算，并计算相应的订正结果，通过各分组比较不同权重系数下订正后的预报值精度来确定该分组的最佳权重系数，分组规则为对该时段每个月份在0到1范围内以步长0.1来调整权重系数进行计算，确定各月最佳权重系数，计算公式为：
其中，t为当前预报时间，P(t)为订正误差值，b(t-1)为前一日主控因子预报值和实测值误差，P(t-1)为前一日订正误差值；
S53、根据步骤S52确定的最佳权重系数对PMT计算模型中的主控因子进行修正，修正公式为：
F(t)＝f(t)-P(t)
其中，F(t)为修正预报值；f(t)为预报值。
S6、将步骤S5中修正完成的PMT模型主控因子预报值和步骤S4收集的其它预测地天气预报值代入PMT模型中进行计算，得到基于主控因子修正的参照作物腾发量预测结果。
接下来以北京大兴为实例来对本发明进行详细说明。
大兴区在北京市南部、华北平原东北部，东临通州区，南临河北省固安县、霸州市等，西与房山区隔永定河为邻，北接丰台、朝阳区。东经116°13′-116°43′，北纬39°26′-39°51′。全境属永定河冲积平原，地势自西向东南缓倾，大部分地区海拔14～52m之间，总面积约1030.57km2。
根据本发明步骤S1到S2进行处理过后，步骤S31相关性分析的结果如图2所示，由图可知，最高温度与ET0的相关系数最大值出现在第178d为0.942，最小值出现在第261d为0.306，且一年中最高温与ET0的相关系数呈显著相关占比达100％；最低温度与ET0的相关系数最大值出现在第60d为0.576，最小值出现在第273d为-0.559，一年中最低温与ET0的相关系数呈不显著相关占比高达77.32％；风速与ET0的相关系数最大值出现在第361d为0.871，最小值出现在第124d为-0.240，一年中风速与ET0的相关系数呈显著相关占比为39.89％呈不显著相关占比达45.63％。由此可知，最高温度对PMT模型的影响最大，但冬季各气象因素对模型的影响均较小且影响程度相近。
步骤S32敏感性分析结果如图3所示，由图可知，自第37-341d，最高温度敏感系数的绝对值均大于最低温度和风速，表明在此期间最高温度对ET0的影响最大，其余时间段最高、最低温度以及风速的敏感系数相差不大，敏感系数变化范围分别为0.044～2.932，-1.681～0.468，0.083～0.316，均值分别为1.488，-0.445，0.208，表明全年时段最高温度对ET0的影响最大，其次是最低温度和风速。
由各气象因子的相关系数以及敏感系数可知，冬季各气象因素对ET0的影响较小，夏季对ET0影响最大，且最高温度发挥的作用要比最低温度和风速大。
步骤S33贡献率分析结果如图4所示，贡献量分析结果如图5所示，由图可知，各气象因素的贡献率主要集中在-1.5mm～1.5mm之间变化，贡献率主要集中在-30％～30％之间变化，最高温度在年内变幅比最低温度和风速大。最高温度贡献率最大为52.5％，贡献量最大为2.47mm，日均贡献量为0.46mm；最低温度贡献率最大为34.2％，贡献量最大为0.31mm，日均贡献量为-0.37mm；风速贡献率最大为13.07％，贡献量最大为0.18mm，日均贡献量为-0.10mm。由此可知最高温度对ET0的贡献最大其次是最低温度和风速，需要指出的是各气象因素的贡献量在夏季较大冬季较小，且冬季各因素贡献量相近。
依据各参数绝对值最大为选择依据，确定最高温度为PMT模型的主控因子。
接下来对修正前PMT模型的预报进行精度评价，图6为修正前PMT预报模型不同预见期的回归系数b、决定系数R2、一致性指数d、相对误差RE、平均绝对误差MAE以及均方根误差RMSE的分析结果。由结果可知，在修正前，模型的b和d在第16d出现了较为明显的下降趋势，回归系数b和一致性指数d分别由15d预报的0.92、0.93下降到0.56、0.76；与此同时模型的各类误差指标出现明显增大，模型的RE、RMSE、MAE值分别由15d预报的0.35mm/day、1.06mm/d、0.80mm/d增大到16d预报的0.57mm/d、1.79mm/d、1.40mm/d，说明PMT模型的精度自第16d起出现显著下降。
然后参照图7的贡献率箱线图，将各个月份贡献率均值按照分类规则进行分类，贡献率在0～10％之间的月份为1、10、11、12月份，贡献率为10％～15％之间的月份为4、6、9月份，贡献率为15％～25％之间的月份为2、5、7、8月份，贡献率大于25％的月份为3月份。将16-30d最高温预报结果按上述分类结果分别代入0～1区间内以0.1为步长的权重系数，将计算结果进行精度评价，找出不同月份的权重系数ωG。
按月份进行分组后各预见期权重系数选择结果见表1。根据表1给出的权重系数选择结果代入基于贡献率的卡尔曼滤波法修正公式得到修正后的预报最高温度并代入PMT模型得到修正后的ET0数值。
表1基于贡献率的卡尔曼滤波法权重系数选择结果
对修正后的ET0数值进行精度评价，修正后的PMT模型(RCD)与为修正前的PMT模型(OR)预测ET0的均方根误差RMSE差异请参照图8，相对误差RE差异请参照图9，平均绝对误差MAE差异请参照图10，从图中可看出在修正后三个值均出现了明显下降，降低范围分别为0.70～0.90mm/d、0.22～0.29mm/d、0.55～0.72mm/d，表明本发明相较于现有技术预测精度有了显著的提升。
综上所述，本发明提供了一种基于主控因子修正的参照作物腾发量预测方法，该方法获取数据方法简单方便，中长期预测结果的精确性相较于现有技术有了明显提升，具有实用性。
本领域的普通技术人员将会意识到，这里所描述的实施例是为了帮助读者理解本发明的原理，应被理解为本发明的保护范围并不局限于这样的特别陈述和实施例。本领域的普通技术人员可以根据本发明公开的这些技术启示做出各种不脱离本发明实质的其它各种具体变形和组合，这些变形和组合仍然在本发明的保护范围内。
Claims (4)
1.一种基于主控因子修正的参照作物腾发量预测方法，其特征在于，包括以下步骤：
S1、收集预测地PMT模型计算所需历史实测气象数据；
S2、将所述步骤S1收集的预测地历史实测气象数据代入PMT模型中计算出历年参照作物腾发量；包括以下分步骤：
S21、根据日最高气温和日最低气温计算出平均温度T，计算公式为：
T＝(Tmax+Tmin)/2
其中，Tmax为日最高气温，Tmin为日最低气温；
S22、根据日最高气温和日最低气温计算出饱和水汽压es，计算公式为：
S23、根据平均温度T计算出饱和水汽压曲线斜率Δ，计算公式为：
S24、根据日最低气温计算出实际水汽压ea，计算公式为：
S25、根据待测地纬度值和当天序数J计算当天天文辐射Ra，计算公式为：
其中，Gsc为太阳常数，δ为太阳倾角，dr为日地相对距离，ωs为日落时角；
根据日最高气温，日最低气温和天文辐射Ra计算出太阳短波辐射Rs，计算公式为：
Rs＝Kr(Tmax-Tmin)0.5Ra
其中，Kr为地区调整系数；
根据日最高气温，日最低气温、天文辐射Ra和太阳短波辐射Rs计算出作物表面上净辐射量Rn，计算公式为：
其中，Z为预测地高程；
S26、根据不同高程的风速值计算出2m高处风速值u2，计算公式为：
其中，uz为z米高处风速值；
S27、根据所述步骤S21到S26的计算结果通过PMT模型计算历年参照作物腾发量，计算公式为：
其中，ET0为参照作物腾发量，G为土壤热通量，γ为湿度计常数；
S3、根据所述步骤S2得到的历年参照作物腾发量采用基于多参数法确定PMT模型中的主控因子；包括以下分步骤：
S31、将日最高温度，日最低温度和风速及对应ET0进行相关性分析，计算出每天各值和ET0的相关系数r，计算公式为：
其中，xi为控制因子，即日最高温度，日最低温度或风速的值，yi为ET0值，为xi的平均值，/>为yi的平均值；
S32、将日最高温度，日最低温度和风速及对应ET0进行敏感性分析，计算出每天各值和ET0的敏感系数Sx，计算公式为：
其中，ΔET0为ET0的变化量，x为各控制因子的值，Δx为x的变化量；
S33、根据步骤S32计算的各控制因子的敏感系数，计算出各控制因子对ET0的贡献率RCx和贡献量Gx，计算公式为：
其中，Trendx为某控制因子x的气候倾向率；|avx|为某控制因子x的多年均值；为ET0的多年日均值；
S34、根据所述步骤S31、步骤S32和步骤S33的结果，按照筛选规则确定PMT模型中的主控因子，筛选规则为：
Main＝Max[abs(r)or abs(Sx)or abs(RCx)or abs(Gx)]
其中，Main为主控因子，Max[]为取最大值对应的气象因子，abs()为取绝对值；
S4、收集预测地PMT模型计算所需天气预报值；
S5、采用基于贡献率的卡尔曼滤波法对所述步骤S4中收集的天气预报值中的PMT模型主控因子预报值进行修正；包括以下分步骤：
S51、对所述步骤S34确定的主控因子进行1至30d预报结果精度计算，确定预报精度明显下降时段；
S52、对所述步骤S51确定的预报精度明显下降时段的主控因子预报值进行修正，依据主控因子在该时段各个月份的贡献率取值范围分成4组，表示为：
其中ωG为权重系数；
以步长0.1在各分组0到1范围内进行权重系数试算，并计算相应的订正结果，通过各分组比较不同权重系数下订正后的预报值精度来确定该分组的最佳权重系数，分组规则为对该时段每个月份在0到1范围内以步长0.1来调整权重系数进行计算，确定各月最佳权重系数，计算公式为：
其中，t为当前预报时间，P(t)为订正误差值，b(t-1)为前一日主控因子预报值和实测值误差，P(t-1)为前一日订正误差值；
S53、根据所述步骤S52确定的最佳权重系数对PMT计算模型中的主控因子进行修正，修正公式为：
F(t)＝f(t)-P(t)
其中，F(t)为修正预报值；f(t)为预报值；
S6、将所述步骤S5中修正完成的PMT模型主控因子预报值和所述步骤S4收集的其它预测地天气预报值代入PMT模型中进行计算，得到基于主控因子修正的参照作物腾发量预测结果。
2.如权利要求1所述的一种基于主控因子修正的参照作物腾发量预测方法，其特征在于：所述步骤S1具体包括收集过去30年以上的长系列预测地日最高气温、日最低气温和风速的实测值。
3.如权利要求1所述的一种基于主控因子修正的参照作物腾发量预测方法，其特征在于，所述步骤S4具体包括收集预测地中长期天气预报风速等级，并解析获得PMT模型计算所需预报风速值，收集预测地中长期天气预报日最高温度预报值，日最低温度预报值。
4.如权利要求1所述的一种基于主控因子修正的参照作物腾发量预测方法，其特征在于，所述步骤S51包括以下分步骤：
S511、计算回归系数b，计算公式为：
其中，oi为第i个实测数据，pi为第i个预测数据，i＝1,2，…，n；n为统计样本数；
S512、计算决定系数R2，计算公式为：
其中，为实测数据组的平均值，/>为预测数据组的平均值；
S513、计算均方根误差RMSE，计算公式为：
S514、计算相对误差RE，计算公式为：
S515、计算平均绝对误差MAE，计算公式为：
S516、计算一致性指数d，计算公式为：
S517、根据所述步骤S511到步骤S516的计算结果，确定预报精度明显下降时段。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011473915.6A CN112507558B (zh) | 2020-12-15 | 2020-12-15 | 一种基于主控因子修正的参照作物腾发量预测方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011473915.6A CN112507558B (zh) | 2020-12-15 | 2020-12-15 | 一种基于主控因子修正的参照作物腾发量预测方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN112507558A CN112507558A (zh) | 2021-03-16 |
| CN112507558B true CN112507558B (zh) | 2024-02-27 |
Family
ID=74973364
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202011473915.6A Active CN112507558B (zh) | 2020-12-15 | 2020-12-15 | 一种基于主控因子修正的参照作物腾发量预测方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN112507558B (zh) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116681158A (zh) * | 2023-05-16 | 2023-09-01 | 西安理工大学 | 基于集成极限学习机的参考作物蒸散量预测方法 |
| CN117521546A (zh) * | 2023-11-07 | 2024-02-06 | 中国水利水电科学研究院 | 一种考虑时间过程的潜在蒸散发变化归因分析方法 |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106355264A (zh) * | 2016-08-11 | 2017-01-25 | 河海大学 | 一种参考作物腾发量组合预测方法 |
| CN111461909A (zh) * | 2020-04-02 | 2020-07-28 | 中国水利水电科学研究院 | 一种农田蒸散量短期预测方法 |
-
2020
- 2020-12-15 CN CN202011473915.6A patent/CN112507558B/zh active Active
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106355264A (zh) * | 2016-08-11 | 2017-01-25 | 河海大学 | 一种参考作物腾发量组合预测方法 |
| CN111461909A (zh) * | 2020-04-02 | 2020-07-28 | 中国水利水电科学研究院 | 一种农田蒸散量短期预测方法 |
Non-Patent Citations (5)
| Title |
|---|
| 基于FAO-Blaney-Criddle方法的河套灌区参考作物蒸散发量估算;张乐昕;丛振涛;;农业工程学报(第16期);103-109 * |
| 基于天气预报的参照作物腾发量中短期预报模型研究;张倩;段爱旺;王广帅;Jha K.Shiva;申孝军;蔡焕杰;;农业机械学报(第05期);113-120 * |
| 基于天气预报的参考作物腾发量LS-SVM预测模型;张展羽;王声锋;段爱旺;王斌;;水科学进展(第01期);65-70 * |
| 基于日光温室内外温湿度的参考作物腾发量计算;郭彬;刘玉春;蔡伟;张芳园;潘秋艳;淮学民;;灌溉排水学报(第S2期);89-93 * |
| 基于气象预报的参考作物蒸发蒸腾量的神经网络预测模型;徐俊增;彭世彰;张瑞美;李道西;;水利学报(第03期);126-129 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN112507558A (zh) | 2021-03-16 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Gneiting et al. | Calibrated probabilistic forecasting at the stateline wind energy center: The regime-switching space–time method | |
| CN112507558B (zh) | 一种基于主控因子修正的参照作物腾发量预测方法 | |
| CN111461909B (zh) | 一种农田蒸散量短期预测方法 | |
| CN111241473B (zh) | 一种提高区域地下水储量估计精度的方法 | |
| Baldi et al. | Hail occurrence in Italy: Towards a national database and climatology | |
| US20220279741A1 (en) | Decision-making Method for Variable Rate Irrigation Management | |
| CN112558188B (zh) | 一种通过同化闪电资料改进强对流预报的方法 | |
| Mhizha et al. | Use of the FAO AquaCrop model in developing sowing guidelines for rainfed maize in Zimbabwe | |
| CN108983324B (zh) | 一种基于卡尔曼滤波的气温预报方法及系统 | |
| CN112800634A (zh) | 耦合干湿状态辨识与多源信息融合的降水估计方法及系统 | |
| CN115018127B (zh) | 一种植被覆盖度变化归因方法 | |
| CN104143043A (zh) | 一种多功能气候数据模型及其应用 | |
| CN109783774B (zh) | 一种温度集合预报方法及系统 | |
| CN116305933B (zh) | 一种基于dem数据的坡面产汇流简易计算方法及装置 | |
| CN117436619A (zh) | 一种基于等效防洪作用的梯级水库防洪库容联合预留方法 | |
| CN115310796A (zh) | 一种不同类型干旱间传播关系的确定方法及系统 | |
| CN110210142B (zh) | 一种南方大型灌区水稻实时需水量测算方法 | |
| Goler et al. | Influence of climate change on river discharge in Austria | |
| Hu et al. | Rainfall erosivity in the Republic of Korea | |
| CN112597698B (zh) | 一种人类活动对流域蒸散发影响的不确定性评估方法 | |
| CN109543911A (zh) | 一种日照辐射预测方法及系统 | |
| CN117236217A (zh) | 一种基于遥感数据的干旱区植被生态需水量定量估算方法 | |
| CN117494419A (zh) | 一种多模型耦合的流域水土流失遥感监测方法 | |
| CN106355264B (zh) | 一种参考作物腾发量组合预测方法 | |
| CN112598266A (zh) | 基于产量构成要素的冬小麦气象灾害产量损失评估方法 |
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 |