基于时间序列数据的南方冬种作物种植面积提取方法

摘要

本发明公开了一种基于时间序列数据的南方冬种作物种植面积提取方法，其包括以下步骤：步骤S1、获取研究地区的遥感影像；步骤S2、对遥感影像预处理；步骤S3、计算不同时相的增强植被指数EVI；步骤S4、得到EVI时间序列影像数据；步骤S5、提取EVI时序特征参考曲线；步骤S6、计算每个像元的EVI时间序列曲线与EVI时序特征参考曲线的光谱夹角值；步骤S7、确定用于提取研究区域作物面积的阈值；步骤S8、将光谱夹角值图像中像元的光谱夹角值与阈值比较；步骤S9、确定监测区冬种作物马铃薯的种植面积。本发明只须在监测区域中测定若干个已知目标作物种植区的经纬度即能够在图像上获得EVI时序特征参考曲线，就能计算机程序自动实现面积提取，省时省力效率高。

权利要求

1.一种基于时间序列数据的冬种马铃薯种植面积提取方法，其特征在于，其包括以下步骤：步骤S1、通过HJ-1A/B CCD传感器获取研究地区的遥感影像；

步骤S2、对所述遥感影像进行预处理，所述预处理包括辐射定标、大气校正、几何校正、镶嵌和剪切；

步骤S3、计算不同时相的增强植被指数EVI，得到各时相的EVI影像：其中，G为背景调节参数，L是土壤修正因子，C1和C2是利用蓝光波段修正红光波段中大气气溶胶散射影响系数，ρNIR、ρred、ρblue分别是近红外波段、红光波段和蓝光波段的反射率；

步骤S4、将各时相的EVI影像叠加在一起，得到EVI时间序列影像数据；

步骤S5、提取EVI时序特征参考曲线，具体方法是：在研究地区的种植区域选取多个样本点，并从该多个样本点中随机抽取部分样本点构成训练样本点，通过GPS实地测量所述训练样本点的经纬度；根据实测的训练样本点经纬度数据，在EVI时间序列影像数据确定这些训练样本点对应在EVI时间序列影像上的位置，提取所述位置对应的EVI时间序列曲线，这些位置对应的EVI时间序列曲线的平均值作为EVI时序特征参考曲线；

步骤S6、使用光谱角度匹配方法，计算每个像元的EVI时间序列曲线与EVI时序特征参考曲线的光谱夹角值，得到每个像元的光谱夹角值图像：即：

其中：θi为第i个像元的光谱夹角值， Xi为第i个像元的EVI时间序列曲线，Y是EVI时序特征参考曲线，xij是第i个像元的EVI时间序列曲线的第j个向量点，yj为EVI时序特征参考曲线第j个向量点，n为向量维数；

步骤S7、确定用于提取研究区域作物面积的阈值：计算训练样本点EVI时间序列曲线与参考曲线的广义夹角的平均值δ和广义夹角的标准差σ：其中，m为训练样本点个数，1≤k≤m，θk′为第k个训练样本点的EVI时间序列曲线与EVI时序特征参考曲线之间的夹角；

利用平均值δ和夹角的标准差σ作为提取研究区域作物面积的阈值ε：ε＝δ+σ   (6)

步骤S8、将光谱夹角值图像中第i个像元的光谱夹角值θi与阈值ε比较，如果光谱夹角值θi小于阈值ε，则把该第i个像元归为目标地物；

步骤S9、根据被归为目标地物的像元数量和遥感影像中每个像元对应的实际面积的乘积确定监测区冬种作物马铃薯的种植面积。

说明书

基于时间序列数据的南方冬种作物种植面积提取方法

技术领域

[0001] 本发明涉及遥感技术领域，具体涉及一种基于物候特征的HJ-1 A/B CCD EVI时间序列数据的南方冬种作物种植面积的遥感监测和估算方法。

背景技术

[0002] 南方冬种作物在我国粮食生产中具有重要地位，确定作物的面积及空间分布是监测农作物种植信息的重要内容。但是，目前南方大部分农作物种植面积和空间分布信息，依赖的还是传统的统计汇报或实地土地测量。这个过程中需要很多的时间和人力物力，并且所得到的农情并不准确、不及时，无法提供农作物种植的空间分布状况。卫星遥感为快速、准确、大范围地获取农作物种植面积信息，监测农作物的空间分布状况提供了重要的技术手段。

[0003] 很多利用时间序列物候特征的土地覆盖制图的研究是借助高时间分辨率传感器的数据，如AVHRR和MODIS；或借助中高空间分辨率的传感器数据，如SPOT和Landsat。然而，在考虑到华南地区多云多雨的气候状况及其相对破碎的耕地地块，这些数据要么空间分辨率太低，要么时间分辨率太低以至于不适用于小范围尺度。我们认为搭载在中国环境与灾害预报卫星上的HJ-1 A/B CCD传感器由于具有2天重复周期的高时间分辨率和30m的高空间分辨率的特点，能够为小尺度的土地覆盖制图研究提供合适的数据源。

