基于多元遥感影像的多云雾区域农作物种植模式提取方法

摘要

本发明属于图像处理技术领域，具体涉及一种基于多元遥感影像的多云雾区域农作物种植模式提取方法；基于多时相的Sentinel‑2数据，结合广汉市物候节律和作物多时相光谱特征确定小麦、油菜等冬季作物阈值函数；基于时序Sentinel‑1数据构建SAR纹理特征数据集，确定水稻、玉米等夏季作物阈值函数；基于不同作物的物候差异，统计农作物的时间‑特征变化曲线；基于农作物时间‑特征变化曲线，确定农作物时间‑特征阈值判定函数；采用面向对象策略,提取农作物类别及时空分布，提取农作物种植模式。

权利要求

1.基于多元遥感影像的多云雾区域农作物种植模式提取方法，其特征在于，包括以下步骤：数据获取：获取多时相多光谱数据以及时序SAR影像数据，并分别对多光谱数据以及SAR影像数据进行预处理，得到光谱特征和后向散射系数；

特征提取：基于光谱特征和后向散射系数进行计算，分别得到光谱指数特征集和SAR纹理特征集；

阈值函数确定：基于光谱指数特征以及训练样本得到冬季作物阈值函数；基于SAR纹理特征以及训练样本得到夏季作物阈值函数；

农作物分类：采用面向对象的分类策略，提取农作物的类型以及空间分布；

所述农作物的分类包括以下步骤：

A.多尺度分割：利用ESP算法计算不同分割尺度下对象同质性局部变化的变化率，用于表示分割效果的最佳参数，计算公式如下：式中：LVn为分割尺度为n时分割对象的局部方差，LVn‑1是分割尺度为n‑1时分割对象的局部方差，当ROC_LVn出现峰值时，该点对应的分割尺度即为最佳分割尺度；

B.面向对象分类：基于实地调查采样所获取的样本建立分类体系；通过多尺度分割提取对象层级的光谱特征，SAR纹理特征，时相物候特征，利用特征向量构建CART决策树分类器自动建立分类规则进行分类；

种植模式提取：基于不同作物的物候差异，统计农作物的时间‑特征变化曲线，将所述光谱指数特征集和SAR纹理特征集进行组合，基于农作物时间‑特征变化曲线，确定农作物时间‑特征阈值判定函数；在多时相影像上提取农作物时间‑空间分布，提取农作物种植模式；种植模式即农作的类型、时空分布及轮作模式。

2.根据权利要求1所述的基于多元遥感影像的多云雾区域农作物种植模式提取方法，其特征在于，所述多光谱数据为Sentinel‑2多光谱数据；

Sentinel‑2多光谱数据的预处理：利用IDL编程语言批量进行波段叠加、影像拼接及裁剪，得到包含7个波段的光谱特征集。

3.根据权利要求1所述的基于多元遥感影像的多云雾区域农作物种植模式提取方法，其特征在于：所述SAR影像数据采用Sentinel‑1卫星数据，SAR影像数据的获取模式为IW模式的GRDH格式产品，空间分辨率为20m，每景SAR影像都具有同极化和交叉极化两种方式；

SAR影像数据的预处理：利用欧空局SNAP软件构建处理模型批量对SAR影像数据进行预处理，包括辐射定标、多普勒地形校正、斑点滤波以及地理编码。

4.根据权利要求1所述的基于多元遥感影像的多云雾区域农作物种植模式提取方法，其特征在于：所述SAR纹理特征包括均值、协同性、方差以及对比度，具体计算公式如下：式中：Mean为均值；Homogeneity为协同性；Variance为方差；Contrast为对比度；n为灰度值的阶数；p(i,j)是n×n的归一化共生矩阵；μ为p(i,j)的均值；i和j表示相距为d的两个灰度级。

5.根据权利要求1所述的基于多元遥感影像的多云雾区域农作物种植模式提取方法，其特征在于，所述光谱指数特征包括归一化植被指数、红边植被指数以及红边叶绿素指数，具体计算公式如下：NDVI＝(ρnir‑ρred)/(ρnir+ρred)；

NDVI705＝(ρ750‑ρ705)/(ρ750+ρ705)；

CI＝(ρ750‑800/ρ690‑725)‑1；

式中：NDVI表示归一化植被指数；NDVI705表示红边植被指数；CI表示红边叶绿素指数；

ρnir为近红外波段反射率ρred为红波段反射率；ρ750为中心波长为750nm波段的反射率；ρ705为中心波长为705nm波段的反射率；ρ750‑800为中心波长区间为750nm～800nm波段的反射率；

ρ690‑725为中心波长区间为690nm～725nm波段的反射率。

说明书

基于多元遥感影像的多云雾区域农作物种植模式提取方法

