一种基于地形和遥感影像融合技术提取冲沟方法

摘要

本发明提供一种基于地形和遥感影像融合技术提取冲沟方法，该方法是基于遥感RGB（红、绿、蓝）色彩空间同HIS（色相、色彩亮度、饱和度）色彩空间之间的变换，以线性标准化地表粗糙度

说明书

技术领域： 

本发明涉及遥感技术领域，具体地讲是一种将地形信息和遥感影像融合技术应用于冲沟沟缘线的提取。该发明可用于基于光学遥感图像，进行冲沟解译、沟蚀研究及制图场合。

背景技术： 

冲沟侵蚀是一种重要的土壤侵蚀方式，它不但是导致土地退化的一个重要过程，还是江河泥沙的一个主要来源。目前国内外学者主要利用地面手工测量法（Casali et al., 2006; Vandekerckhove et al., 2001；）、侵蚀针(桩)监测法（Ionita et al.，2006；Martinez et al.2003）、RTK-GPS测量法（何福红 et al.,2005;YougqiuWu et al.,2005）、航片和数字高程模型（Martinez et al.，2003；Harley et al.,1999）等方法来定量研究沟蚀的发生演变规律。自20世纪60年代初航空摄影技术出现以来，遥感技术已经被广泛应用于土壤侵蚀、冲沟沟蚀研究等领域。利用航天遥感图像解译土壤侵蚀的研究可以追溯到1940年（Smith et al., 1943）。随后的数十年里，诸多学者利用航空照片及摄影测量技术对冲沟侵蚀进行研究（Nachtergaele and Poesen, 1999; Betts and DeRose, 1999；Martínez-Casasnovas et al., 2003）。研究冲沟的方式是对遥感图像，目视解译冲沟沟缘线，由于沟缘线是一条重要的坡度分界线，其上部为坡度较为平缓的沟间地类型，其下部是坡度较为陡峻的沟谷地类型，而遥感二维图像提供了的三维地形信息量不足以准确解译出冲沟沟缘线，因而精度较低。另外，对北半球来说，由于资源卫星大多为太阳同步极轨卫星，其成像时间为地方时十点半。成像时，阳光从东南向射入，山脊两侧的南向坡形成光照面，北向坡为阴影面，而传统的上北下南的影像构图方式，使得阴影面位于光照面上方，这就形成了视觉上的反立体现象（Saraf et al., 1996; Rudnicki, 2000; Patterson, 2004，BoWu et al,2012）。反立体现象，使得遥感图像中的冲沟在视觉上表现为凸起山脊，而山脊则表现为凹陷的冲沟，增加了准确判读冲沟的难度。因而单纯基于遥感图像解译冲沟信息效果有限。

数字高程模型（DEM）是地表三维地形的数字表达方式，是用一组有序数值阵列形式表示地 面高程的一种实体地面模型能，是计算诸如坡度、坡向、坡度变化率等因子的数学模型。但利用DEM数据难以精确解译冲沟沟缘线。目前，基于DEM数据计算机自动识别冲沟沟缘线的方法，但由于数据源、地貌类型的多样性及沟缘线的多级性和层次性均导致自动识别沟缘线算法精度和通用性差。 

因此，将符合人类视觉习惯的三维地形信息融合到遥感图像中，将增加遥感二维图像的三维地形信息，且使得图像空间几何关系及纹理特征符合人类视觉习惯，有助于精确识别冲沟；并辅以由DEM数据计算的冲沟沟缘线坡度阈值数据，将有助于提高解译冲沟沟缘线的精度。目前三维地形信息同遥感图像相融合的方法主要有图像南北向旋转法、像元值逆转法、SRM（Shaded Relief Model）参与的HIS融合3类。像元值逆转法和SRM参与的HIS融合法通过修改图像像素灰度值来实现正立体校正，但均存在很明显的光谱信息损失。图像南北向旋转法尽管不存在光谱信息损失问题，但构图方式采用上南下北构图方式，导致较难推广使用。 

发明内容： 

本发明的目的是克服上述已有技术的不足，而提供一种基于地形和遥感影像融合技术提取冲沟方法，主要解决SRM参与的HIS融合法方法有明显的遥感图像光谱信息损失的问题。

本发明的技术方案是：一种基于地形和遥感影像融合技术提取冲沟的方法，其特殊之处在于，包括以下步骤： 

（１）将图幅相同的光学遥感红、绿、蓝三波段（RGB三波段）图像和DEM（数字高程模型）数据的坐标系统一成高斯-克吕格平面直角坐标系，并进行重采样，使两者具有相同的空间分辨率；

（２）将光学遥感红、绿、蓝三波段图像（RGB色彩空间），变换至彩色亮度（Intensity）、色相（Hue）和饱和度（Saturation）空间（HIS色彩空间），获得彩色亮度 、色相和饱和度三分量数据；

