CN105956599B - 合成孔径雷达图像建模方法及装置 - Google Patents

Abstract

本发明公开了一种合成孔径雷达图像建模方法及装置，方法包括：步骤1，创建灰度等级直方图Y；步骤2，构建广义伽玛混合模型；步骤3，基于HECM‑MML‑GΓMM算法，估计广义伽玛混合模型的混合权重值、各广义伽玛分量的参数以及最优的混合模型分量数，从而确定所构建的广义伽玛混合模型，完成图像建模。本发明能够广泛地应用于幅度和强度、具有单峰或多峰分布特性、具有同纹理或异质纹理特性的SAR图像精确建模；参数具有闭式表达形式，易于实现，计算效率高；保留了图像中的细节，减少了信息丢失，有效地获取图像分类结果。

Description

合成孔径雷达图像建模方法及装置

技术领域

本发明涉及数字图像处理技术领域，具体地，涉及一种合成孔径雷达图像建模方法及装置。

背景技术

作为重要的对地观测技术手段，搭载合成孔径雷达(Synthetic Aperture Radar，SAR)传感器的遥感平台(机载或星载)具有全天候、全时段、全地形条件工作的优势，在军事侦测、资源探测、海洋预报、灾害监测等方面具有广泛的应用。当今，SAR技术的快速发展也极大地扩展了其应用范围。新一代高分辨率搭载于卫星和飞机的SAR系统可以获得较以往更为丰富的图像信息。但是，这也对SAR图像的解译提出了新的挑战，即对于包含不同地物信息，具有不同纹理特性的SAR图像数据进行精确的统计建模和图像分类。

SAR图像的统计建模是一种基于数据的建模方法。根据SAR图像数据计算图像的灰度直方图，建立模型，估计模型中的特征参数，进而实现图像分类，区分图像数据中的地物、场景等信息。该过程主要包括四个步骤：1)图像预处理；2)建立统计模型；3)估计参数；4)后处理及性能评估。本发明主要涉及第2步建立统计模型和第3步估计参数。

截至目前，国内外学者提出了许多SAR图像统计建模方法。按照所采用的方法不同，可以分为三类：

1)非参数方法。其完全根据图像数据信息进行建模，计算量巨大，耗时长。

2)参数化方法。其基本思想是，先验地确定数学模型，将概率密度估计问题转化为参数估计问题。相比非参数方法而言，参数化方法计算简单，易于应用。尽管如此，参数化方法在高分辨率SAR图像处理应用中受到了极大限制，其原因在于高分辨率SAR图像数据中包含了更多的细节信息，其直方图表现出重拖尾或者多峰的特点，使得单个参数化模型难以准确描述原图像的分布情况。

3)半参数化方法。作为一种折衷的方法，半参数化方法可以有效地解决以上问题，将分布模型按照不同的权重进行线性组合，以此估计图像数据的分布情况，进行精确统计建模。半参数化方法相对简单，应用灵活，依托于数据信息，保证了模型的有效性。因而，得到了广泛的研究和应用。

G.Moser等人在文献[G.Moser,J.Zerubia,and S.B.Serpico,“Dictionary-basedstochastic expectation maximization for SAR amplitude probability densityfunction estimation,”IEEE Transactions on Geoscience and Remote Sensing,vol.44,no.1,pp.188-200,2006.]提出了一种基于模型字典的统计建模方法，将随机期望最大化(Stochastic Expectation Maximization，SEM)算法与对数累积量方法(Method ofLog Cumulants，MoLC)结合在一起，同时估计模型参数、混合模型分量数并处理分量模型选择问题。但是，模型字典方法本身具有枚举性质，只有修改字典内待选模型后才可用于新的场景，灵活性有限；字典内的待选模型通常只适用于低分辨率或者中等分辨率的SAR图像，而其他模型不存在解析形式，导致更为庞大的计算量；最后，采用MoLC方法和SEM算法执行参数估计，不能保证对数似然函数单调性，从而难以估计有效混合分量总数。

之后，V.A.Krylov等人在文献[V.A.Krylov,G.Moser,S.B.Serpico,andJ.Zerubia,“Enhanced dictionary-based SAR amplitude distribution estimationand its validation with very high-resolution data,”IEEE Geoscience and RemoteSensing letters,vol.8,no.1,pp.148-152,2011.]中，扩展字典内待选分布模型，将模型字典方法扩展应用于高分辨率SAR图像，同时，改进了确定有效分量总数的方法。这在一定程度上提高了基于字典方法的准确性，但是其灵活性和计算效率仍受到制约。因此，有必要发明一种简单、灵活、计算效率较高的统计建模与参数估计方法，以提高其在实际中的应用价值。

发明内容

本发明的目的是提供一种合成孔径雷达图像建模方法及装置，该方法及装置为基于广义伽玛混合模型(Generalized Gamma Mixture Models，GΓMM)的半参数化的建模方法及装置。

为了实现上述目的，本发明提供了一种合成孔径雷达图像建模方法，包括以下步骤：

