一种窗口遍历的SEVI调节因子自动优化算法的制作方法

本发明涉及一种窗口遍历的sevi调节因子自动优化算法。

背景技术

现有地形阴影消除植被指数tavi中地形调节因子f(△)优化方法主要有3种：“匹配寻优法(国家专利号200910111688x)”、“极值优化法(国家专利号201010180895.3)”和“相关系数法(国家专利号2015108077580)”。

“匹配寻优”算法计算步骤为：(1)影像分类，划分遥感影像中山体的阴坡和阳坡，并选取典型样区；(2)目标识别，借助地面调查资料、实地考察数据、航拍资料或者googleearth的高分辨率影像数据等核实阴坡与阳坡植被的均质性，识别典型样区阴坡与阳坡植被一致或接近的部分；(3)优化匹配，令f(△)从0开始，依次递增，考察tavi在典型样区阴坡与阳坡植被一致部分的植被指数值变化，当二者相等时，即可确定f(△)的最优结果。

“极值优化”算法计算步骤为：(1)影像分类，划分遥感影像中山体的阴坡和阳坡；(2)计算极值，计算阴坡部分tavi的最大值mtavi阴与阳坡部分tavi的最大值mtavi阳；(3)迭代寻优，令f(△)从0开始，依次递增，当满足以下公式条件时，得到f(△)最优值。

|mtavi阴-mtavi阳|≤ε,ε→0,f(△)＝0～∞

“相关系数”算法计算步骤为：(1)选择样区，在复杂地形山区选择地形影响明显、山体阳坡和阴坡对称分布的典型样区；(2)计算样区植被指数tavi、rvi和svi；(3)计算相关系数，包括tavi与cvi的相关系数r1，tavi与svi的相关系数r2；(4)逼近优化，令f(△)从0开始，依次递增，当r1与r2满足以下公式条件时，确定f(△)优化计算结果。

r1-r2≤ε,ε→0,f(△)＝0～∞

上述3种优化算法，tavi无需dem等数据的支持，都能有效削减地形阴影对山区植被信息的影响。但前2种优化算法都需要对遥感影像分类，其中“寻优匹配”算法还需要地面数据等的支持；第3种方法，虽然无需影像分类，但需要人工选择样区，存在较大的不稳定性；此外，这三种方法都容易陷入局部最优解而非全局最优解，这都限制了地形阴影消除植被指数应用的自动化水平，不利于推广应用。

技术实现要素：

有鉴于此，本发明的目的在于提供一种窗口遍历的sevi调节因子自动优化算法。该算法无需dem数据和遥感影像分类，同时不依赖于地面调查数据也不用选择样区，对全景影像计算sevi及其应用推广具有重要的科学意义与经济价值。

为实现上述目的，本发明采用如下技术方案：一种窗口遍历的sevi调节因子自动优化算法，其特征在于包括以下步骤：

步骤s1：观察一景遥感影像上的山体分布，通过阴坡、阳坡判断山体坡长，选择最大坡长确定计算窗口参数k；假设整景遥感影像的尺寸为m行n列；

步骤s2：以整景遥感影像的表观反射率数据计算阴影消除植被指数sevi、比值植被指数rvi和阴影植被指数svi；

步骤s3：计算相关系数，具体如下：

其中：r1为sevi与rvi的相关系数，r2为sevi与svi的相关系数，x、y1、y2分别为遥感影像sevi、rvi及svi计算结果，n为sevi、rvi及svi的影像数据的像元数；

步骤s4：令f(△)从0开始，间隔a，依次递增计算sevi，同时考察sevi与rvi的相关系数r1及sevi与svi的相关系数r2，当r1与r2满足以下条件时，退出内循环，得到单窗口优化解fl：

r1-r2≤ε,ε→0,f(△)＝0～∞；

步骤s5：从遥感影像左上角第一个像元(1，1)开始计算第一个窗口(1:k，1:k)的f(△)优化值；然后先行后列依次递增计算其它窗口的f(△)优化值；最终得到一个行列为(m-k，n-k)的f(△)优化值矩阵；