（３）根据遥感图像成像时刻（年月日时分秒）及图像中心点经纬度坐标，计算遥感图像成像时刻的太阳方位角 和太阳高度角；

（４）基于DEM数据，计算地表粗糙度，获得地表粗糙度数据，并之进行线性标准化处理（式2），（高斯-克吕格平面直角系）获得线性标准化地表粗糙度数据（式3 或式4）；

（５）基于DEM数据，以遥感图像成像时刻太阳方位角加（180°）和遥感图像成像时刻太阳高度角为参数计算研究区太阳西北方向地形阴影图（Shaded Relief Model）；

（６）太阳西北方向地形阴影图和彩色亮度分量数据分别以线性标准化地表粗糙度和作为权重，合并为新的彩色亮度分量数据（式 1）；

   （1）；

（７）将（２）获得的色相和饱和度以及（６）获得的新亮度三分量数据进行HIS->RGB变换，反变换至RGB色彩空间；该步骤实现了将地形信息同遥感图像相融合，形成地形信息同遥感信图像融合图像；

（８）基于DEM数据（高斯-克吕格平面直角系），计算每个像素的坡度，获得地形坡度数据；

（９）野外对冲沟沟缘线坡度值进行实测和统计，确定冲沟沟缘线坡度最小阈值数据；

（１０）根据野外实测冲沟沟缘线坡度最小阈值，对地形坡度数据进行阈值分割，将冲沟沟缘线分割出来，获得冲沟沟缘线图斑数据；

（１１）坡度阈值切割获得的冲沟沟缘线图斑数据，结合地形信息同遥感信图像相融合的图像，基于目视解译方法，精确解译冲沟沟缘线。

 [0007] 进一步的，所述的遥感数据（图像）为航天、航空遥感图像中的可见光波段（红光、绿光和蓝光）。

进一步的，所述的地形坡度数据和地表粗糙度数据是基于DEM数据来计算的；从地形学角度出发，将地面凹凸不平的程度定义为粗糙度，一般定义为地表单元曲面面积与投影面积之比（式2）；地表粗糙度反映了地表地形复杂程度，取值大于等于1；地表坡度为0，即平坦地区的地表粗糙度值为1，而地表坡度大，即区域（山岭）地表粗糙度值高；线性标准化地表粗糙度是由（式3）或（式4）确定， 取值(0,1 ]。： 

                    （2）

其中， 为坡度数据，以弧度为单位；

线性标准化地表粗糙度是由式（3）确定：

             （3）

其中，取值(0,1]；

式（3）亦可以采用式（4）来完成线性标准化：

       （4）。

进一步的，所述的HIS色彩空间向RGB色彩空间变换过程中，新的彩色亮度分量数据取值是太阳西北方向地形阴影图和原始彩色亮度分量数据分别以粗糙度作为权重进行取值；即由式（1）来确定。 

进一步的，所述HIS色彩空间包括HSV色彩空间、HLS色彩空间； 

进一步的，所述RGB色彩空间向HIS色彩空间变换，在有的遥感图像处理软件系统中又被称为HSV变换、HLS变换、USGS Munsell HSV变换；在HIS色彩空间中彩色亮度I分量同HLS色彩空间中的颜色明亮度（Lightness）分量、HSV色彩空间中的色彩明度（Value）分量、USGS Munsell HSV色彩空间中的色彩明度（Value）分量相对应，均表征色彩强弱（明亮），彩色亮度I(Intensity) 分量及公式（3）中的，分别等价为HLS色彩空间中的颜色明亮度（Lightness）分量、HSV色彩空间中的色彩明度（Value）分量。

所述的RGB->HIS色彩变换是遥感图像处理中的彩色空间变换算法，该算法将遥感图像红、绿、蓝色彩空间（即RGB空间）变换至色调（Hue）、彩色亮度（Intensity）、饱和度（Saturation）空间（即HIS空间）；在HIS空间中，亮度（Intensity）反映了照度信息。 

所述彩色亮度、色相、饱和度色彩空间包括HIS色彩空间、HSV色彩空间、HLS色彩空间及USGS Munsell HSV色彩空间；所述红、绿、蓝色彩空间同照度、色相和饱和度色彩空间间的正逆变换，在商业遥感图像处理软件系统中被称为RGB->HIS变换、RGB->HSV变换、RGB->HLS变换以及RGB->USGS Munsell HSV变换、HIS-> RGB变换、HSV-> RGB变换、HLS-> RGB变换以及USGS Munsell HSV-> RGB变换等；在HIS色彩空间中彩色亮度I(Intensity) 分量同HLS色彩空间中的颜色明亮度（Lightness）分量、HSV色彩空间中的色彩明度（Value）分量、USGS Munsell HSV色彩空间中的色彩明度（Value）分量等价，均表征色彩强弱（明亮）。故彩色亮度I(Intensity) 分量及公式（3）中的I，分别等价为HLS色彩空间中的颜色明亮度（Lightness）分量、HSV色彩空间中的色彩明度（Value）分量。 