技术领域

[0001] 本发明属于图像处理技术领域，具体涉及一种基于多元遥感影像的多云雾区域农作物种植模式提取方法。

背景技术

[0002] 快速准确地获取农作物种植结构和时空分布信息对于政府部门制定粮食政策、调整农业结构、保障国家粮食安全十分重要；遥感技术具有覆盖范围大、探测周期短、费用成本低等优势，已成为大尺度范围农作物监测的主要手段之一。

[0003] 目前国内外相关研究多基于高分辨率的光学影像，利用农作物物候特征差异和光谱特征差异提取农作物的时序生长曲线，从而实现农作物分类识别；例如如下现有技术：

[0004] 1.吴静, 吕玉娜, 李纯斌，等. 基于多时相Sentinel‑2A的县域农作物分类 [J]. 农业机械学报, 2019, 50(9): 194‑200；

[0005] 2.杜保佳, 张晶, 王宗明，等. 应用Sentinel‑2A NDVI时间序列和面向对象决策树方法的农作物分类 [J]. 地球信息科学学报, 2019, 21(5): 740‑751；

[0006] 3.潘力, 夏浩铭, 王瑞萌，等. 基于Google Earth Engine的淮河流域越冬作物种植面积制图 [J]. 农业工程学报, 2021, 37(18): 211‑218；

[0007] 现有技术1中基于多基于多时相Sentinel‑2数据对景泰县农作物进行分类，结果表明，Sentinel‑2 特有的红边波段数据及其较高的空间分辨率在县域农作物分类上具有较好的精度。

[0008] 现有技术2中基于Sentinel‑2 NDVI时序数据，利用面向对象决策树分类方法对北安市进行农作物识别，分类精度达96.2%。

[0009] 现有技术3中基于Landsat和Sentinel‑2融合数据，根据不同作物物候期差异，准确提取淮河流域越冬作物空间分布。

[0010] 基于上述可知，光学影像在农作物遥感监测中已经取得了显著成果，但上述研究多集中在地块平坦整齐、光学影像获取质量较高的区域，中国西南区域上空常年多云雾覆盖，获取的高质量光学影像非常有限，因此仅利用光学影像，并不能满足多云雾区域的农作物种植模式提取监测需求；

[0011] 而雷达卫星能够穿透云雾，全天时全天候观测地面，获得各种作物在不同极化方式下的雷达相应特性，从而反映作物的形态、结构和生长等信息，为农作物监测提供数据支持，例如如下现有技术：

[0012] 4.杨沈斌, 李秉柏, 申双和，等. 基于多时相多极化差值图的稻田识别研究 [J]. 遥感学报, 2008, 12(4): 613‑619；

[0013] 5.王晨丞, 王永前, 王利花. 基于SAR纹理信息的农作物识别研究——以农安县为例 [J]. 遥感技术与应用, 2021, 36(2): 372‑380；

[0014] 现有技术4中基于多时相多极化SAR数据对江西省高安地区水稻进行提取，分类精度达到84.92%；

[0015] 现有技术5中以Sentinel‑1双极化数据为数据源，对吉林省农安县典型农作物进行分类制图，结果表明，基于SAR影像数据对农作物分类具有可行性。

[0016] 但目前利用SAR数据进行农作物监测主要集中在空间异质性较低的区域，且多以水稻识别为主，对于地块破碎、种植结构复杂的旱地作物识别研究较少，尚不能满足区域性农作物监测需求。

[0017] 并且若是研究区上空常年覆盖云雾，获取的高质量光学影像有限，难以满足农作物监测需求，仅利用获取的Sentinel‑2光学数据，难以区分水稻、玉米等夏季作物。星载SAR影像获取不受气候状况限制，为多云雾气候区域夏季农作物提取提供了可能性。光学影像和SAR影像在农作物识别提取中各有优势，目前国内外相关研究多基于光学影像或SAR影像进行农作物分类提取，将二者融合协同用于农作物识别的研究很少。

[0018] 并且目前对于农作物的识别依旧停留在农作物种类的识别上，对于如何实现农作物种植模式的提取，目前在国内外还属于空白区域，基于此有必要提供一种农作物种植模式的提取方法，以期为全天候农作物遥感监测提供新的思路。

发明内容

[0019] 本发明公开了一种基于多元遥感影像的多云雾区域农作物种植模式提取方法，拟实现对于目前地块破碎、种植结构复杂、多云雾区域的旱地作物的识别，并实现农作物种植模式的提取。

[0020] 为解决上述技术问题，本发明采用的技术方案如下：

[0021] 基于多元遥感影像的多云雾区域农作物种植模式提取方法，包括以下步骤：

