CN113484245B - 一种丘陵山区水田种植模式遥感快速监测方法、系统及计算机可读存储介质 - Google Patents

Abstract

本发明公开一种丘陵山区水田种植模式遥感快速监测方法、系统及存储介质，包括：步骤1，对耕作初期获取的多光谱光学遥感影像进行预处理，包括进行辐射定标、大气校正和地形辐射校正；步骤2，根据定量化遥感影像提取初始水体，并对初始水体进行非监督分类，然后统计各类别的光谱均值和方差，最后基于光谱均值和方差建立地表水体指数反演模型，计算地表水体指数影像SWI，监测地表是否含有水体；步骤3，叠加土地利用现状数据和地表水体遥感指数影像SWI，统计水田图斑范围内像元的地表水体指数大小，判别水田图斑的种植模式。本发明具有人工参与少、自动化程度高，提取结果准确性高，可快速区分水田的种植模式的优点。

Description

一种丘陵山区水田种植模式遥感快速监测方法、系统及计算 机可读存储介质

技术领域

本发明属于遥感图像处理以及自然资源管理技术领域，具体涉及丘陵山区水田种植模式的监测技术。

背景技术

水田是土地利用现状分类体系中耕地类别下的二级类，是指筑有田埂(坎)，可以蓄水，用于种植水稻、莲藕、荸荠、莼菜等水生作物的土地。丘陵山区的水田大多为梯田，主要种植水稻。丘陵山区地形地貌起伏破碎，水田面积小且分散，难以进行规模化、机械化种植水稻，因此在实际耕作中，因种植经济效益低、农村劳动力不足等因素，较多水田没有种植水稻，而改为旱作(种植旱生作物，如玉米)，或撂荒等，影响国家主粮安全。此外，丘陵山区的水田若长期旱作，可能导致田坎破坏、区域地下水位下降等，从而影响水田的蓄水能力，导致耕地的等级下降。因此，水田的种植模式(水作与旱作)监测，是丘陵山区耕地保护工作中的重要内容。

目前，丘陵山区水田种植模式监测主要有野外调查法和基于高分辨率遥感影像的人工解译法。这两种方法具有各自的不足：

野外实地调查方法是指调查人员携带定位设备到野外实地查看水田的种植模式，并记录相关信息。这种方法需要投入大量人力和财力，难以应用于大范围的监测，一般只能选择少量样本进行抽样监测。

基于高分辨率遥感影像的人工解译法是采用一期或多期高分辨率遥感影像，技术人员在专业软件中叠加遥感影像和土地利用现状数据库的水田范围，根据遥感影像的光谱和纹理信息对水田的种植模式进行遥感人工解译。这种方法的效率虽比外业调查高，可用于丘陵山区大范围水田种植模式监测，但仍存在耗费时间长、效率低的问题，同时遥感解译结果的正确性依赖于解译人员的经验，解译人员的经验和主观性差异必然影响结果精度。

中国专利文献CN112507763A公开了一种基于多源多光谱遥感图像水体提取方法，该方法包括：获取多源遥感图像，并分别对遥感对象进行尺寸预处理；对预处理后的遥感图像进行图像大气校正；对大气校正后的遥感图像进行归一化水体指数NDWI计算；消除归一化水体指数NDWI的异常值，并绘制NDWI图像灰度直方图，根据NDWI图像灰度直方图峰值分布确定水体分割阈值；根据归一化水体指数NDWI与水体分割阈值的预设关系判断遥感图像中的地物是否为水体。但该方法所述的归一化水体指数NDWI仅用了多光谱遥感数据中的两个波段，用于识别河湖库塘等较深的水体较为有效，但难以监测水田中的浅水体。

发明内容

针对现有方法效率低、成本高、主观性大等问题，本发明的目的是提供一种基于多光谱光学遥感数据的丘陵山区水田种植模式定量遥感监测方法、系统及计算机可读存储介质，减少人工参与程度、提高自动化程度，提高提取结果的准确性，快速区分水田的种植模式，为丘陵山区耕地保护工作提供技术支撑。

