CN101595812A - 一种作物收获指数的获取方法 - Google Patents

Abstract

本发明提供了一种作物收获指数的获取方法，包括：获取地面不同观测点的历史实测作物收获指数；获取所述不同观测点不同时间的时序归一化植被指数NDVI数据；根据各个观测点不同时间的NDVI数据结合作物生长过程曲线构建与作物收获指数概念相关参数HI_{NDVI_SUM}；构建HI_{NDVI_SUM}与历史实测作物收获指数之间的函数关系；根据所构建的函数关系和所要计算年份的观测点不同时间的NDVI数据计算该观测点作物收获指数。使用本发明，现利用遥感时序植被指数更加准确地获取空间作物收获指数。

Description

一种作物收获指数的获取方法

技术领域

本发明提供了一种作物收获指数获取方法，特别是指一种基于时序归一化植被指数的作物收获指数遥感获取方法。

背景技术

农作物收获指数作为影响作物单产的重要生物学参数之一，早已引起人们的重视。对粮食作物来说，作物的收获指数即为农作物籽粒产量占农作物地上生物量的百分数。一般而言，正常生长条件下的粮食作物收获指数与作物单产呈正相关关系，同时，作物收获指数与作物光合产物运转、分配及器官的发育建成有密切关系。因此，收获指数是长期以来农学家及育种专家提高作物单产、选育作物新品种和品种改良过程中所需考虑的最重要因素之一。众多研究表明，进几十年来，稻麦等作物单产的不断提高，作物收获指数的不断提高是其中重要原因之一。另外，随着作物生长机理模型的出现和不断发展，收获指数成为作物生长模型通过修正作物地上部生物量得到作物单产所必须输入参数之一。

目前，区域作物收获指数的确定方法主要包括以点代面法、空间内插法或者通过文献查询获取区域收获指数。以点代面法指通过多年定点试验进行作物收获指数测定，然后利用该点多年均值作为该区域作物收获指数的常数。空间内插法主要指通过实际调查获取区域多点作物收获指数，进而通过空间内插的方法得到当年区域空间收获指数分布状况。以上针对作物收获指数的计算或研究的方法，大多基于农学试验在田块尺度进行测量和研究，研究内容主要涉及作物收获指数的数学模拟或作物收获指数对作物生长环境及其管理措施的响应，而对于大范围空间上作物收获指数的信息获取研究，国内外鲜有报道。作物收获指数受育种水平、作物品种、田间管理水平、外界胁迫条件(高温、缺水等)及气候条件影响，虽然在一定时期一定区域具有相对稳定性，但同一作物的不同品种、不同管理水平、不同胁迫条件下均导致收获指数在一定时间内小区域范围内存在较大的空间变异。特别是在我国以农户为基本农业生产单位的特定条件下，作物品种空间分布差异大，随机性强，便增加了作物收获指数在空间上的变异性。

经对现有技术文件的检索发现，在《Agricultural WaterManagement》2003年第58卷的第145-157页，Samarasinghe G.B编写的《Growth and yields of Sri Lanka’s major crop interpreted frompublic domain satellites》一文中，利用统计部门区域平均作物单产数据与遥感数据获取的作物生物量数据获取了区域平均水稻收获指数，并将该参数作为该区作物收获指数常数，这比上述作物收获指数的确定方法更为合理，但该方法不能获取收获指数的空间分布信息，作物收获指数的空间变异性仍无法解决。虽然，在《EuropeanJournal of Agronomy》2007年第2期的第266-274页，Moriondo M，Maselli F和Bindi M编写的《A simple model of regional wheat yieldbased on NDVI data》，提出运用NDVI指数(归一化植被指数，是一个能够通过遥感获取且能够直接、有效地反映作物绿度长势、盖度、生物量和作物单产的最佳参考量之一)进行区域范围提取小麦收获指数的方法，但该方法需要首先确定一个地区的最大收获指数和收获指数可能变幅，上述两个参数在大范围区域内较难准确获取，而且具有一定的不确定性，因此，使得结果具有一定不确定性。

发明内容