步骤1：对输入的SAR图像数据X，创建非零灰度值的非归一化的灰度等级直方图Y，其中，Y＝{h(r)；r＝1，2，...，L-1}，r为所述SAR图像数据X的灰度等级，L为所述灰度等级的个数，h(r)为所述SAR图像数据X中灰度等级为r的数据的个数；

步骤2：根据如下公式，将广义伽玛分布作为广义伽玛混合模型的多个分量，并将所述多个分量进行加权累加，得到所述广义伽玛混合模型：其中，M为所述分量的个数，π_m为第m个分量的混合权重值，满足且π_m≥0，θ_m＝{v_m，κ_m，σ_m；m＝1,2,,M}为第m个分量的参数集合，其中，v_m为幂参数，κ_m为形状参数，σ_m为尺度参数，Θ＝{M，π₁，π₂，…，π_M，θ₁，θ₂，…，θ_M}为所述广义伽玛混合模型的参数集合，GΓD(·)为所述广义伽玛分布的概率密度函数，满足

步骤3：基于HECM-MML-GΓMM算法对所述广义伽玛混合模型的参数集合Θ进行估计，得到分量的个数的估计量各分量的混合权重值的估计量各分量的参数集合的估计量以获得所述广义伽玛混合模型的参数集合的估计量从而完成所述广义伽玛混合模型的建立。

优选地，所述步骤3进一步包括：

步骤3-00：将分量的个数的估计变量以及第一、第二、第三和第四参数子集的估计变量和分别设置为初始值 和并将迭代次数t设置为初始值0，其中，第一参数子集为各分量的所述混合权重值的估计变量的集合，第二参数子集为各分量的所述幂参数的估计变量的集合，第三参数子集为各分量的所述形状参数的估计变量的集合，第四参数子集为各分量的所述尺度参数的估计变量的集合；

步骤3-05：根据如下公式计算灰度等级为r的数据属于第m个分量的后验概率

步骤3-10：对个分量中的每个分量，根据如下公式更新所述混合权重值的估计变量：其中，为所述广义伽玛混合模型的参数集合的估计变量，然后根据如下公式更新所述第一参数子集

步骤3-15：删除所述混合权重值的估计变量小于或等于0的分量，并更新所述广义伽玛混合模型的参数集合的估计变量

步骤3-20：对个分量中的每个分量，根据如下公式更新所述幂参数的估计变量：

其中，然后根据如下公式更新所述第二参数子集

步骤3-30：对个分量中的每个分量，根据如下公式更新所述形状参数的估计变量的具有单调性的双伽玛函数

并根据的单调性计算然后根据如下公式更新所述第三参数子集

步骤3-40：对个分量中的每个分量，根据如下公式更新所述尺度参数的估计变量：然后根据如下公式更新所述第四参数子集

步骤3-50：更新所述迭代次数t＝t+1，记录本次迭代过程得到的所述广义伽玛混合模型的参数集合的估计变量

步骤3-53：对本次迭代过程得到的所述广义伽玛混合模型的参数集合的估计变量根据如下公式计算不完全数据惩罚对数似然函数

其中，

为所述SAR图像数据X中非零灰度值的像素的个数，

为所述混合权重值的估计变量大于零的所有分量的总个数；

步骤3-55：判断所述不完全数据惩罚对数似然函数是否满足收敛条件，若是则继续执行步骤3-56，若否则返回执行步骤3-05；

步骤3-56：将所述不完全数据惩罚对数似然函数记录为有效对数似然函数

步骤3-58：删除所述混合权重值的估计变量在个分量中为最小值的分量，并更新所述广义伽玛混合模型的参数集合的估计变量

步骤3-60：判断所述分量的个数的估计变量是否小于预定的最小分量数M_min，若否则返回执行步骤3-05，若是则执行步骤3-70；

步骤3-70：找出满足如下公式的该为所述广义伽玛混合模型的参数集合的估计量

其中，T为所述有效对数似然函数的个数。

优选地，在所述步骤3-55中所述收敛条件为： 其中∈为预定的极小值正数。

优选地，所述步骤3-00中的初始值和根据如下方法确定：所述分量的个数的估计变量的初始值其中M_max为最大分量数；所述第一参数子集的初始值设置为所述第二参数子集的初始值对于幅度SAR图像数据设置为对于强度SAR图像数据设置为所述第三参数子集的初始值设置为所述第四参数子集的初始值对于幅度SAR图像数据，通过均匀选择M_max个灰度等级分别作为各瑞利分量的众数来估计对于强度SAR图像数据，通过均匀选择M_max个灰度等级分别作为各指数分量的均值来估计其中，

优选地，所述根据的单调性计算具体为采用标准二分法求解

优选地，所述最大分量数M_max的取值范围为20～40，所述预定的最小分量数M_min的取值为1或2。

本发明还提供了一种合成孔径雷达图像建模装置，包括：

第一处理模块，用于对输入的SAR图像数据X，创建非零灰度值的非归一化的灰度等级直方图Y，其中，Y＝{h(r)；r＝1，2，...，L-1}，r为所述SAR图像数据X的灰度等级，L为所述灰度等级的个数，h(r)为所述SAR图像数据X中灰度等级为r的数据的个数；

