CN1790052A - 基于遥感影像和gis数据的面状地物变化检测方法 - Google Patents

Abstract

本发明涉及一种基于遥感影像和GIS数据的面状地物变化检测方法，包括：一、对GIS数据预处理；二、对遥感影像预处理；三、计算面状地物的Label点；四、预测GIS面状地物Label点对应的遥感影像同名点；五、定义一组由粗到细的相似特征量；六、根据相似特征量迭代提取影像上的同名面状地物；七、根据预先定义的相似特征量阈值判断地物是否发生变化，如果发生了变化，则标记该地物，否则计算该地物多边形的重心坐标作为控制点，返回第四步。本发明对感兴趣地物的GIS数据进行上述操作，可以检测出变化地物，对于没有变化的地物，可以获得控制点，利用这些控制点可以实现影像与GIS数据的自动配准。

Description

基于遥感影像和GIS数据的面状地物变化检测方法

技术领域

本发明属于遥感影像与GIS(地理信息系统)相结合应用领域，涉及一种基于遥感影像和GIS数据的面状地物变化检测整体迭代解求方法。

背景技术

自然界的变化和人类的各种活动每天都在改变着地表景观及其利用形式。人口的快速增长及城市化的发展，加速了这种变化。这些变化将对地球资源与环境产生深远的影响，因此及时、有效地监测地表变化，更新相关的地理信息系统，为资源管理与规划和环境保护等职能部门提供科学决策的依据是十分必要的。为了解地球而发展起来的卫星对地观测技术，无疑是监测地表变化的最佳技术手段。现代遥感技术已经进入了一个能够快速、及时提供多种对地观测海量数据的新阶段。要从目前每天接收的数以TB计的遥感影像中，分拣出我们感兴趣的数据，依靠传统的人工解译方法显然不行，这就需要通过某种方法让计算机来理解图像和检测哪儿发生了变化(变化检测)，把我们的注意力引向我们感兴趣的地方。在空间数据框架建成以后，为了实现数据快速更新，保持其现势性，变化检测显得尤为重要，自动变化检测是遥感与GIS领域今后研究的重点之一(李德仁2003)。

现有的变化检测方法主要还是停留在像素级的数据引导上，缺乏知识引导的特征级变化检测方法(李德仁2003)。目前变化检测方法模型主要是对不同时相的遥感影像进行处理，其前提条件是地物的时相变化能够引起图像上像素值的明显变化。这类变化检测方法仅仅利用了图像灰度信息，没有利用检测对象的几何信息和其他知识。由于图像灰度对地物的表达有一定的不确定性，因而使得变化检测的结果出现很多伪变化。尽管已经存在许多变化检测方法，但对于某一特定的应用和研究区域选择一个合适的方法仍然十分困难。GIS数据中包含了丰富的语义和非语义信息，集成遥感图像和GIS数据进行变化检测分析是近年来变化检测方法发展的趋势之一(D.Lu 2004)。

目前一般的变化检测方法步骤可以归纳为：数据预处理、数据配准、特征提取(基于特征级的变化检测方法)、变化检测、结果分析，这些步骤是独立进行的，每一步的误差都会传递到下一步，其中数据配准和特征提取对变化检测有着关键性的影响，使得变化检测的结果精度和可靠性不可控制，而且检测过程也比较盲目。为此本发明提出了一种基于特征级的变化检测方法，把影响变化检测结果的关键步骤：数据配准、特征提取和变化检测同步迭代求解，相互约束、修正。利用GIS数据中的知识进行引导减少特征提取和变化检测的盲目性，提高结果的可靠性和精度。

发明内容

本发明的目的是提供一种基于遥感影像和GIS数据的面状地物变化检测方法，该方法能够实现多源、多传感器、不同时相数据之间的变化检测。

本发明提供的技术方案是，一种基于遥感影像和GIS数据的面状地物变化检测方法，包括以下步骤：

一、对GIS数据预处理；

I.从GIS数据中选取按实际成图比例表达的面状地物，并提取表达其几何形状的多边形；

II.计算I中每一个多边形的以下相似特征量：多边形最小外接矩形的面积；多边形最小外接矩形的宽和高；多边形的面积；多边形的周长；多边形的形状编码；

III.按下述步骤确定I中的多边形的“Label”点：

1)根据每个多边形的坐标进行栅格化，生成多边形外接矩形大小的图像；

2)用腐蚀运算(Erosion)对1)中的图像进行腐蚀处理至凸多边形图像变成一条直线；凹多边形图形变成一条或多条联通的曲线；

