CN110517277B - 基于PCANet与高阶CRF的SAR图像分割方法 - Google Patents

基于PCANet与高阶CRF的SAR图像分割方法 Download PDF

Info

Publication number
CN110517277B
CN110517277B CN201910718215.XA CN201910718215A CN110517277B CN 110517277 B CN110517277 B CN 110517277B CN 201910718215 A CN201910718215 A CN 201910718215A CN 110517277 B CN110517277 B CN 110517277B
Authority
CN
China
Prior art keywords
sar image
pixel point
representing
pixel
stage
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
CN201910718215.XA
Other languages
English (en)
Other versions
CN110517277A (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.)
Xidian University
Original Assignee
Xidian University
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 Xidian University filed Critical Xidian University
Priority to CN201910718215.XA priority Critical patent/CN110517277B/zh
Publication of CN110517277A publication Critical patent/CN110517277A/zh
Application granted granted Critical
Publication of CN110517277B publication Critical patent/CN110517277B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • 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/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
    • G06F18/214—Generating training patterns; Bootstrap methods, e.g. bagging or boosting
    • G06F18/2155—Generating training patterns; Bootstrap methods, e.g. bagging or boosting characterised by the incorporation of unlabelled data, e.g. multiple instance learning [MIL], semi-supervised techniques using expectation-maximisation [EM] or naïve labelling
    • 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
    • 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/12—Edge-based segmentation
    • 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
    • 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/20081—Training; Learning
    • 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/25—Determination of region of interest [ROI] or a volume of interest [VOI]
    • 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/26—Segmentation of patterns in the image field; Cutting or merging of image elements to establish the pattern region, e.g. clustering-based techniques; Detection of occlusion
    • G06V10/267—Segmentation of patterns in the image field; Cutting or merging of image elements to establish the pattern region, e.g. clustering-based techniques; Detection of occlusion by performing operations on regions, e.g. growing, shrinking or watersheds
    • 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