第二处理模块，用于根据如下公式，将广义伽玛分布作为广义伽玛混合模型的多个分量，并将所述多个分量进行加权累加，得到所述广义伽玛混合模型：

其中，

M为所述分量的个数，

π_m为第m个分量的混合权重值，满足且π_m≥0，

θ_m＝{v_m，κ_m，σ_m；m＝1，2，...，M}为第m个分量的参数集合，其中，v_m为幂参数，κ_m为形状参数，σ_m为尺度参数，

Θ＝{M，π₁，π₂，…，π_M，θ₁，θ₂，…，θ_M}为所述广义伽玛混合模型的参数集合，

GΓD(·)为所述广义伽玛分布的概率密度函数，满足

第三处理模块，用于基于HECM-MML-GΓMM算法对所述广义伽玛混合模型的参数集合Θ进行估计，得到分量的个数的估计量各分量的混合权重值的估计量各分量的参数集合的估计量以获得所述广义伽玛混合模型的参数集合的估计量从而完成所述广义伽玛混合模型的建立。

通过上述技术方案，本发明的合成孔径雷达图像建模方法及装置具有如下有益效果：

1)本发明采用广义伽玛分布(GΓD)作为广义伽玛混合模型(GΓMM)的基本分量。对比于现有技术中大多数的已知模型，该经验模型能更加有效地估计SAR图像数据分布，尤其在高分辨率SAR图像处理中能够合理地处理重拖尾和多峰情况，有效区分各种地表图像信息，获得更为清晰、准确的图像分类结果。

2)本发明采用HECM-MML-GΓMM算法高效地进行模型的参数估计，在优选的实施例中，每执行一次基于直方图的条件期望最大化(Expectation-conditional Maximization，ECM)算法，就可以完成参数估计与模型选择的操作。ECM算法将复杂的最大化步分解为多个较易计算的条件最大化步，使算法的最大化步骤更易于实现；同时，各参数的更新方程具有简单紧凑的显示表达，在保证拟合效果的同时，大大缩短了运算时间。

3)本发明采用的ECM算法为基于直方图的ECM算法，该算法构造了一种关于灰度等级的函数，使得算法运算量与图像的大小无关，降低了参数估计的计算复杂度。

4)本发明建立的广义伽玛混合模型(GΓMM)具有高度的灵活性。其灵活性表现在三个方面：对于幅度格式和强度格式的SAR图像建模和分类均有优异的性能表现；可作为通用的建模方法，用于处理8-bit或16-bit量化的SAR图像；对于强度或幅度格式、具有单峰或多峰分布特性和具有不同异质度的SAR图像数据(幅度或强度格式的SAR图像数据)具有良好的拟合效果。

本发明的其他特征和优点将在随后的具体实施方式部分予以详细说明。

附图说明

附图是用来提供对本发明的进一步理解，并且构成说明书的一部分，与下面的具体实施方式一起用于解释本发明，但并不构成对本发明的限制。在附图中：

图1是本发明的合成孔径雷达图像建模方法的流程示意图；

图2是图1中优选的实施例中步骤3的流程示意图；

图3a是仿真实验使用的第一组SAR图像之一，其为由TerraSAR-X传感器获得的某第一地区的卫星图像；

图3b是仿真实验使用的第一组SAR图像之二，其为由TerraSAR-X传感器获得的某第二地区的卫星图像；

图3c是仿真实验使用的第一组SAR图像之三，其为由机载EMISAR系统获得的某第三地区的L波段SAR图像；

图4a是仿真实验使用的第二组SAR图像之一，其为某第四地区的垂直极化CONVAIR-580SAR图像；

图4b是仿真实验使用的第二组SAR图像之二，其为某第五地区的具有16位量化的SAR图像；

图5a、图5b、图5c分别是图3a、图3b、图3c的SAR图像的直方图分布情况以及本发明的建模方法和对比方法建模后拟合得到的分布情况；

图6a和图6c分别是用对比方法处理图3b后的带有分类信息的概率密度曲线以及相应的图像分类的结果；

图6b和图6d分别是用对比方法处理图3c后的带有分类信息的概率密度曲线以及相应的图像分类的结果；

图6e和图6g分别是用本发明的方法处理图3b后的带有分类信息的概率密度曲线以及相应的图像分类的结果；

图6f和图6h分别是用本发明的方法处理图3c后的带有分类信息的概率密度曲线以及相应的图像分类的结果；

图7a是图4a的SAR图像的直方图分布情况以及本发明的建模方法建模后拟合得到的分布情况；

图7b是图4b的SAR图像的直方图分布情况以及本发明的建模方法建模后拟合得到的分布情况；

图8是本发明的合成孔径雷达图像建模装置的框图。

具体实施方式

以下结合附图对本发明的具体实施方式进行详细说明。应当理解的是，此处所描述的具体实施方式仅用于说明和解释本发明，并不用于限制本发明。在阅读这样的描述之后，本领域技术人员很容易构造其他修改或替换，这样的修改或替换应被理解为落入本发明的范围中。

图1给出了本发明的合成孔径雷达图像建模方法的流程示意图，具体实施包含以下步骤：