有鉴于此，本发明的主要目的在于提供一种作物收获指数的获取方法，以实现利用遥感时序植被指数更加准确地获取空间作物收获指数。

本发明提供的一种作物收获指数的获取方法，包括步骤：

A、获取地面不同观测点的历史实测作物收获指数；

B、获取所述不同观测点不同时间的时序归一化植被指数NDVI数据；

C、根据各个观测点不同时间的NDVI数据结合作物生长过程曲线构建与作物收获指数概念相关参数HI_{NDVI_SUM}；

D、构建HI_{NDVI_SUM}与历史实测作物收获指数之间的函数关系；

E、根据所构建的函数关系和所要计算年份的观测点不同时间的NDVI数据计算该观测点作物收获指数。

由此可知，通过利用遥感技术获取作物NDVI数据，根据不同时间的NDVI数据结合作物生长过程曲线构建与作物收获指数概念相关参数HI_{NDVI_SUM}构建HI_{NDVI_SUM}与历史实测作物收获指数之间的函数关系，以计算出作物收获指数的方法，可以更加有效地、准确地获取空间作物收获指数。

上述的方法中，其特征在于，所述步骤A包括子步骤：

A1、在每个观测点内选取至少3个以上的取样点采用作物实割实测获得所述取样点的作物收获指数；

A2、将观测点内各取样点的作物收获指数进行均值处理作为该观测点作物的实测收获指数。

由此可知，在每个观测点内选取至少3个以上的取样点采用作物实割实测获得所述取样点的作物收获指数，并对其进行均值处理作为该观测点作物的实测收获指数，采用最少的工作量获取最接近真实的实测收获指数。

上述的方法中，其特征在于，所述步骤B包括子步骤：

B1、获取观测点的MODIS数据；

B2、根据MODIS数据中的近红外波段反射率和红光波段反射率计算生成观测点的NDVI日数据，其中采用如下计算公式：

$\mathrm{NDVI}=\frac{R_n-R_r}{R_n+R_r}$

其中，R_n为近红外波段的反射率；R_r为红光波段的反射率。

由此可知，MODIS数据因其光谱分辨率高、观测周期短等有点，可以为及时、经济地获取NDVI提供了便捷。

上述的方法中，其特征在于，步骤B2后进一步包括：

B3、采用最大值合成法MVC将NDVI日数据合成NDVI旬数据。

由此可知，通过最大值合成法MVC将NDVI日数据合成NDVI旬数据，可以进一步消除云、大气、太阳高度角等的部分干扰。

上述的方法中，其特征在于，步骤B3后进一步包括：

采用Savitzky-Golay滤波平滑方法对NDVI旬数据时间序列进行平滑去噪处理，所述Savitzky-Golay滤波平滑方法采用如下公式：

$Y_j=\frac{\underset{i=-m}{\overset{i=m}{\mathrm{\Σ}}}{C}_i{X}_{j+i}}{N}$

其中，X是未平滑的NDVI值，C_i是第i个NDVI值的权重系数，N是滤波卷积算子的NDVI值个数，其值大小等于平滑窗口的大小(2m+1)，j是在未平滑数据表中的序数，m值等于平滑窗口宽度的一半。

由此可知，采用Savitzky-Golay滤波平滑方法对NDVI旬数据时间序列进行平滑去噪处理，可以有效去除多时相NDVI遥感数据受云、气溶胶影响造成的噪音，得到质量较高的MODIS-NDVI旬时序数据。

上述的方法中，其特征在于，所述步骤C包括：

利用作物开花前曲线特征和开花后曲线特征，采用如下公式构建与作物收获指数概念相关参数HI_{NDVI_SUM}：

${\mathrm{HI}}_{\mathrm{NDVI}\_\mathrm{SUM}}=\frac{\mathrm{\Σ}{\mathrm{NDVI}}_{\mathrm{post}}}{\mathrm{\Σ}{\mathrm{NDVI}}_{\mathrm{pre}}}$

其中，∑NDVI_post为作物开花后至乳熟初期NDVI累积值；∑NDVI_pre为作物出苗至开花前NDVI累积值。