[0022] 数据获取：获取多时相多光谱数据以及时序SAR影像数据，并分别对多光谱数据以及SAR影像数据进行预处理，得到光谱特征和后向散射系数；

[0023] 特征提取：基于光谱特征和后向散射系数进行计算，分别得到光谱指数特征集和SAR纹理特征集；

[0024] 阈值函数确定：基于光谱指数特征以及训练样本得到冬季作物阈值函数；基于SAR纹理特征以及训练样本得到夏季作物阈值函数；

[0025] 农作物分类：采用面向对象的分类策略，提取农作物的类型以及空间分布；

[0026] 种植模式提取：基于不同作物的物候差异，统计农作物的时间‑特征变化曲线，将所述光谱指数特征集和SAR纹理特征集进行组合，基于农作物时间‑特征变化曲线，确定农作物时间‑特征阈值判定函数；在多时相影像上提取农作物时间‑空间分布，提取农作物种植模式。

[0027] 所述种植模式为农作物的耕种模式，例如：小麦‑水稻、油菜‑水稻、土豆‑大豆、土豆‑玉米4种年内复种模式；粮油作物种植模式多为一年两熟的耕作模式，主要以“小麦‑水稻轮种”、“油菜‑水稻轮种”、“土豆‑玉米轮种”、“土豆‑大豆轮种”为主。

[0028] 优选的，所述多光谱数据为Sentinel‑2多光谱数据；

[0029] Sentinel‑2多光谱数据的预处理：利用IDL编程语言批量进行波段叠加、影像拼接及裁剪，得到包含7个波段的光谱特征集。

[0030] 优选的，所述SAR影像数据采用Sentinel‑1卫星数据，SAR影像数据的获取模式为IW模式的GRDH格式产品，空间分辨率为20m，每景SAR影像都具有同极化和交叉极化两种方式；

[0031] SAR影像数据的预处理：利用欧空局SNAP软件构建处理模型批量对SAR影像数据进行预处理，包括辐射定标、多普勒地形校正、斑点滤波以及地理编码。

[0032] 优选的，所述SAR纹理特征包括均值、协同性、方差以及对比度，具体计算公式如下：

[0033] ；

[0034] ；

[0035] ；

[0036] ；

[0037] 式中：Mean为均值；Homogeneity为协同性；Variance为方差；Contrast为对比度；n为灰度值的阶数；p(i,j)是n×n的归一化共生矩阵；μ为p(i,j)的均值；i和j表示相距为d的两个灰度级。

[0038] 优选的，所述光谱指数特征包括归一化植被指数、红边植被指数以及红边叶绿素指数，具体计算公式如下：

[0039] ；

[0040] ；

[0041] ；

[0042] 式中：NDVI表示归一化植被指数；NDVI1705表示红边植被指数；CI表示红边叶绿素指数； 为近红外波段反射率；  为红波段反射率； 为中心波长为750nm波段的反射率； 为中心波长为705nm波段的反射率； 为中心波长区间为750nm~
800nm波段的反射率； 为中心波长区间为690nm 725nm波段的反射率。
~

[0043] 优选的，所述农作物的分类包括以下步骤：

[0044] A.多尺度分割：利用ESP算法计算不同分割尺度下对象同质性局部变化的变化率，用于表示分割效果的最佳参数，计算公式如下：

[0045] ；

[0046] 式中： 为分割尺度为n时分割对象的局部方差， 是分割尺度为n‑1时分割对象的局部方差，当 出现峰值时，所述峰值对应的分割尺度即为最佳分割尺度；

[0047] B.面向对象分类：基于实地调查采样所获取的样本建立分类体系；通过多尺度分割提取对象层级的光谱特征，SAR纹理特征，时相物候特征，利用特征向量构建CART决策树分类器自动建立分类规则进行分类。

[0048] 综上所述，由于采用了上述技术方案，本发明的有益效果是：

[0049] 针对常年覆盖云雾，光学影像获取受限且空间异质性较高的区域农作物监测的瓶颈，基于多元影像协同，提取夏季作物、冬季作物的光谱特征集和SAR纹理特征集，结合农作物物候特征，实现了农作物种植模式提取，填补了国内多云雾区域农作物种植模式快速提取的空白，为全天候农作物监测提供了新的思路。

附图说明

[0050] 本发明将通过例子并参照附图的方式说明，其中：

[0051] 图1为本发明的整体流程图。

[0052] 图2为研究区位置及作物采样点分布示意图。

[0053] 图3为分割尺度评价图。

[0054] 图4为区域分割效果对比，其中a部分表示裸地水体区域未分割影像，b部分表示裸地水体区域中分割尺度为80的分割影像，c部分表示裸地水体区域中分割尺度为118的分割影像；d部分表示耕地区域未分割影像，e部分表示耕地区域中分割尺度为80的分割影像，f部分表示耕地区域中分割尺度为118的分割影像。