3)对2)生成的结果，从这些直线或者曲线中，选取一个中间点，作为Label点；或者继续利用腐蚀运算，对得到的线图像进行运算，直到线变成一个点，最后得到的这个点作为图像的Label点；

二、对遥感影像预处理；

处理步骤如下：

1)初始化i＝0，输入原始图像f_i(x，y)；

2)用二维卷积算子h(x，y)与图像f_i(x，y)进行卷积，得到f_i+1(x，y)；

其中 $h(x,y)=\left[\begin{array}{ccccc}\frac{1}{256}& \frac{1}{64}& \frac{3}{128}& \frac{1}{64}& \frac{1}{256}\\ \frac{1}{64}& \frac{1}{16}& \frac{3}{32}& \frac{1}{16}& \frac{1}{64}\\ \frac{3}{128}& \frac{3}{32}& \frac{9}{64}& \frac{3}{32}& \frac{3}{128}\\ \frac{1}{64}& \frac{1}{16}& \frac{3}{32}& \frac{1}{16}& \frac{1}{64}\\ \frac{1}{256}& \frac{1}{64}& \frac{3}{128}& \frac{1}{64}& \frac{1}{256}\end{array}\right]$

3)计算第一个小波系数：w_i+1(x，y)＝f_i(x，y)-f_i+1(x，y)；

4)如果i＜n，其中n为给定的分解次数，i＝i+1，返回2)；

5)重复2)，3)，4)直至i＝n；

6)选取1-3个小波面作为遥感影像预处理结果；

对遥感影像边界的处理采用镜像对称的方法，即：

行向：f(-i，j)＝f(i，j)；

      f(i+k，j)＝f(i-k，j)其中i＜＝N，k＝1，2，…，N为图像的总的行数；

列向：f(i，-j)＝f(i，j)；

f(i，j+k)＝f(i，j-k)其中j＜＝N，k＝1，2，…，N为图像的总的列数；

三、建立粗略几何位置关系：根据确定的控制点，通过最小二乘法计算遥感影像与GIS数据之间的几何变换模型参数，建立GIS数据与遥感影像之间的粗略几何位置关系；

四、读取经过GIS数据预处理后的一个多边形，根据三中建立的几何位置关系，预测GIS多边形Label点在遥感影像上对应的同名点的位置，以该点为种子点，采取自适应迭代区域生长方法提取遥感影像上的面状地物多边形特征；灰度一致性的阈值从零开始随着迭代次数增加逐步增大(增大区域提取阈值)，对每次提取的多边形计算最小外接矩形面积(相似量)，计算遥感影像上提取的多边形与GIS多边形最小外接矩形面积之间的差异，直到满足预定阈值条件；记录下每个满足阈值的多边形，作为GIS中多边形的候选同名多边形；

五、面状地物的变化检测：

1.计算候选同名多边形的相似特征量：多边形最小外接矩形的宽和高；多边形的面积；多边形的周长；多边形的形状编码；

2.以多边形最小外接矩形宽和高为相似特征量，从1结果中选择特征相似的多边形；

3.以多边形的面积为相似度量，从2结果中选择特征相似的多边形；

4.以多边形的周长为相似度量，从3结果中选择特征相似的多边形；

5.以多边形的形状编码为相似度量，从4结果中选择特征相似的多边形；选择形状编码的相似性最大的多边形为同名特征；

6.如果5中的结果不满足预先设定的阈值，就认为该面状地物发生了变化，标记该GIS多边形为变化区域；否则(最大相似值大于相似性阈值)该面状地物没有发生变化，取该多边形及遥感影像上的同名多边形的重心为控制点，参与下一个面状地物变化检测，优化粗略几何位置关系；

六、对其骤一I中提取的每个待检测的GIS多边形按照步骤三至五进行迭代处理，直至所有多边形都匹配完成。最后，对检测的结果进行精度评价。

本发明对感兴趣地物的GIS数据进行上述操作，可以检测出变化地物，对于没有变化的地物，可以获得控制点，利用这些控制点可以实现影像与GIS数据的自动配准。本发明利用GIS数据中的知识进行引导减少特征提取和变化检测的盲目性，提高结果的可靠性和精度。本发明能够实现多源、多传感器、不同时相数据之间的变化检测。

附图说明

附图为本发明的流程示意图。

具体实施方式

参见附图，本发明包括以下步骤：

一、对GIS数据的预处理