步骤1：对输入的SAR图像数据X，创建非零灰度值的非归一化的灰度等级直方图Y，其中，Y＝{h(r)；r＝1，2，...，L-1}，r为SAR图像数据X的灰度等级，L为灰度等级的个数，h(r)为SAR图像数据X中灰度等级为r的数据的个数；

步骤2：将M个广义伽玛分布(GΓD)作为广义伽玛混合模型(GΓMM)的多个分量，分别按照相应的混合权重值进行累加，建立GΓMM，即：

其中，M为分量的个数；π_m为第m个分量的混合权重值，满足且π_m≥0；θ_m＝{v_m，κ_m，σ_m；m＝1，2，...，M}为第m个分量的参数集合，其中，v_m为幂参数，κ_m为形状参数，σ_m为尺度参数；Θ＝{M，π₁，π₂，…，π_M，θ₁，θ₂，…，θ_M}为GΓMM的参数集合；GΓD(·)为GΓD的概率密度函数，满足如下公式2：

步骤3：基于HECM-MML-GΓMM算法对GΓMM的参数集合Θ进行估计，得到分量的个数的估计量各分量的混合权重值的估计量各分量的参数集合的估计量以获得GΓMM的参数集合的估计量从而完成广义伽玛混合模型的建立，即完成SAR图像的建模过程。

为了使得本发明的步骤3中对GΓMM的参数集合Θ进行估计的运算量尽可能得小，本发明提供了GΓMM的参数更新的闭式表达式，以下描述GΓMM的参数更新的闭式表达式的推导过程。

该推导过程包括：在GΓMM的基础之上，基于类MML准则，构造具有惩罚函数形式的完全数据的对数似然函数并进一步获得用于HECM-MML-GΓMM算法的函数函数将用于推导GΓMM中各参数的闭式解。该推导过程的具体实施过程可分为如下五个步骤：

(S1)在独立性假设的前提下，待估计的参数的联合先验分布具有如下形式：

p(π₁，π₂，…，π_M)以及p(θ_m)分别表示混合权重参数π_m的联合先验分布和第m个分量的参数集合θ_m的先验分布。

(S2)采用Jeffery无信息先验分布，确定步聚(S1)中的混合权重参数π_m的联合先验分布p(π₁，π₂，…，π_M)和第m个分量的参数集合θ_m的先验分布p(θ_m)，即：

其中，|·|为行列式的运算符号，A为多项式分布的费雪信息矩阵(FisherInformation Matrix)，I⁽¹⁾(θ_m)为与第m个分量相关的单个观测数据的费雪信息矩阵。

(S3)在常规GΓMM的极大似然函数中，引入隐变量Z。隐变量Z＝{z^(r)|r＝1，2，...L-1}是观测数据所对应的二元向量的集合，用于表征灰度等级为r的数据来源。具体而言，对于二元向量有且仅有其余元素均为0表明，灰度等级为r的数据仅来自于混合模型中的第m个分量。从而，得到完全数据的对数似然函数：

(S4)采用类MML准则(MML-like Criterion)，整合(S1)，(S2)，(S3)的结果，得到具有惩罚函数形式的完全数据的对数似然函数以及用于HECM-MML-GΓMM算法的函数。

上述类MML准则具有如下形式：

其中，p(Θ)为参数Θ的先验分布，p(X|Θ)为观测数据的条件分布，I(Θ)为期望费雪信息矩阵，c为GΓMM中自由参数的总数。

为了与极大似然估计法相匹配，对类MML准则的公式7做等效变换：将表达式取反，同时为使目标函数取得极小值，而进行的求解M和Θ的运算(argmin_M,Θ{·})变为使目标函数取极大值，求解M和Θ的运算(argmax_M,Θ{·})的运算。结合(S1)，(S2)，(S3)的结果，具有惩罚函数形式的完全数据的对数似然函数表达式如下：

用于HECM-MML-GΓMM算法的函数如下：

其中，为所述SAR图像数据X中非零灰度值的像素的总数，K_nz为具有非零混合权重值的分量总数，其数学表达为K_nz＝∑_m[π_m＞0]。

(S5)根据函数推导参数的更新方程，得到：

混合权重值的估计变量的更新方程如下：

其中，为灰度等级为r的数据属于第m个分量的后验概率，为所述广义伽玛混合模型的参数集合的估计变量；

幂参数的估计变量的更新方程如下：

其中，

形状参数的估计变量的具有单调性的双伽玛函数(Digamma Function)其定义为常规伽玛函数的对数的导数，如下所示：

尺度参数的估计变量：

下面结合图2描述优选的实施例中步骤3的执行过程，即如何基于HECM-MML-GΓMM算法对所述广义伽玛混合模型的参数集合Θ进行估计，得到分量的个数的估计量各分量的混合权重值的估计量各分量的参数集合的估计量以获得所述广义伽玛混合模型的参数集合的估计量

步骤3-00：将分量的个数的估计变量以及第一、第二、第三和第四参数子集的估计变量和分别设置为初始值以及和并将迭代次数t设置为初始值0，其中，