[0055] 图5为多时相植被指数曲线图，其中图5的a部分表示农作物的归一化植被指数；图5的b部分表示农作物的红边叶绿素指数；图5的c部分表示农作物的红边植被指数。

[0056] 图6为2021‑03‑27光谱特征图。

[0057] 图7为研究区夏季作物后向散射系数时间序列图，其中图7中的a部分表示VH极化后向散射系数时间序列图；图7中的b部分表示VV极化后向散射系数时间系列图。

[0058] 图8为研究区夏季作物纹理特征时间序列图，图8的a部分表示VH极化协同性时间序列图；图8的b部分表示VV极化协同性时间序列图；图8的c部分表示VH对比度时间序列图；图8的d部分表示VV对比度时间序列图；图8的e部分表示VH均值时间序列图；图8的f部分表示VV均值时间序列图；图8的g部分表示VH方差时间序列图；图8的h部分表示VV方差时间序列图。

[0059] 图9为研究区不同分类方案分类结果；图9中的a部分附图为面对对象分类结果；图9中的b部分附图为面对像素的分类结果。

具体实施方式

[0060] 为使本申请实施例的目的、技术方案和优点更加清楚，下面将结合本申请实施例中附图，对本申请实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅是本申请一部分实施例，而不是全部的实施例。通常在此处附图中描述和示出的本申请实施例的组件可以各种不同的配置来布置和设计。因此，以下对在附图中提供的本申请的实施例的详细描述并非旨在限制要求保护的本申请的范围，而是仅仅表示本申请的选定实施例。基于本申请的实施例，本领域技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本申请保护的范围。

[0061] 下面结合附图1对本发明的实施例作详细描述；

[0062] 参见附图1所示，基于多元遥感影像的多云雾区域农作物种植模式提取方法，包括以下步骤：

[0063] 数据获取：获取多时相多光谱数据以及时序SAR影像数据，并分别对多光谱数据以及SAR影像数据进行预处理，得到光谱特征和后向散射系数；

[0064] 所述多光谱数据为Sentinel‑2多光谱数据；

[0065] Sentinel‑2多光谱数据的预处理：利用IDL编程语言批量进行波段叠加、影像拼接及裁剪，得到包含7个波段的光谱特征集。

[0066] 所述SAR影像数据采用Sentinel‑1卫星数据，SAR影像数据的获取模式为IW模式的GRDH格式产品，空间分辨率为20m，每景SAR影像都具有同极化和交叉极化两种方式；

[0067] SAR影像数据的预处理：利用欧空局SNAP软件构建处理模型批量对SAR影像数据进行预处理，包括辐射定标、多普勒地形校正、斑点滤波以及地理编码。

[0068] 特征提取：基于光谱特征和后向散射系数进行计算，分别得到光谱指数特征集和SAR纹理特征集；

[0069] 所述SAR纹理特征包括均值、协同性、方差以及对比度，具体计算公式如下：

[0070] ；

[0071] ；

[0072] ；

[0073] ；

[0074] 式中：Mean为均值；Homogeneity为协同性；Variance为方差；Contrast为对比度； n为灰度值的阶数；p(i,j)是n×n的归一化共生矩阵；μ为p(i,j)的均值；i和j表示相距为d的两个灰度级。

[0075] 所述光谱指数特征包括归一化植被指数、红边植被指数以及红边叶绿素指数，具体计算公式如下：

[0076] ；

[0077] ；

[0078] ；

[0079] 式中：NDVI表示归一化植被指数；NDVI1705表示红边植被指数；CI表示红边叶绿素指数； 为近红外波段反射率；  为红波段反射率； 为中心波长为750nm波段的反射率； 为中心波长为705nm波段的反射率； 为中心波长区间为750nm
~
800nm波段的反射率； 为中心波长区间为690nm 725nm波段的反射率。
~

[0080] 阈值函数确定：基于光谱指数特征以及训练样本得到冬季作物阈值函数；基于SAR纹理特征以及训练样本得到夏季作物阈值函数；

[0081] 农作物分类：采用面向对象的分类策略，提取农作物的类型以及空间分布。

[0082] 所述农作物的分类包括以下步骤：

[0083] A.多尺度分割：利用ESP算法计算不同分割尺度下对象同质性局部变化的变化率，用于表示分割效果的最佳参数，计算公式如下：

[0084] ；

[0085] 式中： 为分割尺度为n时分割对象的局部方差， 是分割尺度为n‑1时分割对象的局部方差，当 出现峰值时，所述峰值对应的分割尺度即为最佳分割尺度；