在HIS色彩空间中彩色亮度I(Intensity) 分量分量取值范围为[0，1]， HLS色彩空间中的颜色明亮度（Lightness）分量和HSV色彩空间中的色彩明度（Value）分量取值范围为[0，1]，而USGS Munsell HSV色彩空间中的色彩明度（Value）分量取值范围为[0，512]。若采用USGS Munsell HSV色彩空间进行地形信息同遥感图像融合时，需将USGS Munsell HSV色彩空间中的色彩明度（Value）分量做归一化处理，使之同HIS色彩空间中的彩色亮度I(Intensity)分量取值范围相同，基于式（1）构建新彩色亮度后，需将新彩色亮度数值范围由[ 0，1]线性变换为[0，512]，然后再由USGS Munsell HSV色彩空间反变换至RGB色彩空间。 

进一步的，所述的太阳西北方向地形阴影图即是由DEM数据模拟太阳位于西北方向时地形阴影图，是由ESRI ENVI4.8软件Topographic Modeling模块中的Shaded Relief程序计算获得。 

进一步的，所述的太阳方位角 和太阳高度角，其单位是弧度。 

进一步的，所述的冲沟沟缘线图斑数据是基于野外实测沟缘线处坡度确定冲沟沟缘线坡度阈值，通过对地形坡度数据进行阈值分割而获得，其形态特征零散、破碎，不便于目视解译冲沟；冲沟沟缘线图斑数据，反映了冲沟沟缘线信息，其中也包含具有较大坡度值的非冲沟地面，因而表现的很凌乱，地形信息同遥感图像相融合后的遥感图像中，冲沟表现出明显的沟谷特征，因而，由地形信息同遥感图像相融合后的遥感图像确定冲沟位置，由冲沟沟缘线图斑数据可以精确确定冲沟沟缘线，两者结合能够高精度的解译出冲沟沟缘线。 

本发明的一种基于地形和遥感影像融合技术提取冲沟方法与已有技术相比具有突出的实质性特点和显著进步：1、将地形三维信息有效的融合到二维的遥感图像中，使得图像中冲沟的纹理特征完全符合人类视觉习惯；2、这种融合技术，使得遥感图像中地形平缓的地区图像基本无光谱信息损失，地形崎岖的地区图像光谱信息损失较少的优点；3、基于地形信息同遥感图像融合图像，结合地形坡度数据阈值数据解译出冲沟沟缘线，其精度更高；在遥感解译、沟蚀分析及制图应用中，具有重要意义。 

附图说明： 

图1是本发明的计算流程图。

具体实施方式： 

为了更好的理解与实施，下面结合附图给出具体实施例详细说明本发明一种基于地形和遥感影像融合技术提取冲沟方法；所举实施例仅用于解释本发明，并非用于限制本发明的范围、

实施例1，参见图1，第一步，首先，将研究区内，光学遥感图像RGB三波段遥感图像和DEM数据坐标系统统一成高斯-克吕格平面直角坐标系，并进行重采样，使两者具有相同的空间分辨率；

第二步，对遥感图像的RGB三色波段进行HIS变换，获得色相 和饱和度 和彩色亮度三分量数据；

第三步，根据DEM数据，结合遥感图像成像时刻（年月日时分秒）图像中心点经纬度坐标，计算成像时刻太阳方位角 和太阳高度角；

第四步，基于DEM数据，计算地表粗糙度，获得地表粗糙度数据，并之进行线性标准化处理（式2），获得线性标准化地表粗糙度数据（式3 或式4）；

第五步，基于DEM数据，以遥感图像成像时刻太阳方位角加（180°）和遥感图像成像时刻太阳高度角为参数计算研究区太阳西北方向地形阴影图（Shaded Relief Model）；

第六步，根据（式 1），太阳西北方向地形阴影图和彩色亮度I 分量数据分别以和作为权重，求和合并为新的彩色亮度分量数据；

第七步，基于第二步获得的色相和饱和度以及第六步获得的新彩色亮度三分量数据，做HIS色彩空间反变换至RGB色彩空间运算，获得三维地形信息同遥感图像相融合图像；

第八步，基于研究区DEM数据，计算每个像素的坡度，获得地形坡度数据；

第九步，通过野外对冲沟沟缘线坡度值的调查、实测，确定冲沟沟缘线坡度最小阈值数据；

第十步，根据野外实测冲沟沟缘线坡度最小阈值，对地形坡度数据进行阈值分割，将冲沟沟缘线分割出来，获得冲沟沟缘线图斑数据；

第十一步，基于地形信息同遥感图像相融合图像可以有效的识别出冲沟的位置，参考坡度阈值切割获得的冲沟沟缘线图斑数据可以精准识别出冲沟的沟缘线，从而完成精确解译冲沟沟缘线。