第一参数子集为各分量的混合权重值的估计变量的集合，

第二参数子集为各分量的幂参数的估计变量的集合，

第三参数子集为各分量的形状参数的估计变量的集合，

第四参数子集为各分量的尺度参数的估计变量的集合。

优选的，初始值和可采用如下方法确定：

a)对于分量的个数的估计变量的初始值可设置为最大分量数M_max，通常M_max的取值范围为20～40。

b)对于第一参数子集的初始值将其中的各个分量的混合权重值设置为等概率分布，即其中，

c)对于第二、第三和第四参数子集和的初始值 和的设置，针对幅度和强度SAR图像数据，分别采用瑞利(Rayleigh)分布和指数分布进行初始化。瑞利分布和指数分布均为广义伽玛分布的特例，瑞利分布为情况下广义伽玛分布的特例，指数分布为情况下广义伽玛分布的特例。

具体为：对于幅度SAR图像数据，设置每个分量的 (此时对应于瑞利分布)，并均匀选择M_max个灰度等级分别作为各瑞利分量的众数进而估计得到对于强度SAR图像数据，设置每个分量的(此时对应于指数分布)，并均匀选择M_max个灰度等级分别作为各指数分量的均值来估计得到

步骤3-05：为了后续步骤的计算，根据如下公式计算灰度等级为r的数据属于第m个分量的后验概率

步骤3-10：对个分量中的每个分量，根据公式10更新混合权重值的估计变量，然后根据如下公式14更新第一参数子集

步骤3-15：删除混合权重值的估计变量小于或等于0的分量，并更新广义伽玛混合模型的参数集合的估计变量包括分量的个数的估计变量以及第一、第二、第三和第四参数子集的估计变量 和即将分量的个数的估计变量的值更新为减去混合权重值小于或等于0的分量的个数，将第一、第二、第三和第四参数子集的估计变量和中混合权重值小于或等于0的分量的参数删除。

步骤3-20：对个分量中的每个分量，根据公式11更新幂参数的估计变量，然后根据如下公式15更新第二参数子集

步骤3-30：对个分量中的每个分量，根据公式12更新形状参数的估计变量的具有单调性的双伽玛函数(Digamma function)并根据的单调特性利用标准二分法计算然后根据如下公式16更新第三参数子集

用标准二分法求解单调函数的零点参数的方法为本领域的公知方法。例如，先选定合适的区间[a，b]，通过不断地把函数f(x)(满足f(a)·f(b)<0)的零点所在的区间一分为二，使区间的两个端点逐步逼近零点，进而得到零点的近似值。

步骤3-40：对个分量中的每个分量，根据公式13更新尺度参数的估计变量，然后根据如下公式17更新第四参数子集

步骤3-50：更新迭代次数t＝t＋1，记录本次迭代过程得到的GΓMM的参数集合的估计变量

步骤3-53：对本次迭代过程得到的广义伽玛混合模型的参数集合的估计变量根据如下公式计算不完全数据惩罚对数似然函数

其中，

为SAR图像数据X中非零灰度值的像素的个数，

为混合权重值的估计变量大于零的所有分量的总个数。

步骤3-55：判断不完全数据惩罚对数似然函数是否满足收敛条件，若是则继续执行步骤3-56，若否则返回执行步骤3-05。

收敛条件为：其中∈为预定的极小值正数。即，当相邻两次的不完全数据惩罚对数似然函数值变化的绝对值小于某一极小值正数∈时，认为算法收敛。

步骤3-56：将不完全数据惩罚对数似然函数记录为有效对数似然函数

步骤3-58：删除混合权重值的估计变量在个分量中为最小值的分量，并更新广义伽玛混合模型的参数集合的估计变量

步骤3-60：判断分量的个数的估计变量是否小于预定的最小分量数M_min(该预定的最小分量数M_min可以取一个较小的值，例如，其取值可以为1或2)，若大于或等于则返回执行步骤3-05继续参数更新的迭代过程，若小于则结束参数更新的迭代过程执行步骤3-70。

步骤3-70：找出满足如下公式的该为所述广义伽玛混合模型的参数集合的估计量

其中，T为有效对数似然函数的个数。

如图2所示，本发明的合成孔径雷达图像建模方法，其包括两个循环，分别为步骤3-55判断完收敛条件然后返回步骤3-05的内循环、和步骤3-60判断完分量的个数然后返回步骤3-05的外循环。

内循环用于估计混合模型参数。内循环收敛表明完成一次建模。收敛时，要记录当前的不完全数据对数似然函数为有效对数似然函数(步骤3-58)，以及记录估计的混合模型参数。

外循环用于剔除具有最小混合权重值的分量，以防止具有最小混合权重值的分量本应被模型中某一分量吸收却被单独划出而产生较高的模型复杂度。剔除具有最小混合权重值的分量后，要再次进入内循环，重新估计模型参数，获得新的模型，记录收敛时的不完全数据对数似然函数及相关参数。如此往复。

直到模型的分量数最终达到预设的最少分量数，然后在全部记录的有效对数似然函数中，选取使最大的那组参数作为最优解(步骤3-70)。