[0086] B.面向对象分类：基于实地调查采样所获取的样本建立分类体系；通过多尺度分割提取对象层级的光谱特征，SAR纹理特征，时相物候特征，利用特征向量构建CART决策树分类器自动建立分类规则进行分类。

[0087] 种植模式提取： 种植模式即农作的类型、时空分布及轮作模式，包括年内复种的耕作模式，如“小麦‑水稻轮种”、“油菜‑水稻轮种”，“瓜菜轮种”或一年一熟单作制，如“莲藕”等。基于不同作物的物候差异，统计农作物的时间‑特征变化曲线，将所述光谱指数特征集和SAR纹理特征集进行组合，基于农作物时间‑特征变化曲线，确定农作物时间‑特征阈值判定函数；在多时相影像上提取农作物时间‑空间分布，提取农作物种植模式。

[0088] 下面以广汉市的农作物分类为例对本发明的实施方式以及本发明所达到的有益效果做详细的说明书；

[0089] 研究区域广汉市位于四川省腹心地带核心区，隶属四川省德阳市，地理位置：104.11 104.50°E，30.89 31.14°N，区域总面积548平方千米（如图2所示）；该区域为沱江冲~ ~
击平原地带，地处龙泉山脉西麓，地势西北高，东南低，地貌以平原、丘陵为主；地处亚热带湿润气候区，年平均气温16.2℃，冬季温暖多雾，夏秋两季多云多雨，降水充沛，年平均降雨量达890mm。广汉市素有“天府粮仓”之称，是成都平原重要的商品粮基地之一，水稻、小麦、油菜是该区域主要粮油作物，采用轮作方式种植，区域内主要农作物及生育周期如表1所示；

[0090]

[0091] 获取研究区域的相关图像数据并进行预处理：

[0092] Sentinel‑2数据介绍及预处理

[0093] Sentinel‑2（哨兵2号）多光谱数据来自欧洲航天局（ESA），产品等级为L2A，即已通过几何校正、辐射校正、大气校正的地表反射率产品；Sentinel‑2多光谱数据覆盖从可见光、近红外到短波红外13个光谱波段，具有不同的空间分辨率（10‑60m），在红边范围含3个波段，为细致地监测农作物生长情况提供了更多可能性；研究区地处西南区域，受多云雾气候影响，高质量光学影像获取困难，分析研究区2021年农作物生育时期（3月‑10月）Sentinel‑2多光谱数据，仅3景影像满足10%以下云量覆盖，影像获取日期见表2；

[0094]

[0095] 考虑到研究区域地块破碎，选取空间分辨率为10 m的可见光、近红外4个波段和空间分辨率为20 m的3个红边波段反射率数据，为保持空间分辨率一致性，采用双线性插值法将20 m分辨率的波段重采样为10 m分辨率，利用IDL编程语言批量进行波段叠加、影像拼接及裁剪，得到包含7个波段的研究区光学影像图（光谱特征）。

[0096] Sentinel‑1 SAR数据介绍及预处理

[0097] 研究区上空常年覆盖云雾，获取的高质量光学影像有限，难以满足农作物监测需求，由表2分析可知，仅利用获取的Sentinel‑2光学数据，难以区分水稻、玉米等夏季作物。星载SAR影像获取不受气候状况限制，为多云雾气候区域夏季农作物提取提供了可能性。

[0098] 本发明的SAR影像采用Sentinel‑1（哨兵1号）卫星数据，Sentinel‑1搭载一个基于C波段的雷达成像系统，总共有4种数据获取模式：条带模式（SM）、干涉宽幅（IW）模式、超宽幅（EW）模式和波模式（WM）；本发明所用数据的获取模式为IW模式的GRDH（地距）格式产品，空间分辨率为20 m，每景影像都具有同极化（VV极化）和交叉极化（VH极化）两种方式。结合地面调查和农作物物候期分析，选取夏季作物生长期的5景Sentinel‑1卫星影像，影像获取时间见表2。利用欧空局SNAP软件构建处理模型批量对数据进行预处理，包括辐射定标、多普勒地形校正、斑点滤波、地理编码等得到地物后向散射系数。

[0099] 地面调查数据

[0100] 地面样本调查数据来源于2021年3月‑9月农作物生长季的4次实地调查工作，在实地采样过程中，为了尽量减小混合像元和样本定位误差对分类精度的影响，采样时均选择大片均一地块，地面调查采样点分布如图2所示。地面样本数据主要记录了所到之处的经纬度信息、农作物类型、物候期和生长状况；利用地面调查数据，结合作物物候期信息，对影像进行目视解译，选取17115个样本点，包括小麦、水稻、油菜等10种类型，如表3所示；采用分层随机抽样方式选择样本集的70%作为训练样本，30%作为验证样本。研究区辣椒、白菜、苋菜、西瓜等蔬菜瓜果作物种植分布破碎，种植模式具有一定的随机性，缺乏规律的物候特征，因此将其合并为蔬菜瓜果一类。研究区粮油作物种植模式多为一年两熟的耕作模式，主要以“小麦‑水稻轮种”、“油菜‑水稻轮种”、“土豆‑玉米轮种”、“土豆‑大豆轮种”为主。