I.从GIS数据中选取按实际成图比例表达的面状地物，并提取表达其几何形状的多边形。

II.根据多边形坐标计算面积，并对所有的多边形按面积由大到小进行排序。

III.计算每一组用于匹配的多边形的相似特征量。相似性特征量的计算见下述步骤四。

IV.确定多边形的“Label”点。对每个多边形按照以下步骤，计算Label点：

1)根据每个多边形的坐标在计算机内存中进行栅格化，生成多边形外接矩形大小的二值图像。

2)用数学形态学中的腐蚀运算(Erosion)对1)中的结果进行多次腐蚀处理。在进行腐蚀的过程中，一般取腐蚀结构元素为[1 1；1 1]。腐蚀结构元素的对称性主要是控制腐蚀的方向。

3)如果图像是一个凸多边形，在经过步骤2)后，会变成一条直线；如果图形是一个凹多边形，腐蚀的结果会是一条或多条联通的曲线。

4)对3)生成的结果，在选取Label点时，有两种方法：一是从这些直线或者曲线中，选取一个中间点，作为Label点；二是继续利用腐蚀运算，对得到的线图像进行运算，直到线变成一个点，最后得到的这个点，就可以作为图像的Label点。

在对3)中的结果进行腐蚀时，腐蚀结构元素一般定义为[1 0；0 1]。

二、基于 小波分解的遥感影像预处理

基于遥感影像与GIS变化检测的整体解求方法中要求从影像上提取面状地物的多边形特征，为此希望影像预处理过程在抑止图像噪声的同时增强边缘信息、平滑同质区域。这里应用小波变换方法研究了一种 小波分解的遥感影像预处理方法。方法原理如下：

设函数ψ(x)∈L²(R)(L²(R)为平方可积空间)，且ψ(x)满足

$C_{\mathrm{\&psi;}}=2\mathrm{\&pi;}\underset{-\&infin;}{\overset{+\&infin;}{\&Integral;}}\frac{{\left|\stackrel{\&OverBar;}{\mathrm{\&psi;}}\left(\mathrm{\&omega;}\right)\right|}^2}{\left|\mathrm{\&omega;}\right|}\mathrm{d\&omega;}<+\&infin;---(2-1)$

或 $\underset{-\&infin;}{\overset{+\&infin;}{\&Integral;}}\mathrm{\&psi;}\left(x\right)\mathrm{dx}=0---(2-2)$

式中： ψ(ω)是ψ(x)的傅立叶变换。

当ψ(x)满足(2-1)或(2-2)式，并且能较快速收敛时，我们称ψ(x)为基本小波。当ψ(x)经过伸缩a和平移b操作后得：

${\mathrm{\&psi;}}_{a,b}\left(x\right)={\left|a\right|}^{-\frac{1}{2}}\mathrm{\&psi;}\left(\frac{x-b}{a}\right)---(2-3)$

式中：a，b∈R且a≠0

此时，ψa，b(x)称为小波。

函数f(x)的小波变换定义为：

$\mathrm{wf}(a,b)=\underset{-\&infin;}{\overset{+\&infin;}{\&Integral;}}f\left(x\right)\&CenterDot;{\left|a\right|}^{\frac{1}{2}}\stackrel{\&OverBar;}{\mathrm{\&psi;}}\left[a^{-1}(x-b)\right]\mathrm{dx}---(2-4)$

式中：f(x)∈L²(R)，b∈R， ψ(x)为ψ(x)的复共轭。

上面介绍了连续小波变换的定义，连续小波不能用来编程计算，必须要离散化。连续小波变换的离散化的方法有多种，但是，对于不同的问题不是所有的离散方法都有效。著名的Mallat算法是利用正交基进行塔式分解，但是经过该方法变换后的图像的大小发生了变化，这种变化在一些图像处理过程中往往是不利的，例如：模式识别，多源影像融合等。

为了使图像经过小波变换后尺寸大小不变，这里采用一种被称为 小波分解的方法把图像分解成不同的小波平面。

小波算法的基本思想是把信号或图像分解为不同频率通道上的近似信号和每一尺度下的细节信号。该细节信号称为小波面，其图像大小与原始图像尺寸相同。

对于一维信号C(x)，假设{C₀(x)}为信号C(x)和尺度函数φ(x)的标量积，尺度函数实际上是一个低通滤波器。信号C(x)经过第一次滤波后得到C₁(x)，w₁(x)＝C₀(x)-C₁(x)包含这两个尺度之间的信息，w₁(x)称为第一小波面，也是对应尺度函数的小波变换的结果。而小波函数ψ(x)与尺度函数φ(x)有如下关系：