以下结合图3a至7b的仿真实验数据来描述本发明的基于HECM-MML-GΓMM的合成孔径雷达图像建模方法的效果。

仿真实验条件：

本发明提出的基于HECM-MML-GΓMM的合成孔径雷达图像建模方法及对比方法均在相同的实验仿真平台上实现，软件平台为MATLAB 7.13。

仿真实验数据：

本发明的仿真实验使用的是真实的幅度SAR图像(图3a、3b和3c所示的第一组SAR图像，包括具有单峰和多峰分布形态特征的数据)和强度SAR图像(图4a)，包括8位(图3a、3b、3c和图4a)和16位(图4b)量化的SAR图像。

结果评价方法：

评价模型的精确性与有效性的方法，包括三类：

(1)基于视觉观察的评价方法。本发明输出的结果可以实现图像数据信息的分类。基于光学遥感图像，判断图像数据分类结果的合理性。

(2)分布比对。将图像数据的直方图与所得模型的分布情况相比较。分布比对方法可以检测数据建模的准确性。

(3)基于指标的评价方法。常用的评价指标有，Kolmogorov-Smirnov(KS)距离和对称Kullback-Leibler(KL)距离。

D_KS：Kolmogorov-Smirnov(KS)距离，表示拟合所得的分布函数F(x)与实际分布函数G(x)之间的差异值，即D_Ks＝max_x∈Ω|F(x)-G(x)|。KS距离表示实际分布与拟合分布的差异情况。

D_KL：对称Kullback-Leibler(KL)距离，表示拟合所得的概率真密度函数f(x)与归一化实际数据直方图h(x)之间的差异程度，即对称KL距离表示拟合概率密度与实际概率密度的不相似程度。

以上两个指标是评价模型准确性与有效性的重要评价依据，能够对微小的性能差异做出明确的判定。

对比仿真实验

为了说明本发明的有效性，本发明的仿真效果与如下两个方法进行比较：

对比方法1，文献[H.C.Li,W.Hong,Y.R.Wu and P.Z.Fan,“On the empirical-statistical modeling of SAR images with generalized Gamma distribution,”IEEEJournal of Selected Topics in Signal Processing,vol.5,no.3,pp.386-397,2011.]采用参数化的统计建模方法，以广义伽玛分布作为数学模型，利用多伽玛函数的一阶和二阶近似表达，获得了具有闭式形式的参数估计量。

对比方法2，文献[V.A.Krylov,G.Moser,S.B.Serpico,and J.Zerubia,“Enhanceddictionary-based SAR amplitude distribution estimation and its validationwith very high-resolution data,”IEEE Geoscience and Remote Sensing letters,vol.8,no.1,pp.148-152,2011.]将待选模型归入模型字典，字典中的模型可作为每个分量的备选分布；基于备选分布，考察归属于第m个有效分量的所有数据的对数似然函数，选择使对数似函数取得最大值的分布作为第m个有效分量的分布；采用SCM算法与MoLC方法执行参数估计，完成图像数据建模，实现图像分类。该方法简称为“EDSEM”。

仿真实验内容：

实验1：SAR图像数据分布估计与分析

本发明的仿真实验采用第一组三幅具有不同复杂分布情况的幅度SAR图像数据集。采用本发明所提出的建模方法，图5a、5b和5c给出了第一组仿真数据所对应的直方图分布情况与建模后所拟合得到的分布情况，同时给出了对比方法1、2和本发明的建模方法的拟合情况。其中，GΓD、EDSEM、GΓMM曲线分别表示对比方法1、对比方法2和本发明建模后拟合的分布情况。

利用D_KS和D_KL两个指标，定量地分析了本发明的准确性。表1列出了仿真实验输出图像的精度评估数据。

表1本发明与对比方法的建模精度评价指标

实验2：SAR图像分类

本实验采用第一组图像数据，基于实验1的结果，采用贝叶斯决策准则，进行数据分类操作，实现图像分类。图6a-6h是带有分类信息的概率密度曲线以及相应的图像分类结果，其中，图6a和图6c是用对比方法2处理图3b后的结果；图6b和图6d是用对比方法2处理图3c后的结果；图6e和图6g是用本发明的方法处理图3b后的结果；图6f和图6h是用本发明的方法处理图3c后的结果。

实验3：强度格式的SAR图像与16位量化的高分辨率SAR图像

本次仿真实验使用的是第二组图像数据，用于证明本发明具有广泛地应用性，能够有效地处理不同格式、不同分辨率的SAR图像数据。图4a和4b分别是强度格式SAR图像和16位量化的高分辨率SAR图像。图7a-7b是本发明对强度格式的SAR图像数据与高分辨率SAR图像数据的拟合结果。其中，图7a是本发明处理强度SAR图像数据分布的拟合结果，图7b是16位量化的高分辨率SAR图像直方图分布与本发明的仿真拟合结果的对比。

仿真结果分析：

定性地分析本发明与对比方法的拟合结果可知，本发明的拟合更为接近原始SAR图像数据的分布直方图；从定量的角度考察，用于评价拟合精度的指标表明，对于不同类型分布的SAR图像数据，本发明的D_KS和D_KL指标均优于对比方法。从而证明，本发明可以准确、广泛用于SAR图像建模。