[0101]

[0102] SAR纹理特征提取

[0103] 灰度共生矩阵（GLCM：Grey‑Level Co‑occurrence Matrix）是一个统计描述影像中的一个局部区域或整个区域相邻像元或一定距离内两像元灰度呈现的某种关系的矩阵，在矩阵中，元素值代表大小为m的基准窗口内方向为θ，相距为d的两个灰度级i和j同时出现的联合条件概率密度。本发明基于农作物生长期内时序Sentinel‑1影像提取均值（Mean）、协同性（Homogeneity）、方差（Variance）、对比度（Contrast）4种纹理特征（计算公式见式1‑4），式中：n为灰度值的阶数， 是 的归一化共生矩阵， 为 的均值；

[0104]   ；

[0105] ；

[0106] ；

[0107] ；

[0108] 光谱指数特征提取

[0109] 归一化植被指数（NDVI）是目前应用最为广泛的植被指数，基于覆盖作物生长期的NDVI时序数据能够反映作物从播种、生长到收割的变化过程，可利用NDVI的差异对不同作物进行识别。Sentinel‑2是光学遥感卫星中唯一一个在红边范围含有三个波段数据的卫星，为农作物精细化提取与监测提供了可能性。本文结合Sentinel‑2影像的波段特点和研究区主要农作物物候特征，选取了2种红边指数进行农作物信息提取，包括红边植被指数（NDVI705）、红边叶绿素指数（CI）；NDVI705、CI指数计算所需用到的中心波长为750nm和708.75nm等波段的反射率在Sentinel‑2中无法获得，因此，采用最邻近波段处的反射率替代；计算公式如下所示：

[0110] ；

[0111] ；

[0112] ；

[0113] 式中：NDVI表示归一化植被指数；NDVI1705表示红边植被指数；CI表示红边叶绿素指数。

[0114] 阈值函数确定：基于光谱指数特征以及训练样本得到冬季作物阈值函数；基于SAR纹理特征以及训练样本得到夏季作物阈值函数；

[0115] 面向对象决策树分类

[0116] 传统的分类方法主要是基于像元的数理统计法，在进行地物分类时考虑的主要是像元的灰度值，对影像空间、纹理信息利用不充分，容易产生“同谱异物”“同物异谱”现象，分类精度较低。面向对象的分类思想，是以对象为分类的基本单元，基于光谱信息、形状信息、纹理信息等特征，将同质性较强的相邻像素合并为一个对象，再对每个对象进行类别划分，这一分类方法明显区别于传统的在像素层级上进行分类的方法，分类精度有明显提高。在面向对象分类的过程中，图像分割和决策树分类是最关键的两个步骤。

[0117] A. 多尺度分割与最优分割尺度选择

[0118] 多尺度分割是面向对象分类的前提，分割尺度决定了分割后对象的细碎程度以及后续分类的精度。本文利用ESP（Estimation of Scale Parameter）算法计算不同分割尺度下对象同质性局部变化（local variance, LV）的变化率（rates of change of LV, ROC_LV），用于表示分割效果的最佳参数，如图3所示。ROC_LV计算公式如式5：

[0119]  ；

[0120] 其中， 为分割尺度为n时分割对象的局部方差， 是分割尺度为n‑1时分割对象的局部方差，当 出现峰值时，所述峰值对应的分割尺度即为最佳分割尺度。

[0121] 由图3可知，变化率存在多个峰值，这是由于不同地物有不同的最优分割尺度。通过调整不同分割尺度进行对比分析，最终确定采用2个尺度进行多尺度分割，多尺度分割效果放大如图4所示；从图中所画的方框处可以看出，对于水体、裸地等均质区域（在图4的a部分视图的基础上，对比图4的b部分视图和c部分视图），在80的分割尺度下，容易出现“过分割”现象，在118的分割尺度下，分割完整性较好；研究区作物种植较为细碎，对于耕地区域（在图4的d部分视图的基础上，对比图4的e部分视图和f部分视图），在118的分割尺度下，出现大量“欠分割”现象，不利于后期农作物信息准确提取，在80的分割尺度下，能识别出地块之间的细微差异，较好的把不同类型的地块分割开，对于后续农作物精细提取有较大帮助。因此本实施例首先在118的分割尺度下提取裸地、水体等大片均质非农作物信息，然后在80的分割尺度下对不同作物进行识别与提取。