发明内容

[0004] 针对现有技术的上述缺陷，本发明的目的在于提供一种基于时间序列数据的南方冬种作物种植面积提取方法，其只须在监测区域中测定若干个已知目标作物种植区的经纬度即能够在图像上获得EVI时序特征参考曲线，就能计算机程序自动实现面积提取，省时省力效率高的同时也避免了传统监督分类通过概率和依赖人为经验的缺点。

[0005] 本发明解决技术问题所采用的技术方案如下：

[0006] 一种基于时间序列数据的南方冬种作物种植面积提取方法，其包括以下步骤：

[0007] 步骤S1、通过HJ-1 A/B CCD传感器获取研究地区的遥感影像；

[0008] 步骤S2、对所述遥感影像进行预处理，所述预处理包括辐射定标、大气校正、几何校正、镶嵌和剪切；

[0009] 步骤S3、计算不同时相的增强植被指数EVI，得到各时相的EVI影像：

[0010]

[0011] 其中，G为背景调节参数，L是土壤修正因子，C1和C2是利用蓝光波段修正红光波段中大气气溶胶散射影响系数，ρNIR、ρred、ρblue分别是近红外波段、红光波段和蓝光波段的反射率；

[0012] 步骤S4、将各时相的EVI影像叠加在一起，得到EVI时间序列影像数据；

[0013] 步骤S5、提取EVI时序特征参考曲线，具体方法是：在研究地区的种植区域选取多个样本点，并从该多个样本点中随机抽取部分样本点构成训练样本点，通过GPS实地测量所述训练样本点的经纬度；根据实测的训练样本点经纬度数据，在EVI时间序列影像数据确定这些训练样本点对应在EVI时间序列影像上的位置，提取所述位置对应的EVI时间序列曲线，这些位置对应的EVI时间序列曲线的平均值作为EVI时序特征参考曲线；

[0014] 步骤S6、使用光谱角度匹配方法，计算每个像元的EVI时间序列曲线与EVI时序特征参考曲线的光谱夹角值，得到每个像元的光谱夹角值图像：

[0015]

[0016] 即：

[0017] 其中：θi为第i个像元的光谱夹角值， Xi为第i个像元的EVI时间序列曲线，Y是EVI时序特征参考曲线，xij是第i个像元的EVI时间序列曲线的第j个向量点，yj为EVI时序特征参考曲线第j个向量点，n为向量维数；

[0018] 步骤S7、确定用于提取研究区域作物面积的阈值：计算训练样本点EVI时间序列曲线与参考曲线的广义夹角的平均值δ和广义夹角的标准差σ：

[0019]

[0020]

[0021] 其中，m为训练样本点个数，1≤k≤m，θ′k为第k个训练样本点的EVI时间序列曲线与EVI时序特征参考曲线之间的夹角；

[0022] 利用平均值δ和夹角的标准差σ作为提取研究区域作物面积的阈值ε：

[0023] ε＝δ+σ    (6)

[0024] 步骤S8、将光谱夹角值图像中第i个像元的光谱夹角值θi与阈值ε比较，如果光谱夹角值θi小于阈值ε，则把该第i个像元归为目标地物；

[0025] 步骤S9、根据被归为目标地物的像元数量和遥感影像中每个像元对应的实际面积的乘积确定监测区冬种作物马铃薯的种植面积。

[0026] 本发明的有益效果是：因为影像数据很确定，EVI的计算有现成软件实现，只须在监测区域中测定若干个已知目标作物种植区的经纬度即能够在图像上获得EVI时序特征参考曲线，就能计算机程序自动实现面积提取，比野外测量汇报的方法省时省力效率高，也避免了传统监督分类通过概率和依赖人为经验的缺点；而且可以通过查找表逐年优化，实现多年重复使用。

附图说明

[0027] 图1研究区稔平半岛地理位置示意图及其真彩色合成遥感影像；

[0028] 图2数据处理与分类方法步骤：从数据获取到提取方法；

[0029] 图3惠州示范区铁涌镇、稔山镇、平海镇和稔平半岛的训练样本点EVI时间序列平均曲线(参考曲线)及其与“早稻-晚稻-冬种马铃薯”物候期的对应关系图；

[0030] 图4研究区稔平半岛冬种马铃薯不同时期的生长情况；

[0031] 图5研究区稔平半岛马铃薯种植面积提取结果；

[0032] 图6研究区稔平半岛马铃薯样本点分布情况与提取结果；