通过分析本发明与对比方法的SAR图像分类结果，并与光学卫星照片相对照可知，分类结果图像显示了地表物体的景观，不同的色块代表不同的景观类型，包括水体、裸露地面、建筑、公路和农作物等。对比方法正确地分类了水体和陆地信息，但是无法对图像中的信息作进一步的分类，比如，无法完全区分地表的建筑、裸露陆地、植被覆盖的陆地。本发明相较于其他方法具备对图像精细分类的优点，保留了更为完整的细节信息，不仅能够有效地对陆地和水体信息进行分类，同时，也能够精致地区分出地表建筑、植被以及农作物。

通过分析本发明对强度SAR图像与16位量化的高分辨率SAR图像处理的结果，可知，本发明不仅能够精确地建模幅度SAR图像数据，而且对于强度SAR图像和高分辨率SAR图像也具有较好的拟合性能。

综上所述，本发明广泛适用于幅度和强度、具有单峰和多峰分布特征、不同纹理程度的SAR图像和高分辨率SAR图像建模和分类。基于广义伽玛分布，建立广义伽玛混合模型，采用基于直方图的HECM-MML-GΓMM算法执行参数估计，完成建模。通过如上分析可知，本方法能够对各类SAR图像数据进行精确建模，保留更多的图像细节信息；有效地区分SAR图像中的水体、地面、建筑、公路、植被、农作物等信息，减少了信息丢失，有效地获取图像分类结果。

本发明还提供一种合成孔径雷达图像建模装置，如图8所示，该装置800可以包括：第一处理模块810，用于对输入的SAR图像数据X，创建非零灰度值的非归一化的灰度等级直方图Y，其中，Y＝{h(r)；r＝1，2，...，L-1}，r为所述SAR图像数据X的灰度等级，L为所述灰度等级的个数，h(r)为所述SAR图像数据X中灰度等级为r的数据的个数；

第二处理模块820，用于根据如下公式，将广义伽玛分布作为广义伽玛混合模型的多个分量，并将所述多个分量进行加权累加，得到所述广义伽玛混合模型：

其中，

M为所述分量的个数，

π_m为第m个分量的混合权重值，满足且π_m≥0，

θ_m＝{v_m，κ_m，σ_m；m＝1，2，...，M}为第m个分量的参数集合，其中，v_m为幂参数，κ_m为形状参数，σ_m为尺度参数，

Θ＝{M，π₁，π₂，…，π_M，θ₁，θ₂，…，θ_M}为所述广义伽玛混合模型的参数集合，

GΓD(·)为所述广义伽玛分布的概率密度函数，满足

第三处理模块830，用于基于HECM-MML-GΓMM算法对所述广义伽玛混合模型的参数集合Θ进行估计，得到分量的个数的估计量各分量的混合权重值的估计量各分量的参数集合的估计量以获得所述广义伽玛混合模型的参数集合的估计量从而完成所述广义伽玛混合模型的建立。

关于上述实施例中的装置，其中各个模块执行操作的具体方式已经在有关该方法的实施例中进行了详细描述，此处将不做详细阐述说明。

以上结合附图详细描述了本发明的优选实施方式，但是，本发明并不限于上述实施方式中的具体细节，在本发明的技术构思范围内，可以对本发明的技术方案进行多种简单变型，这些简单变型均属于本发明的保护范围。

另外需要说明的是，在上述具体实施方式中所描述的各个具体技术特征，在不矛盾的情况下，可以通过任何合适的方式进行组合，为了避免不必要的重复，本发明对各种可能的组合方式不再另行说明。

此外，本发明的各种不同的实施方式之间也可以进行任意组合，只要其不违背本发明的思想，其同样应当视为本发明所公开的内容。

Claims (5)

1.一种合成孔径雷达图像建模方法，其特征在于，包括以下步骤：

步骤1：对输入的SAR图像数据X，创建非零灰度值的非归一化的灰度等级直方图Y，其中，Y＝{h(r)；r＝1，2，...，L-1}，r为所述SAR图像数据X的灰度等级，L为所述灰度等级的个数，h(r)为所述SAR图像数据X中灰度等级为r的数据的个数；

步骤2：根据如下公式，将广义伽玛分布作为广义伽玛混合模型的多个分量，并将所述多个分量进行加权累加，得到所述广义伽玛混合模型：

其中，

M为所述分量的个数，

π_m为第m个分量的混合权重值，满足且π_m≥0，

θ_m＝{v_m，κ_m，σ_m；m＝1，2，...，M}为第m个分量的参数集合，其中，v_m为幂参数，κ_m为形状参数，σ_m为尺度参数，

Θ＝{M，π₁，π₂，...，π_M，θ₁，θ₂，...，θ_M}为所述广义伽玛混合模型的参数集合，

GΓD(·)为所述广义伽玛分布的概率密度函数，满足 其中，Γ(κ_m)是关于参数κ_m的伽玛函数；