$\frac{1}{2}\mathrm{\&psi;}\left(\frac{x}{2}\right)=\mathrm{\&phi;}\left(x\right)-\frac{1}{2}\mathrm{\&phi;}\left(\frac{x}{2}\right)---(2-5)$

相邻的尺度之间相差两倍，经过i次滤波后得到C_i(x)为：

$C_i\left(x\right)=\underset{1}{\mathrm{\&Sigma;}}h\left(1\right)\&CenterDot;C_{i-1}(x+2^{i-1}1)---(2-6)$

离散小波变换小波系数为：

w_i(x)＝C_i-1(x)-C_i(x)               (2-7)

w_i(x)为尺度i下的小波系数(小波面)，C_i(x)为i尺度下的近似信号，h为低通滤波器，它与尺度函数φ(x)满足下列方程：

$\frac{1}{2}\mathrm{\&phi;}\left(\frac{x}{2}\right)=\underset{1}{\mathrm{\&Sigma;}}h\left(1\right)\&CenterDot;\mathrm{\&phi;}(x-1)---(2-8)$

离散小波对信号进行分解，生成一组相邻不同分辨率的小波面{wi}和近似信号之和。

如果选择线性内插的尺度函数即：

φ(x)＝1-|x|   如果x∈[-1，1]

φ(x)＝0       如果 $x\&NotElement;[-\mathrm{1,1}]$

                                             (2-9)

于是可以计算出h(-1)＝1/4，h(0)＝1/2，h(1)＝1/4，得：

$\frac{1}{2}\mathrm{\&phi;}\left(\frac{x}{2}\right)=\frac{1}{4}\mathrm{\&phi;}(x+1)+\frac{1}{2}\mathrm{\&phi;}\left(x\right)+\frac{1}{4}\mathrm{\&phi;}(x-1)$

$C_{i+1}\left(x\right)=\frac{1}{4}{C}_i(x-2^i)+\frac{1}{2}{C}_i\left(x\right)+\frac{1}{4}{C}_i(x+2^i)$

                                             (2-10)

上述的 离散小波分解很容易推广到二维从而得到一个二维的3×3的卷积算子：

$\left[\begin{array}{ccc}\frac{1}{16}& \frac{1}{8}& \frac{1}{16}\\ \frac{1}{8}& \frac{1}{4}& \frac{1}{8}\\ \frac{1}{16}& \frac{1}{8}& \frac{1}{16}\end{array}\right]$

如果采用B3次样条尺度函数，那么二维卷积算子为：

$\left[\begin{array}{ccccc}\frac{1}{256}& \frac{1}{64}& \frac{3}{128}& \frac{1}{64}& \frac{1}{256}\\ \frac{1}{64}& \frac{1}{16}& \frac{3}{32}& \frac{1}{16}& \frac{1}{64}\\ \frac{3}{128}& \frac{3}{32}& \frac{9}{64}& \frac{3}{32}& \frac{3}{128}\\ \frac{1}{64}& \frac{1}{16}& \frac{3}{32}& \frac{1}{16}& \frac{1}{64}\\ \frac{1}{256}& \frac{1}{64}& \frac{3}{128}& \frac{1}{64}& \frac{1}{256}\end{array}\right]$

实现步骤如下：

1)初始化i＝0，输入原始图像f_i(x，y)：

2)用滤波器h(x，y)与图像f_i(x，y)进行卷积，得到f_i+1(x，y)：

f_i+1(x，y)＝f_i(x，y)×h(x，y)；

3)进行第一次小波变化，得到第一个小波系数：

w_i+1(x，y)＝f_i(x，y)-f_i+1(x，y)；

4)如果i＜n(n为给定的分解次数)，i＝i+1，返回2)；

5)重复2)，3)，4)直至i＝n。

对边界的处理采用镜像对称的方法，即：

行向：f(-i，j)＝f(i，j)；

f(i+k，j)＝f(i-k，j)其中i＜＝N，k＝1，2，…，N为图像的总的行数；

列向：f(i，-j)＝f(i，j)；

f(i，j+k)＝f(i，j-k)其中j＜＝N，k＝1，2，…，N为图像的总的列数。

在实际应用中根据需要选用合适的小波面个数。

三、几何变换参数的计算

根据控制点，通过最小二乘法计算遥感影像与GIS数据之间的几何变换模型参数，通过几何变换模型建立GIS数据坐标与遥感影像像素之间的几何位置关系。常用的几何变换模型有：仿射变换，多项式变换和透视投影变换(共线方程)等。