本发明方法是采用耕作初期获取的多光谱光学遥感数据和土地利用现状数据，建立的一种基于多光谱光学遥感数据的丘陵山区水田种植模式定量遥感监测方法，利用光学遥感影像获得多光谱光学遥感数据，其遥感影像是通过航天(卫星)或航空(飞机)成像方式，拍摄的地表的影像(照片)，其实质是获取地面目标的反射、散射或辐射的电磁波信息。光学遥感影像是指获取可见光到红外线这部分电磁波信息进行成像，得到的影像。

本发明的技术方案如下：

一种丘陵山区水田种植模式遥感快速监测方法，所述方法包括如下步骤：

步骤1，多光谱光学遥感影像定量化预处理

对耕作初期获取的多光谱光学遥感影像进行预处理，包括进行辐射定标、大气校正和地形辐射校正，获得地表反射率遥感影像SR，使得遥感影像能定量反映地表真实信息。

步骤2，地表水体指数反演

根据定量化遥感影像提取初始水体，并对初始水体进行非监督分类，然后统计各类别的光谱均值和方差，最后基于光谱均值和方差建立地表水体指数反演模型，然后计算地表水体指数影像SWI，根据地表水体指数监测地表是否含有水体。

步骤3，水田图斑种植模式的自动判别

叠加土地利用现状数据和地表水体遥感指数影像，统计水田图斑范围内像元的地表水体指数大小，利用水体指数影像SWI对土地利用现状数据中的水田图斑进行计算，判别水田图斑的种植模式。

进一步，所述步骤1中，辐射定标是利用多光谱光学遥感影像原始数据的像元亮度值，结合遥感影像元数据中的增益和偏差系数转换为辐亮度数据L，然后利用公式计算星上反射率

式中d为日地距离，L为辐亮度数据，单位是w/m²/sr/μm，E₀为太阳常数，β为成像时太阳天项角。

进一步，所述步骤1中，大气校正是采用辐射传输模型法，公式如下

式中SR为地表反射率，TOA为星上反射率，E₀为太阳常数，L_path为大气的路径辐射，β、v分别为太阳天顶角和卫星观测天顶角，T(β)、T(v)为分别为大气下行辐射总透过率和上行辐射总透过率，S是大气底层向下的半球反照率。

进一步，所述步骤1中，地形辐射校正是采用C辐射校正模型结合DEM数据来消除地形的影响，校正公式为：

c＝a/b (6)

其中，SR_T是倾斜地表(即校正前)的反射率，a和b为拟合方程的系数，i为太阳入射角，β为太阳天顶角。

进一步，所述步骤2的地表水体指数反演具体包括如下步骤：

步骤2.1，初始水体识别：采用归一法或人工交互解译方法提取初始水体。

步骤2.2，初始水体非监督分类：将提取的所有初始水体像元利用k-均值算法、核k-均值算法或迭代自组织数据分析技术(ISODATA)等进行非监督分类。

步骤2.3，初始水体参数统计：结合初始水体非监督分类结果和地表反射率影像，逐波段、逐类别统计初始水体的光谱反射率均值Mean_i，j和均方差STD_i，j，其中i为不同的波段，j为不同的类别。

步骤2.4，计算地表水体指数影像SWI

式中i为不同的波段，j为不同的类别，SR_i为波段i的地表反射率值，m为波段数，n为类别数。

更进一步，所述步骤2.1，采用归一法提取初始水体具体如下：

基于地表反射率影像，逐像元计算归一化水体指数MNDWI，红、绿、近红外三波段最小值Min_RGN，红、绿、蓝、近红外和短波红外最大值Max_Val，公式如下

Min_RGN＝min{SR_R，SR_G，SR_NIR} (8)

Max_Val＝max{SR_R，SR_G，SR_B，SR_NIR，SR_SWIR} (9)

式中SR_R、SR_a、SR_B、SR_NIR、SR_SWIR分别表示红、绿、蓝、近红外和短波红外的反射率，然后提取同时满足以下条件的像元为初始水体

MNDWI＜0.3且Min_RGN＜0.15且Max_Val＜0.2 (10)

进一步，所述步骤3，水田图斑种植模式的自动判别可以采用赋值计算，根据阈值实现水田图斑种植模式的自动判别，或者采用四分位数、均值等判别条件，实现水田图斑种植模式的自动判别。

更进一步，所述步骤3，所述赋值计算方法是叠加地表水体指数影像SWI和土地利用现状数据中的水田图斑，逐图斑统计图斑范围内的SWI像元值，如该图斑范围内90％的SWI像元值都大于1，则该水田图斑耕地模式为水作，否则为旱作。