Landscapes

  • Engineering & Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Evolutionary Biology (AREA)
  • Evolutionary Computation (AREA)
  • Artificial Intelligence (AREA)
  • General Engineering & Computer Science (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明公开一种基于PCANet与高阶CRF的SAR图像分割方法。其方案为：对输入的原始SAR图像进行初始分割；利用PCANet模型，计算SAR图像高阶CRF的一元势能；计算SAR图像的高阶CRF的二元势能和高阶势能，从而推算出SAR图像的后验边缘概率；估计SAR图像高阶CRF的模型参数，最终得到该图像的分割结果。本发明充分利用SAR图像深度特征信息和高阶空域结构信息，增强了SAR图像分割方法的普适性和抗斑点噪声能力，提高了图像分割的区域一致性和边缘定位准确性，可用于SAR图像的场景解析与感兴趣区域提取。

Description

基于PCANet与高阶CRF的SAR图像分割方法
技术领域
本发明属于图像处理技术领域，更进一步涉及一种图像分割方法，可用于合成孔径雷达SAR图像的场景解析与感兴趣区域提取。
背景技术
合成孔径雷达SAR是一种全天候、全天时的高分辨率成像雷达，在地质勘测、灾害控制、参数反演和军事领域等方面应用广泛。图像分割是合成孔径雷达SAR图像分析的一个重要且持续发展的手段，准确的分割对后续地物分类和地表参数的提取具有重要意义。主成分分析网络PCANet是一种用于图像分类的简单深度学习网络，因此可以用来进行图像分割，其优势在于泛化能力强，斑点噪声抑制能力好。高阶条件随机场CRF模型是处理合成孔径雷达SAR图像分割问题的重要手段，其优势在于在合成孔径雷达SAR图像分割过程中图像的区域一致性好，边界定位准确。
中国地质大学在其拥有的专利技术“基于超像素和分形网络演化算法的极化SAR图像分割方法”(申请号：201710411372.7，授权公告号：CN107180434A)中公开了一种利用超像素和分形网络演化算法进行图像分割的方法。该方法对于待分割的极化合成孔径雷达SAR图像，生成超像素作为初始对象，并计算初始对象中相邻对象之间的相似性准则。然后统计每个对象与相邻对象之间相似性准则的最小值，若最小值小于或等于尺度参数，则合并两相邻对象生成新的对象；若最小值大于尺度参数，则不合并，遍历所有对象，完成一次分割，生成新的对象层。最后生成对象多边形，得到最终的分割结果。该方法存在的不足之处是，没有充分利用图像特征信息，使得斑点噪声抑制能力差，分割结果不够好。
西安电子科技大学在其拥有的专利技术“基于子模字典学习的SAR图像分割方法”(申请号：201810031938.8，授权公告号：CN108154511A)中公开了一种利用子模字典学习的方法进行图像分割。该专利技术在以原始SAR图像的像素为中心的邻域内，计算空间金字塔特征。然后从空间金字塔特征中选取训练数据构造图模型，并通过最大化一个子模目标函数，对图模型进行聚类。接着根据聚类结果计算字典，并固定字典，计算训练数据的稀疏编码特征、分类参数矩阵以及所有数据的稀疏编码特征。再根据分类参数矩阵计算类标签向量。最后将类标签向量转换为类标签，得到最终的分割结果。该方法存在的不足之处是，该方法忽略了SAR图像的散射特征，使得分割结果的区域一致性差，边界定位不够准确。
发明内容
本发明的目的在于针对现有技术的不足，提出了一种基于主成分分析网络与高阶条件随机场的SAR图像分割方法，以充分利用SAR图像的深度特征信息和高阶空域结构信息，增强SAR图像分割方法的泛化能力和斑点噪声抑制能力，提高图像分割的区域一致性和边界定位的准确性，从而获得更好的SAR图像分割结果。
本发明的技术思路是：利用滑窗操作，计算SAR图像局部特征向量；利用主成分分析网络PCANet模型，计算该模型中两个阶段的PCA滤波器；利用迭代条件估计ICE法，估计出合成孔径雷达SAR图像的高阶CRF模型参数，进而迭代完成SAR图像分割。其实现步骤包括如下：
(1)输入合成孔径雷达SAR图像；
(2)利用滑窗操作，计算SAR图像局部特征向量；
(3)根据(2)的计算结果，利用模糊C均值FCM算法，计算SAR图像局部类条件概率；
(4)将SAR图像中每个像素点的局部类条件概率的最大值对应的类别，作为该像素点的类别，获得对合成孔径雷达SAR图像的初始分割的结果；
(5)计算原始SAR图像中每个像素点的一元势能：
(5a)利用主成分分析网络PCANet模型，计算该模型中两个阶段的PCA滤波器；
(5b)利用两个阶段的PCA滤波器，分别计算原始SAR图像的基准深度特征向量和测试深度特征向量；
(5c)利用基准深度特征向量和测试深度特征向量的差异，计算原始合成孔径雷达SAR图像中每个像素点的一元势能；
(6)利用波茨Potts模型，计算当前SAR图像中每个像素点的高阶CRF的二元势能；
(7)利用归一化联合直方图法，计算当前SAR图像中每个像素点的高阶CRF的高阶势能；
(8)利用每个像素点的一元势能、二元势能和高阶势能，计算当前SAR图像中每个像素点高阶条件随机场CRF的后验边缘概率ps；
(9)利用迭代条件估计ICE方法，估计出当前SAR图像的模型参数；
(10)分割合成孔径雷达SAR图像：
利用当前的高阶CRF的模型参数，计算SAR图像中每个像素点的后验边缘概率ps的最大值对应的类别，作为该像素点的类别，获得本次对SAR图像的分割结果；
(11)比较本次迭代得到的图像的分割结果与上一次迭代的分割结果，若无明显差异，则结束对SAR图像的分割，得到高阶CRF无监督的SAR图像最终的分割结果；否则，将当前的迭代次数加1后返回步骤(6)。
本发明与现有的技术相比具有以下优点：
第一，由于本发明将主成分分析网络PCANet模型引入高阶条件随机场CRF模型进行图像分割，克服了现有技术泛化能力差，抗噪性能不高的问题，使得本发明利用了主成分分析网络PCANet模型来提取深度层次特征，替代了高阶条件随机场CRF模型中提取简单的图像纹理特征，充分利用了图像的特征信息，提高了图像分割的泛化能力，降低了图像分割的错误率。
第二，由于本发明利用半监督的主成分分析网络PCANet模型计算一元势能，克服了现有技术在分割合成孔径雷达SAR图像时手动提取图像特征以及训练数据不足的问题，使得本发明在只需要提取少量的样本数据训练网络，就可以得到两层PCA卷积滤波器，进而推导计算出高阶条件随机场CRF模型的一元势能，充分利用了图像的先验信息，提高了分类结果的准确性。
第三，由于本发明以高阶条件随机场CRF模型为框架进行图像分割，利用合成孔径雷达SAR图像中每个像素点与其邻域中像素点的联系，计算SAR图像中高阶条件随机场CRF模型的局部特征向量、每个像素点的后验边缘概率以及模型参数，克服了现有技术区域一致性差的问题，增强了图像分割的区域一致性和边界定位的准确性。
附图说明
图1为本发明的实现流程图；
图2为采用本发明方法、现有条件随机场CRF方法和现有卷积神经网络CNN方法对三种同质区域组合的合成孔径雷达SAR图像分割的仿真结果图；
图3为采用本发明方法、现有条件随机场CRF方法和现有卷积神经网络CNN方法对四种同质区域组合的合成孔径雷达SAR图像分割的仿真结果图；
图4为采用本发明方法、现有条件随机场CRF方法和现有卷积神经网络CNN方法对农田合成孔径雷达SAR图像分割的仿真结果图；
图5为采用本发明方法、现有条件随机场CRF方法和现有卷积神经网络CNN方法对AIRSAR图像分割的仿真结果图；
图6为采用本发明方法、现有条件随机场CRF方法和现有卷积神经网络CNN方法对河道合成孔径雷达SAR图像分割的仿真结果图。
具体实施方式
下面结合附图对本发明的实施例和效果做进一步详细描述。
参照图1，对本发明的实施步骤如下：
步骤1，输入待分割的原始合成孔径雷达SAR图像。
步骤2，利用滑窗操作，计算SAR图像局部特征向量。
(2.1)用半径为7个像素点的窗口，在合成孔径雷达SAR图像内，以一个像素点为间隔滑动，对合成孔径雷达SAR图像进行滑窗操作；
(2.2)根据滑窗操作，按照下式计算每个滑窗内的所有像素值的均值：
Figure BDA0002156204880000041
其中，μm表示第m个滑窗内的所有像素值的均值，i表示第m个滑窗内像素点行的序号，j表示第m个滑窗内的像素点列的序号，yi,j表示第m个滑窗内的第i行第j列像素点的像素值，nm表示第m个滑窗内的所有像素点的总数；
(2.3)根据所有像素值的均值，按照下式计算每个滑窗内的所有像素值的方差：
Figure BDA0002156204880000042
其中，
Figure BDA0002156204880000043
表示第m个滑窗内的所有像素值的方差；
(2.4)根据所有像素值的均值和方差，按照下式计算每个滑窗内的所有像素值的偏度：
Figure BDA0002156204880000044
其中，Φm表示第m个滑窗内的所有像素值的偏度；
(2.5)根据所有像素值的均值和方差，按照下式计算每个滑窗内的所有像素值的峰态：
Figure BDA0002156204880000051
其中，Km表示第m个滑窗内的所有像素值的峰态；
(2.6)按照下式，计算每个滑窗内的所有像素值的熵：
Figure BDA0002156204880000052
其中，Tm表示第m个滑窗内的所有像素值的熵；
(2.7)由每个滑窗内的所有像素值的均值μm、方差
Figure BDA0002156204880000053
偏度Φm、峰态Km和熵Tm组成该滑窗的直方图特征；
(2.8)在合成孔径雷达SAR图像内，将每个滑窗中所有像素值的直方图特征组成高阶条件随机场CRF的局部特征向量hs。
步骤3，根据(2)的计算结果，利用模糊C均值FCM法，初始分割合成孔径雷达SAR图像。
初始分割合成孔径雷达SAR图像的现有技术包括模糊C均值FCM法和K均值法。
本步骤采用但不限于模糊C均值法，具体实现如下：
(3.1)按照下式，计算SAR图像中每个像素点的局部类条件概率：
Figure BDA0002156204880000054
其中，Ps表示SAR图像中每个像素点的局部类条件概率；S表示SAR图像中像素点的总数；L表示SAR图像中类别的总数；rl s表示SAR图像中第s个像素点属于第l类的隶属度；φ表示每个隶属度的加权指数；cl表示第l类的聚类中心的位置；
(3.2)将SAR图像中每个像素点的局部类条件概率的最大值对应的类别，作为该像素点的类别，获得对合成孔径雷达SAR图像初始分割的结果。
步骤4，计算原始SAR图像中每个像素点的一元势能。
(4.1)将原始SAR图像中每个像素点的局部类条件概率的最大值对应的类别，作为该像素点的类别，获得对合成孔径雷达SAR图像的初始分割的结果；
(4.2)利用主成分分析网络PCANet模型，计算该模型中两个阶段的PCA滤波器：
主成分分析网络PCANet模型是一种用于图像分类的深度学习网络模型，该模型包括两个阶段的卷积滤波操作，因此需要计算该模型的两个阶段的PCA滤波器，具体实现如下：
(4.2a)在原始SAR图像中，手动选定每一类别对应的任意一块区域，然后分别在各区域中随机提取等量像素点作为当前类别的一个训练样本，再分别以训练样本中的每个像素点为中心构建半径为7个像素点的子图像；
(4.2b)用半径为7个像素点的二维卷积滤波器窗口，在每个子图像内，以一个像素点为间隔滑动进行滑窗操作；
(4.2c)按照下式，计算第一阶段的PCA滤波器：
Figure BDA0002156204880000061
其中，
Figure BDA0002156204880000062
表示第一阶段的第d个主特征向量下的滤波器，该第d个主特征向量表示将特征值按降序排列后第d个特征值对应的特征向量，mat(·)表示可以将列向量转换成方阵的函数，qd(CCT)表示CCT的第d个主特征向量，C表示所有样本的二维卷积滤波器窗去均值后的的组合，k表示二维卷积滤波器窗口的半径大小，设置为7，∈表示属于符号，R表示实数集，Rk×k表示维度为k×k的全部实矩阵的集合，H1表示第一阶段的主特征向量总数；
(4.2d)按照下式，计算第一阶段的基准特征向量：
Figure BDA0002156204880000063
其中，
Figure BDA0002156204880000064
表示用第d个主特征向量提取到的对于每个训练样本中第h个像素点的第一阶段的基准特征向量，Bh表示对应于每个训练样本中第h个像素点的子图像，*表示二维卷积操作，N表示每个训练样本中的像素点总数；
(4.2e)按照下式，计算第二阶段的PCA滤波器：
Wb 2＝mat(qd(JJT))∈Rk×k,1≤b≤H2
其中，Wb 2表示第二阶段的第b个主特征向量下的滤波器，qd(JJT)表示JJT的第d个主特征向量，J表示第一阶段全部特征向量的二维卷积滤波器窗的组合，H2表示第二阶段的主特征向量总数；
(4.3)利用两个阶段的PCA滤波器，分别计算原始SAR图像的基准深度特征向量和测试深度特征向量：
(4.3a)按照下式，计算第二阶段的基准特征向量：
Figure BDA0002156204880000071
其中，
Figure BDA0002156204880000072
表示用第一阶段的第d个主特征向量和第二阶段的第b个主特征向量提取到的对于每个训练样本中第h个像素点的第二阶段基准特征向量；
(4.3b)按照下式，计算基准深度特征向量：
Figure BDA0002156204880000073
其中，
Figure BDA0002156204880000074
表示第l类的第h个像素点的基准深度特征向量，
Figure BDA0002156204880000075
表示用第一阶段的第H1个主特征向量和第二阶段的第H2个主特征向量提取到的对于第l类的第h个像素点的第二阶段基准特征向量，L表示合成孔径雷达SAR图像的总类别；
(4.3c)在原始SAR图像中，以每个像素点为中心构建半径为7个像素点的子图像；
(4.3d)用半径为7个像素点的二维卷积滤波器窗口，在每个子图像内，以一个像素点为间隔进行滑窗操作；
(4.3e)按照下式，计算第一阶段的测试特征向量：
Figure BDA0002156204880000076
其中，
Figure BDA0002156204880000077
表示用第d个主特征向量提取到的合成孔径雷达SAR图像中第s个像素点的第一阶段测试特征向量，Qs表示合成孔径雷达SAR图像中第s个像素点的子图像；
(4.3f)按照下式，计算第二阶段的测试特征向量：
Figure BDA0002156204880000078
其中，Vs d,b表示用第一阶段的第d个主特征向量和第二阶段的第b个主特征向量提取到的合成孔径雷达SAR图像中第s个像素点的第二阶段测试特征向量；
(4.3g)按照下式，计算测试深度特征向量：
Figure BDA0002156204880000079
其中，
Figure BDA0002156204880000081
表示原始SAR图像中第s个像素点属于第l类的基准深度特征向量，
Figure BDA0002156204880000082
表示用第一阶段的第H1个主特征向量和第二阶段的第H2个主特征向量提取到的原始合成孔径雷达SAR图像中第s个像素点属于第l类的第二阶段测试特征向量；
(4.4)利用基准深度特征向量和测试深度特征向量的差异，计算原始合成孔径雷达SAR图像中每个像素点的一元势能：
(4.4a)按照下式，计算基准深度特征向量和测试深度特征向量的差异：
Figure BDA0002156204880000083
其中，Ds,l表示原始合成孔径雷达SAR图像中第s个像素点属于第l类的特征向量差，h表示每个训练样本中像素点的编号，N表示每个训练样本中像素点的总数，||·||2表示求平方和操作；
(4.4b)按照下式，计算每个像素点的一元势能：
Figure BDA0002156204880000084
其中，As,l表示原始合成孔径雷达SAR图像中第s个像素点属于第l类的一元势能，exp(·)表示以自然常数为底的指数操作，
Figure BDA0002156204880000085
表示一个常数，设置为20。
步骤5，计算当前合成孔径雷达SAR图像中每个像素点与其邻域系统中的每个像素点之间的高阶CRF模型的二元势能。
本步骤基于现有的波茨Potts模型实现，其计算步骤如下：，
(5.1)按照下式，计算指数加权平均比算子：
Figure BDA0002156204880000086
其中，μst表示当前合成孔径雷达SAR图像中第s个像素点和与其邻域系统中的第t个像素点间的边缘强度，
Figure BDA0002156204880000087
表示合成孔径雷达SAR图像中第s个像素点的像素值，
Figure BDA0002156204880000088
表示合成孔径雷达SAR图像中第t个像素点的像素值；
(5.2)按照下式，计算当前合成孔径雷达SAR图像中每个像素点与其邻域系统中的每个像素点之间的二元势能：
Figure BDA0002156204880000091
其中，D{s,t}表示当前合成孔径雷达SAR图像中第s个像素点与其邻域系统中第t个像素点之间的二元势能；
αH表示当前SAR图像中第s个像素点与其水平邻域系统中第t个像素点之间的水平相关性参数，该参数初始值为1，之后每次迭代，该参数的取值都取自上一次迭代过程中计算得到的高阶条件随机场CRF的模型参数；
δ(xs,xt)表示一个常数，当合成孔径雷达SAR图像中第s个像素点的类别xs与其邻域系统中第t个像素点的类别xt相同时，δ(xs,xt)＝1，当合成孔径雷达SAR图像中第s个像素点的类别xs与其邻域系统中第t个像素点的类别xt不同时，δ(xs,xt)＝0；
NsH表示当前合成孔径雷达SAR图像中第s个像素点的水平邻域集合；
αV表示当前合成孔径雷达SAR图像中第s个像素点与其垂直邻域集合中的第t个像素点之间的垂直相关性参数，该参数初始值为1，之后每次迭代，该参数的取值都取自上一次迭代过程中计算得到的高阶条件随机场CRF的模型参数；
NsV表示当前合成孔径雷达SAR图像中第s个像素点的垂直邻域集合。
步骤6，计算当前合成孔径雷达SAR图像中每个像素点的高阶条件随机场CRF模型的高阶势能。
计算每个像素点的条件随机场CRF模型的高阶势能的现有技术包括归一化联合直方图法、鲁棒高阶势能法和高阶波茨Potts模型法。
本步骤采用但不限于归一化联合直方图法，具体实现如下：
(6.1)用半径为7个像素点的窗口，在当前合成孔径雷达SAR图像中，以一个像素点为间隔滑动，对当前合成孔径雷达SAR图像进行一次滑窗操作；
(6.2)按照下式，计算每个滑窗内每个类别的概率：
Figure BDA0002156204880000092
其中，
Figure BDA0002156204880000093
表示第
Figure BDA0002156204880000094
个滑窗内第l个类别的概率，
Figure BDA0002156204880000095
表示第
Figure BDA0002156204880000096
个滑窗内属于第l个类别的像素点总数，
Figure BDA0002156204880000101
表示第
Figure BDA0002156204880000102
个滑窗内单个基团C1的总数，该单个基团表示由单个像素点构成的基团；
(6.3)按照下式，计算每个滑窗内每个类别与其他类别相邻的概率：
Figure BDA0002156204880000103
其中，
Figure BDA0002156204880000104
表示第
Figure BDA0002156204880000105
个滑窗内第l个类别与第ε个类别相邻的概率；Π表示求积操作；
τ表示二元基团的种类，该二元基团表示由两个相邻像素点构成的基团，根据τ的不同取值，C2(τ)表示不同的二元基团，当τ＝1时，C2(τ)表示像素点与第l个类别左相邻的二元基团；当τ＝2时，C2(τ)表示像素点与第l个类别右相邻的二元基团；当τ＝3时，C2(τ)表示像素点与第l个类别下相邻的二元基团；当τ＝4时，C2(τ)表示像素点与第l个类别上相邻的二元基团；
Figure BDA0002156204880000106
表示第
Figure BDA0002156204880000107
个滑窗内第l个类别与第ε个类别所构成的C2(τ)的总数；
Figure BDA0002156204880000108
表示第
Figure BDA0002156204880000109
个滑窗内第l个类别的C2(τ)的总数；
(6.4)按照下式，计算每个滑窗中每个像素点的高阶势能：
Figure BDA00021562048800001010
其中，
Figure BDA00021562048800001011
表示第
Figure BDA00021562048800001012
个滑窗中的第s个像素点的高阶势能。
步骤7，计算当前合成孔径雷达SAR图像中每个像素点高阶条件随机场CRF的后验边缘概率。
(7.1)按照下式，计算SAR图像中每个像素点的散射统计：
Figure BDA00021562048800001013
其中，ds表示SAR图像中第s个像素点的散射统计；
βl表示广义伽马Gamma分布模型的形状参数，该参数初始值由初始分割结果得到，之后每次迭代，该参数的取值都取自上一次迭代过程中计算得到的高阶条件随机场CRF的模型参数；
λl表示广义伽马Gamma分布模型的索引形状参数，该参数初始值由初始分割结果得到，之后每次迭代，该参数的取值都取自上一次迭代过程中计算得到的高阶条件随机场CRF的模型参数；
rs表示SAR图像中第s个像素点的强度；
αl表示广义伽马Gamma分布模型的尺度参数，该参数初始值由初始分割结果得到，之后每次迭代，该参数的取值都取自上一次迭代过程中计算得到的高阶条件随机场CRF的模型参数；
Γ(·)表示伽马Gamma函数；
(7.2)按照下式，计算SAR图像中每个像素点高阶条件随机场CRF的后验边缘概率：
Figure BDA0002156204880000111
其中，ps表示当前合成孔径雷达SAR图像中第s个像素点的高阶条件随机场CRF的本次迭代的后验边缘概率。
步骤8，利用迭代条件估计ICE方法，估计出合成孔径雷达SAR图像的模型参数。
(8.1)设置迭代条件估计ICE法的初始迭代次数为1，进入ICE法的迭代过程；
(8.2)用SAR图像的水平相关参数αH、垂直相关参数αV、广义伽马Gamma分布模型的形状参数βl、索引形状参数λl和尺度参数αl，组成高阶条件随机场CRF的模型参数；
(8.3)将当前的高阶条件随机场CRF的模型参数输入到步骤3～步骤7，计算SAR图像中每个像素点的后验边缘概率，将所有像素点的后验边缘概率和原始SAR图像输入到吉布斯Gibbs采样器中，得到2个采样图像；
(8.4)采用相关参数迭代法，计算每个采样图像的相关参数：
(8.4a)设置相关参数迭代法的初始迭代次数为1，进入相关参数迭代法的迭代过程；
(8.4b)按照下式，计算每个采样图像的势能：
Figure BDA0002156204880000112
其中，Yq表示第q个采样图像的势能，q的取值范围为{1,2}，S'表示每个采样图像中像素点的总数，
Figure BDA0002156204880000113
表示第q个采样图像在第
Figure BDA0002156204880000114
个滑窗中第s′个像素点的高阶势能，
Figure BDA0002156204880000115
表示SAR图像中第s′个像素点与其邻域系统中第t'个像素点之间的二元势能，Ns'表示第q个采样图像中第s′个像素点的邻域系统，
Figure BDA0002156204880000121
表示第q个采样图像中第s′个像素点属于第l类的高阶条件随机场CRF的一元势能；
(8.4c)利用当前的水平相关参数αH，垂直相关参数αV、广义伽马Gamma分布模型的形状参数βl、索引形状参数λl和尺度参数αl，组成高阶条件随机场CRF的模型参数；
(8.4d)将当前的高阶条件随机场CRF的模型参数输入到步骤3～步骤7，计算SAR图像中每个像素点的后验边缘概率，将所有像素点的后验边缘概率和原始SAR图像输入到吉布斯Gibbs采样器中，得到本次迭代的采样图像；
(8.4e)按照下式，计算本次迭代中每个采样图像的水平相关参数：
Figure BDA0002156204880000122
其中，
Figure BDA0002156204880000123
表示本次迭代中第q个采样图像的水平相关参数，αH表示水平相关参数，
Figure BDA0002156204880000124
表示本次迭代的上一次迭代中第q个采样图像的水平相关参数，N表示每个采样图像的像素点总数，a表示本次迭代的次数，
Figure BDA0002156204880000125
表示求梯度操作，
Figure BDA0002156204880000126
表示第q个采样图像，
Figure BDA0002156204880000127
表示本次迭代的采样图像；
(8.4f)按照下式，计算本次迭代中每个采样图像的垂直相关参数：
Figure BDA0002156204880000128
其中，
Figure BDA0002156204880000129
表示本次迭代中第q个采样图像的垂直相关参数，αV表示垂直相关参数，
Figure BDA00021562048800001210
表示本次迭代的上一次迭代中第q个采样图像的垂直相关参数；
(8.4g)用本次迭代中得到的每个采样图像的水平相关参数
Figure BDA00021562048800001211
和垂直相关参数
Figure BDA00021562048800001212
组成本次迭代的每个采样图像的相关参数；
(8.4h)比较本次迭代得到的相关参数与上次迭代得到的相关参数的差值与阈值10-3的大小：
若该差值小于阈值10-3，则将本次迭代的每个采样图像的相关参数作为本次ICE法迭代的每个采样图像的相关参数；
否则，将本次迭代得到的每个采样图像的相关参数作为当前的相关参数，并将相关参数迭代法的迭代次数加1后，返回步骤(8.4c)；(8.5)采用特征参数迭代法，计算每个采样图像的特征参数：
(8.5)采用对数累积量MoLC法，计算每个采样图像的广义伽马Gamma参数：
(8.5a)对广义伽马Gamma分布的概率密度函数进行梅林Mellin变换，得到广义伽马Gamma分布的第一个第二类特征函数；
(8.5b)按照下式，计算广义Gamma分布的第二个第二类特征函数：
ξ(ω)＝lnγ(ω)
其中，ξ(ω)是广义伽马Gamma分布的第二个第二类特征函数，ω是广义伽马Gamma分布的概率密度函数的自变量，γ(ω)是广义伽马Gamma分布的第一个第二类特征函数；
(8.5c)按照下式，计算广义伽马Gamma分布的第二个第二类特征函数的一阶对数累积量，二阶对数累积量和三阶对数累积量：
Figure BDA0002156204880000131
其中，cρ表示广义伽马Gamma分布的第二个第二类特征函数的对数累积量；
当ρ＝1时，cρ表示广义伽马Gamma分布的第二个第二类特征函数的一阶对数累积量；
当ρ＝2时，cρ表示广义伽马Gamma分布的第二个第二类特征函数的二阶对数累积量；
当ρ＝3时，cρ表示广义伽马Gamma分布的第二个第二类特征函数的三阶对数累积量；
Figure BDA0002156204880000132
表示y对x求ρ阶偏导；
(8.5d)按照下式，建立起参数估计的一阶矩模型：
Figure BDA0002156204880000133
其中，ψ(0,λl)表示多伽马polygamma函数的参数k′取0，x′取λl时的情况；ns′表示每个采样图像中像素点的总数；rs′表示每个采样的图像在第s′个像素点的强度，S′表示每个采样图像中所有像素点的集合；
(8.5e)按照下式，建立起参数估计的二阶矩模型：
Figure BDA0002156204880000141
其中，ψ(1,λl)表示多伽马polygamma函数的参数k′取1，x′取λl时的情况；c1表示广义伽马Gamma分布的第二个第二类特征函数的一阶对数累积量；
(8.5f)按照下式，建立起参数估计的三阶矩模型：
Figure BDA0002156204880000142
其中，ψ(2,λl)表示多伽马polygamma函数的参数k′取2，x′取λl时的情况；
(8.5g)根据参数估计的一阶矩，二阶矩，三阶矩模型公式，解出广义Gamma分布的参数αl，βl和λl；
(8.6)将本次ICE法迭代得到的每个采样图像的相关参数和广义Gamma分布的参数，组成本次ICE法迭代的每个采样图像的高阶CRF的模型参数；
(8.7)将本次ICE法迭代的两个采样图像的高阶CRF的模型参数取平均，得到本次ICE法迭代的SAR图像的高阶CRF的模型参数；
(8.8)比较本次ICE法迭代得到的模型参数与上次ICE法迭代得到的模型参数的差值与阈值10-3的大小：
若该差值小于阈值10-3，则得到ICE法最终的高阶CRF的模型参数；
否则，将本次ICE法的迭代得到的模型参数作为当前的模型参数，然后将ICE法的迭代次数加1后，返回步骤(8.3)。
步骤9，分割合成孔径雷达SAR图像：
(9.1)将步骤8得到的高阶CRF的模型参数输入到步骤3～步骤7，计算SAR图像中每个像素点的后验边缘概率ps，
(9.2)将每个像素点的后验边缘概率ps中的最大值对应的类别，作为该像素点的类别，获得本次迭代的SAR图像分割结果。
步骤10，获得最终SAR图像分割结果图像。
将本次迭代得到的SAR图像分割结果与上一次迭代的分割结果进行比较：
若两者的分割结果无明显差异，则结束对合成孔径雷达SAR图像的分割，得到高阶条件随机场CRF的最终SAR图像分割结果图像；
否则，将当前的迭代次数加1，并将本次的分割结果图像作为当前的SAR图像，返回步骤5。
下面结合仿真实验对本发明的效果作进一步的说明。
1.仿真实验条件：
本发明仿真实验的硬件测试平台是：处理器为CPU intel Core i7-7700，主频为4.2GHz，内存16GB；软件平台为：Windows 10家庭版，64位操作系统，MATLAB R2018b。
2.仿真内容及仿真结果分析：
仿真实验1：采用本发明、现有条件随机场CRF方法和现有卷积神经网络CNN方法，对三种同质区域组合的合成孔径雷达SAR图像进行仿真分割，结果如图2所示。其中：
图2(a)是仿真实验1输入的原图，该图是从德国慕尼黑附近的帕弗沃伦地区，利用德国航空航天中心ESAR L波段HV数据生成一个模拟图像，对该模拟图像手动选择其中三类的均匀区域，然后将区域放入预设的地面实况中，得到的三种同质区域组合的合成孔径雷达SAR图像。
图2(b)是使用现有条件随机场CRF方法对图2(a)进行图像分割的仿真结果图。
图2(c)是使用现有卷积神经网络CNN方法对图2(a)进行图像分割的仿真结果图。
图2(d)是使用本发明方法对图2(a)进行分割的仿真结果图。
由图2(a)可以看出，该图像的噪声较强，将本发明方法的仿真结果图2(d)与现有现有条件随机场CRF方法的仿真结果图2(b)进行比较，可以看出本发明方法的抗噪能力更强，分割结果更准确；将本发明方法的仿真结果图2(d)与现有卷积神经网络CNN方法的仿真结果图2(c)进行比较可以看出，本发明方法分割出的三种同质区域的区域一致性更好，边界定位较为准确。
仿真实验2：采用本发明、现有条件随机场CRF方法和现有卷积神经网络CNN方法，对四种同质区域组合的合成孔径雷达SAR图像进行仿真分割，结果如图3所示。其中：
图3(a)是仿真实验2输入的原图，该图是由四种同质区域组合的合成孔径雷达SAR图像。
图3(b)是使用现有条件随机场CRF方法对图3(a)进行图像分割的仿真结果图。
图3(c)是使用现有卷积神经网络CNN方法对图3(a)进行图像分割的仿真结果图。
图3(d)是使用本发明方法对图3(a)进行分割的仿真结果图。
由图3(a)可以看出，该图的噪声稍弱，将本发明方法的仿真结果图3(d)与现有条件随机场CRF方法的仿真结果图3(b)进行比较，可以看出本发明方法的泛化能力更强；将本发明方法的仿真结果图3(d)与现有卷积神经网络CNN方法的仿真结果图3(c)进行比较，可以看出本发明方法的仿真结果中四种同质区域的区域边界更为光滑，边缘位置较为准确，分割结果准确性非常高。
仿真实验3：采用本发明方法、现有条件随机场CRF方法和现有卷积神经网络CNN方法，对农田合成孔径雷达SAR图像进行仿真分割，结果如图4所示。其中：
图4(a)是仿真实验3输入的原图，该图是TerraSAR-X波段SAR数据生成的农田合成孔径雷达SAR图像。
图4(b)为使用现有条件随机场CRF方法对图4(a)进行图像分割的仿真结果图。
图4(c)为使用现有卷积神经网络CNN方法对图4(a)进行图像分割的仿真结果图。
图4(d)是使用本发明方法对图4(a)进行分割的仿真结果图。
由图4(a)可以看出，该图像主要包含农田、道路和建筑物，将本发明方法的仿真结果图4(d)与现有条件随机场CRF方法的仿真结果图4(b)进行比较，可以看出本发明方法对斑点噪声的抑制能力更强；将本发明方法的仿真结果图4(d)与现有卷积神经网络CNN方法的仿真结果图4(c)进行比较，可以看出本发明方法的仿真结果中对农田的边界定位明显优于现有的卷积神经网络CNN方法，并且获得了较为平滑的均匀区域。
仿真实验4：采用本发明方法、现有条件随机场CRF方法和现有卷积神经网络CNN方法，对AIRSAR图像进行仿真分割，结果如图5所示。其中：
图5(a)是仿真实验4输入的原图，该图是L波段的AIRSAR数据生成的合成孔径雷达SAR图像。
图5(b)为使用现有条件随机场CRF方法对图5(a)进行图像分割的仿真结果图
图5(c)为使用现有卷积神经网络CNN方法对图5(a)进行图像分割的仿真结果图。
图5(d)是使用本发明方法对图5(a)进行分割的仿真结果图。
由图5(a)可以看出，该图像主要包含城区、农田和道路，将本发明方法的仿真结果图5(d)与现有条件随机场CRF方法的仿真结果图5(b)进行比较，可以看出本发明方法抑制斑点噪声的鲁棒性更强，误分割少；将本发明方法的仿真结果图5(d)与现有卷积神经网络CNN方法的仿真结果图5(c)进行比较，可以看出本发明方法的仿真结果中的区域一致性明显优于现有的卷积神经网络CNN方法，并且对道路的边缘保持较为准确，保存了更多的图像结构信息。
仿真实验5：采用本发明方法、现有条件随机场CRF方法和现有卷积神经网络CNN方法，对河道合成孔径雷达SAR图像进行仿真分割，结果如图6所示。其中：
图6(a)是仿真实验5输入的原图。该图是由美国喷气推进实验室获取的某河道L波段数据生成的河道合成孔径雷达SAR图像。
图6(b)为使用现有条件随机场CRF方法对图6(a)进行图像分割的仿真结果图。
图6(c)为使用现有卷积神经网络CNN方法对图6(a)进行图像分割的仿真结果图。
图6(d)是使用本发明方法对图6(a)进行分割的仿真结果图。
由图6(a)可以看出，该图像主要包含三个区域，农田，城区和河道；
将本发明方法的仿真结果图6(d)与现有条件随机场CRF方法的仿真结果图6(b)进行比较，可以看出本发明方法的仿真结果中在农田区域保持了很好的均匀度，结构信息更为准确；
将本发明方法的仿真结果图6(d)与现有卷积神经网络CNN方法的仿真结果图6(c)进行比较，可以看出本发明方法的仿真结果对河道和农田的边界定位明显优于现有的卷积神经网络CNN方法，并且保留了较多的细节信息。

Claims (10)

1.一种基于PCANet与高阶CRF的SAR图像分割方法，其特征在于，包括如下：
(1)输入合成孔径雷达SAR图像；
(2)利用滑窗操作，计算SAR图像局部特征向量；
(3)根据(2)的计算结果，利用模糊C均值FCM算法，计算SAR图像局部类条件概率；
(4)将SAR图像中每个像素点的局部类条件概率的最大值对应的类别，作为该像素点的类别，获得对合成孔径雷达SAR图像的初始分割的结果；
(5)计算原始SAR图像中每个像素点的一元势能：
(5a)利用主成分分析网络PCANet模型，计算该模型中两个阶段的PCA滤波器；
(5b)利用两个阶段的PCA滤波器，分别计算原始SAR图像的基准深度特征向量和测试深度特征向量；
(5c)利用基准深度特征向量和测试深度特征向量的差异，计算原始合成孔径雷达SAR图像中每个像素点的一元势能；
(6)利用波茨Potts模型，计算当前SAR图像中每个像素点的高阶CRF的二元势能；
(7)利用归一化联合直方图法，计算当前SAR图像中每个像素点的高阶CRF的高阶势能；
(8)利用每个像素点的一元势能、二元势能和高阶势能，计算当前SAR图像中每个像素点高阶条件随机场CRF的后验边缘概率ps；
(9)利用迭代条件估计ICE方法，估计出当前SAR图像的模型参数；
(10)分割合成孔径雷达SAR图像：
利用当前的高阶CRF的模型参数，计算SAR图像中每个像素点的后验边缘概率ps的最大值对应的类别，作为该像素点的类别，获得本次对SAR图像的分割结果；
(11)比较本次迭代得到的图像的分割结果与上一次迭代的分割结果，若无明显差异，则结束对SAR图像的分割，得到高阶CRF无监督的SAR图像最终的分割结果；否则，将当前的迭代次数加1后返回步骤(6)。
2.根据权利要求1所述的方法，其特征在于，步骤(2)中利用滑窗操作，计算SAR图像局部特征向量，其实现如下：
(2a)用半径为7个像素点的窗口，在合成孔径雷达SAR图像内，以一个像素点为间隔滑动，对合成孔径雷达SAR图像进行滑窗操作；
(2b)按照下式，计算每个滑窗内的所有像素值的均值：
Figure FDA0002156204870000021
其中，μm表示第m个滑窗内的所有像素值的均值，i表示第m个滑窗内像素点行的序号，j表示第m个滑窗内的像素点列的序号，yi,j表示第m个滑窗内的第i行第j列像素点的像素值，nm表示第m个滑窗内的所有像素点的总数；
(2c)按照下式，计算每个滑窗内的所有像素值的方差：
Figure FDA0002156204870000022
其中，
Figure FDA0002156204870000023
表示第m个滑窗内的所有像素值的方差；
(2d)按照下式，计算每个滑窗内的所有像素值的偏度：
Figure FDA0002156204870000024
其中，Φm表示第m个滑窗内的所有像素值的偏度；
(2e)按照下式，计算每个滑窗内的所有像素值的峰态：
Figure FDA0002156204870000025
其中，Km表示第m个滑窗内的所有像素值的峰态；
(2f)按照下式，计算每个滑窗内的所有像素值的熵：
Figure FDA0002156204870000026
其中，Tm表示第m个滑窗内的所有像素值的熵；
(2g)由每个滑窗内的所有像素值的均值、方差、偏度、峰态和熵组成该滑窗的直方图特征；
(2h)在合成孔径雷达SAR图像内，将每个滑窗中所有像素值的直方图特征组成高阶条件随机场CRF的局部特征向量。
3.根据权利要求1所述的方法，其特征在于，步骤(3)中利用模糊C均值FCM算法，计算SAR图像局部类条件概率，其公式如下：
Figure FDA0002156204870000031
其中，Ps表示SAR图像中每个像素点的局部类条件概率；S表示SAR图像中像素点的总数；L表示SAR图像中类别的总数；rl s表示SAR图像中第s个像素点属于第l类的隶属度；φ表示每个隶属度的加权指数；cl表示第l类的聚类中心的位置。
4.根据权利要求1所述的方法，其特征在于，步骤(5a)中利用主成分分析网络PCANet模型，计算该模型中两个阶段的PCA滤波器，其实现如下：
(5a1)在合成孔径雷达SAR图像中，手动选定每一类别对应的任意一块区域，然后在每一类别对应的区域中随机提取等数量的像素点作为当前类别的一个训练样本，再分别以每个训练样本中的每个像素点为中心构建半径为7个像素点的子图像；
(5a2)用半径为7个像素点的二维卷积滤波器窗口，在每个子图像内，以一个像素点为间隔滑动，对子图像进行滑窗操作；
(5a3)按照下式，计算第一阶段的PCA滤波器：
Figure FDA0002156204870000032
其中，
Figure FDA0002156204870000033
表示第一阶段的第d个主特征向量下的滤波器，mat(·)表示可以将列向量转换成方阵的函数，qd(CCT)表示CCT的第d个主特征向量，C表示所有样本的二维卷积滤波器窗去均值后的组合，k表示二维卷积滤波器窗口的半径大小，设置为7，R表示实数集，Rk×k表示维度为k×k的全部实矩阵的集合，H1表示第一阶段的主特征向量总数；
(5a4)按照下式，计算第一阶段的基准特征向量：
Figure FDA0002156204870000041
其中，
Figure FDA0002156204870000042
表示用第d个主特征向量提取到的对于每个训练样本中第h个像素点的第一阶段的基准特征向量，Bh表示对应于每个训练样本中第h个像素点的子图像，*表示二维卷积操作，N表示每个训练样本中的像素点总数；
(5a5)按照下式，计算第二阶段的PCA滤波器：
Figure FDA0002156204870000043
其中，
Figure FDA0002156204870000044
表示第二阶段的第b个主特征向量下的滤波器，qd(JJT)表示JJT的第d个主特征向量，J表示第一阶段全部特征向量的二维卷积滤波器窗的组合，H2表示第二阶段的主特征向量总数。
5.根据权利要求1所述的方法，其特征在于，步骤(5b)中利用两个阶段的PCA滤波器，计算原始SAR图像的基准深度特征向量，其实现如下：
(5b1)按照下式，计算第二阶段的基准特征向量：
Figure FDA0002156204870000045
其中，
Figure FDA0002156204870000046
表示用第一阶段的第d个主特征向量和第二阶段的第b个主特征向量提取到的对于每个训练样本中第h个像素点的第二阶段基准特征向量，
Figure FDA0002156204870000047
表示第二阶段的第b个主特征向量下的滤波器，*表示二维卷积操作，
Figure FDA0002156204870000048
表示用第d个主特征向量提取到的对于每个训练样本中第h个像素点的第一阶段的基准特征向量，k表示二维卷积滤波器窗口的半径大小，设置为7，R表示实数集，Rk×k表示维度为k×k的全部实矩阵的集合，H1表示第一阶段的主特征向量总数，H2表示第二阶段的主特征向量总数；
(5b2)按照下式，计算基准深度特征向量：
Figure FDA0002156204870000049
其中，
Figure FDA0002156204870000051
表示第l类的第h个像素点的基准深度特征向量，
Figure FDA0002156204870000052
表示用第一阶段的第H1个主特征向量和第二阶段的第H2个主特征向量提取到的对于第l类的第h个像素点的第二阶段基准特征向量，L表示合成孔径雷达SAR图像的总类别，N表示每个训练样本中的像素点总数。
6.根据权利要求1所述的方法，其特征在于，步骤(5b)中利用两个阶段的PCA滤波器，计算原始SAR图像的测试深度特征向量，其实现如下：
(5b3)在原始合成孔径雷达SAR图像中，以每个像素点为中心构建半径为7个像素点的子图像；
(5b4)用半径为7个像素点的二维卷积滤波器窗口，在每个子图像内，以一个像素点为间隔滑动，对子图像进行滑窗操作；
(5b5)按照下式，计算第一阶段的测试特征向量：
Figure FDA0002156204870000053
其中，
Figure FDA0002156204870000054
表示用第d个主特征向量提取到的合成孔径雷达SAR图像中第s个像素点的第一阶段测试特征向量，Qs表示合成孔径雷达SAR图像中第s个像素点的子图像，*表示二维卷积操作，
Figure FDA0002156204870000055
表示第一阶段的第d个主特征向量下的滤波器；
(5b6)按照下式，计算第二阶段的测试特征向量：
Figure FDA0002156204870000056
其中，Vs d,b表示用第一阶段的第d个主特征向量和第二阶段的第b个主特征向量提取到的合成孔径雷达SAR图像中第s个像素点的第二阶段测试特征向量，
Figure FDA0002156204870000057
表示第二阶段的第b个主特征向量下的滤波器，H1表示第一阶段的主特征向量总数，H2表示第二阶段的主特征向量总数；
(5b7)按照下式，计算测试深度特征向量：
Figure FDA0002156204870000058
其中，
Figure FDA0002156204870000061
表示合成孔径雷达SAR图像中第s个像素点属于第l类的基准深度特征向量，
Figure FDA0002156204870000062
表示用第一阶段的第H1个主特征向量和第二阶段的第H2个主特征向量提取到的合成孔径雷达SAR图像中第s个像素点属于第l类的第二阶段测试特征向量，L表示合成孔径雷达SAR图像的总类别，N表示每个训练样本中的像素点总数。
7.根据权利要求1所述的方法，其特征在于，步骤(5c)中利用基准深度特征向量和测试深度特征向量的差异，计算原始SAR图像中每个像素点的一元势能,其实现如下：
(5c1)按照下式，计算基准深度特征向量和测试深度特征向量的差异：
Figure FDA0002156204870000063
其中，Ds,l表示原始合成孔径雷达SAR图像中第s个像素点属于第l类的特征向量差，h表示每个训练样本中像素点的编号，N表示每个训练样本中像素点的总数；
(5c2)按照下式，计算每个像素点的一元势能：
Figure FDA0002156204870000064
其中，As,l表示原始合成孔径雷达SAR图像中第s个像素点属于第l类的一元势能，exp(·)表示以自然常数为底的指数操作，
Figure FDA0002156204870000065
表示一个常数，设置为20。
8.根据权利要求1所述的方法，其特征在于，步骤6中利用波茨Potts模型，计算当前SAR图像中高阶CRF模型的二元势能，其实现如下：
(6a)按照下式，计算指数加权平均比算子：
Figure FDA0002156204870000066
其中，μst表示合成孔径雷达SAR图像中第s个像素点和与其邻域系统中的第t个像素点间的边缘强度，
Figure FDA0002156204870000067
表示合成孔径雷达SAR图像中第s个像素点的像素值，
Figure FDA0002156204870000068
表示合成孔径雷达SAR图像中第t个像素点的像素值；
(6b)按照下式，计算分割后的合成孔径雷达SAR图像中每个像素点与其邻域系统中的每个像素点之间的二元势能：
Figure FDA0002156204870000071
其中，D{s,t}表示分割后的合成孔径雷达SAR图像中第s个像素点与其邻域系统中第t个像素点之间的二元势能；αH表示分割后的SAR图像中第s个像素点与其水平邻域系统中第t个像素点之间的水平相关性参数，该参数初始值为1，之后每次迭代，该参数的取值都取自上一次迭代过程中计算得到的高阶条件随机场CRF的模型参数；δ(xs,xt)表示一个常数，当分割后的合成孔径雷达SAR图像中第s个像素点的类别xs与其邻域系统中第t个像素点的类别xt相同时，δ(xs,xt)＝1，当分割后的合成孔径雷达SAR图像中第s个像素点的类别xs与其邻域系统中第t个像素点的类别xt不同时，δ(xs,xt)＝0；NsH表示分割后的合成孔径雷达SAR图像中第s个像素点的水平邻域系统；αV表示分割后的合成孔径雷达SAR图像中第s个像素点与其垂直邻域系统中的第t个像素点之间的垂直相关性参数，该参数初始值为1，之后每次迭代，该参数的取值都取自上一次迭代过程中计算得到的高阶条件随机场CRF的模型参数；NsV表示分割后合成孔径雷达SAR图像中第s个像素点的垂直邻域系统。
9.根据权利要求1所述的方法，其特征在于，步骤7中利用归一化联合直方图法，计算当前SAR图像中每个像素点的条件随机场CRF模型的高阶势能，其实现如下：
(7a)用半径为7个像素点的窗口，在分割后的合成孔径雷达SAR图像中，以一个像素点为间隔滑动，对分割后的合成孔径雷达SAR图像进行一次滑窗操作；
(7b)按照下式，计算每个滑窗内每个类别的概率：
Figure FDA0002156204870000072
其中，
Figure FDA0002156204870000073
表示第
Figure FDA0002156204870000074
个滑窗内第l个类别的概率，
Figure FDA0002156204870000075
表示第
Figure FDA0002156204870000076
个滑窗内属于第l个类别的像素点总数，
Figure FDA0002156204870000077
表示第
Figure FDA0002156204870000078
个滑窗内单个基团C1的总数，单个基团表示由单个像素点构成的基团；
(7c)按照下式，计算每个滑窗内每个类别与其他类别相邻的概率：
Figure FDA0002156204870000081
其中，
Figure FDA0002156204870000082
表示第
Figure FDA0002156204870000083
个滑窗内第l个类别与第ε个类别相邻的概率；Π表示求积操作；τ表示二元基团的种类，二元基团表示由两个相邻像素点构成的基团；τ＝1时，C2(τ)表示像素点与第l个类别左相邻的二元基团；τ＝2时，C2(τ)表示像素点与第l个类别右相邻的二元基团；τ＝3时，C2(τ)表示像素点与第l个类别下相邻的二元基团；τ＝4时，C2(τ)表示像素点与第l个类别上相邻的二元基团；
Figure FDA0002156204870000084
表示第
Figure FDA0002156204870000085
个滑窗内第l个类别与第ε个类别所构成的C2(τ)的总数；
Figure FDA0002156204870000086
表示第
Figure FDA0002156204870000087
个滑窗内第l个类别的C2(τ)的总数；
(7d)按照下式，计算每个滑窗中每个像素点的高阶势能：
Figure FDA0002156204870000088
其中，
Figure FDA0002156204870000089
表示第
Figure FDA00021562048700000810
个滑窗中的第s个像素点的高阶势能。
10.根据权利要求1所述的方法，其特征在于，步骤9中利用迭代条件估计ICE方法，估计出合成孔径雷达SAR图像的模型参数，其实现如下：
(9a)利用合成孔径雷达SAR图像的水平相关参数，垂直相关参数，广义伽马Gamma参数，组成高阶条件随机场CRF的模型参数；
(9b)根据当前的高阶条件随机场CRF的模型参数，计算合成孔径雷达SAR中每个像素点的后验边缘概率，将所有像素点的后验边缘概率和原始合成孔径雷达SAR图像输入到吉布斯Gibbs采样器中，得到2个采样图像；
(9c)采用相关参数迭代法，计算每个采样图像的相关参数；
(9d)采用对数累积量MoLC方法，计算每个采样图像的广义伽马Gamma参数；
(9e)将每个采样图像的相关参数和广义Gamma分布的参数，组成每个采样图像的高阶CRF的模型参数；
(9f)将两个采样图像所得的高阶条件随机场CRF的模型参数分别取平均，得到高阶条件随机场CRF的模型参数；
(9g)比较本次迭代得到的参数与上次迭代得到参数的差值，若小于10-3，则得到高阶条件随机场CRF的模型参数；否则，将当前的迭代次数加1后，返回步骤(9b)。
CN201910718215.XA 2019-08-05 2019-08-05 基于PCANet与高阶CRF的SAR图像分割方法 Active CN110517277B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910718215.XA CN110517277B (zh) 2019-08-05 2019-08-05 基于PCANet与高阶CRF的SAR图像分割方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910718215.XA CN110517277B (zh) 2019-08-05 2019-08-05 基于PCANet与高阶CRF的SAR图像分割方法

Publications (2)

Publication Number Publication Date
CN110517277A CN110517277A (zh) 2019-11-29
CN110517277B true CN110517277B (zh) 2022-12-06

Family

ID=68625153

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910718215.XA Active CN110517277B (zh) 2019-08-05 2019-08-05 基于PCANet与高阶CRF的SAR图像分割方法

Country Status (1)

Country Link
CN (1) CN110517277B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111583502B (zh) * 2020-05-08 2022-06-03 辽宁科技大学 基于深度卷积神经网络的人民币冠字号多标签识别方法

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101964112A (zh) * 2010-10-29 2011-02-02 上海交通大学 自适应的基于先验形状的图像分割方法
CN109345537A (zh) * 2018-08-17 2019-02-15 西安电子科技大学 基于高阶多尺度crf半监督的sar图像分割方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106157307B (zh) * 2016-06-27 2018-09-11 浙江工商大学 一种基于多尺度cnn和连续crf的单目图像深度估计方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101964112A (zh) * 2010-10-29 2011-02-02 上海交通大学 自适应的基于先验形状的图像分割方法
CN109345537A (zh) * 2018-08-17 2019-02-15 西安电子科技大学 基于高阶多尺度crf半监督的sar图像分割方法

Also Published As

Publication number Publication date
CN110517277A (zh) 2019-11-29

Similar Documents

Publication Publication Date Title
CN109345508B (zh) 一种基于两阶段神经网络的骨龄评价方法
Lee et al. Unsupervised multistage image classification using hierarchical clustering with a Bayesian similarity measure
Tso et al. A contextual classification scheme based on MRF model with improved parameter estimation and multiscale fuzzy line process
CN108921877B (zh) 一种基于宽度学习的长期目标跟踪方法
US8214177B2 (en) Optimized stochastic resonance signal detection method
US9026404B2 (en) Methods of improving detectors and classifiers using optimized stochastic resonance noise
CN106157330B (zh) 一种基于目标联合外观模型的视觉跟踪方法
CN109801305B (zh) 基于深度胶囊网络的sar图像变化检测方法
CN112446436A (zh) 基于生成对抗网络的抗模糊无人驾驶车辆多目标跟踪方法
CN111008644B (zh) 基于局部动态能量函数fcn-crf模型的生态变化监测方法
Malladi et al. Image denoising using superpixel-based PCA
CN112508963A (zh) 一种基于模糊c均值聚类的sar图像分割方法
CN110517277B (zh) 基于PCANet与高阶CRF的SAR图像分割方法
CN111626380A (zh) 一种基于超像素和卷积网络的极化sar图像分类方法
Sekhar et al. Segmentation of natural images with k-means and hierarchical algorithm based on mixture of Pearson distributions
Bahadure et al. Performance analysis of image segmentation using watershed algorithm, fuzzy c-means of clustering algorithm and simulink design
CN109272515B (zh) 基于高阶多尺度crf无监督的sar图像分割方法
CN108932520B (zh) 结合先验概率估计的sar影像水体概率制图方法
CN108846407B (zh) 基于独立成分高序不确定脑网络的磁共振影像分类方法
Astuti et al. Multilevel thresholding hyperspectral image segmentation based on independent component analysis and swarm optimization methods
Kiwanuka et al. Automatic attribute threshold selection for morphological connected attribute filters
An et al. Unsupervised segmentation of hidden Markov fields corrupted by correlated non-Gaussian noise
CN113838066A (zh) 一种基于改进模糊c均值聚类算法的彩色图像分割方法
CN113160345A (zh) 一种基于ConvLSTM的时间序列影像重建方法
Mirghasemi et al. Severely noisy image segmentation via wavelet shrinkage using PSO and Fuzzy C-Means

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