图像的几何变换模型根据不同的实际情况可以分别采用三种方式：(1)对于地形比较平坦的地区，先采用简单的多项式变换作为遥感影像几何变换粗略模型，待通过匹配确定了足够多的、高精度的控制点后再应用严格的影像纠正模型进行高精度微分纠正；(2)一开始就采用遥感影像几何畸变的严格纠正模型，随着控制点数目和精度的提高不断精化模型参数；

(3)采用(1)和(2)的混合型，先采用简单的多项式变换作为遥感影像几何变换粗略模型，待通过匹配获取到了能够解求严格模型参数最少控制点数目后，采用影像畸变严格改正模型，再随着控制点数目的增加对严格的影像纠正模型参数进行迭代精化。

四、相似特征量

相似特征量的定义十分关键，它是判断同名特征的尺度。相似特征量的定义直接关系到检测结果的可靠性、稳定性和唯一性，也决定了计算量的大小。

对于面状地物的多边形特征这里定义以下五个相似特征量来描述其相似性：

1、多边形最小外接矩形的面积；

2、多边形最小外接矩形的宽和高；

3、多边形的面积；

4、多边形的周长；

5、多边形的形状编码；

这些特征相似量在一起组成了一个由粗到细的相似度量集合，虽然单就每一个测度量而言不能唯一确定同名特征，但这一组相似测定却可以把候选特征缩小到很小的范围内，再加上三中计算的粗略几何位置的约束，就可以唯一确定同名特征。

五、面状地物多边形特征的提取

读取经过GIS数据预处理后的一个多边形，根据三中建立的几何位置关系，预测GIS多边形Label点在遥感影像上对应的同名点的位置，以该点为种子点，采取自适应迭代区域生长方法提取遥感影像上的面状地物多边形特征；灰度一致性的阈值从零开始随着迭代次数增加逐步增大(增大区域提取阈值)，对每次提取的多边形计算最小外接矩形面积(相似量)，计算遥感影像上提取的多边形与GIS多边形最小外接矩形面积之间的差异，直到满足预定阈值条件；记录下每个满足阈值的多边形，作为GIS中多边形的候选同名多边形。这样就自动实现了GIS数据引导下的图像区域多边形特征自适应提取。

六、搜索策略

用GIS多边形Label点的几何变换结果为引导，为了兼顾计算效率和准确性，采用如下的搜索策略：

1、按多边形面积从大到小逐步进行匹配。

面积大的多边形能够允许的几何误差大，按面积由大到小逐步提高几何变换模型的精度，可以保证经过图像几何变换模型计算的GIS多边形Label点落在同名的图像区域内。

2、对同名特征的匹配采用分层匹配，逐步细化的策略。具体方法如下：

1)以GIS多边形的Label点的几何变换结果为图像同名多边形区域提取的种子点，以GIS多边形最小外接矩形面积为相似度量，逐步增大灰度一致性阈值，提取图像多边形区域。最小外接矩形面积相似的度量阈值可以预先设定，阈值的设定可以宽松些。从图像上提取、记录下每个满足阈值的多边形，作为候选的同名多边形特征。

2)以多边形最小外接矩形宽和高为相似度量，从1)结果中选择特征相似的多边形。

3)以多边形的面积为相似度量，从2)结果中选择特征相似的多边形。

4)以多边形的周长为相似度量，从3)结果中选择特征相似的多边形。

5)以多边形的形状编码为相似度量，从4)结果中选择特征相似的多边形。选择形状编码的相似性最大的多边形为同名特征。

七、面状地物的变化检测

如果六5)中的结果不满足预先设定的相似性阈值，就认为该面状地物发生了变化，标记该GIS多边形；否则(即最大相似值大于相似性阈值)该面状地物没有发生变化，取该多边形及遥感影像上的同名多边形的重心为控制点，参与下一个面状地物变化检测，优化粗略几何位置关系，直至步骤一I中提取的待检测的GIS多边形都处理完毕。最后，对检测的结果进行精度评价。

本发明用于不同时相的GIS矢量数据之间做变化检测时，可以省略第二、五步，而直接利用矢量坐标来描述多边形地物。

本发明用于不同时相的栅格数据之间做变化检测时，先对其中一个时相的栅格数据进行地物提取，再按本发明所述的方法处理。

Claims (3)

1.一种基于遥感影像和GIS数据的面状地物变化检测方法，其特征是：包括以下步骤

