CN108257123A - 基于高阶统计特征的多波段雷达图像变化检测方法 - Google Patents
基于高阶统计特征的多波段雷达图像变化检测方法 Download PDFInfo
- Publication number
- CN108257123A CN108257123A CN201810028135.7A CN201810028135A CN108257123A CN 108257123 A CN108257123 A CN 108257123A CN 201810028135 A CN201810028135 A CN 201810028135A CN 108257123 A CN108257123 A CN 108257123A
- Authority
- CN
- China
- Prior art keywords
- image
- formula
- higher order
- sar
- order
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 230000008859 change Effects 0.000 title claims abstract description 74
- 238000001514 detection method Methods 0.000 title claims abstract description 33
- 238000000034 method Methods 0.000 claims abstract description 48
- 230000008569 process Effects 0.000 claims abstract description 19
- 230000004927 fusion Effects 0.000 claims abstract description 8
- 238000012545 processing Methods 0.000 claims description 24
- 238000004364 calculation method Methods 0.000 claims description 18
- 238000000205 computational method Methods 0.000 claims description 12
- 238000013178 mathematical model Methods 0.000 claims description 9
- 230000001427 coherent effect Effects 0.000 claims description 8
- 238000000605 extraction Methods 0.000 claims description 8
- 230000009466 transformation Effects 0.000 claims description 7
- 238000009825 accumulation Methods 0.000 claims description 3
- 230000003044 adaptive effect Effects 0.000 claims description 3
- 238000004422 calculation algorithm Methods 0.000 claims description 3
- 238000012937 correction Methods 0.000 claims description 3
- 230000002708 enhancing effect Effects 0.000 claims description 3
- 238000010606 normalization Methods 0.000 claims description 3
- 238000002203 pretreatment Methods 0.000 claims description 2
- 230000005855 radiation Effects 0.000 claims description 2
- 238000003384 imaging method Methods 0.000 description 15
- 238000012544 monitoring process Methods 0.000 description 4
- 238000013461 design Methods 0.000 description 3
- 230000000694 effects Effects 0.000 description 3
- 238000001914 filtration Methods 0.000 description 3
- 230000002123 temporal effect Effects 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 2
- 239000000284 extract Substances 0.000 description 2
- 230000005764 inhibitory process Effects 0.000 description 2
- 230000007246 mechanism Effects 0.000 description 2
- VMXUWOKSQNHOCA-UKTHLTGXSA-N ranitidine Chemical compound [O-][N+](=O)\C=C(/NC)NCCSCC1=CC=C(CN(C)C)O1 VMXUWOKSQNHOCA-UKTHLTGXSA-N 0.000 description 2
- 238000012360 testing method Methods 0.000 description 2
- 238000007476 Maximum Likelihood Methods 0.000 description 1
- 230000009471 action Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 238000013473 artificial intelligence Methods 0.000 description 1
- 238000013528 artificial neural network Methods 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 238000007621 cluster analysis Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 230000009977 dual effect Effects 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 238000007689 inspection Methods 0.000 description 1
- 238000000513 principal component analysis Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000003595 spectral effect Effects 0.000 description 1
- 238000012706 support-vector machine Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/0002—Inspection of images, e.g. flaw detection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/70—Denoising; Smoothing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/30—Determination of transform parameters for the alignment of images, i.e. image registration
- G06T7/33—Determination of transform parameters for the alignment of images, i.e. image registration using feature-based methods
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10032—Satellite or aerial image; Remote sensing
- G06T2207/10044—Radar image
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Quality & Reliability (AREA)
- Radar Systems Or Details Thereof (AREA)
- Image Processing (AREA)
Abstract
基于高阶统计特征的多波段雷达图像变化检测方法，先输入不同时相的不同波段的SAR遥感图像，再对SAR图像进行预处理，然后对SAR图像进行滤波处理，再对两幅SAR图像进行归一化处理，然后对两幅SAR图像进行配准处理，分别获得各图像的高阶统计特征图，再获得两图像的联合高阶统计特征图；然后产生高阶统计特征差异图，在高阶统计特征差异图中产生检测阈值，再获得由高阶统计特征产生的变化信息，然后对两图像的高阶联合统计特征图进行处理，设置阈值，获取变化信息，对不同特征获得的变化信息进行融合处理，最后输出变化信息图，本发明能够有效地检测到不同时相不同波段SAR图像的变化信息。
Description
技术领域
本发明属于信号处理技术领域，具体涉及基于高阶统计特征的多波段雷达图像变化检测方法。
背景技术
合成孔径雷达(Synthetic Aperture Radar,以下简称SAR)成像遥感是一种非常重要的数据源，用不同时相的SAR图像来获取变化信息是遥感领域的一个重要应用方向，SAR图像变化检测跟普通的光学遥感图像变化检测原理是一样的，首先获取同一区域不同时间的SAR图像，然后通过对两幅或序列的SAR图像进行处理来获取监测区域的动态或变化信息。
目前多时相的SAR图像变化检测方法可以归纳为以下五类：第一类是基于图像像素运算的直接变化检测方法，即在不同时相的遥感图像之间直接进行处理或运算，以获取变化信息，如图像灰度差值法、图像灰度比值法、图像纹理特征差值法、相关系数法、图像回归法和典型相关法等。第二类是基于图像特征域的变化检测方法，首先分别提取不同时相图像的特征，然后对这些特征进行比较和运算处理，以获得变化信息，如最大似然估计法、贝叶斯估计法和统计特征法等。第三类是基于变换域的变化检测方法，这类方法首先对不同时相的图像进行相同的变换，然后在变换域中提取变化信息，如傅里叶变换、小波变换、主成份分析和多尺度几何变换理论等等。第四类是基于图像分类的变化检测方法，这类方法首先对不同时相的遥感图像分别进行分类，然后对不同的类分别进行比较分析，以获得最终的变化信息，如支持向量机、聚类分析和神经网络等。第五类是基于新的理论和交叉学科的变化检测方法，如稀疏表示、模糊理论、人工智能、信息融合和各种群智能等。其中，基于特征的方法主要是利用低阶的统计特性、纹理特征、区域特征、分形特征和结构特征等，几乎没有涉及到图像的高阶统计量特征。
遥感图像存在空间分辨率、辐射分辨率、光谱分辨率和时间分辨率，遥感图像变化检测就是利用不同时间遥感图像来获取变化信息。目前，不同时相的遥感图像进行变化检测，充分利用的就是遥感图像的时间分辨率。同卫星同传感器的遥感图像很难满足实际需求，因为单颗卫星的周期通常比较长，而一些突发事件，例如各种灾害的监测和救援，需要及时了解情况。此时，就需要不同卫星的不同传感器的来弥补单卫星的时间分辨率。不同卫星上安装的传感器一般不同，尤其是SAR成像，每颗卫星的SAR成像传感器发射的微波频率不同，获取的信息有别。因微波成像过程中，频率是影响微波成像的一个重要因素。所以，不同时相的多波段SAR图像变化检测比单波段的变化检测难，而且更复杂。
发明内容
为了克服上述现有技术的缺点，本发明的目的在于提供了基于高阶统计特征的多波段雷达图像变化检测方法，能够有效地检测到监测区域或目标的变化信息。
为了达到上述目的，本发明采取的技术方案为：
基于高阶统计特征的多波段雷达图像变化检测方法，包括以下步骤：
步骤1：输入不同时相的SAR遥感图像；
步骤2：对SAR图像进行预处理；
步骤3：对SAR图像进行滤波处理；
步骤4：对两幅SAR图像进行归一化处理；
步骤5：对两幅SAR图像进行配准处理；
步骤6：分别获得各图像的高阶统计特征图；
步骤7：获得两图像的联合高阶统计特征图；
步骤8：产生高阶统计特征差异图，把步骤6中所获得的高阶统计特征图进行差值或比值运算，获得不同时相SAR图像的特征差异图。
步骤9：在高阶统计特征差异图中产生检测阈值，具体步骤为，
步骤9.1：获取特征差异图的直方图，通过直方图的波峰和波谷确定一个大概的阈值T1；
步骤9.2：采用自适应的期望值最大化(Expectation Maximization,EM)算法对特征差异图进行处理，获得比较精细的阈值T2；
步骤9.3：对步骤9.1和步骤9.2中阈值进行均值处理，获得最后的检测阈值T，即
步骤10：获得由高阶统计特征产生的变化信息；由高阶统计矩特征差异图按步骤9的方法产生的阈值来获取变化信息，具体的判断函数为
式中的C1表示获取的变化信息，(m,n)表示像素的空间坐标，FMD表示高阶矩特征差异图；获得由高阶统计累积量特征差异图产生的变化信息为C2；
步骤11：对两图像的高阶联合统计特征图进行处理，设置阈值，获取变化信息，高阶联合统计特征包括高阶联合矩特征和累积量特征，对联合特征图按步骤9产生阈值，按步骤10进行判断检测，获得的变化信息用C3表示；
步骤12：进行融合处理，变化信息的融合处理按式(3)数学模型进行处理，即
C(m,n)＝α·C1(m,n)+β·C2(m,n)+γ·C3(m,n) (21)
式中的C表示获取的最终的变化信息，α、β和γ表示权系数，它们分别设为0.25，0.25和0.5；
步骤13：输出变化信息图。
所述的步骤2中对输入SAR图像进行预处理步骤为：对输入不同时相的SAR图像对分别进行辐射和几何校正处理。
所述的步骤3中对SAR图像进行滤波处理的具体步骤为；
步骤3.1：选择相干性原理对SAR图像进行滤波处理，去除斑点噪声的同时能够有效保持边缘和几何细节信息；具体实现过程如下，先对SAR图像进行傅里叶变换，在频域内进行错位处理，获得两幅子图像，进行反傅里叶变换，进行相干处理，获得相干处理增强的图像，进行坐标变换，恢复到原坐标下的图像；
步骤3.2：以图像中某点(m，n)为中心，设置一个卷积模板，其窗口大小为K×L，对窗口内的所有像素值进行均值处理，处理结果代替像素点(m，n)的值；所有的像素都按上述规则处理，直到全部像素处理完毕，设置窗口参数K＝L＝5；
步骤3.3：把增强后的图像和均值滤波后的图像合成一幅新的SAR图像。
所述的步骤4中对两幅SAR图像进行归一化处理的具体步骤为：
步骤4.1：先分别计算两幅图像像素灰度值的最大值与最小值；
步骤4.2：按公式(1)计算，实现SAR图像的归一化处理，
式中为归一化后的某像素的灰度值，Ii为归一化前该像素的灰度值，Imax和Imin分别为该图像中像素灰度值的最大值和最小值。
所述的步骤5中按Hu-SIFT方法对两幅图像实施配准处理。
所述的步骤6中获得各SAR图像的高阶统计特征图的具体步骤为：
步骤6.1：获取SAR图像的高阶距统计特性，对于一个随机变量x，其k阶原点矩mk和k阶中心矩μk的计算方法分别为式(2)和式(3)，其中式中E[·]表示求数学期望值；
mk＝E[xk] (2)
其中表示均值，一幅图像的大小为N×N，即像素个数为N2，图像的k阶原点矩mk和k阶中心矩μk的计算方法分别为式(4)和式(5)；
式中的为均值，n表示统计特性的阶次；
当n＝1时，获得的是均值特征图；当n＝2时，获得的是方差特征图；
步骤6.2：邻域窗口大小选择，在计算SAR图像的高阶统计特征时，在式(4)和式(5)中，N表示图像的所有像素个数，对于计算某个像素(m，n)的统计特性时，以该像素为中心，选取一个邻域范围内的像素来计算该像素的统计特性，选择计算高阶统计特性的邻域窗口大小设置为5×5；
步骤6.3：获取高阶矩特征图，对提取窗口大小为5的图像三阶和四阶中心矩特征图的计算公式如下，
式中的阶数为3，即n＝3，用来获取三阶高阶矩统计特征图；当n＝4时，得到计算四阶统计矩特征的公式，即
步骤6.4：获取高阶累积量特征图，随机变量的二阶和三阶累积量与对应的二阶和三阶中心矩的计算方法是一样的，即图像的二阶和三阶累积量特征图就是图像的二阶和三阶矩特征图；随机变量的四阶累积量计算方法为
同样，取每个像素点的5×5邻域块来获取图像的四阶累积量，得到图像的四阶累积量特征图，计算公式如下：
所述的步骤7中获取两图像的联合高阶统计特征图的具体步骤为：
步骤7.1：两幅SAR图像的高阶联合矩特征图的获取，对于任意两个随机变量x和y，它们的k1+k2联合矩计算公式为
则三阶和四阶联合矩的计算公式分别如下：
由式(10)-(12)推导配准后的不同波段的两幅SAR图像中对应像素点的三阶和四阶联合矩特征提取的数学模型分别为式(13)和式(14)，
步骤7.2：获取两幅SAR图像的高阶联合累积量特征图，由统计理论知随机变量的三阶和四阶累积量计算公式分别为
c111＝cum(x,y,z)＝E[xyz] (15)
c1111＝cum(x,y,z,t)＝E[xyzt]-E[xy]E[zt]-E[xz]E[yt]-E[xt]E[yz] (16)
从式(15)和式(16)式中可知，在求三阶和四阶联合累积量特征时，需要三个或四个变量参数参与求解，选择三个或四个不同波段的SAR图像来计算它们的联合高阶累积量，针对两个临近波段SAR图像来提取变化信息，把某一波段图像视为第三个波段的图像，就求得了它们的三阶联合累积量特征，c111＝E[xyz]＝E[xy2]＝E[x2y]，这样得到的三阶联合累积量与三阶联合矩是一样的；
同理，得到四阶联合累积量统计特征的数学模型为
在两幅配准后的SAR图像中，取移到窗口大小为5的区域，计算两幅图像对应像素的联合四阶累积量，具体计算公式如式(18)所示，
本发明与现有技术相比的优点在于：
(1)方法设计的优点：由于SAR成像是相干性成像，SAR图像中不可避免地会产生大量的斑点噪声，为了有效去除斑点噪声或者尽量降低斑点对变化信息提取的影响，因此，在设计上采用双重去除斑点噪声的策略。首先，根据SAR相干性成像机理，对SAR图像进行相干性滤波处理，这是第一次对噪声的滤除。其次，在获取SAR图像的高阶统计特征时，采用邻域窗的方法，即以某像素为中心，对周围像素进行处理，获得的结果代替该像素的值，进一步对噪声进行抑制处理，即第二次滤波。采用的是异源数据来获取变化信息，难度比同源数据的处理大多了。因此，设计上不是直接对SAR图像处理，而是通过高阶特征来处理。
(2)数据处理的优点：处理的不是原始SAR图像数据，而是通过原始图像获取其高阶统计特征数据，有利于变化信息的获取。另外，处理的数据不是同波段的数据，而是不同波段的数据，即不同时间获取的数据来源于多个波段，即异源遥感数据，可以弥补实际的时间分辨率不足的问题，与一般的变化检测有所不同，尤其是数据描述的内容。
(3)变化信息获取的优点：用本方法对不同波段不时相的SAR图像进行处理，以获取变化或动态监测信息。本方法通过不同层次和不同角度来获取变化信息，不是单纯地利用某种高阶统计特征来获取变化信息，因此最后得到的变化信息更加准确。
附图说明
图1为基于高阶统计特征的多波段雷达图像变化检测方法流程图。
图2为多波段SAR图像配准的Hu-SIFT方法流程图。
具体实施方式
下面结合附图和实施例对本发明作详细描述。
参照图1，基于高阶统计特征的多波段雷达图像变化检测方法，包括以下步骤：
步骤1：输入不同时相的SAR遥感图像；
步骤2：对SAR图像进行预处理：对输入不同时相的SAR图像对分别进行辐射和几何校正处理，提高图像整体质量，有利于后续处理和信息提取；
步骤3：对SAR图像进行滤波处理，具体步骤为；
步骤3.1：因为SAR成像是相干成像，所以选择相干性原理对SAR图像进行滤波处理，目的是去除斑点噪声的同时能够有效保持边缘和几何细节信息；具体实现过程如下，先对SAR图像进行傅里叶变换，在频域内进行错位处理，获得两幅子图像，进行反傅里叶变换，进行相干处理，获得相干处理增强的图像，进行坐标变换，恢复到原坐标下的图像；
步骤3.2：以图像中某点(m，n)为中心，设置一个卷积模板，其窗口大小为K×L，对窗口内的所有像素值进行均值处理，处理结果代替像素点(m，n)的值；所有的像素都按上述规则处理，直到全部像素处理完毕，在本方法中设置窗口参数K＝L＝5；
步骤3.3：把增强后的图像和均值滤波后的图像合成一幅新的SAR图像；
步骤4：对两幅SAR图像进行归一化处理，具体步骤为：
步骤4.1：先分别计算两幅图像像素灰度值的最大值与最小值；
步骤4.2：按公式(1)计算，实现SAR图像的归一化处理，
式中为归一化后的某像素的灰度值，Ii为归一化前该像素的灰度值，Imax和Imin分别为该图像中像素灰度值的最大值和最小值；
步骤5：对两幅SAR图像进行配准处理，如图2所示，具体步骤为：
步骤5.1：在参考的SAR图像中提取尺度不变特征转换(Scale Invariant FeatureTransform，简称SIFT)特征点，并确定SIFT特征点的方向；
步骤5.2：获取特征点邻域图像块的Hu矩，构成一个8维的描述子；
步骤5.3：用同样的方法获得待配准SAR图像的Hu矩和8维描述子；
步骤5.4：用描述算子对两幅SAR图像进行配准处理，并输出配准后的SAR图像；
步骤6：分别获得各图像的高阶统计特征图，具体步骤为：
步骤6.1：获取SAR图像的高阶距统计特性，对于一个随机变量x，其k阶原点矩mk和k阶中心矩μk的计算方法分别为式(2)和式(3)，其中式中E[·]表示求数学期望值；
mk＝E[xk] (2)
其中表示均值，一幅图像的大小为N×N，即像素个数为N2，图像的k阶原点矩mk和k阶中心矩μk的计算方法分别为式(4)和式(5)；
式中的为均值，n表示统计特性的阶次；
当n＝1时，获得的是均值特征图；当n＝2时，获得的是方差特征图；
步骤6.2：邻域窗口大小选择，在计算SAR图像的高阶统计特征时，在式(4)和式(5)中，N表示图像的所有像素个数，对于计算某个像素(m，n)的统计特性时，以该像素为中心，选取一个邻域范围内的像素来计算该像素的统计特性，例如3×3、5×5和7×7等；统计特征的阶数越高，其对灰度值变化较大的区域越敏感，所以，窗口越小，获取的边缘或轮廓的高阶特征越清晰，且定位越准确；但是窗口太小，又不利于噪声的抑制和滤除；因此，在大量实验的基础上和实验结果的比较后，选择计算高阶统计特性的邻域窗口大小设置为5×5，这样既能获得好的边缘现象，又能对噪声起到滤波作用；
步骤6.3：获取高阶矩特征图，选择窗口大小为5时，得到的统计值不会受到远离中心点的像素点的影响，所以，对提取窗口大小为5的图像三阶和四阶中心矩特征图的计算公式如下，
式中的阶数为n＝3，即用来获取三阶高阶矩统计特征图；当n＝4时，得到计算四阶统计矩特征的公式，即
步骤6.4：获取高阶累积量特征图，由高阶累积量理论可知，随机变量的二阶和三阶累积量与对应的二阶和三阶中心矩的计算方法是一样的，即图像的二阶和三阶累积量特征图就是图像的二阶和三阶矩特征图；下面主要研究图像的四阶累积量的获取，随机变量的四阶累积量计算方法为
c4＝m4-3m2 2-4m1m3+12m1 2m2-6m1 4 (8)
同样，取每个像素点的5×5邻域块来获取图像的四阶累积量，得到图像的四阶累积量特征图，计算公式如下：
步骤7：获得两图像的联合高阶统计特征图，具体步骤为：
步骤7.1：两幅SAR图像的高阶联合矩特征图的获取，对于任意两个随机变量x和y，它们的k1+k2联合矩计算公式为
则三阶和四阶联合矩的计算公式分别如下：
由式(10)-(12)推导配准后的不同波段的两幅SAR图像中对应像素点的三阶和四阶联合矩特征提取的数学模型分别如式(13)和式(14)表示，
步骤7.2：获取两幅SAR图像的高阶联合累积量特征图，由统计理论知随机变化的三阶和四阶累积量计算公式分别为
c111＝cum(x,y,z)＝E[xyz] (15)
c1111＝cum(x,y,z,t)＝E[xyzt]-E[xy]E[zt]-E[xz]E[yt]-E[xt]E[yz] (16)
从式(15)和式(16)式可知，在求三阶和四阶联合累积量特征时，需要三个或四个变量参数参与求解，本发明中，可以选择三个或四个不同波段的SAR图像来计算它们的联合高阶累积量，但是SAR成像受波长(即频率)的影响，特别是当两个波长的频率相差比较大时，SAR图像中包含的信息差别也比较大，所以，本发明针对两个临近波段SAR图像来提取变化信息，例如X波段和C波段，把某一波段图像视为第三个波段的图像，就求得了它们的三阶联合累积量特征，c111＝E[xyz]＝E[xy2]＝E[x2y]，这样得到的三阶联合累积量与三阶联合矩是一样的；
同理，得到四阶联合累积量统计特征的数学模型为
在两幅配准后的SAR图像中，取移到窗口大小为5的区域，计算两幅图像对应像素的联合四阶累积量，具体计算公式如式(18)所示，
步骤8：产生高阶统计特征差异图，把步骤6中所获得的高阶统计特征图进行差值或比值运算，获得不同时相SAR图像的特征差异图；
步骤9：在高阶统计特征差异图中产生检测阈值，具体步骤为：
步骤9.1：获取特征差异图的直方图，通过直方图的波峰和波谷确定一个大概的阈值T1；
步骤9.2：采用自适应的期望值最大化(Expectation Maximization,EM)算法对特征差异图进行处理，获得比较精细的阈值T2；
步骤9.3：对步骤9.1和步骤9.2中阈值进行均值处理，获得最后的检测阈值T，即
步骤10：获得由高阶统计特征产生的变化信息，具体步骤为：
由高阶统计矩特征差异图按步骤9的方法产生的阈值来获取变化信息，具体的判断函数为
式中的C1表示获取的变化信息，(m,n)表示空间坐标，FMD表示高阶矩特征差异图；获得由高阶统计累积量特征差异图产生的变化信息为C2；
步骤11：对两图像的高阶联合统计特征图进行处理，设置阈值，获取变化信息，具体步骤为：高阶联合统计特征包括高阶联合矩特征和累积量特征，由于它们是直接进行的联合特征运算，所以获得的是联合特征图，而不是特征差异图，对联合特征图按步骤9产生阈值，按步骤10进行判断检测，获得的变化信息用C3表示；
步骤12：进行融合处理，变化信息的融合处理按式(21)数学模型进行处理，即
C(m,n)＝α·C1(m,n)+β·C2(m,n)+γ·C3(m,n) (21)
式中的C表示获取的最终的变化信息，α、β和γ表示权系数，它们分别设为0.25，0.25和0.5；
步骤13：输出变化信息图。
SAR图像变化检测的难点在于SAR复杂的成像机理，相干性成像使SAR图像中包含着大量的斑点噪声，同时成像***参数和地物参数也制约SAR成像。为了尽可能降低斑点噪声对变化信息提取的影响，本发明提出了两次滤波的概念，即在预处理阶段用相干滤波方法对SAR图像进行滤波处理，在提取SAR图像高阶统计特征时采用领域的方法进一步降低斑点噪声的影响。因为是不同的波段，意味着合成孔径雷达发射的电磁波的波长不同，即频率不同，同一地物在不同频率的SAR图像中体现是有差别的，为了降低这种差别的影响，因此对不同时相不同波段的SAR图像进行归一化处理。在本发明中，核心内容是SAR图像高阶统计特征图的获取，考虑到像素的空间结构和关系，采用邻域运算的方式产生SAR图像的高阶统计特征图。高阶统计特征图包括单幅SAR的三阶和四阶矩及累积量，两幅图像的联合高阶距和累积量。高阶统计特征的最大优势是能够扑捉差异信息。为了使提取的信息更加准确，除了用不同时相的SAR图像的统计特征产生变化信息外，还用联合高阶统计特征，然后把它们的结果进行融合，得到最终的变化信息。实验表明该方法是一种非常有效的变化信息获取方法，对多源遥感图像处理，特别是遥感大数据的处理和协同利用，有着重要借鉴和启发作用，具有极大的应用潜力。
Claims (7)
1.基于高阶统计特征的多波段雷达图像变化检测方法，其特征在于，包括以下步骤：
步骤1：输入不同时相的SAR遥感图像；
步骤2：对SAR图像进行预处理；
步骤3：对SAR图像进行滤波处理；
步骤4：对两幅SAR图像进行归一化处理；
步骤5：对两幅SAR图像进行配准处理；
步骤6：分别获得各图像的高阶统计特征图；
步骤7：获得两图像的联合高阶统计特征图；
步骤8：产生高阶统计特征差异图，把步骤6中所获得的高阶统计特征图进行差值或比值运算，获得不同时相SAR图像的特征差异图。
步骤9：在高阶统计特征差异图中产生检测阈值，具体步骤为，
步骤9.1：获取特征差异图的直方图，通过直方图的波峰和波谷确定一个大概的阈值T1；
步骤9.2：采用自适应的期望值最大化(Expectation Maximization,EM)算法对特征差异图进行处理，获得比较精细的阈值T2；
步骤9.3：对步骤9.1和步骤9.2中阈值进行均值处理，获得最后的检测阈值T，即
步骤10：获得由高阶统计特征产生的变化信息；由高阶统计矩特征差异图按步骤9的方法产生的阈值来获取变化信息，具体的判断函数为
式中的C1表示获取的变化信息，(m,n)表示像素的空间坐标，FMD表示高阶矩特征差异图；获得由高阶统计累积量特征差异图产生的变化信息为C2；
步骤11：对两图像的高阶联合统计特征图进行处理，设置阈值，获取变化信息，高阶联合统计特征包括高阶联合矩特征和累积量特征，对联合特征图按步骤9产生阈值，按步骤10进行判断检测，获得的变化信息用C3表示；
步骤12：进行融合处理，变化信息的融合处理按式(21)数学模型进行处理，即
C(m,n)＝α·C1(m,n)+β·C2(m,n)+γ·C3(m,n) (21)
式中的C表示获取的最终的变化信息，α、β和γ表示权系数，它们分别设为0.25，0.25和0.5；
步骤13：输出变化信息图。
2.根据权利要求1所述的基于高阶统计特征的多波段雷达图像变化检测方法，其特征在于，所述的步骤2中对输入SAR图像进行预处理步骤为：对输入不同时相的SAR图像对分别进行辐射和几何校正处理。
3.根据权利要求1所述的基于高阶统计特征的多波段雷达图像变化检测方法，其特征在于，所述的步骤3中对SAR图像进行滤波处理的具体步骤为；
步骤3.1：选择相干性原理对SAR图像进行滤波处理，去除斑点噪声的同时能够有效保持边缘和几何细节信息；具体实现过程如下，先对SAR图像进行傅里叶变换，在频域内进行错位处理，获得两幅子图像，进行反傅里叶变换，进行相干处理，获得相干处理增强的图像，进行坐标变换，恢复到原坐标下的图像；
步骤3.2：以图像中某点(m，n)为中心，设置一个卷积模板，其窗口大小为K×L，对窗口内的所有像素值进行均值处理，处理结果代替像素点(m，n)的值；所有的像素都按上述规则处理，直到全部像素处理完毕，设置窗口参数K＝L＝5；
步骤3.3：把增强后的图像和均值滤波后的图像合成一幅新的SAR图像。
4.根据权利要求1所述的基于高阶统计特征的多波段雷达图像变化检测方法，其特征在于，所述的步骤4中对两幅SAR图像进行归一化处理的具体步骤为：
步骤4.1：先分别计算两幅图像像素灰度值的最大值与最小值；
步骤4.2：按公式(1)计算，实现SAR图像的归一化处理，
式中为归一化后的某像素的灰度值，Ii为归一化前该像素的灰度值，Imax和Imin分别为该图像中像素灰度值的最大值和最小值。
5.根据权利要求1所述的基于高阶统计特征的多波段雷达图像变化检测方法，其特征在于，所述的步骤5中按Hu-SIFT方法对两幅图像实施配准处理。
6.根据权利要求1所述的基于高阶统计特征的多波段雷达图像变化检测方法，其特征在于，所述的步骤6中获得各SAR图像的高阶统计特征图的具体步骤为：
步骤6.1：获取SAR图像的高阶距统计特性，对于一个随机变量x，其k阶原点矩mk和k阶中心矩μk的计算方法分别为式(2)和式(3)，其中式中E[·]表示求数学期望值；
mk＝E[xk] (2)
其中表示均值，一幅图像的大小为N×N，即像素个数为N2，图像的k阶原点矩mk和k阶中心矩μk的计算方法分别为式(4)和式(5)；
式中的为均值，n表示统计特性的阶次；
当n＝1时，获得的是均值特征图；当n＝2时，获得的是方差特征图；
步骤6.2：邻域窗口大小选择，在计算SAR图像的高阶统计特征时，在式(4)和式(5)中，N表示图像的所有像素个数，对于计算某个像素(m，n)的统计特性时，以该像素为中心，选取一个邻域范围内的像素来计算该像素的统计特性，选择计算高阶统计特性的邻域窗口大小设置为5×5；
步骤6.3：获取高阶矩特征图，对提取窗口大小为5的图像三阶和四阶中心矩特征图的计算公式如下，
式中的阶数为n＝3，即用来获取三阶高阶矩统计特征图；当n＝4时，得到计算四阶统计矩特征的公式，即
步骤6.4：获取高阶累积量特征图，随机变量的二阶和三阶累积量与对应的二阶和三阶中心矩的计算方法是一样的，即图像的二阶和三阶累积量特征图就是图像的二阶和三阶矩特征图；随机变量的四阶累积量计算方法为
同样，取每个像素点的5×5邻域块来获取图像的四阶累积量，得到图像的四阶累积量特征图，计算公式如下：
。
7.根据权利要求1所述的基于高阶统计特征的多波段雷达图像变化检测方法，其特征在于，所述的步骤7中获取两图像的联合高阶统计特征图的具体步骤为：
步骤7.1：两幅SAR图像的高阶联合矩特征图的获取，对于任意两个随机变量x和y，它们的k1+k2联合矩计算公式为
则三阶和四阶联合矩的计算公式分别如下：
由式(10)-(12)推导配准后的不同波段的两幅SAR图像中对应像素点的三阶和四阶联合矩特征提取的数学模型分别如式(13)和式(14)表示，
步骤7.2：获取两幅SAR图像的高阶联合累积量特征图，由统计理论知随机变化的三阶和四阶累积量计算公式分别为
c111＝cum(x,y,z)＝E[xyz] (15)
c1111＝cum(x,y,z,t)＝E[xyzt]-E[xy]E[zt]-E[xz]E[yt]-E[xt]E[yz] (16)
从式(15)和式(16)式可知，在求三阶和四阶联合累积量特征时，需要三个或四个变量参数参与求解，选择三个或四个不同波段的SAR图像来计算它们的联合高阶累积量，针对两个临近波段SAR图像来提取变化信息，把某一波段图像视为第三个波段的图像，就求得了它们的三阶联合累积量特征，c111＝E[xyz]＝E[xy2]＝E[x2y]，这样得到的三阶联合累积量与三阶联合矩是一样的；
同理，得到四阶联合累积量统计特征的数学模型为
在两幅配准后的SAR图像中，取移到窗口大小为5的区域，计算两幅图像对应像素的联合四阶累积量，具体计算公式如式(18)所示，
。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810028135.7A CN108257123B (zh) | 2018-01-11 | 2018-01-11 | 基于高阶统计特征的多波段雷达图像变化检测方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810028135.7A CN108257123B (zh) | 2018-01-11 | 2018-01-11 | 基于高阶统计特征的多波段雷达图像变化检测方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN108257123A true CN108257123A (zh) | 2018-07-06 |
| CN108257123B CN108257123B (zh) | 2022-01-18 |
Family
ID=62726349
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201810028135.7A Active CN108257123B (zh) | 2018-01-11 | 2018-01-11 | 基于高阶统计特征的多波段雷达图像变化检测方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN108257123B (zh) |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110121109A (zh) * | 2019-03-22 | 2019-08-13 | 西安电子科技大学 | 面向监控***数字视频实时溯源方法、城市视频监控*** |
| CN112818966A (zh) * | 2021-04-16 | 2021-05-18 | 武汉光谷信息技术股份有限公司 | 多模态遥感影像数据检测方法及*** |
| CN116051398A (zh) * | 2022-11-23 | 2023-05-02 | 广东省国土资源测绘院 | 面向多源多模态遥感数据调查监测特征库构建方法及装置 |
| CN117452367A (zh) * | 2023-12-21 | 2024-01-26 | 西安电子科技大学 | 基于宽带成像雷达的sar载荷辐射信号提取方法及装置 |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101738607A (zh) * | 2009-12-07 | 2010-06-16 | 西安电子科技大学 | 基于聚类的高阶累量交叉熵的sar图像变化检测方法 |
| CN102750705A (zh) * | 2012-07-08 | 2012-10-24 | 西安电子科技大学 | 基于图像融合的光学遥感图像变化检测 |
| CN103927737A (zh) * | 2013-10-31 | 2014-07-16 | 王浩然 | 基于非局部均值的sar图像变化检测方法 |
| CN104680549A (zh) * | 2015-03-24 | 2015-06-03 | 西安电子科技大学 | 基于高阶邻域tmf模型的sar图像变化检测方法 |
-
2018
- 2018-01-11 CN CN201810028135.7A patent/CN108257123B/zh active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101738607A (zh) * | 2009-12-07 | 2010-06-16 | 西安电子科技大学 | 基于聚类的高阶累量交叉熵的sar图像变化检测方法 |
| CN102750705A (zh) * | 2012-07-08 | 2012-10-24 | 西安电子科技大学 | 基于图像融合的光学遥感图像变化检测 |
| CN103927737A (zh) * | 2013-10-31 | 2014-07-16 | 王浩然 | 基于非局部均值的sar图像变化检测方法 |
| CN104680549A (zh) * | 2015-03-24 | 2015-06-03 | 西安电子科技大学 | 基于高阶邻域tmf模型的sar图像变化检测方法 |
Non-Patent Citations (6)
| Title |
|---|
| LIU JIA等: "Change detection in synthetic aperture radar images based on unsupervised artificial immune systems", 《APPLIED SOFT COMPUTING》 * |
| 刘永春等: "基于二维概率密度函数比较的SAR图像变化检测方法", 《电子与信息学报》 * |
| 吕雁: "基于高阶累积量的红外图像时域检测", 《激光与红外》 * |
| 盛钊等: "多波段遥感影像变化检测中差异影像构造研究", 《电子测量技术》 * |
| 范宜强等: "一种基于高阶统计量的图像混和加权滤波方法", 《计算机工程与应用》 * |
| 黄世奇等: "基于小波变换的多时相SAR图像变化检测技术", 《测绘学报》 * |
Cited By (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110121109A (zh) * | 2019-03-22 | 2019-08-13 | 西安电子科技大学 | 面向监控***数字视频实时溯源方法、城市视频监控*** |
| CN112818966A (zh) * | 2021-04-16 | 2021-05-18 | 武汉光谷信息技术股份有限公司 | 多模态遥感影像数据检测方法及*** |
| CN116051398A (zh) * | 2022-11-23 | 2023-05-02 | 广东省国土资源测绘院 | 面向多源多模态遥感数据调查监测特征库构建方法及装置 |
| CN116051398B (zh) * | 2022-11-23 | 2023-09-22 | 广东省国土资源测绘院 | 面向多源多模态遥感数据调查监测特征库构建方法及装置 |
| CN117452367A (zh) * | 2023-12-21 | 2024-01-26 | 西安电子科技大学 | 基于宽带成像雷达的sar载荷辐射信号提取方法及装置 |
| CN117452367B (zh) * | 2023-12-21 | 2024-03-26 | 西安电子科技大学 | 基于宽带成像雷达的sar载荷辐射信号提取方法及装置 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN108257123B (zh) | 2022-01-18 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Noyel et al. | Morphological segmentation of hyperspectral images | |
| Zhang et al. | A categorization of multiscale-decomposition-based image fusion schemes with a performance study for a digital camera application | |
| CN108257123A (zh) | 基于高阶统计特征的多波段雷达图像变化检测方法 | |
| Hariharan et al. | Image fusion and enhancement via empirical mode decomposition | |
| Gagnon et al. | Speckle noise reduction of airborne SAR images with symmetric daubechies wavelets | |
| Ranjbarzadeh et al. | LNPSS: SAR image despeckling based on local and non-local features using patch shape selection and edges linking | |
| Mustafa et al. | Image enhancement technique on contrast variation: a comprehensive review | |
| Baraiya et al. | Comparative study of different methods for brain tumor extraction from MRI images using image processing | |
| Yadav et al. | A review on image fusion methodologies and applications | |
| CN106815819A (zh) | 多策略粮虫视觉检测方法 | |
| Ouma et al. | Multiscale remote sensing data segmentation and post-segmentation change detection based on logical modeling: Theoretical exposition and experimental results for forestland cover change analysis | |
| Lu et al. | Human body flexibility fitness test based on image edge detection and feature point extraction | |
| Chen et al. | Attention-based hierarchical fusion of visible and infrared images | |
| CN114066816B (zh) | 基于静态小波变换提取的sar图像无监督变化检测方法 | |
| Sharma et al. | Change Detection from Remotely Sensed Images Based on Stationary Wavelet Transform. | |
| Akbarpour et al. | Medical image fusion based on nonsubsampled shearlet transform and principal component averaging | |
| Alsulami et al. | Detection and tracking of dim objects in infrared (IR) images using Support Vector Machine | |
| Zhang et al. | Infrared polarization and intensity image fusion algorithm based on the feature transfer | |
| Liu et al. | Infrared aerial small target detection with NSCT and two-dimensional property histogram | |
| Bhujle | Weighted-average fusion method for multiband images | |
| Khan et al. | Efficient target detection in cluttered FLIR imagery | |
| Rahman et al. | Multisensor multichannel image fusion based on fuzzy logic | |
| Newey et al. | Self-supervised speckle reduction GAN for synthetic aperture radar | |
| Wu et al. | Fusion of remote sensing images based on nonsubsampled contourlet transform and region segmentation | |
| Cao et al. | Image Recognition Based on Denoising and Edge Detection |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |