本发明属于光学遥感影像预处理的技术领域，尤其涉及一种通用的基于遥感影像分割单元的地形校正优化方法。
背景技术：
地表物体的表面反射率是定量遥感进行计算的重要依据，但是遥感数据在成像过程中容易受诸多因素的干扰，导致山区不同位置接收到的有效太阳辐射差别明显，即使坡度相同的同一地物，由于坡向的影响会显示不同灰度值；不同的地物由于阴影的存在可能表现出相似灰度值，即“同物异谱”和“异物同谱”的现象。因此对山区遥感影像进行地形校正具有重要的现实意义，特别是在我国山区面积约占国土陆域总面积69％的具体国情下，地形校正是必不可少的前提。
为消除地形效应对遥感影像的影响，国内外学者从上世纪80年代开始致力于各种地形校正模型的研究，概括起来主要有三类：经验模型、物理模型和半经验模型。经验模型具有简单、易实现、适用性强的优点，但其理论基础不完备；物理模型侧重于机理研究，具有明确的物理意义，但其结构复杂、输入参数较多；半经验模型综合了前两者的优点，广泛应用于山区遥感影像的地形校正，然而其只适用于地形起伏较小的区域，地形起伏较大的区域精度仍然较低。随着地形校正研究的深入，一些改进的模型陆续用于复杂山区的地形校正，如scs c校正、旋转校正，坡度匹配模型。虽然这些模型均能取得较好的校正效果，但它们的性能往往过度依赖由像元反射率和太阳入射角余弦拟合的经验参数。
综上所述，现有地形校正模型均过度依赖研究区，且相同的地形校正模型在不同的试验区具有明显不同的校正效果，没有一种普适的地形校正模型；校正效果好的模型往往来自一些面积相对较小的研究区，且整幅影像所有像元使用一个全局参数，存在地域性强、通用性差的缺陷，无法满足大区域遥感影像的需要。
解决上述技术问题的难度：分割参数的设置直接影响到地形校正模型经验参数的估计精度。对于异质性高的遥感影像不存在普适的分割参数，需要构建一套完整的理论技术进行分割参数的优选，参数的敏感性与稳定性需要做系统的评价。
解决上述技术问题的意义：本发明提出的基于遥感影像分割单元的地形校正优化方法能够更好地适应复杂山区地形的起伏变化，恢复山区遥感影像的真实光谱信息，避免全局参数估计的盲目性，进而为传统的地形校正模型优化提供一种参考，服务于大区域森林监测及相近领域的光学遥感影像预处理。
技术实现要素：
针对现有技术存在的问题，本发明提供了一种通用的基于遥感影像分割单元的地形校正优化方法。
本发明是这样实现的，一种通用的基于遥感影像分割单元的地形校正优化方法，所述通用的基于遥感影像分割单元的地形校正优化方法包括：
步骤一，输入光学遥感影像、默认参数、dem数据以及其他测试数据；
步骤二，处理遥感数据，利用dem数据生成相应的坡度图和坡向图；计算给定光学遥感数据的光照图；
步骤三，利用meanshift进行遥感影像的自动分割，随后采用k-mean算法进行聚类融合，基于聚类融合结果，结合gps实地测量结果评价分割精度；
步骤四，基于步骤二计算得到的光照图，利用分割后的矢量图层逐一裁剪，获得基于分割单元的icmap；
步骤五，基于步骤三、步骤四，按照逐单元从左到右自上而下的顺序对每个分割单元的经验参数进行估算，保持整幅影像每个分割单元均有一个经验参数；
步骤六，以旋转校正模型为例，构建基于分割单元的地形校正优化模型；
步骤七，从视觉评价、相关性分析、参数统计分析方面，进行校正效果的综合评价。
进一步，步骤一中，所述默认参数为太阳参数，所述太阳参数包括太阳方位角和太阳高度角；所述太阳参数从遥感影像元文件中读取得到。
进一步，步骤二中，所述遥感数据初处理包括正射校正、几何精校正、彩色合成、数据融合。
进一步，步骤二中，所述给定光学遥感数据的光照图计算方法包括：
式中，z表示太阳天顶角，
ic值的变化范围介于-1到1之间，相对于水平面而言ic＝cos(z)。
进一步，步骤三中，所述利用meanshift进行遥感影像的自动分割包括：
利用三个具有一定物理意义的参数从下至上控制分割过程，进行光学遥感影像的自动分割；
所述三个具有一定物理意义的参数为空间尺度hs,色彩尺度hr,最小区域面积尺度m。
进一步，步骤五中，所述经验参数估算公式如下：
li(λ)i＝α(λ)i*ici b(λ)i；
式中，i代表遥感影像分割单元；λ表示波段号；li(λ)i表示校正前波段λ中分割单元i的反射率；a(λ)i,b(λ)i分别表示波段λ中分割单元i线性回归的斜率和截距；ici表示分割单元i的光照系数。
进一步，步骤六中，所述以旋转校正模型为例，构建基于分割单元的地形校正优化模型包括：
lh(λ)i＝li(λ)i-a(λ)i*(ic-cos(z))。
式中，i代表遥感影像分割单元；λ表示波段号；lh(λ)i,li(λ)i分别表示波段λ中分割单元i校正前和校正后的反射率；a(λ)i表示波段λ中分割单元i线性回归的斜率；ici表示分割单元i的光照系数，z表示太阳天顶角。
本发明的另一目的在于提供一种存储在计算机可读介质上的计算机程序产品，包括计算机可读程序，供于电子装置上执行时，提供用户输入接口以实施通用的基于遥感影像分割单元的地形校正优化方法。
本发明的另一目的在于提供一种计算机可读存储介质，储存有指令，当所述指令在计算机上运行时，使得计算机执行通用的基于遥感影像分割单元的地形校正优化方法。
本发明的另一目的在于提供一种实施通用的基于遥感影像分割单元的地形校正优化方法的山区森林资源的遥感监测仪。
综上所述，本发明的优点及积极效果为：本发明以中高分辨率光学遥感影像为主要数据源，采用基于遥感影像分割单元的经验参数估计替代传统的全局经验参数估计，发明一套通用的地形校正优化方法。
本发明具有普遍适用性，可为传统的地形校正模型优化提供一种参考，提出的基于遥感影像分割单元的地形校正优化方法能够更好地适应复杂山区地形的起伏变化，避免全局参数估计的盲目性。
本发明思路清晰，可移植性强，可用于大区域森林监测及相近领域的光学遥感影像预处理。
附图说明
图1是本发明实施例提供的通用的基于遥感影像分割单元的地形校正优化方法流程图。
图2是本发明实施例提供的通用的基于遥感影像分割单元的地形校正优化方法原理图。
图3是本发明实施例提供的研究区地理位置图。
图4是本发明实施例基于meanshift算法的影像分割结果图。
图5是本发明实施例提供的例地形校正前后的放大效果图(对应图3中的区域a)。
图6是本发明实施例提供的例地形校正前后的放大效果图(对应图3中的区域b)。
图7是本发明实施例提供的不同模型地形校正前后光照系数ic和近红外nir波段的相关系数图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。
现有校正方法地域性强、通用性差的缺陷，无法适应大区域山区光学遥感影像的地形校正。
针对现有技术存在的问题，本发明提供了一种通用的基于遥感影像分割单元的地形校正优化方法，下面结合附图对本发明作详细的描述。
如图1所示，本发明实施例提供的通用的基于遥感影像分割单元的地形校正优化方法包括：
s101，输入光学遥感影像、默认参数、dem数据以及其他测试数据。
s102，处理遥感数据，利用dem数据生成相应的坡度图和坡向图；计算给定光学遥感数据的光照图。
s103，利用meanshift进行遥感影像的自动分割，基于高分辨率遥感影像数据的多参数分割结果，结合gps实地测量结果评价分割精度。
s104，基于步骤s102计算得到的光照图，利用分割后的矢量图层逐一裁剪，获得基于分割单元的icmap。
s105，基于步骤s103、步骤s104，按照逐单元从左到右自上而下的顺序对每个分割单元的经验参数进行估算，保持整幅影像每个分割单元均有一个经验参数。
s106，以旋转校正模型为例，构建基于分割单元的地形校正优化模型。
s107，从视觉评价、相关性分析、参数统计分析方面，进行校正效果的综合评价。
步骤s101中，本发明实施例提供的默认参数为太阳参数，所述太阳参数包括太阳方位角和太阳高度角；所述太阳参数从遥感影像元文件中读取得到。
步骤s102中，本发明实施例提供的遥感数据初处理包括正射校正、几何精校正、彩色合成、数据融合。
步骤s102中，本发明实施例提供的给定光学遥感数据的光照图计算方法包括：
式中，z表示太阳天顶角，
ic值的变化范围介于-1到1之间，相对于水平面而言ic＝cos(z)。
步骤s103中，本发明实施例提供的利用meanshift进行遥感影像的自动分割包括：
利用三个具有一定物理意义的参数从下至上控制分割过程，进行光学遥感影像的自动分割；
所述三个具有一定物理意义的参数为空间尺度hs,色彩尺度hr,最小区域面积尺度m。
步骤s105中，本发明实施例提供的经验参数估算公式如下：
li(λ)i＝α(λ)i*ici b(λ)i；
式中，i代表遥感影像分割单元；λ表示波段号；li(λ)i表示校正前波段λ中分割单元i的反射率；a(λ)i,b(λ)i分别表示波段λ中分割单元i线性回归的斜率和截距；ici表示分割单元i的光照系数。
步骤s106中，本发明实施例提供的以旋转校正模型为例，构建基于分割单元的地形校正优化模型包括：
lh(λ)i＝li(λ)i-a(λ)i*(ic-cos(z))。
式中，i代表遥感影像分割单元；λ表示波段号；lh(λ)i,li(λ)i分别表示波段λ中分割单元i校正前和校正后的反射率；a(λ)i表示波段λ中分割单元i线性回归的斜率；ici表示分割单元i的光照系数，z表示太阳天顶角。
图2是本发明实施例提供的通用的基于遥感影像分割单元的地形校正优化方法原理图。
本发明实施例提供的通用的基于遥感影像分割单元的地形校正优化方法测试方法包括：
通过更换地形校正模型，采用常用的半经验地形校正模型如minnaert校正、c校正、scs c校正和sec校正模型等，测试所述通用的基于遥感影像分割单元的地形校正优化方法的通用性和稳定性。
下面结合具体实施例对本发明的技术方案作进一步说明。
实施例：
本发明实施例提供的通用的基于遥感影像分割单元的地形校正优化方法，适用于大部分半经验参数的地形校正模型，具体步骤如下：
步骤1：测试数据输入
光学遥感影像；太阳方位角(sun_azimuth)和太阳高度角(sun_elevation)作为默认参数；dem数据。
步骤2：数据初处理
遥感数据初处理包括正射校正、几何精校正、彩色合成、数据融合等常规初处理过程；利用dem数据生成相应的坡度图(slope)和坡向图(aspect)。
步骤3：计算光照图
利用公式(1)计算给定光学遥感数据的光照图。
式中，z表示太阳天顶角，
步骤4：影像分割
本发明引入一种非参数密度估计——meanshift实现遥感影像的分割，分割过程由三个具有一定物理意义的参数从下至上控制(hs,hr,m分别表示空间尺度、色彩尺度和最小区域面积尺度)，进而实现光学遥感影像的自动分割。利用高分辨率遥感影像数据的多参数分割结果，结合gps实地测量结果评价分割精度。
步骤5：矢量裁剪
结合步骤3的光照图，利用分割后的矢量图层逐一裁剪，获得基于分割单元的icmap。
步骤6：经验参数估计
基于步骤4和步骤5，对每个分割单元的经验参数α(λ)i进行估算，按照逐单元从左到右自上而下的顺序移动，最后整幅影像的每个分割单元均有一个经验参数。具体公式如下：
li(λ)i＝α(λ)i*ici b(λ)i
式中，i代表遥感影像分割单元；λ表示波段号；li(λ)i表示校正前波段λ中分割单元i的反射率；a(λ)i,b(λ)i分别表示波段λ中分割单元i线性回归的斜率和截距；ici表示分割单元i的光照系数。
步骤7：地形校正
以旋转校正模型为例，构建基于分割单元的地形校正优化模型。
lh(λ)i＝li(λ)i-a(λ)i*(ic-cos(z))
式中，i代表遥感影像分割单元；λ表示波段号；lh(λ)i,li(λ)i分别表示波段λ中分割单元i校正前和校正后的反射率；a(λ)i表示波段λ中分割单元i线性回归的斜率；ici表示分割单元i的光照系数，z表示太阳天顶角。
步骤8：性能评价
从视觉评价、相关性分析、参数统计分析(如最大值、最小值、平均值、标准差和变异系数等)方面，综合评价本发明方法的校正效果。
步骤9：适用性分析
更换地形校正模型，采用常用的半经验地形校正模型如minnaert校正、c校正、scs c校正和sec校正模型等，测试基于遥感影像分割单元的地形校正优化方法的通用性和稳定性。
下面结合具体实验及仿真对本发明作进一步描述。
如图3所示，是选择的研究区地理位置图。选择的研究区位于湖南省攸县黄丰桥国有林场，研究区主要部分均为陡峭的斜坡和崎岖的地形，海拔落差大，坡度0～59度，具有典型的中国山地特征。图3中的遥感数据为2013年9月17日接收的landsat8影像，dem数据为30米空间分辨率的srtm数据。
图3本发明试验例的研究区地理位置图，区域a和区域b分别对应后续的放大效果图：如图5和图6。
图4是本发明实施例的影像分割结果图。采用的是基于meanshift的分割算法，该算法通过控制核带宽参数h＝(hs,hr,m)来决定分割精度。具体分割过程：首先确定分割参数，对遥感影像进行初始分割；然后确定合并参数，将光谱属性相同或相似的斑块合并为一类，减少初始分割斑块数量。经过反复试验，分割尺度为(15,25,200)、合并参数为(15,90,200)时，可以获得较好的分割效果。
图5、图6分别是本发明实施例旋转地形校正前后图3中区域a和区域b的放大效果图。观察可知原始影像校正后，图像的山体立体效果趋于水平，地形起伏与山体阴影现象得到明显消除。全局旋转校正模型在阴影区域出现了一定程度的过度校正，而且局部区域的地形较为破碎，发明者提出的基于分割单元的局部旋转校正模型，局部区域过度校正现象明显改善。
图5中(a)原始遥感影像(b)基于全局参数的地形校正(c)基于分割单元的地形校正，以旋转校正为例。图6中(a)原始遥感影像(b)基于全局参数的地形校正(c)基于分割单元的地形校正，以旋转校正为例。
如图7本发明实施例不同模型地形校正前后光照系数ic和近红外nir波段的相关系数所示，以近红外波段为例，计算原始影像和校正后影像的反射率与太阳入射角余弦值的pearson相关系数。观察可知，校正前影像的反射率与太阳入射角余弦值相关性高(r2为0.167，相关系数达0.447)，整体红黄色偏多；地形校正后影像反射率与太阳入射角余弦值的相关性明显减低，r2均小于0.05，相关系数从0.447分别减小到0.065和0.028，基于分割单元的旋转校正模型减小幅度最大，达85.41％。
图7中(a)原始影像；(b)基于全局参数的地形校正；(c)基于分割单元的地形校正，以旋转校正为例。图例0.1～0.9的渐变色系为可见光的蓝色到红色，波长由短逐渐变长。
如表1，为了进一步检验基于分割单元的旋转模型校正效果，研究以近红外波段为例，分析遥感影像近红外波段校正前后的参数统计结果。分析可知，地形校正前后影像的最小值、最大值基本没变，分别为0.000和1.000；均值有轻微增大，从校正前的0.805变为校正后的0.807和0.806；标准差从原始影像到全局旋转校正和局部旋转校正，呈现逐渐减小的趋势，说明基于分割单元的旋转校正模型在消除地形和背景因素干扰的同时，最大限度地保留了原始影像的光谱信息。变异系数的逐渐减小进一步佐证了该结论。
表1遥感影像近红外波段校正前后的参数统计结果
在上述实施例中，可以全部或部分地通过软件、硬件、固件或者其任意组合来实现。当使用全部或部分地以计算机程序产品的形式实现，所述计算机程序产品包括一个或多个计算机指令。在计算机上加载或执行所述计算机程序指令时，全部或部分地产生按照本发明实施例所述的流程或功能。所述计算机可以是通用计算机、专用计算机、计算机网络、或者其他可编程装置。所述计算机指令可以存储在计算机可读存储介质中，或者从一个计算机可读存储介质向另一个计算机可读存储介质传输，例如，所述计算机指令可以从一个网站站点、计算机、服务器或数据中心通过有线(例如同轴电缆、光纤、数字用户线(dsl)或无线(例如红外、无线、微波等)方式向另一个网站站点、计算机、服务器或数据中心进行传输)。所述计算机可读取存储介质可以是计算机能够存取的任何可用介质或者是包含一个或多个可用介质集成的服务器、数据中心等数据存储设备。所述可用介质可以是磁性介质，(例如，软盘、硬盘、磁带)、光介质(例如，dvd)、或者半导体介质(例如固态硬盘solidstatedisk(ssd))等。
以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。
1.一种通用的基于遥感影像分割单元的地形校正优化方法，其特征在于，所述通用的基于遥感影像分割单元的地形校正优化方法包括：
步骤一，输入光学遥感影像、默认参数、dem数据以及测试数据；
步骤二，处理遥感数据，利用dem数据生成相应的坡度图和坡向图；计算给定光学遥感数据的光照图；
步骤三，利用meanshift进行遥感影像的自动分割，基于高分辨率遥感影像数据的多参数分割结果，结合gps实地测量结果评价分割精度；
步骤四，基于步骤二计算得到的光照图，利用分割后的矢量图层逐一裁剪，获得基于分割单元的icmap；
步骤五，基于步骤三、步骤四，按照逐单元从左到右自上而下的顺序对每个分割单元的经验参数进行估算，使整幅影像每个分割单元均有经验参数；
步骤六，构建基于分割单元的地形校正优化模型；
步骤七，进行校正效果的综合评价。
2.如权利要求1所述通用的基于遥感影像分割单元的地形校正优化方法，其特征在于，步骤一中，所述默认参数为太阳参数，所述太阳参数包括太阳方位角和太阳高度角；所述太阳参数从遥感影像元文件中读取得到。
3.如权利要求1所述通用的基于遥感影像分割单元的地形校正优化方法，其特征在于，步骤二中，所述遥感数据初处理包括正射校正、几何精校正、彩色合成、数据融合。
4.如权利要求1所述通用的基于遥感影像分割单元的地形校正优化方法，其特征在于，步骤二中，所述给定光学遥感数据的光照图计算方法包括：
式中，z表示太阳天顶角，
ic值的变化范围介于-1到1之间，相对于水平面而言ic＝cos(z)。
5.如权利要求1所述通用的基于遥感影像分割单元的地形校正优化方法，其特征在于，步骤三中，所述利用meanshift进行遥感影像的自动分割包括：
利用三个具有一定物理意义的参数从下至上控制分割过程，进行光学遥感影像的自动分割；
所述三个具有一定物理意义的参数为空间尺度hs,色彩尺度hr,最小区域面积尺度m。
6.如权利要求1所述通用的基于遥感影像分割单元的地形校正优化方法，其特征在于，步骤五中，所述经验参数估算公式如下：
li(λ)i＝α(λ)i*ici b(λ)i；
式中，i代表遥感影像分割单元；λ表示波段号；li(λ)i表示校正前波段λ中分割单元i的反射率；a(λ)i,b(λ)i分别表示波段λ中分割单元i线性回归的斜率和截距；ici表示分割单元i的光照系数。
7.如权利要求1所述通用的基于遥感影像分割单元的地形校正优化方法，其特征在于，步骤六中，构建基于分割单元的地形校正优化模型为：
lh(λ)i＝li(λ)i-a(λ)i*(ici-cos(z))；
式中，i代表遥感影像分割单元；λ表示波段号；lh(λ)i,li(λ)i分别表示波段λ中分割单元i校正前和校正后的反射率；a(λ)i表示波段λ中分割单元i线性回归的斜率；ici表示分割单元i的光照系数，z表示太阳天顶角。
8.一种存储在计算机可读介质上的计算机程序产品，包括计算机可读程序，供于电子装置上执行时，提供用户输入接口以实施如权利要求1～7任意一项通用的基于遥感影像分割单元的地形校正优化方法。
9.一种计算机可读存储介质，储存有指令，当所述指令在计算机上运行时，使得计算机执行如权利要求1～7任意一项通用的基于遥感影像分割单元的地形校正优化方法。
10.一种实施权利要求1～7任意一项通用的基于遥感影像分割单元的地形校正优化方法的山区森林资源的遥感监测仪。
技术总结