[0122] 本实施例基于实地调查采集的样本建立分类体系，分类样本主要包含农作物及非农作物两大类，其中农作物样本包括水稻、玉米、大豆等夏季作物，冬小麦、冬油菜、土豆等冬季作物，以及蔬菜瓜果等作物；非农作物样本包括果树、裸地、水体。通过面向对象多尺度分割提取对象层级的光谱特征，SAR纹理特征，时相物候特征，利用特征向量构建CART决策树分类器自动建立分类规则进行分类。

[0123] 农作物种植模式提取

[0124] 基于多时相Sentinel‑2多光谱影像提取植被指数特征，通过训练样本获取各类别NDVI、NDVI705、CI三种植被指数的均值，并绘制时间‑植被指数特征曲线；参见图5的a、b和c三部分可知，各类农作物在3月27日植被指数差异较为明显，且该时期多光谱影像质量较好，因此实施例选择3月27日Sentinel‑2多光谱影像作为典型时相多光谱数据，进行光谱特征分析，如图6所示。

[0125] 由图5的a、b和c三部分可知，裸地、水体的植被指数变化较小，裸地植被指数趋近于0，水体小于0,利用三个时相的植被指数特征能将水体、裸地与耕地等区分开。研究区果树以桃树、李树、柚子树为主，3月至8月是三类果树生长期，植被指数逐渐升高，在8月达到最大值，而双熟制耕地植被指数均为先下降后上升，两者上升趋势差异明显，基于多时相植被指数特征可以将果树与双熟制耕地进行有效区分。

[0126]

[0127] 基于表4提取农作物的种植模式，参见图5的a、b和c三部分所示，双熟制耕地植被指数特征均为先下降后上升，在5月出现波谷，由研究区农作物物候期可知，4月中下旬至5月上旬是冬季作物成熟收获与夏季作物播种发芽的交替时期，随着夏季作物的生长，5月以后植被指数快速上升，在8月达到峰值。对于水稻、大豆、玉米三类夏季作物（生长期4月‑9月），水稻的植被指数上升速度和数值明显不同，8月份水稻植被指数高于大豆、玉米，该时间点3种植被指数均可作为水稻识别的重要特征；玉米、大豆生育期接近，5月和8月的各类植被指数（3月尚未播种，不做分析）差异较小，仅依靠植被指数特征难以将二者区分开。对于小麦、油菜、土豆三类冬季作物（生长期10月‑次年5月），区分的关键时间节点为3月下旬，此时，小麦处于拔节期，生长较快，植被指数明显高于其他类别，容易识别；油菜处于盛花期，土豆处于结薯期，两者植被指数相近，但近红外及红边波段反射率差异明显（如图6），可利用光谱特征将其有效区分；蔬菜瓜果作物3月下旬气温较低，为提高土壤温度，蔬菜瓜果种植时多覆有地膜或大棚，植被指数较低，与小麦、油菜、土豆差异明显，识别度较高。

[0128] 综上所述，基于获取的多时相多光谱数据可对研究区非植被地表、冬季作物、果树进行有效识别，但由于缺乏夏季作物生长季关键时期（6月、7月、9月）的光谱影像，大豆、玉米不能被有效区分与识别。

[0129] 基于覆盖夏季作物（水稻、大豆、玉米）完整生长季的SAR影像，分析不同极化条件下，水稻、玉米、大豆的后向散射系数差异。如图7中的a部分和b部分所示，在VV极化和VH极化下，水稻都表现出较低的后向散射信号，分析其原因为水稻田需要经常灌水，土壤表面含水量较多，易发生镜面反射，使得后向散射信号较弱，后向散射系数小，VV极化下，水稻后向散射系数表现出更明显的差异，区分度更高，对水稻提取更有效，而在两种极化条件下，玉米、大豆的后向散射系数差异均较小，难以区分。

[0130] 本实施例基于灰度共生矩阵的统计分析方法获取SAR纹理特征参数（均值、协同性、方差、对比度）来扩大玉米、大豆在SAR影像上的差异性，以期提高玉米、大豆的可分离性，结果如图8所示。

[0131] 参见图8中的a部分和b部分视图，协同性：三种农作物VV极化和VH极化下的协同性曲线差异很小，农作物图像纹理均匀，难以均分玉米、大豆。

[0132] 参见图8中的c部分和d部分视图，对比度：VH极化水稻和大豆对比度曲线基本重合与玉米差异明显，有利于玉米提取，最佳提取时相为5月5日至7月16日。

