CN115375036A - 遥感、光能利用率模型与气象融合的作物成熟期预测方法 - Google Patents
遥感、光能利用率模型与气象融合的作物成熟期预测方法 Download PDFInfo
- Publication number
- CN115375036A CN115375036A CN202211116114.3A CN202211116114A CN115375036A CN 115375036 A CN115375036 A CN 115375036A CN 202211116114 A CN202211116114 A CN 202211116114A CN 115375036 A CN115375036 A CN 115375036A
- Authority
- CN
- China
- Prior art keywords
- lai
- crop
- model
- data
- maturity
- 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.)
- Withdrawn
Links
- 238000000034 method Methods 0.000 title claims abstract description 75
- 230000004927 fusion Effects 0.000 title claims description 5
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 claims abstract description 56
- 238000011160 research Methods 0.000 claims abstract description 55
- 230000010287 polarization Effects 0.000 claims abstract description 26
- 238000001914 filtration Methods 0.000 claims abstract description 13
- 239000002689 soil Substances 0.000 claims description 44
- 239000011159 matrix material Substances 0.000 claims description 17
- 230000008569 process Effects 0.000 claims description 7
- 238000004088 simulation Methods 0.000 claims description 6
- 230000000694 effects Effects 0.000 claims description 4
- 238000012897 Levenberg–Marquardt algorithm Methods 0.000 claims description 3
- 238000012937 correction Methods 0.000 claims description 3
- 238000000605 extraction Methods 0.000 claims description 3
- 238000005070 sampling Methods 0.000 claims description 3
- 230000035945 sensitivity Effects 0.000 claims description 3
- 101150077190 sinI gene Proteins 0.000 claims description 3
- 238000004422 calculation algorithm Methods 0.000 claims description 2
- 238000005457 optimization Methods 0.000 claims description 2
- 230000007704 transition Effects 0.000 claims description 2
- 238000003672 processing method Methods 0.000 claims 1
- 238000012544 monitoring process Methods 0.000 abstract description 5
- 230000035800 maturation Effects 0.000 abstract description 4
- 238000005094 computer simulation Methods 0.000 abstract 1
- 230000003287 optical effect Effects 0.000 description 5
- 238000012986 modification Methods 0.000 description 4
- 230000004048 modification Effects 0.000 description 4
- 238000009825 accumulation Methods 0.000 description 2
- 238000011161 development Methods 0.000 description 2
- 230000003993 interaction Effects 0.000 description 2
- 238000007781 pre-processing Methods 0.000 description 2
- 239000002028 Biomass Substances 0.000 description 1
- 208000036855 Left sided atrial isomerism Diseases 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 239000003337 fertilizer Substances 0.000 description 1
- 230000002452 interceptive effect Effects 0.000 description 1
- 230000000737 periodic effect Effects 0.000 description 1
- 230000005855 radiation Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q10/00—Administration; Management
- G06Q10/04—Forecasting or optimisation specially adapted for administrative or management purposes, e.g. linear programming or "cutting stock problem"
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
- G01S13/88—Radar or analogous systems specially adapted for specific applications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
- G01S13/88—Radar or analogous systems specially adapted for specific applications
- G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
- G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01W—METEOROLOGY
- G01W1/00—Meteorology
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q50/00—Systems or methods specially adapted for specific business sectors, e.g. utilities or tourism
- G06Q50/02—Agriculture; Fishing; Mining
-
- 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
Abstract
本发明实施例公开了一种遥感、光能利用率模型与气象融合的作物成熟期预测方法,包括建立SAR多极化模式下的水云模型;选择拟合最好的极化模式下的所述水云模型,建立LAI与雷达后向散射系数之间的关系;根据研究区、研究作物历史多年LAI积分曲线确定成熟期阈值;标定LINTUL作物模型;用气象数据和气象预报数据驱动所述LINTUL模型,用经雷达后向散射系数反演的LAI,采用集合卡尔曼滤波法对LINTUL模型进行数据同化,来模拟截至未来时间段逐天LAI;计算截至未来时间的模拟LAI的积分曲线面积比值,判定是否到达成熟期;利用气象监测和预测数据驱动LINTUL动态模拟叶面积指数,通过计算模拟的叶面积指数积分面积比值来判断是否达到成熟期阈值,从而实现对作物成熟期的预测。
Description
技术领域
本发明实施例涉及农业遥感监测技术领域,具体涉及一种遥感、光能利用率模型与气象融合的作物成熟期预测方法。
背景技术
叶面积指数是反演作物物候期的关键,随着卫星遥感技术的发展,基于光学遥感的作物叶面积指数常被用作作物成熟期预测的关键指标,但是光学影像受天气影响大,而作物生长发育期间多云雨,严重阻碍了卫星影像的使用;单纯用气象监测及预测要素(如积温、积雨、太阳辐射)驱动的叶面积指数增长曲线并不能可靠地反映出作物与包括气象要素在内的多种要素的交互关系;利用作物生长模型的成熟期模拟由于模型要素多,误差的逐步积累严重影响模型的准确度,引入数据同化技术可以提高模型准确度,实现区域尺度的模拟,但是目前的数据同化多基于前述光学影像获得的叶面积指数,同样严重受云雨天气制约,虽然作物模型给定了众多的参数,但是模型本身总有不确定性和特定方向上的倾向性,并且其内置的物候模块由温度驱动。
现有的方法中,虽然有通过遥感技术结合气象信息对植物的成熟期预测的方案,但此种方法在实施的过程中,由于作物生长季多云雨,严重阻碍了基于光学遥感影像的植被指数的获取时间和范围;单纯从LAI 变化轨迹考虑的模型如冠层结构动力学模型等,以温度等气象变量为单一驱动,难以反映出作物特性及其与大气、土壤等的复杂交互过程,不利于在精准农业方向的应用;直接利用作物生长模型进行预测的方法,一是需要大量的参数标定,二是其物候模块依赖于积温,并且作物成熟所需的积温是人为给定的,这导致了单纯使用作物模型的物候模块进行预测的局限性;这些都不利于大面积及多种作物的推广应用。
发明内容
为此,本发明实施例提供一种遥感、光能利用率模型与气象融合的作物成熟期预测方法,解决传统利用光学影像同化作物模型时影像获取易受云雨影响、模型对作物生长的地域性考虑欠缺、作物模型内置的物候模块过于简化等问题。
为了实现上述目的,本发明的实施方式提供如下技术方案:
在本发明的实施方式的第一方面中,提供了一种遥感、光能利用率模型与气象融合的作物成熟期预测方法,包含以下步骤:
S1、实测研究区、研究作物的LAI与土壤含水量,获取对应时间、对应位置的合成孔径雷达(SAR)影像,通过水云模型建立卫星VV、VH两种极化模式下各自的后向散射系数与研究区、研究作物的LAI、土壤含水量之间的关系,选择拟合好的极化模式下的水云模型;
S2、通过拟合好的所述水云模型,生成LAI与雷达后向散射系数的查找表,并用支持向量回归方法建立LAI与雷达后向散射系数之间的关系;
S3、获取卫星后向散射系数,反演研究区、研究作物历史三年LAI,根据LAI积分曲线确定成熟期阈值;
S4、根据农业气象站点的气象数据标定作物LINTUL模型;
S5、用气象数据和未来15天气象预报数据驱动所述LINTUL模型,用经雷达后向散射系数反演的LAI,采用集合卡尔曼滤波法对LINTUL模型进行数据同化,来模拟截至未来15天的逐天LAI;
S6、计算截至未来15天的模拟LAI的积分曲线面积比值,判定是否到达成熟期,从而实现提前15天对作物成熟期的预测;
S7、获取新的经过卫星反演的LAI,更新未来15天气象预报数据,重复上述步骤S5到S6。
进一步地,步骤S1中建立水云模型的步骤如下:
提取卫星后向散射系数,根据研究区、研究作物生育期制作卫星影像数据;
在与上述卫星影像对应的样方实地采样,建立水云模型,并对建立的水云模型优化。
进一步地,所述卫星影像的制作方法包括
下载研究区、研究作物生育期内的卫星干涉宽带模式数据,所述数据包括VV和VH两种极化模式;
使用SNAP、ENVI软件对数据进行边缘噪声移除、热噪声移除、斑点滤波、辐射定标、地形校正、地理编码和影像裁剪预处理;
其中,所述预处理过程包括如下步骤
提取雷达亮度值β0,公式如下:
β0 j=20*lg(DNj/A2j) (1)
式中,DNj是第j个像元的灰度值,A2j是系统增益,β0 j指的是像元j的亮度值;
提取雷达后向散射系数σ0,公式如下:
σ0 j=β0 j+10lg(sinIj) (2)
式中,Ij是第j个像元的入射角,单位为(°),σ0 j则为像元j的雷达后向散射系数。
进一步地,所述水云模型的建立方法包括
在研究区均匀布设样方,用GPS记录仪记录样方的边界,测量样方内土壤湿度和作物LAI;
建立水云模型(WCM),公式如下:
σ0=σ0 veg+τ2σ0 soil (3)
σ0 veg=A*V1*cosθ*(1–exp(-2B*V2/cosθ) (4)
τ2=exp(-2B*V2/cosθ) (5)
式(3)–(5)中,σ0 veg为作物后向散射,σ0 soil为土壤后向散射,τ2为衰减因子,V1、V2是对植被冠层结构的定量描述,θ为入射角(°),σ0 soil为土壤后向散射系数(m2/m2)参数A、B由作物种类决定。
进一步地,所述水云模型的优化方法包括
用作物LAI作为植被结构的代表、土壤体积含水量SM作为土壤含水量的代表,公式如下:
V1=LAIE,V2=LAIF (6)
σ0 soil=CMv+D (7)
式中,E、F为拟合参数,与拟合参数A、B一起反映因植被结构相关的散射和衰减作用,参数C可以理解为SAR对土壤湿度的敏感性,D反映土壤粗糙表面造成的后向散射;
由以上公式可以得到下式:
σ0=A*LAIE*cosθ*(1–exp(-2B LAIF/cosθ))+(CMv+D)*exp(-2B LAIF/cos θ) (8)
利用实测数据对水云模型进行参数拟合
利用Levenberg-Marquardt算法对VV和VH极化模式下的上述水云模型分别进行非线性最小二乘拟合,求出两种极化模式下各自的最优拟合参数A、B、C、D、E、F;
选择两种极化模式下决定系数最高的水云模型。
进一步地,步骤S2中建立LAI与雷达后向散射系数之间的关系的步骤如下:
生成查找表:改变拟合好的水云模型中的LAI数值,生成(LAI,σ0)查找表;
支持向量回归:利用上述步骤生成的查找表,用支持向量回归法方法进行回归,公式如下:
f(xi)=wφ(xi)+b (9)
式中,φ(xi)为高维特征空间,w和b为模型参数,利用ε -insensitive loss算法,引入稀疏矩阵ξi,ξi *,i=1,…,n,则支持向量回归问题可转变为:
min(0.5*||w||2)+C∑ni=1(ξi+ξi *) (10)
满足:
yi–f(xi,w)≤ε+ξi * (11)
f(xi,w)-yi≤ε+ξi (12)
ξi,ξ≥0,i=1,2,…,n (13)
将LUT中得出的LAI与σ0数据代入上述公式中,进行支持向量回归,从而得到LAI与σ0之间的关系。
进一步地,步骤S3中成熟期阈值的判定方法如下
根据经SVR后得到的LAI与σ0之间的关系,调取研究区、研究作物生育期历史三年卫星数据,反演LAI时间序列,并计算作物LAI达到最大值时到成熟时LAI积分面积占作物出苗时到LAI达到最大值时的积分面积的比值作为该作物的成熟期阈值,记为Mthreshold,公式如下:
式中,t1为作物LAI达到最大值的时间,t2为作物成熟时间,t0为作物出苗时间。
进一步地,步骤S5中模拟LAI的方法如下:
对作物LINTUL模型LINTUL进行参数标定,针对研究作物进行参数标定,包括作物参数、土壤参数,利用气象数据和未来15天气象预报数据驱动模型;
采用步骤S2中的方法,反演研究区不同时期作物LAI,并用集合卡尔曼滤波方法对LINTUL模型进行数据同化,实现逐天、多尺度的LAI预测, EnKF公式为:
Aa=Af+K*(Dt-HA) (15)
K=Ac*HT*(HAc*HT+Dc)-1,K∈Rn*N (16)
其中,Aa为分析矩阵,Af是预报矩阵,K是集合卡尔曼增益系数,Dt是观测变量矩阵,Dc是观测变量的协方差矩阵,指的是由卫星水云模型得到的LAI观测值,H是非线性算子,A来自于预测方程,指LINTUL模型对LAI的模拟,HA为观测变量的误差协方差矩阵。
进一步地,步骤S6中,成熟期判定的方法如下:
计算作物出现LAI最大值到预测期,模型模拟出的LAI值的积分值占比,公式如下:
式中tp为预测截时间点,作物成熟的判定条件为:LAIr>=Mthreshold。
进一步地,步骤S7中,成熟期预测动态更新的方法如下:
逐日更新气象数据和气象预测数据,当新的卫星数据可获得时,同时更新根据卫星数据反演的LAI,对LINTUL模型进行同化更新,进而更新模拟LAIr,从而实现作物成熟期提前15天的预测。
与现有技术相比,根据本发明的实施方式,该方法具有如下优点:
1、该方法包括建立合成孔径雷达多极化模式下的水云模型;选择拟合最好的极化模式下的所述水云模型,建立LAI与雷达后向散射系数之间的关系;根据研究区、研究作物历史多年LAI积分曲线确定成熟期阈值;标定LINTUL作物模型;用气象数据和气象预报数据驱动所述LINTUL 模型,用经雷达后向散射系数反演的LAI,采用集合卡尔曼滤波法对LINTUL模型进行数据同化,来模拟截至未来时间段逐天LAI;计算截至未来时间的模拟LAI的积分曲线面积比值,判定是否到达成熟期;利用气象监测和预测数据驱动基于光能利用率的LINTUL作物生长模型动态模拟叶面积指数,通过计算模拟的叶面积指数积分面积比值来判断是否达到成熟期阈值,从而实现对作物成熟期的预测;
2、本方法用的合成孔径雷达数据(不受云雨影响),是通过作物模型模拟LAI,通过LAI积分面积比值,结合集合卡尔曼滤波将雷达LAI 与作物模型结合的系统化方法综合判断成熟与否,相对于现有的光学遥感(极易受云雨影响)反演LAI、通过作物模型的物候模块以及通过作物模型直接模拟和仅利用单一变量进行预测等方法的实施,能够充分考虑作物与环境要素的交互关系易于实施遥感同化,有更精确的时间和周期性的动态预测过程,易于推广应用;
3、本方法综合了光能利用率模型对作物干物质积累模拟、叶面积指数物候判定、SAR的全天候、全天时观测的优势,以研究区、研究作物成熟期叶面积指数曲线积分面积比值的多年平均作为阈值,利用基于光能利用率的LINTUL作物生长模型动态模拟叶面积指数,利用SAR卫星后向散射系数反演的叶面积指数作为同化模型的观测值,利用气象监测和未来15天预测数据驱动模型模拟叶面积指数,通过计算模拟的叶面积指数积分面积比值来判断是否达到阈值,从而实现对作物成熟期的预测。
附图说明
为了更清楚地说明本发明的实施方式或现有技术中的技术方案,下面将对实施方式或现有技术描述中所需要使用的附图作简单地介绍。显而易见地,下面描述中的附图仅仅是示例性的,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据提供的附图引申获得其它的实施附图。
本说明书所绘示的结构、比例、大小等,均仅用以配合说明书所揭示的内容,以供熟悉此技术的人士了解与阅读,并非用以限定本发明可实施的限定条件,故不具技术上的实质意义,任何结构的修饰、比例关系的改变或大小的调整,在不影响本发明所能产生的功效及所能达成的目的下,均应仍落在本发明所揭示的技术内容得能涵盖的范围内。
图1为本发明实施例提供的遥感、光能利用率模型与气象融合的作物成熟期预测方法的流程图;
图2为本发明实施例提供的预测方法中通过SAR数据反演LAI并提取成熟期阈值方法的流程图;
图3为本发明实施例提供的利用SAR同化作物生长模型动态模拟叶面积指数积分曲线的成熟期预测流程。
具体实施方式
以下由特定的具体实施例说明本发明的实施方式,熟悉此技术的人士可由本说明书所揭露的内容轻易地了解本发明的其他优点及功效,显然,所描述的实施例是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
本说明书中所引用的如“上”、“下”、“左”、“右”、“中间”等的用语,亦仅为便于叙述的明了,而非用以限定本发明可实施的范围,其相对关系的改变或调整,在无实质变更技术内容下,当亦视为本发明可实施的范畴。
如图1所示,其示出了本发明实施例提供的遥感、光能利用率模型与气象融合的作物成熟期预测方法,包含以下步骤:
S1、实测研究区、研究作物的LAI与土壤含水量,获取对应时间、对应位置的合成孔径雷达影像,通过水云模型建立卫星VV (vertical-vertical polarization)、VH(vertical-horizontal)两种极化模式下各自的后向散射系数与研究区、研究作物的叶面积指数 (LAI)、土壤含水量(SM)之间的关系,选择拟合好的极化模式下的水云模型(WCM);本实施例中的合成孔径雷达卫星采用哨兵一号卫星,其中,建立水云模型的步骤如下:
1.1提取哨兵一号后向散射系数,包括
(1)下载研究区作物生育期内的哨兵一号干涉宽带模式 (Interferometric WideSwath Mode,IW)数据,包括VV (vertical-vertical polarization)和VH(vertical-horizontal) 两种极化模式;
(2)使用SNAP、ENVI软件对数据进行边缘噪声移除、热噪声移除、斑点滤波、辐射定标、地形校正、地理编码和影像裁剪等预处理;
(3)提取雷达亮度值β0,公式如下:
β0 j=20*lg(DNj/A2j) (1)
式中,DNj是第J个像元的灰度值,A2j是系统增益,β0 j指的是像元j的灰度值;
(4)提取雷达后向散射系数σ0,公式如下:
σ0 j=β0 j+10lg(sinIj) (2)
式中,Ij是第j个像元的入射角,单位为(°),σ0 j则为像元j的雷达后向散射系数。
1.2在与上述哨兵一号影像对应的样方实地采样,建立水云模型,步骤如下,
(1)在研究区均匀布设样方,用GPS记录仪记录样方的边界,测量样方内土壤湿度(SM)和
作物叶面积指数(LAI);
(2)建立水云模型(WCM),公式如下:
σ0=σ0 veg+τ2σ0 soil (3)
σ0 veg=A*V1*cosθ*(1–exp(-2B*V2/cosθ) (4)
τ2=exp(-2B*V2/cosθ) (5)
式(3)–(5)中,σ0 veg为作物后向散射,σ0 soil为土壤后向散射,τ2为衰减因子,V1、V2是对植被冠层结构的定量描述,θ为入射角(°),σ0 soil为土壤后向散射系数(m2/m2)参数A、B由作物种类决定。用作物叶面积指数LAI作为植被结构的代表、土壤体积含水量SM作为土壤含水量的代表,公式如下:
V1=LAIE,V2=LAIF (6)
σ0 soil=CMv+D (7)
式中,E、F为拟合参数,与拟合参数A、B一起反映因植被结构相关的散射和衰减作用,参数C可以理解为SAR对土壤湿度的敏感性,D反映土壤粗糙表面造成的后向散射。
由以上公式可以得到下式:
σ0=A*LAIE*cosθ*(1–exp(-2B LAIF/cosθ))+(CMv+D)*exp(-2B LAIF/cos θ) (8)
(3)利用实测数据对水云模型进行拟合
利用Levenberg-Marquardt算法对VV和VH极化模式下的上述水云模型分别进行非线性最小二乘拟合,求出两种极化模式下最优拟合参数 A、B、C、D、E、F;
选择两种极化模式下决定系数最高的水云模型。
S2、通过拟合好的所述水云模型,生成LAI与雷达后向散射系数的查找表,并用支持向量回归方法建立LAI与雷达后向散射系数之间的关系;
通过上述拟合好的水云模型,生成叶面积指数与雷达后向散射系数与对应的叶面积指数的查找表(LUT,look up table),并用支持向量回归(SVR)方法建立叶面积指数与雷达后向散射系数之间的关系。
S3、获取哨兵一号卫星后向散射系数,反演研究区、研究作物历史三年叶面积指数,根据叶面积指数积分曲线确定成熟期阈值;
如图2所示,成熟期阈值的判定方法如下:
根据经SVR后得到的LAI与σ0之间的关系,调取研究区、研究作物生育期历史三年哨兵一号卫星数据,反演LAI时间序列,并计算作物LAI 达到最大值时到成熟时LAI积分面积占作物出苗时到LAI达到最大值时的积分面积的比值作为该作物的成熟期阈值,记为Mthreshold,公式如下:
式中,t1为作物LAI达到最大值的时间,t2为作物成熟时间,t0为作物出苗时间。
S4、根据农业气象站点的数据标定作物光能利用率模型(LINTUL),其中的数据信息包括:气象数据、施用水、肥记录数据等,标定方法为常规的方法,一般作物生长模型如WOFOST、光能利用率模型LINTUL,标定一般指该作物从生长到成熟所需要的积温、生物量在各部位的积累比例、作物生长过程中施水、肥的量和施用时间,所处土壤的性质等,可以根据农气站点记录或相关地区文献记录等进行标定;
S5、用气象数据和未来15天气象预报数据驱动光能利用率模型,用经雷达后向散射系数反演的叶面积指数,采用集合卡尔曼滤波法对 LINTUL模型进行数据同化,来模拟截至未来15天的逐天的叶面积指数;
其中模拟叶面积指数的方法如下:
5.1对作物光能利用率模型LINTUL进行参数标定,针对研究作物进行参数标定,包括作物参数、土壤参数等,利用气象数据和未来15天气象预报数据驱动模型;
5.2采用步骤S2中的方法,反演研究区不同时期作物LAI,并用集合卡尔曼滤波方法对LINTUL模型进行数据同化,实现逐天、多尺度的LAI预测,EnKF公式为:
Aa=Af+K*(Dt-HA) (15)
K=Ac*HT*(HAc*HT+Dc)-1,K∈Rn*N (16)
其中,Aa为分析矩阵,Af是预报矩阵,K是集合卡尔曼增益系数, Dt是观测变量矩阵,Dc是观测变量的协方差矩阵,本方法中指的是由卫星水云模型得到的LAI观测值,H是非线性算子,A来自于预测方程,指 LINTUL模型对LAI的模拟,HA为观测变量的误差协方差矩阵。
S6、计算截至未来15天的模拟LAI的积分曲线面积比值,判定是否到达成熟期,从而实现提前15天对作物成熟期的预测;
其中,成熟期判定的方法如下:
计算作物出现LAI最大值到预测期,模型模拟出的LAI值的积分值占比,公式如下:
式中tp为预测截止时间点,作物成熟的判定条件为:LAIr>=Mthreshold。
S7、获取新的经过卫星反演的LAI,更新未来15天气象预报数据,重复上述步骤S5到S6。
成熟期预测动态更新的方法如下:
逐日更新气象数据和气象预测数据,当新的SAR卫星数据可获得时,同时更新根据SAR卫星数据反演的LAI,对LINTUL模型进行同化更新,进而更新模拟LAIr,从而实现作物成熟期提前15天的预测。
具体的实施例中,如图3所示,获取研究区、研究作物的合成孔径雷达(SAR)卫星影像,进行预处理,提取哨兵1后向散射系数,反演研究区、研究作物历史三年叶面积指数,根据叶面积指数积分曲线确定成熟期阈值,根据获取到的气象信息以及研究区、研究作物的数据信息标定作物光能利用率模型,反演研究区不同时期作物LAI,并用集合卡尔曼滤波方法对LINTUL模型进行数据同化,实现逐天、多尺度的LAI预测,计算作物出现LAI最大值到预测期,模型模拟出的LAI值的积分值在整个LAI积分曲线的占比,判定是否到达成熟期,从而实现作物成熟期预测,获取新的经过SAR卫星反演的叶面积指数,对LINTUL模型进行同化更新,进而更新模拟,实现动态循环对作物成熟期的预测。
虽然,上文中已经用一般性说明及具体实施例对本发明作了详尽的描述,但在本发明基础上,可以对之作一些修改或改进,这对本领域技术人员而言是显而易见的。因此,在不偏离本发明精神的基础上所做的这些修改或改进,均属于本发明要求保护的范围。
Claims (10)
1.一种遥感、光能利用率模型与气象融合的作物成熟期预测方法,其特征在于,包含以下步骤:
S1、实测研究区、研究作物的LAI与土壤含水量,获取对应时间、对应位置的合成孔径雷达影像,通过水云模型建立卫星VV、VH两种极化模式下各自的后向散射系数与研究区、研究作物的LAI、土壤含水量之间的关系,选择拟合好的极化模式下的水云模型;
S2、通过拟合好的所述水云模型,生成LAI与雷达后向散射系数的查找表,并用支持向量回归方法建立LAI与雷达后向散射系数之间的关系;
S3、获取卫星后向散射系数,反演研究区、研究作物历史三年LAI,根据LAI积分曲线确定成熟期阈值;
S4、根据农业气象站点的数据标定作物LINTUL模型;
S5、用气象数据和未来15天气象预报数据驱动所述LINTUL模型,用经雷达后向散射系数反演的LAI,采用集合卡尔曼滤波法对LINTUL模型进行数据同化,来模拟截至未来15天的逐天LAI;
S6、计算截至未来15天的模拟LAI的积分曲线面积比值,判定是否到达成熟期,从而实现提前15天对作物成熟期的预测;
S7、获取新的经过卫星反演的LAI,更新未来15天气象预报数据,重复上述步骤S5到S6。
2.如权利要求1所述的遥感、LINTUL模型与气象融合的作物成熟期预测方法,其特征在于,步骤S1中建立水云模型的步骤如下:
提取卫星后向散射系数,根据研究区、研究作物生育期准备卫星影像数据;
在与上述卫星影像对应的样方实地采样,建立水云模型,并对建立的水云模型优化。
3.如权利要求2所述的遥感、光能利用率模型与气象融合的作物成熟期预测方法,其特征在于,所述卫星影像的处理方法包括
下载研究区、研究作物生育期内的卫星干涉宽带模式数据,所述数据包括VV和VH两种极化模式;
使用SNAP、ENVI软件对数据进行边缘噪声移除、热噪声移除、斑点滤波、辐射定标、地形校正、地理编码和影像裁剪预处理;
其中,所述预处理过程包括如下步骤
提取雷达亮度值β0,公式如下:
β0 j=20*lg(DNj/A2j) (1)
式中,DNj是第j个像元的灰度值,A2j是系统增益,β0 j指的是像元j的亮度值;
提取雷达后向散射系数σ0,公式如下:
σ0 j=β0 j+10lg(sinIj) (2)
式中,Ij是第j个像元的入射角,单位为(°),σ0 j则为像元j的雷达后向散射系数。
4.如权利要求3所述的遥感、光能利用率模型与气象融合的作物成熟期预测方法,其特征在于,所述水云模型的建立方法包括
在研究区均匀布设样方,用GPS记录仪记录样方的边界,测量样方内土壤湿度和作物LAI;
建立水云模型(WCM),公式如下:
σ0=σ0 veg+τ2σ0 soil (3)
σ0 veg=A*V1*cosθ*(1–exp(-2B*V2/cosθ) (4)
τ2=exp(-2B*V2/cosθ) (5)
式(3)–(5)中,σ0 veg为作物后向散射,σ0 soil为土壤后向散射,τ2为衰减因子,V1、V2是对植被冠层结构的定量描述,θ为入射角,σ0 soil为土壤后向散射系数(m2/m2)参数A、B由作物种类决定。
5.如权利要求4所述的遥感、光能利用率模型与气象融合的作物成熟期预测方法,其特征在于,所述水云模型的优化方法包括
用作物LAI作为植被结构的代表、土壤体积含水量SM作为土壤含水量的代表,公式如下:
V1=LAIE,V2=LAIF (6)
σ0 soil=CMv+D (7)
式中,E、F为拟合参数,与拟合参数A、B一起反映因植被结构相关的散射和衰减作用,参数C可以理解为SAR对土壤湿度的敏感性,D反映土壤粗糙表面造成的后向散射;
由以上公式可以得到下式:
σ0=A*LAIE*cosθ*(1–exp(-2B LAIF/cosθ))+(CMv+D)*exp(-2B LAIF/cosθ) (8)
利用实测数据对水云模型进行参数拟合
利用Levenberg-Marquardt算法对VV和VH极化模式下的上述水云模型分别进行非线性最小二乘拟合,求出两种极化模式下各自的最优拟合参数A、B、C、D、E、F;
选择两种极化模式下决定系数最高的水云模型。
6.如权利要求2所述的遥感、光能利用率模型与气象融合的作物成熟期预测方法,其特征在于,步骤S2中建立LAI与雷达后向散射系数之间的关系的步骤如下:
生成查找表:改变拟合好的水云模型中的LAI数值,生成(LAI,σ0)查找表;
支持向量回归:利用上述步骤生成的查找表,用支持向量回归法方法进行回归,公式如下:
f(xi)=wφ(xi)+b (9)
式中,φ(xi)为高维特征空间,w和b为模型参数,利用ε-insensitive loss算法,引入稀疏矩阵ξi,ξi *,i=1,…,n,则支持向量回归问题可转变为:
min(0.5*||w||2)+C∑ni=1(ξi+ξi *) (10)
满足:
yi–f(xi,w)≤ε+ξi * (11)
f(xi,w)-yi≤ε+ξi (12)
ξi,ξ≥0,i=1,2,…,n (13)
将LUT中得出的LAI与σ0数据代入上述公式中,进行支持向量回归,从而得到LAI与σ0之间的关系。
8.如权利要求1所述的遥感、光能利用率模型与气象融合的作物成熟期预测方法,其特征在于,步骤S5中模拟LAI的方法如下:
对作物LINTUL模型LINTUL进行参数标定,针对研究作物进行参数标定,包括作物参数、土壤参数,利用气象数据和未来15天气象预报数据驱动模型;
采用步骤S2中的方法,反演研究区不同时期作物LAI,并用集合卡尔曼滤波方法对LINTUL模型进行数据同化,实现逐天、多尺度的LAI预测,EnKF公式为:
Aa=Af+K*(Dt-HA) (15)
K=Ac*HT*(HAc*HT+Dc)-1,K∈Rn*N (16)
其中,Aa为分析矩阵,Af是预报矩阵,K是集合卡尔曼增益系数,Dt是观测变量矩阵,Dc是观测变量的协方差矩阵,指的是由卫星水云模型得到的LAI观测值,H是非线性算子,A来自于预测方程,指LINTUL模型对LAI的模拟,HA为观测变量的误差协方差矩阵。
10.如权利要求1所述的遥感、光能利用率模型与气象融合的作物成熟期预测方法,其特征在于,步骤S7中,成熟期预测动态更新的方法如下:
逐日更新气象数据和气象预测数据,当新的卫星数据可获得时,同时更新根据卫星数据反演的LAI,对LINTUL模型进行同化更新,进而更新模拟LAIr,从而实现作物成熟期提前15天的预测。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202211116114.3A CN115375036A (zh) | 2022-09-14 | 2022-09-14 | 遥感、光能利用率模型与气象融合的作物成熟期预测方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202211116114.3A CN115375036A (zh) | 2022-09-14 | 2022-09-14 | 遥感、光能利用率模型与气象融合的作物成熟期预测方法 |
Publications (1)
Publication Number | Publication Date |
---|---|
CN115375036A true CN115375036A (zh) | 2022-11-22 |
Family
ID=84071501
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202211116114.3A Withdrawn CN115375036A (zh) | 2022-09-14 | 2022-09-14 | 遥感、光能利用率模型与气象融合的作物成熟期预测方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN115375036A (zh) |
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116071644A (zh) * | 2022-12-20 | 2023-05-05 | 中化现代农业有限公司 | 逐日叶面积指数数据反演方法、装置、设备及存储介质 |
CN116305875A (zh) * | 2023-02-17 | 2023-06-23 | 中国科学院地理科学与资源研究所 | 数值天气预报模式中叶面积指数的数据处理方法及装置 |
-
2022
- 2022-09-14 CN CN202211116114.3A patent/CN115375036A/zh not_active Withdrawn
Cited By (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116071644A (zh) * | 2022-12-20 | 2023-05-05 | 中化现代农业有限公司 | 逐日叶面积指数数据反演方法、装置、设备及存储介质 |
CN116071644B (zh) * | 2022-12-20 | 2023-08-08 | 中化现代农业有限公司 | 逐日叶面积指数数据反演方法、装置、设备及存储介质 |
CN116305875A (zh) * | 2023-02-17 | 2023-06-23 | 中国科学院地理科学与资源研究所 | 数值天气预报模式中叶面积指数的数据处理方法及装置 |
CN116305875B (zh) * | 2023-02-17 | 2023-08-29 | 中国科学院地理科学与资源研究所 | 数值天气预报模式中叶面积指数的数据处理方法及装置 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Monteiro et al. | Assessment of NASA/POWER satellite‐based weather system for Brazilian conditions and its impact on sugarcane yield simulation | |
Skakun et al. | Early season large-area winter crop mapping using MODIS NDVI data, growing degree days information and a Gaussian mixture model | |
Liang et al. | A long-term Global LAnd Surface Satellite (GLASS) data-set for environmental studies | |
Reichle et al. | Connecting satellite observations with water cycle variables through land data assimilation: Examples using the NASA GEOS-5 LDAS | |
CN115375036A (zh) | 遥感、光能利用率模型与气象融合的作物成熟期预测方法 | |
Tian et al. | Estimating zero-plane displacement height and aerodynamic roughness length using synthesis of LiDAR and SPOT-5 data | |
Holmes et al. | Microwave implementation of two-source energy balance approach for estimating evapotranspiration | |
Houborg et al. | Adapting a regularized canopy reflectance model (REGFLEC) for the retrieval challenges of dryland agricultural systems | |
Yuan et al. | Vegetation-specific model parameters are not required for estimating gross primary production | |
WO2018107245A1 (en) | Detection of environmental conditions | |
Dehkordi et al. | Yield gap analysis using remote sensing and modelling approaches: Wheat in the northwest of Iran | |
CN116450700B (zh) | 极轨卫星地表温度时间归一化方法、装置及电子设备 | |
Ma et al. | Continuous evaluation of the spatial representativeness of land surface temperature validation sites | |
Hashimoto et al. | High‐resolution mapping of daily climate variables by aggregating multiple spatial data sets with the random forest algorithm over the conterminous United States | |
Chen et al. | Potential of remote sensing data-crop model assimilation and seasonal weather forecasts for early-season crop yield forecasting over a large area | |
Raoult et al. | Evaluating and optimizing surface soil moisture drydowns in the ORCHIDEE land surface model at in situ locations | |
CN112785035A (zh) | 融合多元信息的中短期水文预报方法及系统 | |
Chen et al. | Importance of shaded leaf contribution to the total GPP of Canadian terrestrial ecosystems: evaluation of MODIS GPP | |
Pereira et al. | Solar irradiance modelling using an offline coupling procedure for the Weather Research and Forecasting (WRF) model | |
CN110008621B (zh) | 基于双重流依赖集合平方根滤波同化算法的作物模型遥感同化估产方法 | |
Khesali et al. | A method in near-surface estimation of air temperature (NEAT) in times following the satellite passing time using MODIS images | |
Huang et al. | Topographic effects on estimating net primary productivity of green coniferous forest in complex terrain using Landsat data: a case study of Yoshino Mountain, Japan | |
CN115759524B (zh) | 一种基于遥感影像植被指数的土壤生产力等级识别方法 | |
Kjaersgaard et al. | Comparison of the performance of net radiation calculation models | |
CN114819737B (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 | ||
WW01 | Invention patent application withdrawn after publication |
Application publication date: 20221122 |
|
WW01 | Invention patent application withdrawn after publication |