CN113129292B - 基于迭代马尔科夫的合成孔径雷达影像变化检测方法 - Google Patents
基于迭代马尔科夫的合成孔径雷达影像变化检测方法 Download PDFInfo
- Publication number
- CN113129292B CN113129292B CN202110459965.7A CN202110459965A CN113129292B CN 113129292 B CN113129292 B CN 113129292B CN 202110459965 A CN202110459965 A CN 202110459965A CN 113129292 B CN113129292 B CN 113129292B
- Authority
- CN
- China
- Prior art keywords
- image
- binary image
- binary
- class
- iterative
- 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.)
- Active
Links
- 230000008859 change Effects 0.000 title claims abstract description 74
- 238000001514 detection method Methods 0.000 title claims abstract description 46
- 238000000034 method Methods 0.000 claims abstract description 67
- 230000004927 fusion Effects 0.000 claims abstract description 15
- 238000013139 quantization Methods 0.000 claims abstract description 10
- 230000011218 segmentation Effects 0.000 claims description 7
- 239000011159 matrix material Substances 0.000 claims description 6
- 238000010586 diagram Methods 0.000 claims description 5
- 238000010276 construction Methods 0.000 claims description 4
- 238000009826 distribution Methods 0.000 claims description 4
- 238000000354 decomposition reaction Methods 0.000 claims description 3
- 238000003064 k means clustering Methods 0.000 claims description 3
- 239000013598 vector Substances 0.000 claims description 3
- 238000012804 iterative process Methods 0.000 claims 1
- 230000008569 process Effects 0.000 abstract description 3
- 238000004088 simulation Methods 0.000 abstract description 3
- 238000012545 processing Methods 0.000 description 4
- 238000012544 monitoring process Methods 0.000 description 3
- 230000009286 beneficial effect Effects 0.000 description 2
- 230000001427 coherent effect Effects 0.000 description 2
- 230000000052 comparative effect Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000004364 calculation method Methods 0.000 description 1
- 238000010835 comparative analysis Methods 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 238000003384 imaging method Methods 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
Images
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
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/25—Fusion techniques
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/11—Region-based segmentation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/136—Segmentation; Edge detection involving thresholding
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/143—Segmentation; Edge detection involving probabilistic approaches, e.g. Markov random field [MRF] modelling
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
- G06V10/00—Arrangements for image or video recognition or understanding
- G06V10/20—Image preprocessing
- G06V10/28—Quantising the image, e.g. histogram thresholding for discrimination between background and foreground patterns
-
- 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
-
- 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/20—Special algorithmic details
- G06T2207/20076—Probabilistic image processing
-
- 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/20—Special algorithmic details
- G06T2207/20212—Image combination
- G06T2207/20221—Image fusion; Image merging
Landscapes
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Data Mining & Analysis (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Evolutionary Biology (AREA)
- Evolutionary Computation (AREA)
- Bioinformatics & Computational Biology (AREA)
- General Engineering & Computer Science (AREA)
- Artificial Intelligence (AREA)
- Life Sciences & Earth Sciences (AREA)
- Multimedia (AREA)
- Quality & Reliability (AREA)
- Probability & Statistics with Applications (AREA)
- Software Systems (AREA)
- Image Analysis (AREA)
Abstract
一种基于迭代马尔科夫的合成孔径雷达影像变化检测方法,由构建对数比值差异图D、确定先验概率和先验二值图、构建特征空间、确定特征数值的量化级、确定似然函数、确定变化类BC的后验概率、确定后验概率对应的二值图、马尔科夫迭代融合、输出最终结果步骤组成。本发明采用对同一地区、不同时刻的两幅合成孔径雷达影像进行变化检测,将一次贝叶斯过程中的先验概率图二值图和后验概率二值图进行迭代融合,更新,用后验概率代替先验概率进行马尔科夫融合,更新变化信息。融合了先验概率图与后验概率图,能检测到变化信息。本发明经对比仿真实验结果表明,本发明具有变化检测结果噪声更少,检测到的变化部分更接近真实结果优点。
Description
技术领域
本发明属于遥感影像应用技术领域,具体地涉及到同一地区不同时刻两幅合成孔径雷达影像的变化检测方法。
背景技术
用多时相遥感影像探测来监测地球表面发生的变化是遥感技术最重要的应用之一。其中合成孔径雷达(SAR)因为具有全天候、全天时的特点,可方便地获得同一地区不同时刻的图像,被广泛应用于地震区域的定位和灾害评估、对农作物生长状况的监测、对土地使用的监测等等。SAR图像变化检测也成了它的一个重要的应用领域,通过对不同时期图像的比较分析,根据图像之间的差异得到所需要的变化信息。
由于SAR系统采用合成相干方式成像,乘性相干斑总是伴随SAR影像而存在。因此,如何有效地设计基于SAR影像的变化检测方法得到了越来越多学者的关注。
现有的SAR图像处理存在的主要技术问题是:SAR影像的乘性噪声对后续的构造差异图及差异图分类都有较大影像,会导致变化检测的精度降低,造成错分、误分的现象。
在图像处理技术领域,进行合成孔径雷达影像处理时,要充分考虑相干斑的影响,如何做到既能消除相干斑的影响又能保留真实变化信息不受到损失,是图像处理技术领域当前需迫切解决的一个技术问题。
发明内容
本发明所要解决的技术问题在于提供一种鲁棒性好、正确率高、分类精度高的基于迭代马尔科夫的合成孔径雷达影像变化检测方法。
解决上述技术问题所采用的技术方案由下述步骤组成:
(1)构建对数比值差异图D
用对数比值方法将同一地区不同时刻的合成孔径雷达影像I1和合成孔径雷达影像I2构建对数比值差异图D:
(2)确定先验概率和先验二值图
其中,l表示当前迭代次数,δ是参数,δ≠0,Df为归一差异图,T为将先验概率经过最大类间方差法进行阈值分割得到的分割阈值,将先验概率中大于阈值T的像素点置为1,小于阈值T的像素点置为0,得到先验二值图Tp。
(3)构建特征空间
采用3×3不重叠分块方法将差异图Df分成n个图像块,将图像块拉成列向量构建成自相关矩阵,并采用奇异值分解方法构建成K维特征空间:
其中,λj是自相关矩阵的第j个特征值,k为像素点个数,k∈[1,2,…,8],n为有限的正整数。
(4)确定特征数值的量化级
将特征空间中的像素点用K-means聚类方法分为变化类C、非变化类UC以及混淆类MC三类,并通过下式将变化类C与非变化类UC分布间的K-L距离
DL(C,UC)=K(C|UC)+K(UC|C)
其中,K(C|UC)表示从变化类C到非变化类UC的KL散度,K(UC|C)表示从非变化类UC到变化类C的KL散度,DL(C,UC)表示变化类C与非变化UC之间的K-L距离;将K-L距离最大化,获得变化类C与非变化类UC分布间相应特征值的量化级。
(5)确定似然函数
将量化级采用直方图方法得到先验二值图Tp中变化类BC的似然函数Pl(Fi|BC)和非变化类BUC的似然函数Pl(Fi|BUC):
其中表示像素i处的第k个特征属于变化类BC的条件概率,表示像素i处的第k个特征属于非变化类BUC的条件概率,是像素i处的第k个特征,ωk表示第k维特征的权重,和是第k维特征的变化部分和非变化部分的均值。
(6)确定变化类BC的后验概率
按下式确定差异图Df中各像素属于变化类BC的后验概率Pl(BC|Fo):
(7)确定后验概率对应的二值图
通过最大类间方差阈值法对后验概率Pl(BC|Fi)进行二值化,得到后验概率对应的后验二值图Ta。
(8)马尔科夫迭代融合
利用能量最小化原则,采用迭代方式融合先验二值图Tp和后验二值图Ta,得到融合后二值图Y。
(9)输出最终结果
Tp=Y
重复步骤(5)-(8),直至到最大迭代次数lmax,停止迭代,输出变化检测结果图Y。
在本发明的步骤(8)中,所述的迭代方式融合先验二值图Tp和后验二值图Ta的方法为:
1)用先验二值图Tp和后验二值图Ta,将能量函数U1最小化,得到初始化二值图B0,能量函数U1如下:
其中,bmn表示图像B0中像素点(m,n)的类标,AS(m,n)表示二值图A中像素点(m,n)的邻域,A(p,q)表示当前位置处属于AS(m,n)中的像素点,T表示计算二值图A所使用的阈值,xpq表示二值图A所对应差异图在(p,q)位置处的强度值,γ是一个调节参数,γ>0。
2)将能量函数U2(Bt,Bt-1,Tp,Ta)最小化得到二值图Bt,能量函数U2如下:
3)重复执行步骤2),当二值图Bt-1和二值图Bt中不同标记的像素个数小于给定的阈值Td或者达到最大迭代次数tmax时,停止迭代,Td的取值小于总体像素点个数,二值图Bt为融合后二值图Y。
在本发明的(8)步骤的1)步骤中,所述的γ取值为(0,10];在(8)步骤的2)步骤中,所述的β取值为(0,1]。
在本发明的(8)步骤的1)步骤中,所述的γ取值最佳为5。
在本发明的(8)步骤的3)步骤中,所述的tmax取值为500。
在本发明的步骤(9)中,所述的lmax取值最佳为3。
本发明的有益效果如下:
本发明采用基于迭代马尔科夫融合方法对同一地区、不同时刻的两幅合成孔径雷达影像进行变化检测,马尔科夫模型融合空间信息,可以将一次贝叶斯过程中的先验概率图二值图和后验概率二值图进行迭代融合,再通过迭代更新,利用后验概率代替先验概率继续进行马尔科夫融合,不断更新变化信息。本发明融合了先验概率图与后验概率图,能够最大限度检测到变化信息。本发明方法与PCA+K-means方法、TLFCM方法、MRF方法进行了对比仿真实验,实验结果表明,本发明具有变化检测结果噪声更少,检测到的变化部分更接近真实结果优点。
附图说明
图1是本发明实施例1的流程图。
图2是实施例1变化前的影像。
图3是实施例1变化后的影像。
图4是本发明对图2、图3进行变化检测的结果图。
图5是PCA+K-means方法对图2、图3进行变化检测的结果图。
图6是TLFCM方法对图2、图3进行变化检测的结果图。
图7是MRF方法对图2、图3进行变化检测的结果图。
图8是实施例1真实的变化图。
具体实施方式
下面结合附图和实施例对本发明进进一步详细说明,但本发明不限于下述实施方式。
实施例1
以1997年5月和1997年8月加拿大Ottawa地区的合成孔径雷达影像为例,本实施例的基于迭代马尔科夫的合成孔径雷达影像变化检测方法步骤如下(参见图1):
(1)构建对数比值差异图D
用对数比值方法将同一地区不同时刻的合成孔径雷达影像I1和合成孔径雷达影像I2构建对数比值差异图D:
(2)确定先验概率和先验二值图
其中,l表示当前迭代次数,δ是参数,本实施例的δ取值为4×10-3,Df为归一差异图,T为将先验概率经过最大类间方差法进行阈值分割得到的分割阈值,将先验概率中大于阈值T的像素点置为1,小于阈值T的像素点置为0,得到先验二值图Tp。
(3)构建特征空间
采用3×3不重叠分块方法将差异图Df分成n个图像块,将图像块拉成列向量构建成自相关矩阵,并采用奇异值分解方法构建成K维特征空间:
其中λj是自相关矩阵的第j个特征值,k为像素点个数,k∈[1,2,…,8],n取值为12168。
(4)确定特征数值的量化级
将特征空间中的像素点用K-means聚类方法分为变化类C、非变化类UC以及混淆类MC三类,并通过下式将变化类C与非变化类UC分布间的K-L距离
DL(C,UC)=K(C|UC)+K(UC|C)
其中,K(C|UC)表示从变化类C到非变化类UC的KL散度,K(UC|C)表示从非变化类UC到变化类C的KL散度,DL(C,UC)表示变化类C与非变化UC之间的K-L距离;将K-L距离最大化,获得变化类C与非变化类UC分布间相应特征值的量化级。
(5)确定似然函数
将量化级采用直方图方法得到先验二值图Tp中变化类BC的似然函数Pl(Fi|BC)和非变化类BUC的似然函数Pl(Fi|BUC):
其中表示像素i处的第k个特征属于变化类BC的条件概率,表示像素i处的第k个特征属于非变化类BUC的条件概率,是像素i处的第k个特征,ωk表示第k维特征的权重,和是第k维特征的变化部分和非变化部分的均值。
(6)确定变化类BC的后验概率
按下式确定差异图Df中各像素属于变化类BC的后验概率Pl(BC|Fi)::
(7)确定后验概率对应的二值图
通过最大类间方差阈值法对后验概率Pl(BC|Fi)进行二值化,得到后验概率对应的后验二值图Ta。
(8)马尔科夫迭代融合
利用能量最小化原则,采用迭代方式融合先验二值图Tp和后验二值图Ta,得到融合后二值图Y。
本实施例迭代方式融合先验二值图Tp和后验二值图Ta的方法为:
1)用先验二值图Tp和后验二值图Ta,将能量函数U1最小化,得到初始化二值图B0,能量函数U1如下:
其中,bmn表示图像B0中像素点(m,n)的类标,AS(m,m)表示二值图A中像素点(m,n)的邻域,A(p,q)表示当前位置处属于AS(m,n)中的像素点,T表示计算二值图A所使用的阈值,xpq表示二值图A所对应差异图在(p,q)位置处的强度值,γ是一个调节参数,γ>0,本实施例的M取值为313,N取值为352,γ取值为5。
2)将能量函数U2(Bt,Bt-1,Tp,Ta)最小化得到二值图Bt,能量函数U2如下:
3)重复执行步骤(2),二值图Bt-1和二值图Bt中不同标记的像素个数小于给定的阈值Td或者达到最大迭代次数tmax时,停止迭代,Td的取值小于总体像素点个数,二值图Bt为融合后二值图Y,本实施例的Td取值为2,tmax取值为500。
(9)输出最终结果
Tp=Y
重复步骤(5)-(8),直至到最大迭代次数lmax,停止迭代,本实施例的lmax取值为3,输出变化检测结果图Y。
实施例2
以1997年5月和1997年8月加拿大Ottawa地区的合成孔径雷达影像为例,本实施例的基于迭代马尔科夫的合成孔径雷达影像变化检测方法步骤如下:
(1)-(7)步骤与实施例1相同。
(8)马尔科夫迭代融合
利用能量最小化原则,采用迭代方式融合先验二值图Tp和后验二值图Ta,得到融合后二值图Y。
本实施例迭代方式融合先验二值图Tp和后验二值图Ta的方法为:
1)用先验二值图Tp和后验二值图Ta,将能量函数U1最小化,得到初始化二值图B0,能量函数U1如下:
其中,bmn表示图像B0中像素点(m,n)的类标,AS(m,n)表示二值图A中像素点(m,n)的邻域,A(p,q)表示当前位置处属于AS(m,n)中的像素点,T表示计算二值图A所使用的阈值,xpq表示二值图A所对应差异图在(p,q)位置处的强度值,γ是一个调节参数,γ>0,本实施例的M取值为313,N取值为352,γ取值为0.5。
2)将能量函数U2(Bt,Bt-1,Tp,Ta)最小化得到二值图Bt,能量函数U2如下:
3)重复执行步骤(2),二值图Bt-1和二值图Bt中不同标记的像素个数小于给定的阈值Td或者达到最大迭代次数tmax时,停止迭代,Td的取值小于总体像素点个数,二值图Bt为融合后二值图Y,本实施例的Td取值为2,tmax取值为500。
其它步骤与实施例1相同,输出变化检测结果图Y。
实施例3
以1997年5月和1997年8月加拿大Ottawa地区的合成孔径雷达影像为例,本实施例的基于迭代马尔科夫的合成孔径雷达影像变化检测方法步骤如下:
(1)-(7)步骤与实施例1相同。
(8)马尔科夫迭代融合
利用能量最小化原则,采用迭代方式融合先验二值图Tp和后验二值图Ta,得到融合后二值图Y。
本实施例迭代方式融合先验二值图Tp和后验二值图Ta的方法为:
1)用先验二值图Tp和后验二值图Ta,将能量函数U1最小化,得到初始化二值图B0,能量函数U1如下:
其中,bmn表示图像B0中像素点(m,n)的类标,AS(m,n)表示二值图A中像素点(m,n)的邻域,A(p,q)表示当前位置处属于AS(m,n)中的像素点,T表示计算二值图A所使用的阈值,xpq表示二值图A所对应差异图在(p,q)位置处的强度值,γ是一个调节参数,γ>0,本实施例的M取值为313,N取值为352,γ取值为10。
2)将能量函数U2(Bt,Bt-1,Tp,Ta)最小化得到二值图Bt,能量函数U2如下:
3)重复执行步骤(2),二值图Bt-1和二值图Bt中不同标记的像素个数小于给定的阈值Td或者达到最大迭代次数tmax时,停止迭代,Td的取值小于总体像素点个数,二值图Bt为融合后二值图Y,本实施例的Td取值为2,tmax取值为500。
其它步骤与实施例1相同,输出变化检测结果图Y。
为了验证本发明的有益效果,发明人采用本发明实施例1的基于迭代马尔科夫的合成孔径雷达影像变化检测方法与PCA+K-means方法、TLFCM方法、MRF方法进行了对比仿真实验,实验时,选取1997年5月和1997年8月加拿大Ottawa地区的合成孔径雷达影像,图像大小为352×313,如图2、图3。本发明方法对图2、图3进行变化检测的结果图见图4,PCA+K-means方法对图2、图3进行变化检测的结果见图5,TLFCM方法对图2、图3进行变化检测的结果见图6,MRF方法对图2、图3进行变化检测的结果见图7。
4种方法的接测结果与真实的变化图像图8进行比较,按下式计算误检像素总数(OE)、检测准确率(Pcc)、Kappa系数(Kappa):
OE=FP+FN
其中,FP是非变化类被误检为变化类的像素个数,FN是变化类被误检为非变化类的像素个数,TP是正确检测为变化类像素个数,TN是正确检测为非变化类像素个数,N为像素总数,NC为真实的变化像素个数,NU为真实的非变化像素个数。
计算结果见表1、图4-7。
表1本发明方法和三种已知方法对两种图像的检测结果
检测方法 | FP(个) | FN(个) | OE(个) | Pcc(%) | Kappa |
PCA+K-means | 956 | 1941 | 2897 | 97.371 | 0.89083 |
TLFCM | 684 | 2638 | 3322 | 96.985 | 0.87147 |
MRF | 1376 | 1499 | 2875 | 97.391 | 0.8941 |
本发明 | 1426 | 1270 | 2696 | 97.553 | 0.90142 |
由表1可见,本发明获得了OE、Pcc和kappa的最佳值,OE是误检像素个数总数,本发明OE数值较其他三种方法最小,说明本发明误检的像素个数最少。Pcc可以体现分类的正确性,本发明的Pcc较其他三种方法取得了最优值。Kappa系数是能更准确衡量分类准确度的参数,本发明具有较高的Kappa系数,说明在实际的多时相SAR图像的变化检测中具有较好的效果。先用马尔科夫迭代融合先验二值图和后验二值图,再通过贝叶斯继续更新先验二值图,在迭代的过程中融合空间信息具有较好的变化检测结果。本发明与现有的技术相比,考虑了空间信息,具有方法简单、鲁棒性高、准确度高等优点。
在图4中,白色部分为本发明检测出的发生变化区域,黑色部分为本发明检测出的未发生变化区域。与真实的变化部分参考图图8相比,本发明能够更加完整地检测到变化区域。
在图5中,白色部分为PCA+K-means方法检测出的发生变化区域,黑色部分为PCA+K-means方法检测出的未发生变化区域。与真实的变化部分参考图图8相比,PCA+K-means方法能检测出大部分变化区域,但仍存在漏检的情况。
在图6中,白色部分为TLFCM方法检测出的发生变化区域,黑色部分为TLFCM方法检测出的未发生变化区域。与真实的变化部分参考图图8相比,TLFCM方法检测的结果虽然去除了噪声,但也缺失了部分变化信息。
在图7中,白色部分为MRF方法检测出的发生变化区域,黑色部分为MRF方法检测出的未发生变化区域。与真实的变化部分参考图图8相比,MRF方法也存在检测变化区域不完整的情况,另外检测的结果中存在一些噪声,影响变化检测结果的准确性。
以上实验结果表明,实施例1的方法与三种已知方法相比,本发明具有变化检测结果噪声更少,检测到的变化部分更接近真实结果的优点。
Claims (6)
1.一种基于迭代马尔科夫的合成孔径雷达影像变化检测方法,其特征在于由下述步骤组成:
(1)构建对数比值差异图D
用对数比值方法将同一地区不同时刻的合成孔径雷达影像I1和合成孔径雷达影像I2构建对数比值差异图D:
(2)确定先验概率和先验二值图
其中,l表示当前迭代次数,δ是参数,δ≠0,Df为归一差异图,T为将先验概率经过最大类间方差法进行阈值分割得到的分割阈值,将初始化的先验概率中大于阈值T的像素点置为1,小于阈值T的像素点置为0,得到先验二值图Tp;
(3)构建特征空间
采用3×3不重叠分块方法将差异图Df分成n个图像块,将图像块拉成列向量构建成自相关矩阵,并采用奇异值分解方法构建成K维特征空间:
其中,λj是自相关矩阵的第j个特征值,k为像素点个数,k∈[1,2,…,8],n为有限的正整数;
(4)确定特征数值的量化级
将特征空间中的像素点用K-means聚类方法分为变化类C、非变化类UC以及混淆类MC三类,并通过下式计算变化类C与非变化类UC分布间的K-L距离:
DL(C,UC)=K(C|UC)+K(UC|C)
其中,K(C|UC)表示从变化类C到非变化类UC的KL散度,K(UC|C)表示从非变化类UC到变化类C的KL散度,DL(C,UC)表示变化类C与非变化UC之间的K-L距离;将K-L距离最大化,获得变化类C与非变化类UC分布间相应特征值的量化级;
(5)确定似然函数
将量化级采用直方图方法得到先验二值图Tp中变化类BC的似然函数Pl(Fi|BC)和非变化类BUC的似然函数Pl(Fi|BUC):
其中表示像素i处的第k个特征属于变化类BC的条件概率,表示像素i处的第k个特征属于非变化类BUC的条件概率,是像素i处的第k个特征,ωk表示第k维特征的权重,和是第k维特征的变化部分和非变化部分的均值;
(6)确定变化类BC的后验概率
按下式确定差异图Df中各像素属于变化类BC的后验概率Pl(BC|fi)::
(7)确定后验概率对应的二值图
通过最大类间方差阈值法对后验概率Pl(BC|Fi)进行二值化,得到后验概率对应的后验二值图Ta;
(8)马尔科夫迭代融合
利用能量最小化原则,采用迭代方式融合先验二值图Tp和后验二值图Ta,得到融合后二值图Y;
(9)输出最终结果
Tp=Y
重复步骤(5)-(8),直至到最大迭代次数lmax,停止迭代,输出变化检测结果图Y。
2.根据权利要求1所述的基于迭代马尔科夫的合成孔径雷达影像变化检测方法,其特征在于在步骤(8)中,所述的迭代方式融合先验二值图Tp和后验二值图Ta的方法为:
1)用先验二值图Tp和后验二值图Ta,将能量函数U1最小化,得到初始化二值图B0,能量函数U1如下:
其中,bmn表示图像B0中像素点(m,n)的类标,As(m,n)表示二值图A中像素点(m,n)的邻域,A(p,q)表示当前位置处属于AS(m,n)中的像素点,T表示计算二值图A所使用的阈值,xpq表示二值图A所对应差异图在(p,q)位置处的强度值,γ是一个调节参数,γ>0;
2)将能量函数U2(Bt,Bt-1,Tp,Ta)最小化得到二值图Bt,能量函数U2如下:
3)重复执行步骤2),当二值图Bt-1和二值图Bt中不同标记的像素个数小于给定的阈值Td或者达到最大迭代次数tmax时,停止迭代,Td的取值小于总体像素点个数,二值图Bt为融合后二值图Y。
3.根据权利要求2所述的基于迭代马尔科夫的合成孔径雷达影像变化检测方法,其特征在于:在(8)步骤的1)步骤中,所述的γ取值为(0,10];在(8)步骤的2)步骤中,所述的β取值为(0,1]。
4.根据权利要求2或3所述的基于迭代马尔科夫的合成孔径雷达影像变化检测方法,其特征在于:在(8)步骤的1)步骤中,所述的γ取值为5。
5.根据权利要求2所述的基于迭代马尔科夫的合成孔径雷达影像变化检测方法,其特征在于:在(8)步骤的3)步骤中,所述的tmax取值为500。
6.根据权利要求1所述的基于迭代马尔科夫的合成孔径雷达影像变化检测方法,其特征在于:在步骤(9)中,所述的lmax取值为3。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110459965.7A CN113129292B (zh) | 2021-04-27 | 2021-04-27 | 基于迭代马尔科夫的合成孔径雷达影像变化检测方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110459965.7A CN113129292B (zh) | 2021-04-27 | 2021-04-27 | 基于迭代马尔科夫的合成孔径雷达影像变化检测方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN113129292A CN113129292A (zh) | 2021-07-16 |
CN113129292B true CN113129292B (zh) | 2023-04-07 |
Family
ID=76780224
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202110459965.7A Active CN113129292B (zh) | 2021-04-27 | 2021-04-27 | 基于迭代马尔科夫的合成孔径雷达影像变化检测方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN113129292B (zh) |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106203521A (zh) * | 2016-07-15 | 2016-12-07 | 西安电子科技大学 | 基于差异图自步学习的sar图像变化检测方法 |
CN107644413A (zh) * | 2017-08-25 | 2018-01-30 | 西安电子科技大学 | 基于邻域比值和自步学习的sar图像变化区域检测方法 |
Family Cites Families (8)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101694719B (zh) * | 2009-10-13 | 2011-12-07 | 西安电子科技大学 | 基于非参数密度估计的遥感图像变化检测方法 |
CN104751453A (zh) * | 2015-03-11 | 2015-07-01 | 西安电子科技大学 | 基于cuda和平稳小波变换的sar图像变化检测方法 |
CN104778717A (zh) * | 2015-05-05 | 2015-07-15 | 西安电子科技大学 | 基于导向差异图的sar图像变化检测方法 |
CN106713935B (zh) * | 2017-01-09 | 2019-06-11 | 杭州电子科技大学 | 一种基于贝叶斯决策的hevc块划分快速方法 |
CN107368781B (zh) * | 2017-06-09 | 2019-08-20 | 陕西师范大学 | 基于子空间划分的合成孔径雷达影像变化检测方法 |
CN108764119B (zh) * | 2018-05-24 | 2022-03-18 | 西安电子科技大学 | 基于迭代最大类间方差的sar图像变化检测方法 |
US10984532B2 (en) * | 2018-08-24 | 2021-04-20 | Ordnance Survey Limited | Joint deep learning for land cover and land use classification |
CN109697474B (zh) * | 2018-12-30 | 2023-05-23 | 陕西师范大学 | 基于迭代贝叶斯的合成孔径雷达影像变化检测方法 |
-
2021
- 2021-04-27 CN CN202110459965.7A patent/CN113129292B/zh active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106203521A (zh) * | 2016-07-15 | 2016-12-07 | 西安电子科技大学 | 基于差异图自步学习的sar图像变化检测方法 |
CN107644413A (zh) * | 2017-08-25 | 2018-01-30 | 西安电子科技大学 | 基于邻域比值和自步学习的sar图像变化区域检测方法 |
Also Published As
Publication number | Publication date |
---|---|
CN113129292A (zh) | 2021-07-16 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN103353989B (zh) | 基于先验和融合灰度与纹理特征的sar图像变化检测方法 | |
CN105787935B (zh) | 一种基于Gamma分布的模糊聚类SAR图像分割方法 | |
CN108717694B (zh) | 基于模糊c均值聚类的电阻抗层析成像图像质量评价方法 | |
US20080040083A1 (en) | System and Method for Solid Component Evaluation in Mixed Ground Glass Nodules | |
CN111695468B (zh) | 一种基于K-shot学习的极化SAR地物分类方法 | |
CN111445484B (zh) | 一种基于图像级标注的工业图像异常区域像素级分割方法 | |
CN104732552B (zh) | 基于非平稳条件场的sar图像分割方法 | |
CN104680536B (zh) | 利用改进的非局部均值算法对sar图像变化的检测方法 | |
CN106023229B (zh) | 结合半高斯模型与高斯模型的sar影像变化检测方法 | |
CN106600607A (zh) | 一种基于水平集分割极化sar影像的水体精确提取方法 | |
CN104680549B (zh) | 基于高阶邻域tmf模型的sar图像变化检测方法 | |
CN116258877A (zh) | 土地利用场景相似度变化检测方法、装置、介质及设备 | |
CN112508963A (zh) | 一种基于模糊c均值聚类的sar图像分割方法 | |
CN113822361B (zh) | 一种基于汉明距离的sar图像相似程度度量方法和系统 | |
Luo et al. | Unsupervised multiregion partitioning of fully polarimetric SAR images with advanced fuzzy active contours | |
CN109697474B (zh) | 基于迭代贝叶斯的合成孔径雷达影像变化检测方法 | |
CN113129292B (zh) | 基于迭代马尔科夫的合成孔径雷达影像变化检测方法 | |
CN112348750B (zh) | 基于阈值融合和邻域投票的sar图像变化检测方法 | |
Irving et al. | Classification of targets in synthetic aperture radar imagery via quantized grayscale matching | |
Li et al. | Crack Detection and Recognition Model of Parts Based on Machine Vision. | |
Lv et al. | Fuzzy active contour model using fractional-order diffusion based edge indicator and fuzzy local fitted image | |
CN115829996A (zh) | 基于深度特征图的无监督合成孔径雷达图像变化检测方法 | |
Ali et al. | Object based image segmentation using fuzzy clustering | |
CN106709921B (zh) | 一种基于空间Dirichlet混合模型的彩色图像分割方法 | |
McCarter | The Kernel Density Integral Transformation |
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 |