由此可知，利用作物开花前曲线特征和开花后曲线特征获取作物收获指数相关参数HI_{NDVI_SUM}，避免了受与粮食作物籽粒干物质积累过程和作物茎、叶等干物质积累过程关系较少的作物生长阶段的影响，提高了参数HI_{NDVI_SUM}计算的准确性。

上述的方法中，其特征在于，所述作物为越冬作物时，所述∑NDVI_pre为作物返青至开花前NDVI累积值。

由此可知，对于越冬作物而言，返青至开花前NDVI累积值与越冬作物地上生物量相关性较好。

上述的方法中，其特征在于，步骤D所述构建的方法包括：

采用数据的直线拟合方法构建。

由此可知，数据的直线拟合方法具有重复精度高、稳定性好的优点，能有效抑制噪音的影响。

附图说明

图1为本发明作物收获指数的获取方法的流程图；

图2为本发明的实施区域位置和地面被观测点示意图；

图3为本发明的NDVI旬数据平滑滤波图；

图4为本发明的作物NDVI旬数据-时间曲线示意图；

图5为本发明的参数与作物收获指数关系图；

图6为本发明的2008年研究区域作物收获指数空间分布图；

图7为本发明的作物收获指数反演精度验证图。

具体实施方式

作物收获指数是指农作物籽粒产量占农作物地上生物量的百分数，从作物生理机制角度看，收获指数即为碳素从源分到籽粒库的比例，因此，作物抽穗开花前作物累积生物量和抽穗开花后光合产物向穗部转移水平决定了作物收获指数的大小。

由上，作物开花前阶段营养生长主要与作物地上生物量积累有密切关系，而开花后阶段的生殖生长主要与作物籽粒产量有密切关系。又由于作物生长关键生育期NDVI与作物的绿度、生物量和作物单产等均具有较好的关系(下文论证)，故，本发明通过分析冬小麦开花前和开花后的NDVI数据累积值与冬小麦地上生物量、冬小麦单产之间关系提供了一种作物收获指数获取的方法。

其中，本发明进行作物收获指数信息提取时，将作物的生长阶段划为开花前和开花后两个阶段考虑。其中，开花前阶段主要指作物出苗至开花前。开花后阶段主要指开花期至乳熟初期，这是因为乳熟中后期开始作物叶绿素迅速减少，光合作用变弱或停止，作物开始呈现黄色，通过遥感技术获取的反映作物绿度的植被指数已经不能敏感反映作物生长状况。

下面对本发明作物收获指数获取的方法举例进行详细介绍。

本发明实施例的实施区域为(E115.19°～116.53°，N37.09°～38.36°)区域，位于中国北方粮食生产基地黄淮平原区内的河北省衡水市11个县(市)，覆盖面积为8815km²。如图2示出的实施区域位置和地面被观测点示意图，该区域属于温带半湿润季风气候，大于0℃积温4200～5500℃，年累积辐射量约为5.0×10⁶～5.2×10⁶kJ/m²，无霜期为170～220天，年降水量平均为500～600mm，降水主要集中在夏季的7～9月，该区主要粮食作物为冬小麦、夏玉米，一年两熟轮作制度。其中，冬小麦种植时间为上年的9月底～10月初，返青时间为下一年3月上中旬，收获期为当年的6月上中旬。本实施例以冬小麦收获指数为例，2004年、2007年和2008年地面实测冬小麦收获指数观测点共117个，其中2004年观测点为29个，2007年的观测点为42个，2008年的观测点为46个。2004和2007两年的数据用来建立模型进行收获指数的提取，2008年的数据主要是用来验证本发明的方法的实施及精度。

其中，为了使得本发明具有普通作物的适用性，选取观测点时不仅考虑了冬小麦长势和产量的代表性，同时考虑观测点在实施区中分布的均匀性。

具体为，每个观测点面积不小于500m×500m，观测点内种植结构较为单一，观测点的定位采用差分GPS进行精确定位，并且为减小误差，对地面被观测点做500m缓冲区，在对相应遥感参数进行提取时，得到与地面被观测点相对应的500m范围内遥感参数并进行均值处理。