步骤s6：计算f(△)优化值由高往低前m％数量的分界值f，从而得到全景影像的f(△)全局最优解fg。

在本发明一实施例中，窗口参数k为50、100、150或200。

在本发明一实施例中，a为0.001。

在本发明一实施例中，m为3。

在本发明一实施例中，阴影消除植被指数sevi、比值植被指数rvi和阴影植被指数svi的计算公式为：

其中：f(△)为调节因子；br为遥感影像红光波段数据，bnir为遥感影像近红外波段数据。

本发明与现有技术相比具有以下有益效果：

1、可用于整景影像计算：本发明的自动优化算法可以计算整景影像的f(△)最优解，而非样区的局部优化解，对于sevi的实际应用和工程化推广具有重要的科学价值与经济效益。

2、地形校正效果明显：本发明确定的f(△)全局(全景)最优解，保证sevi能有效消除地形本影和落影对山区植被信息的干扰。

3、数据需求少，成本低：本发明只需要遥感影像自身携带的波段数据，无需地面调查数据或实地考察数据等的支持，数据成本与时间成本实现最小化。

4、流程简单，可操作性强：本发明主要由“窗口选择、植被指数计算、相关系数计算、逼近计算、窗口遍历、全局优化值确定”等步骤就可确定f(△)全局(全景)最优解，流程简单，操作容易，无需遥感影像分类、样区选择等环节，大大提升了sevi应用自动化水平。

附图说明

图1是本发明技术流程示意图。

具体实施方式

下面结合附图及实施例对本发明做进一步说明。

请参照图1，本发明提供一种窗口遍历的sevi调节因子自动优化算法，其包括以下步骤：

步骤s1：窗口选择：观察一景遥感影像上的山体分布，通过阴坡、阳坡判断山体坡长，选择最大坡长确定计算窗口参数k。结合实用性与计算效率，参照landsat影像30米空间分辨率，窗口参数k可以选择50、100、150、200等来表征影像上窗口大小。

步骤s2：植被指数计算：以整景遥感影像(m行n列)表观反射率数据计算阴影消除植被指数、比值植被指数和阴影植被指数；

其中：sevi为阴影消除植被指数；rvi为比值植被指数；svi为阴影植被指数；f(△)为调节因子；br为遥感影像红光波段数据，bnir为遥感影像近红外波段数据。

步骤s3：相关系数计算，具体如下：

其中：r1为sevi与cvi的相关系数，r2为sevi与svi的相关系数，x、y1、y2分别为遥感影像sevi、cvi及svi计算结果，n为sevi、cvi及svi的影像数据的像元数。

步骤s4：逼近优化，令f(△)从0开始，以a(a可以为0.001)为间隔，依次递增进行植被指数的计算，同时考察sevi与cvi的相关系数r1及sevi与svi的相关系数r2，当r1与r2满足以下条件时，退出内循环，得到单窗口优化解fl：

r1-r2≤ε,ε→0,f(△)＝0～∞。

步骤s5：窗口遍历，从遥感影像左上角第一个像元(1，1)开始计算第一个窗口(1:k,1:k)的f(△)优化值；然后先行后列依次递增计算其它窗口的f(△)优化值。最终得到一个行列为(m-k，n-k)的f(△)优化值矩阵。

步骤s6：全局最优值确定，计算f(△)优化值由高往低前m％(m可以取3)数量的分界值f，从而得到全景影像的f(△)全局最优值fg。

主要流程示意图参见图1。

进一步的，所述植被指数为阴影消除植被指数sevi、比值植被指数rvi和阴影植被指数svi，相应计算公式如下：

其中：br为遥感影像红光波段数据，bnir为遥感影像近红外波段数据。

通过对该方法在landsat8oli全景影像的验证，表明本发明计算的sevi在地形本影和落影的误差都低于5％，明显优于c地形校正和6s+c大气和地形复合校正的效果如表1所示。

表1

注：toa为表观反射率数据；c为c地形校正数据；6s+c为6s大气校正和c地形校正数据。以上所述仅为本发明的较佳实施例，凡依本发明申请专利范围所做的均等变化与修饰，皆应属本发明的涵盖范围。