[0033] 图7(a)研究区稔平半岛马铃薯提取结果与马铃薯专业镇及专业村分布情况；(b)铁涌镇样本点及马铃薯网点分布情况。

具体实施方式

[0034] 下面结合附图和具体实施方式对本发明的内容做进一步详细说明。

[0035] 实施例：

[0036] 请参照图2所示，基于时间序列数据的南方冬种作物种植面积提取方法，其包括以下步骤：

[0037] 1.监测区域：如图1所示，本发明以示范区广东省惠州市惠东县地处南海之滨的马铃薯种植为例，对本发明的保护范围进行解释和说明，该地区属亚热带海洋性气候，全年无霜，冬季少雨、阳光充足、气候暖和是冬种马铃薯的理想之地，已成为全国冬种面积最大、产量最高的县区。全县0.867万公顷马铃薯基本集中在沿海，即本文的监测区稔平半岛。稔平半岛包括稔山镇、铁涌镇、平海镇和巽寮镇四个镇，其中稔山镇、铁涌镇和平海镇被评为广东省马铃薯专业镇。

[0038] 2.监测区域物候特征

[0039] 惠东县稔平半岛采用早稻—晚稻—冬种马铃薯的轮作模式。总体而言，早稻是清明后(4月上旬，约第95天)抛秧，7月底(约第206天)收割；晚稻是8月上中旬(约第214至227天)抛秧，沿海地区10月底(约第298天)收割，其他地区11月初(约第312天)收割；冬种马铃薯播种期，一般在11月初(约第312天)至11月下旬(约第331天)播种，次年3月底(约第85天)前全部完成收割。

[0040] 相关影像选择的物候依据：

[0041] 1)1月8日(第8天)：研究区马铃薯均未收获，处于茂盛阶段，EVI高；

[0042] 2)1月31日(第31天)和2月8日(第39天)：研究区马铃薯即将收获，马铃薯叶逐渐变黄，EVI稍有下降，参照图4(c)；

[0043] 3)3月26日(第85天)：研究区马铃薯均已收获，并且早稻均未开始抛秧，耕地处于闲置状态，EVI很低；

[0044] 4)4月11日和4月19日(约105天)：研究区马铃薯均已收获，并且早稻均未开始抛秧，耕地处于闲置状态，EVI很低；

[0045] 5)5月24日(约144天)：研究区早稻均已经抛秧完毕，并生长了一段时间，水稻处于茂盛时期，EVI较高；

[0046] 6)8月2日(第214天)：研究区早稻均已收获，并且晚稻均未开始抛秧，耕地处于闲置状态，EVI很低；

[0047] 7)9月18日(第261天)：研究区晚稻均已经抛秧完毕，并生长了一段时间，水稻处于茂盛时期，EVI较高；

[0048] 8)10月16日(第289天)：研究区晚稻已成熟，准备收割，EVI较低。

[0049] 9)11月8日(第312天)：研究区晚稻均已收获，马铃薯播种前，耕地处于闲置状态，EVI很低。

[0050] 10)11月27日(第331天)：研究区马铃薯均已播种完毕，马铃薯逐渐进入出苗期生长阶段，EVI有所升高。

[0051] 11)12月10日(第344天)、12日18日(第352天)和12月25日(第359天)：马铃薯生长了一段时间，处于茂盛阶段，EVI偏高，参照图4(a)和(b)。

[0052] 为了避免马铃薯集中上市而导致的马铃薯滞销，当地种植人员有意识的错开马铃薯的播种和收获时间。稔平半岛冬种马铃薯播种时间在11月初(约第312天)开始，铁涌镇最早，稔山镇在铁涌镇之后，平海镇最晚播种，大约11月下旬(约第329天)播种工作基本完成。经过70-80天的种植后，铁涌镇的收获时间大约是在2月底3月初(约第60天)，之后是稔山镇，最迟收获的平海镇大约在3月中下旬(约第85天)完成全部收割。示范区稔平半岛各地区冬种马铃薯物候期见表1。

[0053] 表1 稔平半岛冬种马铃薯物候期

[0054]

[0055] 3.获取对应时期的遥感影像

[0056] 本文采用影像数据是由中国资源卫星应用中心(China Centre For Resources Satellite Data and Application)提供的HJ-1 A/B CCD传感器数据。为了有效体现监测区一年的种植模式，本研究选择2011年监测区无云或者云少的数据为主要数据，用2010年或2012年的数据补充缺少数据的月份。所使用的HJ-1 A/B CCD传感器数据列表见表2：

[0057] 表2 HJ-1 A/B CCD传感器数据列表

[0058]