当对作物收获指数进行提取时，参见图1示出的本发明作物收获指数提取方法的流程图，包括以下步骤：

步骤101：针对所选取的各个观测点，获取2004、2007和2008年的地面冬小麦的实测收获指数。

每个观测点内冬小麦收获指数实测点取样面积为1m²，且实测样点不少于3个。各个取样点收获指数获取采用冬小麦实割实测获得，并将观测点内各样点的收获指数进行均值处理作为该观测点冬小麦的实测收获指数。

步骤102：对应所述各个观测点，利用遥感技术，获取各个观测点的不同时间的NDVI数据。

本实施例中，即获取2004年和2007年实施区11个县(市)冬小麦地面被观测点的NDVI日数据，并据此生成NDVI旬数据(即10日数据)。具体方法为：

针对被观测点，通过遥感技术获取其地面分辨率为250米的MODIS数据(MODIS数据即卫星遥感数据，由terra和aqua卫星向全世界免费实时广播发送)，根据MODIS数据中的近红外波段反射率和红光波段反射率计算生成被观测点的NDVI指数。其计算方法为：

$\mathrm{NDVI}=\frac{R_n-R_r}{R_n+R_r}$

其中，R_n为近红外波段的反射率；R_r为红光波段的反射率。

然后，以10日为单位，采用通用的最大值合成法(MVC)将NDVI日数据合成NDVI旬数据(即10日数据)，以进一步减少云、大气、太阳高度角等的干扰。

需要说明的是，在本实施例中2004年、2007年和2008年各年3月至6月上旬MODIS原始日数据源于中国农业科学院农业资源与农业区划研究所卫星接收系统接收并存档的MODIS数据。该数据的预处理由数据接收处理系统完成：首先对接收的MODIS原始数据进行辐射校正和定位校正得到MODIS 1B数据；然后对MODIS 1B数据进行BOWTIE处理、大气校正和几何精校正；最后，利用每天MODIS第二近红外波段反射率和第一红光波段反射率计算得到日MODIS-NDVI。

合成的NDVI旬数据仍可能存在云的干扰或其它原因造成的数据丢失，故，可以进一步采用Savitzky-Golay滤波平滑方法，对NDVI旬数据时间序列进行平滑去噪处理，从而有效去除多时相NDVI遥感数据受云、气溶胶影响造成的噪音。最终，本实施例得到了较高质量的三年冬小麦MODIS-NDVI旬时序数据。如图3示出了冬小麦1月上旬至成熟期间NDVI平滑效果。

其中，Savitzky-Golay滤波平滑方法可以采用下述公式：

$Y_j=\frac{\underset{i=-m}{\overset{i=m}{\mathrm{\Σ}}}{C}_i{X}_{j+i}}{N}$

其中，X是未平滑的NDVI值，C_i是第i个NDVI值的权重系数，N是滤波卷积算子的NDVI值个数，其值大小等于平滑窗口的大小(2m+1)，j是在未平滑数据表中的序数，m值等于平滑窗口宽度的一半。

步骤103：根据步骤102获取的各个观测点的NDVI旬数据与时间变量构建NDVI旬数据-时间的曲线作为各个观测点的作物生长过程曲线，并据此构建与作物收获指数概念相关的参数HI_{NDVI_SUM}。具体如下：

首先，构建出如图4所示的一个观测点的NDVI旬数据-时间的作物生长过程曲线，图中并注明了该冬小麦的生育期，如图4中示出的返青期、开花期、乳熟初期、成熟期等。

然后，利用作物开花前曲线特征和开花后曲线特征，采用如下公式构建与作物收获指数概念相关的参数HI_{NDVI_SUM}：

${\mathrm{HI}}_{\mathrm{NDVI}\_\mathrm{SUM}}=\frac{\mathrm{\Σ}{\mathrm{NDVI}}_{\mathrm{post}}}{\mathrm{\Σ}{\mathrm{NDVI}}_{\mathrm{pre}}}$