[0133] 参见图8中的e部分和f部分视图，均值：VH极化下，水稻和大豆均值曲线相似，玉米纹理均值在7月16日、8月21日与水稻、大豆差异明显，根据研究区农作物物候节律可知玉米在7月下旬逐步成熟，先于水稻和大豆开始收获，此时土壤的散射逐渐占主导作用，与水稻、大豆冠层散射占主导作用的纹理均值差异较大，至8月下旬水稻、大豆成熟后差异逐渐缩小，相对于VV极化，VH极化纹理均值更有利于识别玉米，最佳识别时相为7月16日至8月21日。

[0134] 参见图8中的g部分和h部分视图，方差：VV极化和VH极化下，水稻和大豆纹理方差差异很小，玉米的方差整体大于水稻和大豆，说明玉米纹理强度表现更为相对分散，变化较大，最佳识别时相为7月16日至9月2日。

[0135] 综上所述，本实施例选择VV极化后向散射系数，结合VH极化对比度、均值、方差时序纹理特征，构建决策树完成夏季作物的提取。

[0136] 结果分析

[0137] 在本实施例中分别采用面向对象和面向像素两种尺度的决策树分类方法对研究区农作物进行提取，以验证不同尺度的分类对象对分类结果的影响；从图9的a和b部分可以看到：

[0138] （1）裸地主要分布在中部平原地带，少量零星分布在田间及东部山谷间，方案一提取的裸地整体性更好，精度更高（如表5所示），方案二部分裸地（阴影区域）错分为水体。

[0139]

[0140] （2）果树主要分布于东部丘陵区域，方案一果树分类精度为86.95%远高于方案二的提取精度，方案二虽然能够提取大片果树区域，但“孤岛”现象较为普遍。

[0141] （3）研究区蔬菜瓜果一类包括黄瓜、辣椒、白菜、草莓等多种作物，种植计划随气温、降雨量等变化会随时调整，作物物候特征变化存在一定的不确定性，光谱特征变化存在随机性，且部分作物覆有地膜或大棚进行保温，光谱特征与裸地相似，因此两种方案中蔬菜瓜果分类精度均比较低，存在大量混分、错分现象。

[0142] （4）四种双熟制耕地中，小麦‑水稻、油菜‑水稻提取精度最高，这与其种植模式稳定，光谱特征明显有关；研究区玉米、大豆存在大量带状间作套种情况，两个作物带间距小于2米，Sentinel‑1空间分辨率为20米，虽重采样至10米，但仍然存在大量混合像元，对区分玉米、大豆影响较大，因此土豆‑大豆、土豆‑玉米的提取精度较低。

[0143] 综上所述，面向对象的分类对不同作物的分类精度均优于面向像素的分类方法，且保持了大面积作物的完整性，减少破碎，避免了分类结果的“椒盐现象”。

[0144] 基于上述可知多云雾区域难以获取覆盖作物完整生长季的全时序光学影像，农作物遥感监测受到制约。针对上述问题，本文基于多元遥感影像，构建典型农作物光谱特征、SAR时序特征，并以广汉市主要农作物种植模式遥感提取为实例，探讨光学与雷达遥感数据协作提取多云雾地区农作物种植模式的可行性，得到如下结论：

[0145] （1）广汉市农作物以水稻、小麦、油菜等粮油作物为主，其中水稻种植面积最大，广泛分布在平原各地；农作物种植模式主要为小麦‑水稻、油菜‑水稻、土豆‑大豆、土豆‑玉米4种稳定的年内复种模式，其中大豆、玉米多为带状间作套种模式；蔬菜瓜果随气候、降水等外部环境变化呈现无固定规律年内复种模式。

[0146] （2）SAR影像包含丰富的空间纹理信息，且不受天气条件限制，覆盖作物完整生育期的SAR时序后向散射特征及纹理特征，能够反映出广汉市夏季作物生长发育过程，满足夏季农作物分类要求，验证了Sentinel‑1 SAR时序影像辅助用于农作物提取的技术可行性。

[0147] （3）基于典型时相多光谱数据和SAR时序影像，结合面向对象的决策树分类方法，能更为准确地提取广汉市农作物种植模式信息，总体精度达85.49%，与传统的面向像素分类结果相比，明显改善了作物的分类精度；面向对象的分类方法极大程度减少了“同谱异物”和“同物异谱”现象，同时避免了分类结果出现“椒盐”现象，具有面向像素分类不具备的优势。

[0148] 基于本实施例，基于光学与雷达遥感数据协作的多云雾区域农作物种植模式提取在国内尚属首次，为大区域农作物监测及非农非粮化监测提供了新的思路。

[0149] 以上所述实施例仅表达了本申请的具体实施方式，其描述较为具体和详细，但并不能因此而理解为对本申请保护范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本申请技术方案构思的前提下，还可以做出若干变形和改进，这些都属于本申请的保护范围。