[0059] 其中，2011年4月的数据均有云层覆盖部分监测区，然而4月是一个表现监测区物候规律的重要时间节点。因此我们使用2011年4月11日和2011年4月19日两景时间相近、物候情况差异较小的影像中无云部分进行镶嵌，得到2011年4月监测区的影像数据。因此在分析过程中，EVI时间序列曲线中有14个时间节点。

[0060] 4.影像EVI时序数据的获取

[0061] 对获取的影像进行预处理，包括辐射定标和大气校正，通过下面的公式计算不同时相影像的EVI值：

[0062]

[0063] 式中，G为背景调节参数，L是土壤修正因子，C1和C2是利用蓝光波段修正红光波段中大气气溶胶散射影响的系数，ρNIR、ρred、ρblue分别是近红外波段(0.7-1.1μm)、红光波段(0.6-0.7μm)和蓝光波段(0.45-0.52μm)的反射率。一般而言，G＝2.5，C1＝6.0，C2＝7.5，and L＝1。将各时相的EVI影像叠加在一起，得到EVI时间序列影像数据。

[0064] 5.EVI时间序列参考曲线的获取

[0065] 我们从130个种马铃薯种植区样本点中抽取26个(20％)的样本点作为训练样本点，在EVI的时间序列数据中取得每个训练样本点对应的EVI时间序列曲线。在每个时间节点，求得各个镇的训练样本点除去最大值和最小值的平均值，所有时间节点的平均值组合在一起作为该镇的EVI时间序列参考曲线，如图3所示。

[0066] 6.计算光谱夹角

[0067] 使用光谱角度匹配方法，通过以下公式计算像元的EVI时间序列与参考曲线的广义夹角，得到每个像元的光谱夹角值图像：

[0068]

[0069] 即：

[0070] 式中，θi为第i个像元的光谱夹角值， Xi为第i个像元的EVI时间序列曲线，Y是EVI时序特征参考曲线，xij是第i个像元的EVI时间序列曲线的第j个向量点，yj为EVI时序特征参考曲线第j个向量点，n为向量维数。

[0071] 7.像元的光谱夹角值与光谱角阈值的比较

[0072] 计算训练样本点EVI时间序列曲线(向量)弧度的平均值δ和弧度的标准差σ：

[0073]

[0074]

[0075] m为训练样本点个数，m＝26，k＝1，2，……26，θ′k为第k个训练样本点的EVI时间序列曲线与EVI时序特征参考曲线之间的夹角，训练样本点EVI时间序列曲线(向量)弧度即为训练样本点的EVI时间序列曲线与EVI时序特征参考曲线之间的夹角。

[0076] 利用平均值δ和夹角的标准差σ作为提取研究区域作物面积的阈值ε：

[0077] ε＝δ+σ    (12)

[0078] 将光谱夹角值图像中第i个像元的光谱夹角值θi与阈值ε比较，如果光谱夹角值θi小于阈值ε，则把该第i个像元归为目标地物。

[0079] 8.获得监测区域马铃薯分布图和计算估算面积：

[0080] 监测区域马铃薯分布图见附图5；

[0081] 根据HJ-1 A/B CCD影像的空间分辨率(即每个像元对应的实际面积)和被归为目标地物的像元数目统计，获得估算的面积大小：

[0082] 估算面积大小＝30m×30m×目标地物像元数目。

[0083] 9.分类结果精度评定

[0084] 我们监测区调研一共获取130个实地测量的马铃薯样本点用以评价冬种马铃薯种植分布图的精度，其中110个被正确识别(正确率为84.6％)；130个样本点的分布情况和每个镇的正确率分别如图6和表3所示。

[0085] 表3 基于冬种马铃薯种植区样本点的提取精度评价表

[0086]

[0087] 另外，也可将提取结果与2011年的马铃薯种植面积数据(来源于2011年惠东县统计资料)进行比较，稔平半岛各镇的提取精度见表4。提取结果中，稔山镇冬种马铃薯种植面积最大(1727.19公顷)，其次是铁涌镇(1594.26公顷)和平海镇(1008.72公顷)，总体次序与统计数据一致。稔平半岛总体的相对误差为10.82％。稔山镇和铁涌镇冬种马铃薯面积之和为3321.45公顷，占稔平半岛冬种马铃薯面积的73.96％，它们的相对误差分别是-4.09％和-4.39％。

[0088] 表4 基于统计数据的提取精度评价表

[0089]

[0090] 从图7(a)和(b)中可见，提取结果中集中连片的马铃薯种植区与马铃薯“专业镇/村”的位置一致；铁涌镇的提取结果与该镇马铃薯种植网点的空间分布情况也一致。

[0091] 应当理解的是，本发明的应用不限于上述的举例，对本领域普通技术人员来说，可以根据上述说明加以改进或变换，所有这些改进和变换都应属于本发明所附权利要求的保护范围。