其中，∑NDVI_post指作物NDVI旬数据-时间曲线峰后累积值，即作物开花后至乳熟初期NDVI累积值，该指标反映作物籽粒干物质积累过程；∑NDVI_pre指作物NDVI旬数据-时间曲线峰前累积值，即作物出苗至开花前NDVI累积值，该指标反映了作物茎、叶等干物质积累过程，本例中对冬小麦而言，由于冬小麦越冬期间地上部分干枯，NDVI已经很小，此时的NDVI值更多地表达了裸地NDVI，而非冬小麦NDVI。而返青后，冬小麦开始返青和分蘖，冬小麦NDVI逐渐增大。因此，对冬小麦而言，仅考虑返青-开花期NDVI累积值。但对其他不越冬作物，应从出苗开始计算。故，∑NDVI_post/∑NDVI_pre恰好反映了作物收获指数含义，即得到了籽粒产量与茎、叶等地上生物量的比值。

步骤104：分析步骤103所构建的作物收获指数相关参数HI_{NDVI_SUM}与步骤101历史实测作物收获指数之间的函数关系，并将此函数关系和HI_{NDVI_SUM}作为计算作物收获指数的依据。

本实施例在提取2004年和2007年与地面观测点相对应的HI_{NDVI_SUM}参数值基础上，分别建立HI_{NDVI_SUM}与地面冬小麦实测收获指数的关系，主要利用SPSS统计软件曲线拟合模块进行参数HI_{NDVI_SUM}与历史实测作物收获指数间的直线统计关系拟合。关于数据的直线拟合问题，将在下文介绍。如图5所示，所构建的参数与冬小麦收获指数具有很好的正相关关系，y＝0.4943x+0.2532，R²＝0.4598。其中，x为参数HI_{NDVI_SUM}，y为冬小麦实测收获指数，R为相关系数。当HI_{NDVI_SUM}为最小值0时，收获指数最小值为0.25。

步骤105：根据步骤104所建立的参数HI_{NDVI_SUM}与所述历史实测作物收获指数之间的函数关系，以及要计算年份的所述观测点的NDVI数据计算出该观测点作物收获指数。

在建立HI_{NDVI_SUM}参数与步骤101历史实测作物收获指数定量关系的基础上，本实施例利用平滑后的2008年3月至6月上旬NDVI旬数据，提取了2008年冬小麦HI_{NDVI_SUM}参数。然后，分别代入HI_{NDVI_SUM}与冬小麦实测收获指数的定量关系模型y＝0.4943x+0.2532，R²＝0.4598，其中x为各观测点2008年冬小麦HI_{NDVI_SUM}参数，y为需要计算得出的各观测点的冬小麦收获指数，R为相关系数，最终得到利用HI_{NDVI_SUM}参数预测的如图6所示的2008年本发明实施区域冬小麦收获指数空间分布图。

最后，本实施例对由参数HI_{NDVI_SUM}求得的收获指数进行精度验证。首先从2008年求得的收获指数结果中，提取与2008年实测相对应的500m缓冲区内冬小麦收获指数均值。然后将参数HI_{NDVI_SUM}求得的冬小麦收获指数与2008年的地面实测冬小麦进行对比(采用数据的直线拟合，将在下文介绍)，可以得出如图7所示的参数HI_{NDVI_SUM}预测冬小麦收获指数与实测冬小麦收获指数间的较好的相关关系，其中，参数HI_{NDVI_SUM}求得冬小麦收获指数的平均相对误差(相对误差/样本数)为2.40％(相对误差＝(预测值-真实值)/真实值*100％)，参数HI_{NDVI_SUM}求得冬小麦收获指数均方根误差

为0.02。可见，参数H_{INDVI_SUM}在发明实施区域内冬小麦收获指数均达到了较好的预测效果。

