专利名称：基于序列非线性滤波的遥感影像水体专题信息提取方法
技术领域：
本发明属于遥感技术领域，具体涉及一种基于序列非线性滤波的遥感影像水体专题信 息提取方法。
背景技术：
水体是重要的专题信息，它对于地学分析有重要的价值，是一类基础的专题信息数据。 对于各地理信息系统的管理与建设有重要意义。利用遥感影像进行水体信息提取是一种快 速、且现势性好的方式。目前利用遥感影像提取村镇专题信息的方法有人工方法、半自 动跟踪提取方法，以及自动模式分类识别方法；这些方法中人工方法、半自动获取方法需 要大量的人工操作效率低；利用多光谱信息分类识别方法是一种重要的方法，但现有的分 类识别方法不是专门针对遥感影像上的水体信息提取而研究的，而且这种识别分类方法对 于一般分辨率(空间分辨率20m，或更低)的多波段遥感影像才比较适用；对于高分辨率 (分辨率达到5m及以上)且波段较少(4个波段或更少)的遥感影像，其对水体的分类识 别提取结果不能满足后续成图及应用要求。目前，在高分辨率遥感影像分类识别方面，德 国的eCognition系统的分类方法是比较好的一种，其基于多尺度分割，分类时需要有分 类协议，每一次分类的分类协议用于其它分类时需经过修改，最终实现分类过程的自动化， 这里的分类协议类似于知识支持。其特点是分类的细化，但对分类协议有较大依赖(有协 议库支持较好)，并且也不是专门针对水体这类信息而进行的，因此应用于水体信息提取 时效率不高。为提高水体提取效率，避免提取过程中的人工干扰，实现一种全自动的水体 提取方法，本方法针对上述的情况产生了利用灰度数学形态学运算，构成一种序列非线性 滤波模型，直接进行非线性滤波，有效提取遥感影像中的水体区域；在此基础上，根据水 体的相关知识，对上面提取出的水体区域，计算不同标定区域的各项特征值(面积，区域 灰度均值，区域的灰度直方图)，根据水体成像特性和水体区域相关知识，建立水体目标 区域选择准则，以此进行水体区域自动选择提取，得到水体专题信息。
发明内容
本发明的目的在于，针对现有的利用高分辨率遥感影像进行水体专题信息提取方法存 在的受各种细节信息干扰严重、提取效率低、效果差等缺点，提出一种基于序列非线性滤 波的遥感影像水体专题信息提取方法，本发明在保留利用遥感影像进行水体专题信息提取 现势性强、重复性高等优点的同时，对高分辨率多波段遥感影像的水体信息提取更加适用，其效率更高。克服了一般提取方法需要先验知识，对噪音及细节干扰信息敏感的缺点，利 用灰度数学形态学运算来构成一种序列非线性滤波模型，直接对高分辨率遥感影像进行非 线性滤波，统计标定目标区域的特征值，根据水体所具有的特征，来进行遥感影像水体专 题信息自动识别。
本发明提出的基于序列非线性滤波的遥感影像水体专题信息提取方法，通过从灰度数 学形态学基本运算入手，构成一种序列非线性滤波模型，直接进行非线性滤波，提取遥感 影像中的水体区域；在此基础上，根据水体的相关知识，对序列非线性滤波提取出的水体 区域，计算不同区域的各项特征值(面积，区域灰度均值，区域的灰度直方图)，根据水 体成像特性和水体区域相关知识，建立水体目标区域选择准则，从而实现遥感影像上水体 专题信息的自动快速提取，具体包括遥感影像的预处理与序列非线性滤波，滤波结果影 像区域标定和各标定区域各项特征值的计算，水体目标区域选择准则的建立，水体区域的 自动提取及基于形态学的水体专题提取信息后处理四个步骤，具体如下
(1) 遥感多波段影像的预处理与序列非线性滤波
根据水体光谱特性和遥感影像波段设置，选择红外或红波段、蓝或绿波段作为备选波 段，根据公式(l)进行预处理，并对预处理结果做直方图均衡化，生成待处理数据图层， BW0 = BL + k x BL / (BH + w)(1)
其中，BL表示红外或红波段亮度值，BH表示蓝或绿波段亮度值，BW。表示预处理后 的像素亮度值，k为比例变换系数，w是为了避免除零问题，而対Bh加上的一个很小的 正数。
采用序列非线性滤波方法，对生成待处理数据图层首先进行最大值滤波处理，消除影
像中非连续面状分布的暗的像素，同时呈面状连续分布的水体区域也被削弱；然后进行中 值滤波进行降噪处理；进而采用最小值滤波处理恢复最大值滤波处理中被削弱的水体区 域。
(2) 对滤波结果影像进行区域标定、计算各个标定区域的各项特征值 采用区域生长算法，对步骤(1)中的结果影像先标定候选水体像素组成连续区域，
然后计算各个区域的各项统计特征区域面积、区域灰度均值、区域灰度直方图， 区域灰度直方图为
邵]二X!1i = 0,l，.."255 (2) 叩]
区域面积A为,=o
区域灰度均值 1 255
(4) 直方图峰值灰度级
OP^尸/， H[PI]=MAX{H[i]，i = 0,l,...,255} (5) 其中，R[i]为区域中灰度值为i的像素集合，H[i]为直方图第i级灰度的像素个数，A 为区域面积，GM为区域灰度均值，GP为直方图峰值灰度级； (3)水体目标区域选择准则的建立
①区域面积A < A0的区域为非水体区域；②区域灰度均值GM > GM0的区域为非水 体；③区域峰值灰度GP >区域均值灰度GM的区域为非水体区域。其中A0及GM0取 值均有较大的冗余度，可取八0=100， GM0 = 64。
(4)水体区域的自动提取，基于形态学的水体专题提取信息后处理
根据步骤(3)建立的水体区域准则，及各个区域的特征值，对标记的各个区域进行 判定，从而自动剔除杂波及其他目标，得到水体目标区域；通过形态开、闭运算，对生成 的专题图像进行滤波，除去比结构元素小的特定图像细节，在形态学滤波的基础上对提取 的矢量进行长度及面积统计，对提取出的矢量长度及所围面积小于阈值的进行剔除。阈值 可调整，对于高分辨率航空遥感影像一般长度阈值为25，面积阈值为500，最后提取的专 题信息分层表示为ArcGIS软件的形文件，即shp文件格式。
以下对本发明方法作进一步的限定，具体如下
1、遥感多波段影像的预处理与序列非线性滤波
根据水体光谱特性和遥感影像波段设置，选择红外(IR)及蓝(B)两个波段作为备选
波段；针对不同遥感影像，如果没有波段设置或缺少上述波段，也可选择红(R)代替红
外(IR)或绿(G)代替蓝(B)波段。则可根据下式进行预处理并对预处理结果作直方图
均衡化，生成待处理数据图层。
BWo = BL + k x BL / (BH + w)
其中，Bl表示紅外(IR)或红(R)波段亮度值，BH表示蓝(B)或绿(B)波段亮
度值，BWo表示预处理后的像素亮度值，k为比例变换系数，w是为了避免除零问题，而
対Bh加上的一个很小的正数。
应用灰度形态学算子与非线性统计滤波算子，在灰度形态学运算中定义一个特殊的结构元素，将灰度形态学运算转化为非线性滤波处理，通过序列非线性滤波，将形态学灰度 膨胀运算转化为最大值滤波处理，灰度腐蚀运算转化为最小值滤波处理。采用组合处理步 骤，对遥感影像首先进行最大值滤波处理，消除影像中非连续面状分布的暗的像素，同时 呈面状连续分布的水体区域也被削弱；然后进行中值滤波进行降噪处理；进而采用最小值 滤波处理恢复最大值滤波处理中被削弱的水体区域。处理过程中邻域大小可以灵活选取， 也可以连续进行多次最大值滤波，最小值滤波时也连续进行相同次数的处理，则对于连续 面状分布的水体目标可以得到保留，而消除掉大部分非水体区域，提取遥感影像中水体目 标区域。
2、对滤波结果影像进行区域标定、计算各个标定区域的各项特征值 采用区域生长算法，对(1)中的结果影像先标定候选水体像素组成连续区域，然后计 算各个区域的各项统计特征区域面积、区域灰度均值、区域灰度直方图。以此为自动建 立各个区域的判断规则做准备，避免了影像全局阈值设定困难的问题。
3、 水体目标区域选择准则的建立
根据水体成像特性和水体区域相关知识，建立水体目标区域选择准则，根据水体为连 续面状分布，建立区域面积判定规则，对面积小于额定像素的区域判定为非水体区域；根 据水体区域亮度比较统一，其直方图应为单峰形状，建立灰度分布判定规则，判定出非水 体区域。这里亮度判断可以设定比较宽松的标准，结合直方图单峰特征，先排除掉明显非 水体区域。再根据区域直方图峰值灰度和区域均值灰度的关系，自动判断水体区域。自动 确定以下水体选择准则1)区域面积A〈AO的区域为非水体区域；2)区域灰度均值GM 〉GMO的区域为非水体；3)区域峰值灰度GP〉区域均值灰度GM的区域为非水体区域。 其中AO及GMO取值均有较大的冗余度，可取AO-IOO， GM0 = 64。
4、 水体区域的自动提取，基于形态学的水体专题提取信息后处理 根据以上建立的水体区域准则，及各个区域的特征值，对标记的各个区域进行判定，
从而自动剔除杂波及其它目标，得到水体目标区域；在此基础上，为了消除结果上因细小 地物以及噪声影响，而产生的细小地物图斑，通过形态开、闭运算，对生成的专题图像进 行滤波(结构元素的大小可以选择有3X3、 5X5,7X7, 15X15)，除去比结构元素小的特 定图像细节，同时保证不产生全局的集合失真，在形态学滤波的基础上对提取的矢量进行 长度及面积统计，对提取出的矢量长度及所围面积小于阈值的进行剔除。阈值可调整(对 于高分辨率航空遥感影像一般长度阈值为25，面积阈值为500)，最后提取的专题信息分 层表示为ArcGIS软件的形文件，即shp文件格式。
与现有技术相比，本发明采用了序列非线性滤波方法，避免了传统方法二值化中的阈值计算确定，速度更快，而且对高分辨率卫星影像更加适用，效率更高，采用区域直方图 形状及灰度峰值均值相对关系判断方法，有效的避免了水体判断时的阈值问题，避免提取 过程中作业者的人工干预，具有更好的鲁棒性和自适应特性。本发明可以应用于高分辨率 多波段遥感影像的水体专题信息快速自动提取，在水体专题信息获取与变化检测方面发挥 作用。
图1为本发明实施例1中TM(专题成像仪)图像数据除热红外波段以外的其它谱段1、 2、 3、 4、 5、 7，其中(a)为1波段谱段(蓝)，(b)为2波段谱段，(c)为3波段谱段，(d) 为4波段谱段，(e)为5波段谱段，(f)为7波段谱段(红外)。
具体实施例方式
以下结合一个实施例子对具体实现方法进行说明，即要对一幅高分辨率多波段遥感影 像进行水体专题信息识别，获得满足ArcGIS软件要求的水体专题信息矢量数据(.shp文 件)。
实施例1:
这里取某一TM (专题成像仪)数据——图象大小为2516X2395像素，含7个波段， 其己具有地理编码信息——图象左上角点坐标为(285950， 3475350)，像素分辨率为25 米，图l显示的是除热红外波段以外的其它谱段l、 2、 3、 4、 5、 7。由于数据量极大，例 子中U)、 (2)、 (3)中的显示数据是该图象局部自图象像素行列坐标(753， 723)至(768， 738)范围的16X16大小图像块的处理结果。
(1)遥感多波段影像的预处理和序列非线性滤波
根据水体光谱特性和遥感影像波段设置，选择红外(IR)及蓝(B)两个波段作为备 选波段；针对不同遥感影像，如果没有波段设置缺少上述波段，也可选择红(R)代替红 外(IR)或绿(G)代替蓝(B)波段。定义BL表示上述红外(IR)或红(R)波段亮度 值，BH表示上述蓝(B)或绿(B)波段亮度值，BWo表示预处理后的像素亮度值。则可 以下式进行预处理，并对预处理结果作直方图均衡化生成待处理数据图层。 BW0 = BL + k* BL/ (BH + w)
其中，参数k及w取值都有很大的冗余度，这里取k-20,w-0.1。如果BWo〉255，则 限制BWo-255。
然后，采用序列非线性滤波方法，针对上述预处理后图层，可根据公式(6) (8) 依次进行具有非线性特性的最大值滤波，中值滤波及最小值滤波。
必^ = AWJV{i/8cg^;)} (6)5『2 = (fi『,)} (7)
5= Ar/7V{J7S (丑)} (8)
其中，"—)表示像素的8邻域，M^(')， M￡"(')， M^(')分别表示取邻域内像
素的最大值、中值及最小值，BWQ， BW,， BW2， BW3分别表示预处理后图层、最大值滤 波后、中值滤波后及最小值滤波后的像素亮度值。
这里，选用TM图像的第1波段(蓝)和第7波段(红外)作为BH和B^波段灰度 数据。表l及表2为选取的一个16xl6大小区域影像波段数据示例。表3及表4为线性变 换及直方图均衡化后结果数据。表5为根据公式(6) ~ (8)进行序列非线性滤波后结果 数据，其中已将大于128的像素灰度置为背景灰度255。
表1蓝色波段
102104108110107105106103105110109103102100100102
100104107107106106105102104107107103104102100103
10010511011010710610811010711010910410210199102
97103107105105106107106104106107104101101102103
10110210198103105104104103101103105104102103102
101100100102104103100100101101101104104103104102
9799102103103103100999899103103100102104102
97100100999910210310210110110110299101103101
99991009897100103104101102102102101102102100
100100100989910010210110210310210110410510299
100999998100999910010210110210210310310199
9999101101999999100100101賜101100104102101
1001011011009910110010010110099101102105103104
99100999910010010099100102100102103102101101
99100101101100101100100100102102100102102102102
9910010310110010210410310210210299101101102103
表2红外波段
22 335360615346515957503325273026
33 395968646454424949604320332921
32 426168574853524851603815251512
29 3851605449525451565733131598
31 34273041494134303036372517149
24 262413213027252319253933161422
10 1318141291217202421212111817
537161915101825251233434
3304152019262614500210
24237162126164252211
233105131361240230
1340133323101112
2132142022223312
0021001121221220
0330121340111221
3300351132203413
表3线性转换后结果
926396270726354607067593929323531
39467080767664505858715123393425
38497280675762615660704517291814
3445607164586164606667391517109
37403236485848403535424429201610
28312815253532292722294639191626
1215211614101420242825252513920
638192217112129291433434
3304182322303116500210
24238192531194252211
233106151571240230
1340133323101112
2132142022223312
0021001121221220
0330121340111221
3300351132203413
表4线性转换并经直:方图均衡化后结果
601101822032071851571762031951731107285恥81
110132203221214214187144170170205147501109256
10614120722119516718217916417620312941724238
9212917620518717017918717619319511039413433
1031138510013817013811396恥12012672454034
6781673956968572644972132110444060
35394740383438455367565656373345
29273244494134477272382727282728
2727028425049768140280025220
25282527324456814428252825252222
25272722029393930222528025270
222728022272727252722022222225
2522272522282502525252527272225
002522002222252225252225250
02727022252227280222222252522
27270027282222272525027282227
表5序列非线性滤波后结果
255255255255255255255255255255255255255255255255
255255255255255255255255255255255255255255255255
255255255255255255255255255255255255255255255255
25525525525525525525525525525525525525525525533
25525525525525525525525525525525525525525525534
2552552552552552552552552552552552552552554060
2552552552552552552552552552552552552552553345
292732444925525525525525525525527282728
27270284250255255255255280025220
25282527324425525525528252825252222
2527272202925525530222528025270
222728022272727252722022222225
2522272522282502525252527272225
002522002222252225252225250
02727022252227280222222252522
27270027282222272525027282227(2)对滤波结果影像进行区域标定、计算各个标定区域的各项特征值 对非线性滤波结果影像采用区域生长算法进行区域标定，然后统计各个标定区域的灰 度直方图，根据下面4式，分别计算各个区域的面积、灰度均值和直方图峰值灰度。
则二》p0，1,…,255 对']
255
1 255
力,.=i
OP = P/ ， H[PI] =MAX{H[i],i = 0,l,...,255}
其中，R[i]为区域中灰度值为i的像素集合，H[i]为直方图第i级灰度的像素个数，A为区 域面积，GM为区域灰度均值，GP为直方图峰值灰度级。
(3) 水体目标区域选择准则的建立
根据第(2)歩计算的特征值，可以自动确定以下水体选择准则1)区域面积A〈AO 的区域为非水体区域；2)区域灰度均值GM > GMO的区域为非水体；3)区域峰值灰度 GP>区域均值灰度GM的区域为非水体区域。其中AO及GMO取值均有较大的冗余度， 这里取AO-IOO， GMO = 64。
(4) 水体区域的自动提取，基于形态学的水体专题提取信息后处理
根据(3)建立的水体区域准则以及各个区域的特征值，自动对标记的各个区域进行
判定，从而自动剔除杂波及其它目标，从备选区域中选择出水体区域。在此基础上，通过
形态开、闭运算，对选出的水体区域进行形态学滤波(结构元素的大小可以选择有3X3、
5X5,7X7, 15X15),去掉微小矢量分支(除去比结构元素小的特定图像细节)，同时保证
不产生全局的集合失真，在形态学滤波的基础上，提取区域边缘矢量，对提取的矢量进行
长度及面积统计，矢量长度及所围面积小于阈值的进行剔除。阈值可调整(对于高分辨率
航空遥感影像一般长度阈值为25，面积阈值为500)，最后提取的专题信息分层表示为
ArcGIS软件的形文件，即shp文件格式。表6为生成的某一地物专题矢量数据示例。
表6生成的某一地物专题矢量数据
点号 专题矢量结点坐标(X， Y)
1 ~ 4 300600,3452250 300550,3452200 300475,3452200 300475， 3452075 5~8 300450,3452050300450,3451950300425,3451925 300400,3451925<image>image see original document page 12</image>
权利要求
1、一种基于序列非线性滤波的遥感影像水体专题信息提取方法，其特征在于具体步骤如下(1)遥感多波段影像的预处理与序列非线性滤波根据水体光谱特性和遥感影像波段设置，选择红外或红波段、蓝或绿波段作为备选波段，根据公式(1)进行预处理生成待处理数据图层，BW0＝BL+k×BL/(BH+w) (1)其中，BL表示红外或红波段亮度值，BH表示蓝或绿波段亮度值，BW0表示预处理后的像素亮度值，k为比例变换系数，w是为了避免除零问题，而对BH加上的一个很小的正数；采用序列非线性滤波方法，根据预处理生成待处理数据图层公式(1)，对遥感影像首先进行最大值滤波处理，消除影像中非连续面状分布的暗的像素，同时呈面状连续分布的水体区域也被削弱；然后进行中值滤波进行降噪处理；进而采用最小值滤波处理恢复最大值滤波处理中被削弱的水体区域；(2)对滤波结果影像进行区域标定、计算各个标定区域的各项特征值采用区域生长算法，对步骤(1)中的结果影像先标定候选水体像素组成连续区域，然后利用公式(2)～(5)，计算各个区域的各项统计特征区域面积、区域灰度均值、区域灰度直方图，区域灰度直方图为区域面积A为区域灰度均值直方图峰值灰度级GP＝PI，H[PI]＝MAX{H[i]，i＝0，1，...，255}(5)其中，R[i]为区域中灰度值为i的像素集合，H[i]为直方图第i级灰度的像素个数，A为区域面积，GM为区域灰度均值，GP为直方图峰值灰度级；(3)水体目标区域选择准则的建立①区域面积A<A0的区域为非水体区域；②区域灰度均值GM>GMO的区域为非水体；③区域峰值灰度GP>区域均值灰度GM的区域为非水体区域；其中A0及GM0取值均有较大的冗余度；(4)水体区域的自动提取，基于形态学的水体专题提取信息后处理根据步骤(3)建立的水体区域准则，及各个区域的特征值，对标记的各个区域进行判定，从而自动剔除杂波及其他目标，得到水体目标区域；通过形态开、闭运算，对生成的专题图像进行滤波，除去比结构元素小的特定图像细节，在形态学滤波的基础上对提取的矢量进行长度及面积统计，对提取出的矢量长度及所围面积小于阈值的进行剔除，最后提取的专题信息分层表示为ArcGIS软件的形文件，即shp文件格式。
全文摘要
本发明属于遥感技术领域，具体涉及一种基于序列非线性滤波的遥感影像水体专题信息提取方法。具体步骤为遥感影像的预处理与序列非线性滤波，滤波结果影像区域标定和各标定区域各项特征值的计算，水体目标区域选择准则的建立，水体区域的自动提取及基于形态学的水体专题提取信息后处理四个步骤，最后生成满足通用GIS软件要求的专题信息矢量数据。本发明从灰度数学形态学基本运算入手，构成一种序列非线性滤波模型，直接进行非线性滤波，提取遥感影像中的水体区域；根据水体的相关知识，对序列非线性滤波提取出的水体区域，计算不同区域的各项特征值根据水体成像特性和水体区域相关知识，建立水体目标区域选择准则，从而实现遥感影像上水体专题信息的自动快速提取，本发明受噪声和误差影响较少，鲁棒性强，处理提取水体过程全自动进行，无须人工干预，能处理的数据广泛，耗时少。
文档编号G06T7/00GK101546431SQ200910050718
公开日2009年9月30日 申请日期2009年5月7日 优先权日2009年5月7日
发明者勤 叶, 张绍明, 栋 王, 锋 谢, 邵永社 申请人:同济大学