步骤3：基于HECM-MML-GΓMM算法对所述广义伽玛混合模型的参数集合Θ进行估计，得到分量的个数的估计量各分量的混合权重值的估计量各分量的参数集合的估计量以获得所述广义伽玛混合模型的参数集合的估计量从而完成所述广义伽玛混合模型的建立，所述HECM-MML-GΓMM算法是指针对广义伽玛混合模型GΓMM，采用基于直方图条件期望最大化HECM和最小信息长度准则MML对参数集合Θ进行估计的无监督学习算法。

2.根据权利要求1所述的方法，其特征在于，所述步骤3进一步包括：

步骤3-00：将分量的个数的估计变量以及第一、第二、第三和第四参数子集的估计变量和分别设置为初始值 和并将迭代次数t设置为初始值0，其中，

第一参数子集为各分量的所述混合权重值的估计变量的集合，

第二参数子集为各分量的所述幂参数的估计变量的集合，

第三参数子集为各分量的所述形状参数的估计变量的集合，

第四参数子集为各分量的所述尺度参数的估计变量的集合；

步骤3-05：根据如下公式计算灰度等级为r的数据属于第m个分量的后验概率

步骤3-10：对个分量中的每个分量，根据如下公式更新所述混合权重值的估计变量：

其中，为所述广义伽玛混合模型的参数集合的估计变量，

然后根据如下公式更新所述第一参数子集

步骤3-15：删除所述混合权重值的估计变量小于或等于0的分量，并更新所述广义伽玛混合模型的参数集合的估计变量

步骤3-20：对个分量中的每个分量，根据如下公式更新所述幂参数的估计变量：

其中，

然后根据如下公式更新所述第二参数子集

步骤3-30：对个分量中的每个分量，根据如下公式更新所述形状参数的估计变量的具有单调性的双伽玛函数

并根据的单调性计算

然后根据如下公式更新所述第三参数子集

步骤3-40：对个分量中的每个分量，根据如下公式更新所述尺度参数的估计变量：

然后根据如下公式更新所述第四参数子集

步骤3-50：更新所述迭代次数t＝t+1，记录本次迭代过程得到的所述广义伽玛混合模型的参数集合的估计变量

步骤3-53：对本次迭代过程得到的所述广义伽玛混合模型的参数集合的估计变量根据如下公式计算不完全数据惩罚对数似然函数

其中，

为所述SAR图像数据X中非零灰度值的像素的个数，

为所述混合权重值的估计变量大于零的所有分量的总个数；

步骤3-55：判断所述不完全数据惩罚对数似然函数是否满足收敛条件，若是则继续执行步骤3-56，若否则返回执行步骤3-05；

步骤3-56：将所述不完全数据惩罚对数似然函数记录为有效对数似然函数

步骤3-58：删除所述混合权重值的估计变量在个分量中为最小值的分量，并更新所述广义伽玛混合模型的参数集合的估计变量

步骤3-60：判断所述分量的个数的估计变量是否小于预定的最小分量数M_min，若否则返回执行步骤3-05，若是则执行步骤3-70；

步骤3-70：找出满足如下公式的该为所述广义伽玛混合模型的参数集合的估计量

其中，T 为所述有效对数似然函数的个数。

3.根据权利要求2所述的方法，其特征在于，在所述步骤3-55中所述收敛条件为：其中∈为预定的极小值正数。

4.根据权利要求2所述的方法，其特征在于，所述根据的单调性计算具体为采用标准二分法求解

5.一种合成孔径雷达图像建模装置，其特征在于，包括：

第一处理模块，用于对输入的SAR图像数据X，创建非零灰度值的非归一化的灰度等级直方图Y，其中，Y＝{h(r)；r＝1，2，...，L-1}，r为所述SAR图像数据X的灰度等级，L为所述灰度等级的个数，h(r)为所述SAR图像数据X中灰度等级为r的数据的个数；

第二处理模块，用于根据如下公式，将广义伽玛分布作为广义伽玛混合模型的多个分量，并将所述多个分量进行加权累加，得到所述广义伽玛混合模型：

其中，

M为所述分量的个数，

π_m为第m个分量的混合权重值，满足且π_m≥0，

θ_m＝{v_m，κ_m，σ_m；m＝1，2，...，M}为第m个分量的参数集合，其中，v_m为幂参数，κ_m为形状参数，σ_m为尺度参数，

Θ＝{M，π₁，π₂，...，π_M，θ₁，θ₂，...，θ_M}为所述广义伽玛混合模型的参数集合，

GΓD(·)为所述广义伽玛分布的概率密度函数，满足 其中，Γ(κ_m)是关于参数κ_m的伽玛函数；

第三处理模块，用于基于HECM-MML-GΓMM算法对所述广义伽玛混合模型的参数集合Θ进行估计，得到分量的个数的估计量各分量的混合权重值的估计量各分量的参数集合的估计量以获得所述广义伽玛混合模型的参数集合的估计量从而完成所述广义伽玛混合模型的建立，所述HECM-MML-GΓMM算法是指针对广义伽玛混合模型GΓMM，采用基于直方图条件期望最大化HECM和最小信息长度准则MML对参数集合Θ进行估计的无监督学习算法。