CN114898224B - 一种基于物理散射机制的变化检测方法 - Google Patents

一种基于物理散射机制的变化检测方法 Download PDF

Info

Publication number
CN114898224B
CN114898224B CN202210520453.1A CN202210520453A CN114898224B CN 114898224 B CN114898224 B CN 114898224B CN 202210520453 A CN202210520453 A CN 202210520453A CN 114898224 B CN114898224 B CN 114898224B
Authority
CN
China
Prior art keywords
sar image
span
detected
polarization sar
change
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
Application number
CN202210520453.1A
Other languages
English (en)
Other versions
CN114898224A (zh
Inventor
殷君君
张秀婷
杨健
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Tsinghua University
University of Science and Technology Beijing USTB
Original Assignee
Tsinghua University
University of Science and Technology Beijing USTB
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Tsinghua University, University of Science and Technology Beijing USTB filed Critical Tsinghua University
Priority to CN202210520453.1A priority Critical patent/CN114898224B/zh
Publication of CN114898224A publication Critical patent/CN114898224A/zh
Application granted granted Critical
Publication of CN114898224B publication Critical patent/CN114898224B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V20/00—Scenes; Scene-specific elements
    • G06V20/10—Terrestrial scenes
    • G06V20/13—Satellite images
    • 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/40—Extraction of image or video features
    • G06V10/44—Local feature extraction by analysis of parts of the pattern, e.g. by detecting edges, contours, loops, corners, strokes or intersections; Connectivity analysis, e.g. of connected components
    • 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/70—Arrangements for image or video recognition or understanding using pattern recognition or machine learning
    • G06V10/762—Arrangements for image or video recognition or understanding using pattern recognition or machine learning using clustering, e.g. of similar faces in social networks
    • G06V10/763—Non-hierarchical techniques, e.g. based on statistics of modelling distributions
    • 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/70—Arrangements for image or video recognition or understanding using pattern recognition or machine learning
    • G06V10/764—Arrangements for image or video recognition or understanding using pattern recognition or machine learning using classification, e.g. of video objects
    • 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/70—Arrangements for image or video recognition or understanding using pattern recognition or machine learning
    • G06V10/77—Processing image or video features in feature spaces; using data integration or data reduction, e.g. principal component analysis [PCA] or independent component analysis [ICA] or self-organising maps [SOM]; Blind source separation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
    • G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
    • G06F18/2135—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods based on approximation criteria, e.g. principal component analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/23—Clustering techniques
    • G06F18/232—Non-hierarchical techniques
    • G06F18/2321—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions
    • G06F18/23213—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions with fixed number of clusters, e.g. K-means clustering
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/24—Classification techniques

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Multimedia (AREA)
  • Software Systems (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Computing Systems (AREA)
  • Medical Informatics (AREA)
  • General Health & Medical Sciences (AREA)
  • Evolutionary Computation (AREA)
  • Databases & Information Systems (AREA)
  • Health & Medical Sciences (AREA)
  • Artificial Intelligence (AREA)
  • Probability & Statistics with Applications (AREA)
  • Astronomy & Astrophysics (AREA)
  • Remote Sensing (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明公开了一种基于物理散射机制的变化检测方法，包括：获取待检测的全极化SAR图像；对所述待检测的全极化SAR图像进行预处理；基于预处理后的数据，得到所述待检测的全极化SAR图像对应的基于物理散射机制的参数αB和极化总功率span；基于待检测的全极化SAR图像的αB和span，获得差异图算子，得到差异图；对所述差异图进行PCA和k‑Means聚类，得到变化图，实现变化检测任务。本发明可提高检测率、降低虚警，且可使变化目标保持较好边缘形状。并可检测更加精细的变化细节信息，同时受散射随机性的影响更小，检测性能鲁棒性更强。

Description

一种基于物理散射机制的变化检测方法
技术领域
本发明涉及极化合成孔径雷达图像处理技术领域，特别涉及一种基于物理散射机制的变化检测方法。
背景技术
极化合成孔径雷达(synthetic aperture radar,SAR)是对地观测的一个新兴而有效的工具，近年来对极化SAR图像的目标检测与分类、特征提取等方面研究的不断深入，使得极化SAR图像在更多的应用场景下发挥了巨大的作用。而极化SAR变化检测便是一个重要的研究领域。
极化SAR变化检测是一种针对同一区域多时相极化SAR图像的数据融合技术，目的在于从多时相图像中检测出地表特征发生变化的区域，在灾害监测、城市化监测以及森林采伐等方面具有重要的应用。
通常来说，传统极化SAR图像变化检测可以分为以下几类：
1)基于统计距离和检验统计量的变化检测方法
在众多变化检测方法中，基于统计距离和检验统计量的方法是最早发展起来的，其中最具有代表性的是似然比检验统计量(the Likelyhood Ratio test statistic,LRT)和HLT(the Hotelling-Lawley trace)。Conradsen等人于2003年提出LRT方法，这是在变化检测检验统计距离方面里程碑式的成就。该方法利用LRT对两个复协方差矩阵的相等性进行检验，通过对比假设不变化的零假设和假设变化的备择假设，从而确定是否变化(Conradsen K,Nielsen AA,Schou J,et al.A test statistic in the complex Wishartdistribution and its application to change detection in polarimetric SAR data[J].IEEE Transactions on Geoscience and Remote Sensing,2003,41(1):4-19.)。后来，Akbari等人提出了一种多场景下更简单的变化检测统计距离HLT。该检验假设协方差矩阵服从复Wishart分布，HLT用来衡量两个协方差矩阵之间的相似度(Akbari V,Anfinsen SN,Doulgeris A P,et al.Polarimetric SAR change detection with the complexHotelling–Lawley trace statistic[J].IEEE Transactions on Geoscience andRemote Sensing,2016,54(7):3953-3966.)。由于HLT卓越的性能，Bouhlel等人在2020年根据HLT的思路，提出了行列式比值检验统计量(Determinant Ratio Test,DRT)。其中心思想在于将HLT中统计距离求矩阵的迹转变为求矩阵的行列式，原理虽然简单，但是在UAVSAR数据集上取得了不错的效果(Bouhlel N,Akbari V,Méric S.Change detection inmultilook polarimetric SAR imagery with determinant ratio test statistic[J].IEEE Transactions on Geoscience and Remote Sensing,2020,60:1-15.)。
2)基于极化散射信息的变化检测方法
随着变化检测和极化目标分解技术的快速发展，学者们开始探究利用目标分解或者目标增强等方法得到的极化散射信息来作为变化检测过程中的重要信息。2012年，Tang等人提出了一种基于Yamaguchi四成分分解的方法，主要利用表面体散射成分和二次散射成分，完成定性变化检测。(Tang Y,Zhang H,Wang C,et al.Change detection based onpolarization decomposition using RADARSAT-2quad-pol data[C].//2012IEEEInternational Geoscience and Remote Sensing Symposium.IEEE,2012:6275-6278.)。考虑到前人仅用单一的极化特征来做变化检测，Najafi等人在2019年提出了一种基于组合极化分解的变化检测方法，该方法通过Yamaguchi四成分分解和H/α分解，获得六个表示极化散射信息的分量，再进行后面的操作。(Najafi A,Hasanlou M and Akbari V.Changedetection using distance-based algorithms between synthetic aperture radarpolarimetric decompositions[J].International Journal of Remote Sensing,2019,40(15):6084-6097.)。同年，Pirrone等人提出了一种结合α和span信息的变化检测方法，通过实验验证了该方法的有效性(Pirrone D,De S,Bhattacharya A,et al.AnUnsupervised Approach to Change Detection in Built-Up Areas by MultitemporalPolSAR Images[J].IEEE Geoscience and Remote Sensing Letters,2019,17(11):1914-1918.)。此外，一些学者也借鉴了目标增强的一些方法。2016年，Yin等人根据雷达接收能量的最优极化对比增强技术(the Optimization of Polarimetric ContrastEnhancement)，提出了一种新的变化检测算子，并在地物较为复杂的岛屿数据上取得了不错的效果(Yin J and Yang J.A new change detector in PolSAR imagery[C].//2016IEEE International Geoscience and Remote Sensing Symposium.IEEE,2016.)。因为雷达的接受能量总是为正，所以最小化两个目标A和B之间的能量差相当于最大化目标B和A之间的能量差，也就是说前者的最小化问题可以应用OPCE的流程来解决。因此，优化前者算子的形式，新的实验结果显示该算子能够很好检测出因为雨水产生的建筑物、海湾等自然变化。(Yin J and Yang J.A change detector based on the optimization ofpolarimetric contrast[C].//2017IEEE International Geoscience and RemoteSensing Symposiu.IEEE,2017:5334-5337.)
基于统计距离和统计分布的经典变化检测方法是较早被提出来的，且由于其效果较好，一直在变化检测领域处于主导地位。但是当前变化检测算子普遍存在漏检率高、虚警高以及变化细节很难被检测等问题。
发明内容
本发明提供了一种基于物理散射机制的变化检测方法，以解决变化检测结果检测率较低、虚警较高、很多变化细节不能被检测出来的技术问题，并达到降低受到散射随机性影响的目的。
为解决上述技术问题，本发明提供了如下技术方案：
一方面，本发明提供了一种基于物理散射机制的变化检测方法，包括：
获取待检测的全极化SAR图像；其中，所述待检测的全极化SAR图像包括预设区域的第一时相全极化SAR图像和预设区域的第二时相全极化SAR图像；
对所述待检测的全极化SAR图像进行预处理；
基于预处理后的数据，得到所述待检测的全极化SAR图像对应的基于物理散射机制的参数αB和极化总功率span；
基于待检测的全极化SAR图像的αB和span，获得差异图算子，得到差异图；
对所述差异图进行PCA和k-Means聚类，得到变化图，实现变化检测任务。
进一步地，所述第一时相全极化SAR图像和所述第二时相全极化SAR图像分别包括HH，VV，HV，VH四种极化通道数据；其中，H表示水平极化，V表示垂直极化。
进一步地，对所述待检测的全极化SAR图像进行预处理，包括：
从待处理图像中依次读取每个像素点所对应的S矩阵；其中，所述待处理图像为所述第一时相全极化SAR图像或所述第二时相全极化SAR图像；
利用Pauli基进行共轭相乘，将读取的S矩阵转化为需要的相干矩阵T；
对转化出的T矩阵进行精细Lee滤波，得到滤波处理之后的数据；
所述基于预处理后的数据，得到所述待检测的全极化SAR图像对应的基于物理散射机制的参数αB和极化总功率span，为：
基于精细Lee滤波后的T矩阵，得到所述待处理图像对应的αB和span。
进一步地，S矩阵的表达式为：
Figure BDA0003643139930000045
其中，SHH表示水平方式发射、水平方式接收的极化通道数据，SHV表示水平方式发射、垂直方式接收的极化通道数据，SVH表示垂直方式发射、水平方式接收的极化通道数据，SVV表示垂直方式发射、垂直方式接收的极化通道数据。
进一步地，所述相干矩阵T的表达式为：
Figure BDA0003643139930000041
Figure BDA0003643139930000042
Figure BDA0003643139930000043
其中，k为Pauli基，上标H表示共轭转置，*表示复共轭，<·>表示求均值。
进一步地，对转化出的T矩阵进行精细Lee滤波时，滤波窗口大小为3×3。
进一步地，所述基于预处理后的数据，得到所述待检测的全极化SAR图像对应的基于物理散射机制的参数αB和极化总功率span，包括：
通过下式计算αB：
Figure BDA0003643139930000044
其中，T11、T22、T33为精细Lee滤波后的T矩阵中的对角元素；
通过下式计算span：
span＝T11+T22+T33＝λ1+λ2+λ3
其中，λ1、λ2、λ3为精细Lee滤波后的T矩阵的特征值，并且满足λ1＞λ2＞λ3≥0。
进一步地，基于待检测的全极化SAR图像的αB和span，获得差异图算子，包括：
通过下式基于待检测的全极化SAR图像的αB和span，计算差异图算子Δ(1)：
Figure BDA0003643139930000051
其中，span(1)表示所述第一时相全极化SAR图像的span值，span(2)表示所述第二时相全极化SAR图像的span值，
Figure BDA0003643139930000052
表示所述第一时相全极化SAR图像的αB值，
Figure BDA0003643139930000053
表示所述第二时相全极化SAR图像的αB值。
进一步地，对所述差异图进行PCA和k-Means聚类，得到变化图，包括：
将所述差异图分割为h×h个不重叠的图像块；其中，h取值为4；
在分割出的h×h个不重叠的图像块上使用主成分分析技术PCA提取M个正交特征向量，从而产生特征向量空间；其中，M≤h2，M取值为3；
将h×h的差异图投影到生成的特征向量空间；
利用k-Means算法，将特征向量空间聚为变化和未变化两类；
根据每个像素的特征向量与簇均值特征向量的最小欧氏距离，将差异图的每个像素分到最近的簇内，最终产生变化图。
进一步地，在对所述差异图进行PCA和k-Means聚类，得到变化图，实现变化检测任务之后，所述基于物理散射机制的变化检测方法还包括：
采用预设的参数指标对变化图进行量化评价；其中，所述预设的参数指标包括：准确率、虚警、检测率和Kappa系数中的任意一种或多种的组合。
再一方面，本发明还提供了一种电子设备，其包括处理器和存储器；其中，存储器中存储有至少一条指令，所述指令由处理器加载并执行以实现上述方法。
又一方面，本发明还提供了一种计算机可读存储介质，所述存储介质中存储有至少一条指令，所述指令由处理器加载并执行以实现上述方法。
本发明提供的技术方案带来的有益效果至少包括：
本发明的基于物理散射机制的变化检测方法，使用精细Lee滤波进行数据处理；然后使用αB和span参数作为变化检测算子，获得差异图；最终根据PCA和k-Means进行聚类获得变化图。本发明不仅能够提高检测率、降低虚警，而且能够使变化目标保持较好的边缘形状。并且可以检测更加精细的变化细节信息，同时受到散射随机性的影响更小，所提算法提高了检测性能的鲁棒性。
附图说明
为了更清楚地说明本发明实施例中的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1是本发明实施例提供的基于物理散射机制的变化检测方法的流程图；
图2是本发明实施例提供的两组UAVSAR变化检测实验数据集；其中，(a)、(b)以及(d)、(e)分别是两组数据的伪图，(c)、(f)为两组数据真值图；
图3是本发明实施例提供的基于物理散射机制的变化检测方法在第一组数据集上的差异图结果；其中，(a)是Δ(2)在第一组数据上的差异图，(b)是Δ(1)在第一组数据上的差异图；
图4是本发明实施例提供的基于物理散射机制的变化检测方法在第一组数据集上的变化图结果；其中，(a)是Δ(2)在第一组数据上的变化图，(b)是Δ(1)在第一组数据上的变化图；
图5是本发明实施例提供的基于物理散射机制的变化检测方法在第二组数据集上的差异图结果；其中，(a)是Δ(2)在第二组数据上的差异图，(b)是Δ(1)在第二组数据上的差异图；
图6是本发明实施例提供的基于物理散射机制的变化检测方法在第二组数据集上的变化图结果；其中，(a)是Δ(2)在第二组数据上的变化图，(b)是Δ(1)在第二组数据上的变化图；
图7是本发明方法和对比方法在两组数据集所获变化图的细节对比图；其中，(a)、(b)是两种方法在第一组数据集上的细节对比图；(c)、(d)是两种方法在第二组数据集上的细节对比图。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明实施方式作进一步地详细描述。
第一实施例
本实施例提供了一种基于物理散射机制的变化检测方法，该方法可以由电子设备实现，该方法的执行流程如图1所示，包括以下步骤：
S1，获取待检测的全极化SAR图像；其中，待检测的全极化SAR图像包括预设区域的第一时相全极化SAR图像和预设区域的第二时相全极化SAR图像；
具体地，本实施例中，如图2所示，待检测图像为UAVSAR在预设区域的两组已配准的全极化影像数据，该数据包括HH,VV,HV,VH四种极化通道数据。
S2，对所述待检测的全极化SAR图像进行预处理，接着得到所述待检测的全极化SAR图像对应的基于物理散射机制的参数αB和极化总功率span；
需要说明的是，考虑到极化散射信息对变化检测也是非常重要的，该信息能够反映出极化SAR图像的一些内在属性，所以合理利用极化物理散射信息有助于得到更准确的变化检测结果。对于极化散射信息获取的方法，一般通过极化目标分解或者目标增强。在极化目标分解中，Cloude-Pottier分解(Cloude S Rand Pottier E.An entropy basedclassification scheme for land applications of polarimetric SAR[J].IEEETransactions on Geoscience and Remote Sensing,1997,35(1):68-78.)，即基于特征值-特征向量的方法，它通过几个参数来表征主要平均散射机制。极化熵H、极化角α和极化各向异性A是从散射相干矩阵的特征分解中提取的三个重要参数。在Cloude-Pottier分解中，极化α角仅利用了每个特征向量的第一个元素。但实际中，目标的物理散射特性与特征向量的三个元素都相关。所以本实施例考虑利用一种新的目标特征描述方法αB/ΔαB(YinJ,Moon W M and Yang J.Novel Model-Based Method for Identification ofScattering Mechanisms in Polarimetric SAR Data[J].IEEE Transactions onGeoscience and Remote Sensing,2015.)，该方法能够同时用到此特征向量的三个元素，以更好区分不同的散射机制。具体地，获取参数αB和极化总功率span的过程如下：
S21，从待处理图像中依次读取每个像素点所对应的S矩阵；其中，所述待处理图像为所述第一时相全极化SAR图像或所述第二时相全极化SAR图像；
Figure BDA0003643139930000071
其中，H表示水平极化，V表示垂直极化，四个分量分别表示用两种极化天线接收和发送电磁波的目标后向散射系数。SHH表示水平方式发射、水平方式接收的极化通道数据，SHV表示水平方式发射、垂直方式接收的极化通道数据，SVH表示垂直方式发射、水平方式接收的极化通道数据，SVV表示垂直方式发射、垂直方式接收的极化通道数据。在单站SAR系统中，发射天线和接收天线在同一位置，即对于所有满足互易性的媒质，S矩阵变为对称矩阵，满足SHV＝SVH。
S22，利用Pauli基进行共轭相乘，将读取的S矩阵转化为需要的相干矩阵T；
Figure BDA0003643139930000081
Figure BDA0003643139930000082
Figure BDA0003643139930000083
其中，k为Pauli基，上标H表示共轭转置，*表示复共轭，<·>表示求均值。
S23，对转化出的T矩阵进行精细Lee滤波，得到滤波处理后的数据；其中，对数据使用精细Lee滤波进行预处理需要适当调整滤波窗口的大小；滤波窗口大小的选择方法为：精细Lee滤波能够在对相干斑做较充分滤波的同时保留良好的图像边缘，是一种经典的极化SAR相干斑滤波器，此处滤波窗口选择3×3。
S24，基于精细Lee滤波后的T矩阵，得到对应的αB和span。具体如下：
根据相干矩阵的形式，有：
Figure BDA0003643139930000084
其中，
Figure BDA0003643139930000085
其中，θ是方位角。t11、t22、t33和t12是相干矩阵在0°方位角下的元素。该相干矩阵T秩为3，与共极化通道到交叉极化通道的能量衰减相关。
可以定义一个新参数αB，它是一个旋转不变量：
Figure BDA0003643139930000091
这相当于运用了αB/ΔαB方法中的αB参数。与之相似的是Cloude-Pottier分解方法中的极化α角，但该参数仅利用了每个特征向量的第一个元素。但实际中，目标的物理散射特性与特征向量的三个元素都相关，故我们选择αB作为参数，从而可更好地区分不同的散射机制，更有利于产生贴合实际变化的差异图。
另一个参数是极化总功率span，由于T矩阵是半正定厄米特矩阵，可以通过特征值和特征矢量生成其对角化形式
Figure BDA0003643139930000092
其中，
Figure BDA0003643139930000093
矩阵Σ的对角元素是T矩阵的特征值，且满足λ1＞λ2＞λ3≥0。U是具有正交特征向量列的酉特征向量矩阵。每个特征向量可用4个角来参数化，如下：
Figure BDA0003643139930000099
相干矩阵T可以分解成3个一阶子相干矩阵T1，T2，T3的和，
Figure BDA0003643139930000094
表示单个散射机制。因此，我们可以得到：
span＝T11+T22+T33＝λ1+λ2+λ3
S3，基于待检测的全极化SAR图像的αB和span获得差异图算子，得差异图；
具体地，通过组合αB和span，可得差异图算子Δ(1)为
Figure BDA0003643139930000098
其中，span(1)表示所述第一时相全极化SAR图像的span值，span(2)表示所述第二时相全极化SAR图像的span值，
Figure BDA0003643139930000095
表示所述第一时相全极化SAR图像的αB值，
Figure BDA0003643139930000096
表示所述第二时相全极化SAR图像的αB值。
此外，为了更近一步证明本发明方法的有效性，本实施例利用Cloude-Pottier分解方法中的极化α角来进行对比实验，对比实验的差异图算子Δ(2)表示为：
Figure BDA0003643139930000097
S4，对差异图进行PCA和k-Means聚类，得到变化图，实现变化检测任务。
具体地，根据PCA和k-Means获得变化图的过程主要分为以下五个步骤：
S41，将所述差异图分割为h×h个不重叠的图像块；其中，h取值为4；
S42，在分割出的h×h个不重叠的图像块上使用主成分分析技术PCA提取M个正交特征向量，从而产生特征向量空间；其中，M≤h2，M取值为3；
S43，通过将每个像素周围的数据块投影到特征向量空间，在整个差异图上产生一个特征向量空间；差异图中的每个像素由M维的特征向量表示，即将h×h的差异图投影到生成的特征向量空间；
S44，利用k-Means算法，设置分类个数K＝2，将特征向量空间聚为变化和未变化两类，与ωu和ωc相照应；
S45，根据每个像素的特征向量与簇均值特征向量的最小欧氏距离，将差异图的每个像素分到最近的簇内，最终产生变化图。“1”表示相应像素位置发生了改变，而“0”表示未变化，这个过程可以视为无监督阈值分类，即
Figure BDA0003643139930000101
其中，||||2是欧式距离。
执行本实施例的基于物理散射机制的变化检测方法后，第一组数据集上的差异图结果如图3所示，第一组数据集上的变化图结果如图4所示，第二组数据集上的差异图结果如图5所示，第二组数据集上的变化图结果如图6所示。
本发明方法和对比方法在两组数据集所获变化图的细节对比图如图7所示。
进一步地，在对所述差异图进行PCA和k-Means聚类，得到变化图，实现变化检测任务之后，所述基于物理散射机制的变化检测方法还包括：
S5，采用预设的参数指标对变化图进行量化评价；
其中，所述预设的参数指标包括：准确率(Accuracy)、虚警(False Alarm,FA)、检测率(True Positive Rate,TPR)和Kappa系数。具体地，对变化图进行评价的方法为：对二分类问题，通常可用2×2的混淆矩阵来衡量分类器精确度。
表1分类结果混淆矩阵
Figure BDA0003643139930000102
其中，P表示Positive，N表示Negative，TP、FN、FP、TN分别表示将正样本分为正样本、将正样本分为负样本、将负样本分为正样本和将负样本分为负样本的数目。在变化检测中，通常采用以下的性能评价指标：
Accuracy表示总体的分类正确率：
Figure BDA0003643139930000111
FA表示非变化像素中错分为变化像素的百分比：
Figure BDA0003643139930000112
TPR表示变化像素被正确分为变化像素的比例：
Figure BDA0003643139930000113
Kappa系数是综合考虑了检测率和漏检率的一种衡量检测精度的关键性指标，在0-1范围内，值越大代表变化检测精度越高。
变化检测结果量化评价如表2所示。
表2变化检测结果量化评价
Figure BDA0003643139930000114
综上，本实施例的基于物理散射机制的变化检测方法，使用精细Lee滤波进行数据处理，然后使用αB和span参数作为变化检测算子，获得差异图；最终根据PCA和k-Means进行聚类获得变化图。本方法不仅可提高检测率、降低虚警，且能够使变化目标保持较好的边缘形状。并且可以检测更加精细的变化细节信息，同时受到散射随机性的影响更小，所提算法提高了检测性能的鲁棒性。
第二实施例
本实施例提供一种电子设备，其包括处理器和存储器；其中，存储器中存储有至少一条指令，所述指令由处理器加载并执行，以实现第一实施例的方法。
该电子设备可因配置或性能不同而产生比较大的差异，可以包括一个或一个以上处理器(central processing units，CPU)和一个或一个以上的存储器，其中，存储器中存储有至少一条指令，所述指令由处理器加载并执行上述方法。
第三实施例
本实施例提供一种计算机可读存储介质，该存储介质中存储有至少一条指令，所述指令由处理器加载并执行，以实现上述第一实施例的方法。其中，该计算机可读存储介质可以是ROM、随机存取存储器、CD-ROM、磁带、软盘和光数据存储设备等。其内存储的指令可由终端中的处理器加载并执行上述方法。
此外，需要说明的是，本发明可提供为方法、装置或计算机程序产品。因此，本发明实施例可采用完全硬件实施例、完全软件实施例或结合软件和硬件方面的实施例的形式。而且，本发明实施例可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质上实施的计算机程序产品的形式。
本发明实施例是参照根据本发明实施例的方法、终端设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、嵌入式处理机或其他可编程数据处理终端设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理终端设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理终端设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。这些计算机程序指令也可装载到计算机或其他可编程数据处理终端设备上，使得在计算机或其他可编程终端设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程终端设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
还需要说明的是，在本文中，诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者终端设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者终端设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所述要素的过程、方法、物品或者终端设备中还存在另外的相同要素。
最后需要说明的是，以上所述是本发明优选实施方式，应当指出，尽管已描述了本发明优选实施例，但对于本技术领域的技术人员来说，一旦得知了本发明的基本创造性概念，在不脱离本发明所述原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。所以，所附权利要求意欲解释为包括优选实施例以及落入本发明实施例范围的所有变更和修改。

Claims (4)

1.一种基于物理散射机制的变化检测方法，其特征在于，包括：
获取待检测的全极化SAR图像；其中，所述待检测的全极化SAR图像包括预设区域的第一时相全极化SAR图像和预设区域的第二时相全极化SAR图像；
对所述待检测的全极化SAR图像进行预处理；
基于预处理后的数据，得到所述待检测的全极化SAR图像对应的基于物理散射机制的参数αB和极化总功率span；
基于待检测的全极化SAR图像的αB和span，获得差异图算子，得到差异图；
对所述差异图进行PCA和k-Means聚类，得到变化图，实现变化检测任务；
所述第一时相全极化SAR图像和所述第二时相全极化SAR图像分别包括HH，VV，HV，VH四种极化通道数据；其中，H表示水平极化，V表示垂直极化；
对所述待检测的全极化SAR图像进行预处理，包括：
从待处理图像中依次读取每个像素点所对应的S矩阵；其中，所述待处理图像为所述第一时相全极化SAR图像或所述第二时相全极化SAR图像；
利用Pauli基进行共轭相乘，将读取的S矩阵转化为需要的相干矩阵T；
对转化出的T矩阵进行精细Lee滤波，得到滤波处理之后的数据；
所述基于预处理后的数据，得到所述待检测的全极化SAR图像对应的基于物理散射机制的参数αB和极化总功率span，为：
基于精细Lee滤波后的T矩阵，得到所述待处理图像对应的αB和span；
S矩阵的表达式为：
Figure FDA0003959777440000011
其中，SHH表示水平方式发射、水平方式接收的极化通道数据，SHV表示水平方式发射、垂直方式接收的极化通道数据，SVH表示垂直方式发射、水平方式接收的极化通道数据，SVV表示垂直方式发射、垂直方式接收的极化通道数据；
所述相干矩阵T的表达式为：
Figure FDA0003959777440000012
Figure FDA0003959777440000021
Figure FDA0003959777440000022
其中，k为Pauli基，上标H表示共轭转置，*表示复共轭，<·>表示求均值；
所述基于预处理后的数据，得到所述待检测的全极化SAR图像对应的基于物理散射机制的参数αB和极化总功率span，包括：
通过下式计算αB：
Figure FDA0003959777440000023
其中，T11、T22、T33为精细Lee滤波后的T矩阵中的对角元素；
通过下式计算span：
span＝T11+T22+T33＝λ1+λ2+λ3
其中，λ1、λ2、λ3为精细Lee滤波后的T矩阵的特征值，并且满足λ1＞λ2＞λ3≥0；
基于待检测的全极化SAR图像的αB和span，获得差异图算子，包括：
通过下式基于待检测的全极化SAR图像的αB和span，计算差异图算子Δ(1)：
Figure FDA0003959777440000024
其中，span(1)表示所述第一时相全极化SAR图像的span值，span(2)表示所述第二时相全极化SAR图像的span值，
Figure FDA0003959777440000025
表示所述第一时相全极化SAR图像的αB值，
Figure FDA0003959777440000026
表示所述第二时相全极化SAR图像的αB值。
2.如权利要求1所述的基于物理散射机制的变化检测方法，其特征在于，所述对转化出的T矩阵进行精细Lee滤波时，滤波窗口大小为3×3。
3.如权利要求1所述的基于物理散射机制的变化检测方法，其特征在于，对所述差异图进行PCA和k-Means聚类，得到变化图，包括：
将所述差异图分割为h×h个不重叠的图像块；其中，h取值为4；
在分割出的h×h个不重叠的图像块上使用主成分分析技术PCA提取M个正交特征向量，从而产生特征向量空间；其中，M≤h2，M取值为3；
将h×h的差异图投影到生成的特征向量空间；
利用k-Means算法，将特征向量空间聚为变化和未变化两类；
根据每个像素的特征向量与簇均值特征向量的最小欧氏距离，将差异图的每个像素分到最近的簇内，最终产生变化图。
4.如权利要求1～3任一项所述的基于物理散射机制的变化检测方法，其特征在于，在对所述差异图进行PCA和k-Means聚类，得到变化图，实现变化检测任务之后，所述基于物理散射机制的变化检测方法还包括：
采用预设的参数指标对变化图进行量化评价；其中，所述预设的参数指标包括：准确率、虚警、检测率和Kappa系数中的任意一种或多种的组合。
CN202210520453.1A 2022-05-13 2022-05-13 一种基于物理散射机制的变化检测方法 Active CN114898224B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210520453.1A CN114898224B (zh) 2022-05-13 2022-05-13 一种基于物理散射机制的变化检测方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210520453.1A CN114898224B (zh) 2022-05-13 2022-05-13 一种基于物理散射机制的变化检测方法

Publications (2)

Publication Number Publication Date
CN114898224A CN114898224A (zh) 2022-08-12
CN114898224B true CN114898224B (zh) 2023-04-21

Family

ID=82721946

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210520453.1A Active CN114898224B (zh) 2022-05-13 2022-05-13 一种基于物理散射机制的变化检测方法

Country Status (1)

Country Link
CN (1) CN114898224B (zh)

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10317349B2 (en) * 2015-11-30 2019-06-11 The Boeing Company X-ray scatter systems and methods for detecting structural variations
CN105988113A (zh) * 2016-07-06 2016-10-05 天津大学 极化合成孔径雷达图像变化检测方法
CN107358162B (zh) * 2017-06-12 2020-04-07 西安电子科技大学 基于深度堆栈网络的极化sar影像变化检测方法
CN112580687A (zh) * 2020-11-23 2021-03-30 北京化工大学 基于毫米波全极化sar影像的无监督分类方法
CN113628234B (zh) * 2021-08-16 2024-04-19 西安电子科技大学 基于综合邻域信息的显著性极化sar图像变化检测方法

Also Published As

Publication number Publication date
CN114898224A (zh) 2022-08-12

Similar Documents

Publication Publication Date Title
Matteoli et al. A tutorial overview of anomaly detection in hyperspectral images
He et al. A novel automatic PolSAR ship detection method based on superpixel-level local information measurement
Alonso-Gonzalez et al. Processing multidimensional SAR and hyperspectral images with binary partition tree
Xiang et al. Visual attention and background subtraction with adaptive weight for hyperspectral anomaly detection
CN105989597B (zh) 基于像素选择过程的高光谱图像异常目标检测方法
Sun et al. SAR image change detection based on nonlocal low-rank model and two-level clustering
CN106096506A (zh) 基于子类类间判别双字典的sar目标识别方法
Ma et al. Hyperspectral anomaly detection based on low-rank representation with data-driven projection and dictionary construction
Gui et al. Eigenvalue statistical components-based PU-learning for PolSAR built-up areas extraction and cross-domain analysis
Ren et al. Unsupervised classification of polarimetirc SAR image via improved manifold regularized low-rank representation with multiple features
Zhang et al. Vector phase symmetry for stable hyperbola detection in ground-penetrating radar images
CN114898224B (zh) 一种基于物理散射机制的变化检测方法
Kekre et al. SAR image segmentation using vector quantization technique on entropy images
CN113869119A (zh) 一种多时相sar船只目标追踪方法、系统、设备、介质
Li et al. Block sparse Bayesian learning over local dictionary for robust SAR target recognition
Yan et al. Floating small target detection in sea clutter based on the singular value decomposition of low rank perturbed random matrices
Kang et al. Comparative study on different dimension reduction methods in remote sensing ground object recognition
CN109344767B (zh) 一种基于多方位多特征协同表示的sar目标识别方法
CN111458689A (zh) 基于极化散射中心的多径散射特征分类方法
Ren et al. Unsupervised Classification for Polarimetric SAR Data Using Variational Bayesian Wishart Mixture Model with Inverse Gamma-Gamma Prior
Zhou et al. Research on modified SVM for classification of SAR images
Hu et al. Evaluation of polsar similarity measures with spectral clustering
CN117333468B (zh) 面向多模式时序PolSAR影像的洪涝灾害监测方法
Theoharatos et al. Detection of ship targets in polarimetric SAR data using 2D-PCA data fusion
Yang et al. Unsupervised Classification of Polarimetric SAR Images Based on Scattering Mechanism

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