本实施例中，为了验证上述作物生长关键生育期NDVI与作物的绿度、生物量和作物单产等均具有较好的关系，首先利用2004年和2007年实施区域内11个县(市)冬小麦地面实测地上生物量和冬小麦单产数据以及相应旬NDVI数据，建立了实施区域小麦返青-开花前旬NDVI累积值与冬小麦收获期地上生物量以及开花期-乳熟初期旬NDVI累积值与冬小麦单产间的关系(采用数据的直线拟合，将在下文介绍)，并进行了精度验证，最终得到结果如下：

y_{ab_bio}＝212.53x₁+996.38   (n＝71，R²＝0.7581，SigF＝0.000)；

y_yield＝20.988x₂+448.95    (n＝71，R²＝0.6174，SigF＝0.048)。

其中，x₁为返青期至开花前旬NDVI累积值；x₂为开花期至乳熟初期旬NDVI累积值；y_{ab_bio}为冬小麦收获期地上生物量(g·m^-2)；y_yield为收获期冬小麦单产(kg·ha^-1)；n为样本数；R为相关系数；SigF为显著性标志。上述参数作用仅是说明自变量和因变量相关性显著。否则，说明两个变量没有线性关系。

同时，通过利用2008年46个地面调查点冬小麦实测生物量数据和单产数据及相应旬NDVI数据对实施区域上述模型验证可知，实施区域上述统计关系冬小麦生物量预测平均相对误差(相对误差/样本数)＝-3.00％，相对误差((预测值-真实值)/真实值*100％)范围-11.16％～8.7％，均方根误差

冬小麦单产预测平均相对误差＝-0.29％，相对误差范围-11.90％～9.39％，均方根误差为284.11kg·ha^-1。可见，冬小麦返青-开花前旬NDVI累积值、开花-乳熟初期旬NDVI累积值与冬小麦地上生物量和冬小麦收获期籽粒产量具有较好的相关关系。

上述HI_{NDVI_SUM}与地面冬小麦实测收获指数关系的建立与实施区域小麦返青-开花前旬NDVI累积值与冬小麦收获期地上生物量以及开花期-乳熟初期旬NDVI累积值与冬小麦单产间的关系，均采用了数据的直线拟合--常用一种以最小二乘法为基础的实验数据处理方法。最小二乘法的最佳经验公式y＝a+bx中a、b的求解：

通过实验，等精度地测得一组实验数据(x_i，y_i，i＝1，2...n)，设定此两物理量x、y满足线性关系，且假定实验误差主要出现在y_i上，设拟合直线公式为y＝f(x)＝a+bx，当所测各y_i值与拟合直线上各估计值f(x_i)＝a+bx_i之间偏差的平方和最小，即s＝∑[y_i-f(x_i)]²＝∑[y_i-(a+bx_i)]²→min时，所得拟合公式即为最佳经验公式。据此由 $\frac{\∂s}{\∂a}=-2\mathrm{\Σ}(y_i-a-{\mathrm{bx}}_i)=0,$ $\frac{\∂s}{\∂b}=-2\mathrm{\Σ}(y_i-a-{\mathrm{bx}}_i)x_i=0$ 解得 $a=\frac{\mathrm{\Σ}{x}_i{y}_i\mathrm{\Σ}{x}_i-\mathrm{\Σ}{y}_i\mathrm{\Σ}{x}_i^2}{{\left(\mathrm{\Σ}{x}_i\right)}^2-\mathrm{n\Σ}{x}_i^2},$ $b=\frac{\mathrm{\Σ}{x}_i\mathrm{\Σ}{y}_i-\mathrm{n\Σ}{x}_i{y}_i}{{\left(\mathrm{\Σ}{x}_i\right)}^2-\mathrm{n\Σ}{x}_i^2}.$ 最小二乘法处理数据除给出a、b外，还应给出相关系数r，r定义为 $r=\frac{\mathrm{\Σ}(x_i-\stackrel{\‾}{x})\mathrm{\Σ}(y_i-\stackrel{\‾}{y})}{\sqrt{\mathrm{\Σ}{(x_i-\stackrel{\‾}{x})}^2}\sqrt{\mathrm{\Σ}{(y_i-\stackrel{\‾}{y})}^2}},$ 其中 $\stackrel{\‾}{x}=\frac{\mathrm{\Σ}{x}_i}{n},$ $\stackrel{\‾}{y}=\frac{\mathrm{\Σ}{y}_i}{n}.$ r表示两变量之间的函数关系与线性的符合程度，r∈[-1，1]。|r|→1，x、y间线性关系好，|r|→0，x、y间无线性关系，拟合无意义。同理，也可证假定实验误差主要出现在x_i上。

