发明名称 一种对象级高分辨率SAR影像洪水灾害变化检测方法
摘要 本发明公开一种对象级高分辨率SAR影像洪水灾害变化检测方法,针对SAR影像中可能存在的与水体区域具有相似几何特征及波谱特征的虚假目标,造成大量的“伪变化”给洪水灾害变化检测造成的困难与干扰。因此,对每个时相影像进行轮廓波变换。在保持图像边缘特性的同时压制图像噪声,同时通过简单的样本训练选择最佳分解尺度,并在该尺度上通过分块直方图统计快速获取水体可能存在区域的标记点位置。进而采用基于标记点的分水岭分割以及区域合并策略获得水体的轮廓信息,并利用基于多种特征的判别规则进一步消除虚假目标的干扰。最后,根据配准结果直接比较多时相影像中所提取的水体轮廓,获得水体发生变化的区域。
申请公布号 CN104361582A 申请公布日期 2015.02.18
申请号 CN201410573002.X 申请日期 2014.10.23
申请人 河海大学 发明人 王超;石爱业;高红民;徐立中;黄凤辰
分类号 G06T7/00(2006.01)I 主分类号 G06T7/00(2006.01)I
代理机构 南京苏高专利商标事务所(普通合伙) 32204 代理人 李玉平
主权项 一种对象级高分辨率SAR影像洪水灾害变化检测方法,其特征在于,主要由几何配准、基于轮廓波变换的噪声抑制及标记点提取、基于标记点的分水岭分割及区域合并、虚假目标消除及变化检测组成;其中图像配准采用美国ITT视觉解决方案的ENVI软件;轮廓波变换主要包括两个步骤:Step1:采用多尺度拉普拉斯金字塔分解影像,捕捉影像中的奇异点信息;每次分解过程都首先采用低通滤波器对上一尺度进行滤波,进而进行下采样获得当前尺度的低频影像;对当前尺度低频影像做上采样操作,并进一步利用低通滤波进行平滑,最后与上一尺度低频影像相减可获得当前尺度下的高频影像;Step2:采用方向滤波器组对高频影像进行方向滤波;针对SAR影像中固有的相干斑噪声,首先采用轮廓波变换对配准后的多时相SAR影像进行多尺度去噪;在轮廓波分解后的影像中,噪声对应较小的轮廓波系数,采用硬阈值方法对噪声加以区分;采用如下多尺度策略:Step1:若原始影像尺寸为N×N像素,定义影像中像素的灰度值为f<sub>i,j</sub>(i,j=1,2,3....N)。定义分解层数M,进行方向为2<sup>M</sup>的轮廓波变换;Step2:利用公式(1)对分解后第m(m=1,2...M)层的系数的方差δ<sup>2</sup>(m)进行估计,如下所示:<maths num="0001" id="cmaths0001"><math><![CDATA[<mrow><msup><msub><mi>&delta;</mi><msub><mi>f</mi><mrow><mi>i</mi><mo>,</mo><mi>j</mi></mrow></msub></msub><mn>2</mn></msup><mrow><mo>(</mo><mi>m</mi><mo>)</mo></mrow><mo>=</mo><mfrac><mrow><mi>Med</mi><mrow><mo>(</mo><mo>|</mo><msub><mi>W</mi><msub><mi>f</mi><mrow><mi>i</mi><mo>,</mo><mi>j</mi></mrow></msub></msub><mrow><mo>(</mo><mi>m</mi><mo>)</mo></mrow><mo>|</mo><mo>)</mo></mrow></mrow><mi>C</mi></mfrac><mo>-</mo><mo>-</mo><mo>-</mo><mrow><mo>(</mo><mn>1</mn><mo>)</mo></mrow></mrow>]]></math><img file="FDA0000592676020000011.GIF" wi="1186" he="178" /></maths>其中,C为常数;<img file="FDA0000592676020000012.GIF" wi="194" he="116" />为第m层的轮廓波变换系数。Step3:利用公式(2)计算第m层各高频子带的阈值,如下所示:<maths num="0002" id="cmaths0002"><math><![CDATA[<mrow><mi>T</mi><mrow><mo>(</mo><mi>m</mi><mo>)</mo></mrow><mo>=</mo><msup><msub><mi>&delta;</mi><msub><mi>f</mi><mrow><mi>i</mi><mo>,</mo><mi>j</mi></mrow></msub></msub><mn>2</mn></msup><mrow><mo>(</mo><mi>m</mi><mo>)</mo></mrow><msqrt><mn>21</mn><mi>gL</mi><mrow><mo>(</mo><mi>m</mi><mo>)</mo></mrow></msqrt><mo>-</mo><mo>-</mo><mo>-</mo><mrow><mo>(</mo><mn>2</mn><mo>)</mo></mrow></mrow>]]></math><img file="FDA0000592676020000013.GIF" wi="1126" he="113" /></maths>其中,L(m)为第m层分解系数的数量;Step4:采用硬阈值方法对不同方向、不同子带的系数进行处理;Step5:进行轮廓波反变换,获得去噪影像;最佳轮廓波分解尺度选择:通过计算所选取水体样本的局部方差均值,选择均值最小的尺度为最佳分解尺度,并认为在此尺度下水体内部的均质度最高,因此仅在此尺度下搜索可能存在水体的区域,并提取标记点;最佳轮廓波分解尺度选择过程如下:在图像去噪中,进行M层轮廓波变换后,每层取K=2<sup>M</sup>个方向,则对图像分解得到M×K个高频子带的轮廓波系数W<sub>m,k</sub>,k=1,2,3,...K表示其对应的高频子带;分别利用公式(3)计算各尺度下各个波段的相位一致值PC<sub>m,k</sub>:<maths num="0003" id="cmaths0003"><math><![CDATA[<mrow><msub><mi>PC</mi><mrow><mi>m</mi><mo>,</mo><mi>k</mi></mrow></msub><mo>=</mo><mfrac><mrow><msub><mi>&Sigma;</mi><mi>n</mi></msub><msub><mi>A</mi><mi>n</mi></msub><mrow><mo>(</mo><mi>x</mi><mo>)</mo></mrow><mi>cos</mi><mrow><mo>(</mo><msub><mi>&phi;</mi><mi>n</mi></msub><mrow><mo>(</mo><mi>x</mi><mo>)</mo></mrow><mo>-</mo><mover><mi>&phi;</mi><mo>&OverBar;</mo></mover><mrow><mo>(</mo><mi>x</mi><mo>)</mo></mrow><mo>)</mo></mrow></mrow><mrow><munder><mi>&Sigma;</mi><mi>n</mi></munder><msub><mi>A</mi><mi>n</mi></msub><mrow><mo>(</mo><mi>x</mi><mo>)</mo></mrow></mrow></mfrac><mo>-</mo><mo>-</mo><mo>-</mo><mrow><mo>(</mo><mn>3</mn><mo>)</mo></mrow></mrow>]]></math><img file="FDA0000592676020000021.GIF" wi="1307" he="201" /></maths>其中,A<sub>n</sub>为W<sub>m,k</sub>傅里叶分解的n次谐波分量的幅度,φ<sub>n</sub>(x)表示该相位偏移量分量在x处的局部相位,<img file="FDA0000592676020000022.GIF" wi="107" he="80" />为所在点的傅里叶分量的加权平均相位;若所有傅里叶分量都具有一致的相位则该比值为1;反之该比值最小值为0。对各尺度下各个波段的相位一致值PC<sub>m,k</sub>,求梯度G<sub>m,k</sub>;采用公式(4)融合个尺度下不同波段相位一致梯度值:<maths num="0004" id="cmaths0004"><math><![CDATA[<mrow><msub><mover><mi>G</mi><mo>~</mo></mover><mi>m</mi></msub><mo>=</mo><mi>max</mi><mrow><mo>(</mo><msub><mi>G</mi><mrow><mi>m</mi><mn>1</mn></mrow></msub><mo>,</mo><msub><mi>G</mi><mrow><mi>m</mi><mn>2</mn></mrow></msub><mo>,</mo><mo>.</mo><mo>.</mo><mo>.</mo><mo>,</mo><msub><mi>G</mi><mi>mK</mi></msub><mo>)</mo></mrow><mo>-</mo><mo>-</mo><mo>-</mo><mrow><mo>(</mo><mn>4</mn><mo>)</mo></mrow></mrow>]]></math><img file="FDA0000592676020000023.GIF" wi="1156" he="91" /></maths>在原始图像中选取c个大小为d×d像素的水体区域作为样本,为了使个尺度下的样本具有相同的地物面积,在各个尺度所对应相位一致梯度图的相同位置用不同尺寸的窗口对水体样本进行采样,用公式(5)计算个尺度下水体的局部方差均值:<maths num="0005" id="cmaths0005"><math><![CDATA[<mrow><msub><msup><mi>&sigma;</mi><mn>2</mn></msup><mi>m</mi></msub><mo>=</mo><mfrac><mn>1</mn><mi>c</mi></mfrac><mi>&Sigma;</mi><mfrac><mn>1</mn><mrow><mi>g</mi><mo>&times;</mo><mi>g</mi></mrow></mfrac><munderover><mi>&Sigma;</mi><mrow><mi>x</mi><mo>=</mo><mn>1</mn></mrow><mi>g</mi></munderover><munderover><mi>&Sigma;</mi><mrow><mi>y</mi><mo>=</mo><mn>1</mn></mrow><mi>g</mi></munderover><msup><mrow><mo>[</mo><msub><mover><mi>G</mi><mo>~</mo></mover><mi>m</mi></msub><mrow><mo>(</mo><mi>x</mi><mo>,</mo><mi>y</mi><mo>)</mo></mrow><mo>-</mo><msub><mover><mi>G</mi><mover><mo>~</mo><mo>&OverBar;</mo></mover></mover><mi>m</mi></msub><mo>]</mo></mrow><mn>2</mn></msup><mo>-</mo><mo>-</mo><mo>-</mo><mrow><mo>(</mo><mn>5</mn><mo>)</mo></mrow></mrow>]]></math><img file="FDA0000592676020000024.GIF" wi="1204" he="167" /></maths>式中,σ<sup>2</sup><sub>m</sub>为尺度m下水体地物的局部方差均值;g×g像素为原始影像中样本区域在当前尺度下的对应区域的尺寸;<img file="FDA0000592676020000025.GIF" wi="172" he="79" />表示m尺度下坐标为(x,y)的像元的相位一致梯度值;<img file="FDA0000592676020000026.GIF" wi="76" he="90" />为采样样本的相位一致梯度均值;取σ<sub>m</sub><sup>2</sup>值最小的尺度作为最佳分解尺度,并在此尺度下确定标记点位置。
地址 211100 江苏省南京市江宁区佛城西路8号