与现有技术相比，本发明方法建立了一种水体指数的定量遥感模型，然后通过阈值自动判别水田图斑的种植模式，大大提高了效率，降低了成本，该方法人工参与少、自动化程度高，提取结果准确性高，可快速区分水田的种植模式，为丘陵山区耕地保护工作提供技术支撑。

附图说明

图1是本方法的流程框图。

具体实施方式

以下结合附图进一步详细说明本发明的技术细节：

本发明所述方法从流程上主要分为三大部分：

第一部分是对水稻耕作初期获取的多光谱(包含可见光近红外和短波红外波段)光学遥感影像进行预处理，获得地表反射率遥感影像SR，使得遥感影像能定量反映地表真实信息。在中国来说，耕作初期指的是4-5月份，而其它国家和地区时间会有所不同。这里的的多光谱包括可见光近红外和短波红外光。

第二部分是根据定量化遥感影像提取初始水体，并根据初始水体光谱统计参数建立水体指数反演模型，然后计算地表水体指数影像SWI。

第三部分是利用水体指数影像SWI对土地利用现状数据中的水田图斑进行赋值计算，并根据阈值实现水田图斑种植模式的自动判别。

参见图1，对本方法的具体步骤描述如下：

1多光谱光学遥感影像定量化预处理

光学遥感成像时，由于大气层的存在，大气层中所含大气分子、尘埃及水汽等对太阳辐射的吸收、散射等作用，以及丘陵山区地形起伏引起阴坡、阳坡的太阳光照强度不同，使得传感器的测量值与目标的光谱反射率或光谱辐射亮度等物理量不一致，所获得的影像数据不能真实反映地物的反射或发射光谱特征。因此，通过光学遥感影像获得地表的定量信息，必须消除上述影响，对数据进行定量化处理得到地表反射率。

定量化预处理按以下步骤进行：

①辐射定标

多光谱光学遥感影像原始数据为像元亮度值，结合遥感影像元数据中的增益和偏差系数即可转换为辐亮度数据L。然后，利用公式计算星上反射率

式中d为日地距离，L为辐亮度数据(单位：w/m²/sr/μm)，E₀为太阳常数，β为成像时太阳天顶角。

②大气校正

大气校正，就是采用一定的算法，消除由于大气效应所引起的遥感像元值偏差。本实施例中，大气校正采用辐射传输模型法，即利用电磁波在大气中辐射传输的原理建立辐射传输方程，输入方程参数，对方程求解，从而消除大气影响，公式如下

式中SR为地表反射率，TOA为星上反射率，E₀为太阳常数，L_path为大气的路径辐射，β、v分别为太阳天顶角和卫星观测天顶角，T(β)、T(v)为分别为大气下行辐射总透过率和上行辐射总透过率，S是大气底层向下的半球反照率。大气校正可以采用通用遥感影像处理软件或专门的大气校正软件进行处理。

③地形辐射校正

地形辐射校正的主要方法是结合DEM数据来消除地形的影响，本实施例采用C辐射校正模型。对于山区地表，影像像元对应的太阳入射角为太阳直射光与坡面法线的夹角，其反射率和太阳入射角余弦值之间的关系为：

SR_T＝a+b×cosi (3)

cosi＝cosθcosβ+sinθsinβcos(λ-ω) (4)

式中a和b为拟合方程的系数，i为太阳入射角，β和λ分别为太阳天顶角，λ为太阳方位角，θ和ω分别为地表坡度角和坡向角。利用不同地形条件下的像元进行回归分析，可得到系数a和b，然后将倾斜地表的反射率投影到对应的水平地表反射率，校正公式为：

c＝a/b (6)。

2地表水体指数反演

光学遥感影像经过定量化预处理得到地表反射率影像，其像元值可以定量反演地表参数特征，包括植被覆盖、土壤水分等。本发明建立了一种地表水体指数反演方法，步骤如下

①初始水体识别

可以采用采用归一法或人工交互解译方法等来提取初始水体，本实施例采用了归一法。

基于地表反射率影像，逐像元计算归一化水体指数MNDWI，红、绿、近红外三波段最小值Min_RGN，红、绿、蓝、近红外和短波红外最大值Max_Val，公式如下