由上述实施例可知，构建参数HI_{NDVI_SUM}在区域范围内反演冬小麦收获指数取得了很好的效果，证明本发明利用构建参数HI_{NDVI_SUM}反演区域冬小麦收获指数信息方法的可行性。同时，本实施例仅以冬小麦为例，阐述和实施了如何利用遥感时序植被指数提取区域冬小麦收获指数的方法，但该方法对于一般的农作物收获指数遥感反演具有普遍的适用性。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (8)

1.一种作物收获指数的获取方法，其特征在于，所述方法包括：

A、获取地面不同观测点的历史实测作物收获指数；

B、获取所述不同观测点不同时间的时序归一化植被指数NDVI数据；

C、根据各个观测点不同时间的NDVI数据结合作物生长过程曲线构建与作物收获指数概念相关参数HI_{NDVI_SUM}；

D、构建HI_{NDVI_SUM}与历史实测作物收获指数之间的函数关系；

E、根据所构建的函数关系和所要计算年份的观测点不同时间的NDVI数据计算该观测点作物收获指数。

2.根据权利要求1所述的方法，其特征在于，所述步骤A包括子步骤：

A1、在每个观测点内选取至少3个以上的取样点采用作物实割实测获得所述取样点的作物收获指数；

A2、将观测点内各取样点的作物收获指数进行均值处理作为该观测点作物的实测收获指数。

3.根据权利要求1或2所述的方法，其特征在于，所述步骤B包括子步骤：

B1、获取观测点的MODIS数据；

B2、根据MODIS数据中的近红外波段反射率和红光波段反射率计算生成观测点的NDVI日数据，其中采用如下计算公式：

$\mathrm{NDVI}=\frac{R_n-R_r}{R_n+R_r}$

其中，R_n为近红外波段的反射率；R_r为红光波段的反射率。

4.根据权利要求3所述的方法，其特征在于，步骤B2后进一步包括：

B3、采用最大值合成法MVC将NDVI日数据合成NDVI旬数据。

5.根据权利要求4所述的方法，其特征在于，步骤B3后进一步包括：

采用Savitzky-Golay滤波平滑方法对NDVI旬数据时间序列进行平滑去噪处理，所述Savitzky-Golay滤波平滑方法采用如下公式：

$Y_j=\frac{\underset{i=-m}{\overset{i=m}{\mathrm{\Σ}}}{C}_i{X}_{j+i}}{N}$

其中，X是未平滑的NDVI值，C_i是第i个NDVI值的权重系数，N是滤波卷积算子的NDVI值个数，其值大小等于平滑窗口的大小(2m+1)，j是在未平滑数据表中的序数，m值等于平滑窗口宽度的一半。

6.根据权利要求1、4或5所述的方法，其特征在于，所述步骤C包括：

利用作物开花前曲线特征和开花后曲线特征，采用如下公式构建与作物收获指数概念相关参数HI_{NDVI_SUM}：

${\mathrm{HI}}_{\mathrm{NDVI}\_\mathrm{SUM}}=\frac{\mathrm{\Σ}{\mathrm{NDVI}}_{\mathrm{post}}}{\mathrm{\Σ}{\mathrm{NDVI}}_{\mathrm{pre}}}$

其中，∑NDVI_post为作物开花后至乳熟初期NDVI累积值；∑NDVI_pre为作物出苗至开花前NDVI累积值。

7.根据权利要求6所述的方法，其特征在于，所述作物为越冬作物时，

所述∑NDVI_pre为作物返青至开花前NDVI累积值。

8.根据权利要求1所述的方法，其特征在于，步骤D所述构建的方法包括：

采用数据的直线拟合方法构建。

Images