一、对GIS数据预处理；

I.从GIS数据中选取按实际成图比例表达的面状地物，并提取表达其几何形状的多边形；

II.计算I中每一个多边形的以下相似特征量：多边形最小外接矩形的面积；多边形最小外接矩形的宽和高；多边形的面积；多边形的周长；多边形的形状编码；

III.按下述步骤确定I中的多边形的“Label”点：

1)根据每个多边形的坐标进行栅格化，生成多边形外接矩形大小的图像；

2)用腐蚀运算对1)中的图像进行腐蚀处理至凸多边形图像变成一条直线；凹多边形图形变成一条或多条联通的曲线；

3)对2)生成的结果，从这些直线或者曲线中，选取一个中间点，作为Label点；或者继续利用腐蚀运算，对得到的线图像进行运算，直到线变成一个点，最后得到的这个点作为图像的Label点；

二、对遥感影像预处理；

处理步骤如下：

1)初始化i＝0，输入原始图像f_i(x，y)；

2)用二维卷积算子h(x，y)与图像f_i(x，y)进行卷积，得到f_i+1(x，y)：

其中 $h(x,y)=\left[\begin{array}{ccccc}\frac{1}{256}& \frac{1}{64}& \frac{1}{128}& \frac{1}{64}& \frac{1}{256}\\ \frac{1}{64}& \frac{1}{16}& \frac{3}{32}& \frac{1}{16}& \frac{1}{64}\\ \frac{3}{128}& \frac{3}{32}& \frac{9}{64}& \frac{3}{32}& \frac{3}{128}\\ \frac{1}{64}& \frac{1}{16}& \frac{3}{32}& \frac{1}{16}& \frac{1}{64}\\ \frac{1}{256}& \frac{1}{64}& \frac{1}{128}& \frac{1}{64}& \frac{1}{256}\end{array}\right]$

3)计算第一个小波系数：w_i+1(x，y)＝f_i(x，y)-f_i+1(x，y)；

4)如果i＜n，其中n为给定的分解次数，i＝i+1，返回2)；

5)重复2)，3)，4)直至i＝n；

6)选取1-3个小波面作为遥感影像预处理结果；

对遥感影像边界的处理采用镜像对称的方法，即：

行向：f(-i，j)＝f(i，j)；

      f(i+k，j)＝f(i-k，j)其中i＜＝N，k＝1，2，…，N为图像的总的行数；

列向：f(i，-j)  ＝f(i，j)；

      f(i，j+k)＝f(i，j-k)其中j＜＝N，k＝1，2，…，N为图像的总的列数；

三、根据确定的控制点，通过最小二乘法计算遥感影像与GIS数据之间的几何变换模型参数；

四、采取自适应迭代区域生长方法提取遥感影像上的面状地物多边形特征；遥感影像上区域生长的种子点由GIS多边形的Label点对应的几何变换结果给出，灰度一致性的阈值从零开始随着迭代次数增加逐步增大，对每次提取的多边形计算最小外接矩形面积，计算遥感影像上提取的多边形与GIS多边形最小外接矩形面积之间的差异，直到满足预定阈值条件；记录下每个满足阈值的多边形，作为GIS中多边形的候选同名多边形；

五、面状地物的变化检测：

1)计算候选同名多边形的相似特征量：多边形最小外接矩形的宽和高；多边形的面积；多边形的周长；多边形的形状编码；

2)以多边形最小外接矩形宽和高为相似特征量，从1)结果中选择特征相似的多边形；

3)以多边形的面积为相似度特征量，从2)结果中选择特征相似的多边形；

4)以多边形的周长为相似度特征量，从3)结果中选择特征相似的多边形；

5)以多边形的形状编码为相似度特征量，从4)结果中选择特征相似的多边形；选择形状编码的相似性最大的多边形为同名特征；

6)如果5)中的结果不满足预先设定的阈值，就认为该面状地物发生了变化；否则该面状地物没有发生变化，取该多边形及遥感影像上的同名多边形的重心为控制点，参与下一个面状地物变化检测；

六、对步骤一I中提取的每个待检测的GIS多边形按照步骤三至五进行迭代处理。

2.根据权利要求1所述的方法，其特征是：在步骤一.III.2)进行腐蚀的过程中，取腐蚀结构元素为[11；11]。

3.根据权利要求1或2所述的方法，其特征是：在步骤一.III.3)进行腐蚀的过程中，取腐蚀结构元素为[10；01]。