Min_RGN＝min{SR_R，SR_G，SR_NIR} (8)

Max_Val＝max{SR_R，SR_G，SR_B，SR_NIR，SR_SWIR} (9)

式中SR_R、SR_G、SR_B、SR_NIR、SR_SWIR分别表示红、绿、蓝、近红外和短波红外的反射率，然后提取同时满足以下条件的像元为初始水体

MNDWI＜0.3且Min_RGN＜0.15且Max_Val＜0.2 (10)。

②初始水体非监督分类

将上一步提取的所有初始水体像元进行k-均值非监督分类，类别数设置为10，最大迭代次数为20，变化阈值为0.05。

当然，也可以采用其它方法如核k-均值聚类算法或迭代自组织数据分析技术(ISODATA)进行非监督分类。

③初始水体参数统计

结合初始水体非监督分类结果和地表反射率影像，逐波段、逐类别统计初始水体的光谱反射率均值Mean_i，j和均方差STD_i，j，其中i为不同的波段，j为不同的类别。

④计算地表水体指数影像SWI

式中i为不同的波段，j为不同的类别，SR_i为波段i的地表反射率值，m为波段数，n为类别数。

3水田图斑种植模式的自动判别

结合土地利用现状数据和地表水体定量遥感监测指数，叠加地表水体指数影像SWI和土地利用现状数据中的水田图斑，逐图斑统计图斑范围内的SWI像元值，根据阈值自动判别水田种植模式，如该图斑范围内90％的SWI像元值都大于1，则该水田图斑耕地模式为水作，否则为旱作。

另一方面，本发明实施例还提供一种丘陵山区水田种植模式遥感快速监测系统，其包括：

预处理单元，用于对耕作初期获取的多光谱光学遥感影像进行定量化预处理，包括进行辐射定标、大气校正和地形辐射校正，获得地表反射率遥感影像SR，使得遥感影像能定量反映地表真实信息。

反演单元，用于根据定量化遥感影像提取初始水体，并对初始水体进行非监督分类，然后统计各类别的光谱均值和方差，最后基于光谱均值和方差建立地表水体指数反演模型，计算地表水体指数影像SWI，根据地表水体指数监测地表是否含有水体。

自动判别单元，用于叠加土地利用现状数据和地表水体遥感指数影像，统计水田图斑范围内像元的地表水体指数大小，利用水体指数影像SWI对土地利用现状数据中的水田图斑进行计算，判别水田图斑的种植模式。

本发明进一步的实施例还提供一种计算机可读存储介质，所述计算机可读存储介质中包括丘陵山区水田种植模式遥感快速监测方法程序，所述丘陵山区水田种植模式遥感快速监测方法程序被处理器执行时，实现以上实施例所述的丘陵山区水田种植模式遥感快速监测方法的步骤。

以上实施例仅仅是示意性的，仅为说明本发明的技术思想，不能以此限定本发明的保护范围，本领域技术人员应理解，按本发明提出的技术思想，在技术方案基础上所做的任何改动，均应当落入本发明保护范围之内。

Claims (6)

1.一种丘陵山区水田种植模式遥感快速监测方法，其特征在于，所述方法包括如下步骤：

步骤1，多光谱光学遥感影像定量化预处理

对耕作初期获取的多光谱光学遥感影像进行预处理，包括进行辐射定标、大气校正和地形辐射校正，获得地表反射率遥感影像SR，使得遥感影像能定量反映地表真实信息；

步骤2，地表水体指数反演

根据定量化遥感影像提取初始水体，并对初始水体进行非监督分类，然后统计各类别的光谱均值和方差，最后基于光谱均值和方差建立地表水体指数反演模型，并应用于地表反射率影像SR计算地表水体指数影像SWI，根据地表水体指数监测地表是否含有水体；

具体包括如下步骤：

步骤2.1，初始水体识别：采用归一法或人工交互解译方法提取初始水体；

基于地表反射率影像，逐像元计算归一化水体指数MNDWI，红、绿、近红外三波段最小值Min_RGN，红、绿、蓝、近红外和短波红外最大值Max_Val，公式如下

Min_RGN＝min{SR_R,SR_G,SR_NIR} (8)

Max_Val＝max{SR_R,SR_G,SR_B,SR_NIR,SR_SWIR} (9)

式中SR_R、SR_G、SR_B、SR_NIR、SR_SWIR分别表示红、绿、蓝、近红外和短波红外的反射率，然后提取同时满足以下条件的像元为初始水体

MNDWI<0.3且Min_RGN<0.15且Max_Val<0.2 (10)；

步骤2.2，初始水体非监督分类：将提取的所有初始水体像元利用k—均值算法、核k—均值算法或迭代自组织数据分析技术进行非监督分类；

步骤2.3，初始水体参数统计：结合初始水体非监督分类结果和地表反射率影像，逐波段、逐类别统计初始水体的光谱反射率均值Mean_i,j和均方差STD_i,j，其中i为不同的波段，j为不同的类别；

步骤2.4，计算地表水体指数影像SWI

式中i为不同的波段，j为不同的类别，SR_i为波段i的地表反射率值，m为波段数，n为类别数；

步骤3，水田图斑种植模式的自动判别

叠加土地利用现状数据和地表水体遥感指数影像，统计水田图斑范围内像元的地表水体指数大小，利用水体指数影像SWI对土地利用现状数据中的水田图斑进行计算，判别水田图斑种植模式；

水田图斑种植模式的自动判别采用赋值计算，根据阈值实现水田图斑种植模式的自动判别，或者采用四分位数、均值判别条件，实现水田图斑种植模式的自动判别；

所述赋值计算是叠加地表水体指数影像SWI和土地利用现状数据中的水田图斑，逐图斑统计图斑范围内的SWI像元值，如该图斑范围内90％的SWI像元值都大于1，则该水田图斑耕地模式为水作，否则为旱作。

2.根据权利要求1所述的丘陵山区水田种植模式遥感快速监测方法，其特征在于，所述步骤1中，辐射定标是利用多光谱光学遥感影像原始数据的像元亮度值，结合遥感影像元数据中的增益和偏差系数转换为辐亮度数据L，然后利用公式计算星上反射率

式中d为日地距离，L为辐亮度数据，单位是w/m²/sr/μm，E₀为太阳常数，β为成像时太阳天顶角。

3.根据权利要求1所述的丘陵山区水田种植模式遥感快速监测方法，其特征在于，所述步骤1中，大气校正是采用辐射传输模型法，公式如下

式中SR为地表反射率，TOA为星上反射率，E₀为太阳常数，L_path为大气的路径辐射，β、ν分别为太阳天顶角和卫星观测天顶角，T(β)、T(ν)为分别为大气下行辐射总透过率和上行辐射总透过率，S是大气底层向下的半球反照率。

4.根据权利要求1所述的丘陵山区水田种植模式遥感快速监测方法，其特征在于，所述步骤1中，地形辐射校正是采用C辐射校正模型结合DEM数据来消除地形的影响，校正公式为：

C＝a/b (6)

其中，SR_T是倾斜地表即校正前的反射率，a和b为拟合方程的系数，i为太阳入射角，β为太阳天顶角。

5.一种实现权利要求1-4任一项所述丘陵山区水田种植模式遥感快速监测方法的系统，其特征在于，包括：

预处理单元，对耕作初期获取的多光谱光学遥感影像进行定量化预处理，包括进行辐射定标、大气校正和地形辐射校正，获得地表反射率遥感影像SR，使得遥感影像能定量反映地表真实信息；

反演单元，根据定量化遥感影像提取初始水体，并对初始水体进行非监督分类，然后统计各类别的光谱均值和方差，最后基于光谱均值和方差建立地表水体指数反演模型，计算地表水体指数影像SWI，根据地表水体指数监测地表是否含有水体；

自动判别单元，叠加土地利用现状数据和地表水体遥感指数影像，统计水田图斑范围内像元的地表水体指数大小，利用水体指数影像SWI对土地利用现状数据中的水田图斑进行计算，判别水田图斑的种植模式。

6.一种计算机可读存储介质，其特征在于，所述计算机可读存储介质中包括丘陵山区水田种植模式遥感快速监测方法程序，所述丘陵山区水田种植模式遥感快速监测方法程序被处理器执行时，实现如权利要求1至4中任一项所述的丘陵山区水田种植模式遥感快速监测方法的步骤。

Images