Y)
[0064] 其中，Sp为Pw的变化测度，X为主SAR图像内的子区域，MLLp (X)表示区域X内 Pw的最大对数似然函数，Y为辅SAR图像内的子区域，MLLp(Y)表示区域Y内Pw的最大对 数似然函数，区域X与区域Y为主、辅SAR图像内的对应区域。其中，在X与Y均未变化的 特殊情况下，Sp = 0。该Sp的推导方法在下文中具体描述，该Sp中同时考虑了极不均匀纹 理和斑噪的变化，区别于现有技术中将纹理视为均匀的情况，当然，现有技术中纹理视为均 匀的情况为本发明所考虑的极不均匀纹理的一种特殊情况，因而本发明也可适用于纹理均 匀的图像处理情况。
[0065] 此外，在步骤S6之前，根据本发明实施例的方法还可以包括：步骤S5,即通过对差 异图像进行阀值分割来确定其区域和边界，具体地，步骤S5包括：
[0066] S51、对差异图像进行自适应KI阀值分割得到区别变化区域与不变区域的作为二 值化图像的二值化图像（也可以称为二值化掩膜图像），在实际应用中，变化区域与不变区 域通常在二值化图像中以两种不同的颜色或其他有视觉差异的方式进行表达，；
[0067] S52、通过连通区搜索算法去除所述二值化图像中的变化面积（实际单块面积）小 于预定值（优选的，预定值为100平方米）的变化区域，在实际应用根据本发明的方法时， SAR通常为以米级为精度单位，因而根据SAR图像与实际面积的比例设置，可得到大小不同 的图像，本文中的"实际单块面积"指SAR采集的真实地表的面积。在SAR图像中，二值化 图像中100平方米的变化区域通常仅在肉眼中示出为一点，消除这些点使得变化区域变得 连续平滑是为了后续操作的需求服务；
[0068] S53、通过形态学的闭运算中的膨胀操作对变化区域进行膨胀得到膨胀区域，并对 该膨胀区域的周界进行记录；
[0069] S54、通过形态学的闭运算中的腐蚀操作，采用与所述膨胀操作相对应的比例对所 述膨胀区域进行腐蚀，得到腐蚀后区域；
[0070] S55、将腐蚀后区域的周界与膨胀区域的周界进行作差得到差异图像中所述变化 区域的边界。
[0071] 此外，在步骤S6之后可以进一步包括步骤S7 :将变化区域与主、辅SAR图像进行 叠合，确定主、辅SAR图像内的变化部分，即通过两个图像的叠合可以在像素点的精度里判 定哪些范围有变化，而哪些范围没有变化。并且，可以将得到的最终的变化区域与原来的主 辅图像进行叠合，输出土地利用图。
[0072] 根据本发明的另一个实施例，还提供一种检测全极化SAR图像的变化区域的装 置，全极化SAR图像基于极不均匀纹理，检测装置包括：
[0073] 参数确定模块，其配置为使用基于极不均匀纹理的多视极化SAR乘性模型对先后 接收的全极化主、辅SAR图像进行拟合，以分别得到主、辅SAR图像的纹理图像、斑噪协方差 图像，并分别根据主、辅SAR图像的纹理图像确定其对应的概率密度函数的分布参数，主、 辅SAR图像为包括同一对象区域的高分辨率全极化SAR图像；
[0074] 差异信息确定模块，其配置为根据主、辅SAR图像各自的纹理图像、斑噪协方差图 像和概率密度函数的分布参数，并通过基于极不均匀纹理的多视极化SAR乘性模型的相似 性度量对比主、辅SAR图像得到差异信息，基于差异信息确定相对于主、辅SAR图像的差异 图像；
[0075] 比对模块，其配置为将差异图像与主SAR图像和/或辅SAR图像进行比对，以检测 出变化区域。
[0076] 作为优选，上述装置进一步包括：
[0077] 区域匹配模块，其配置为对主、辅SAR图像进行同一区域的对应匹配。
[0078] 如图2所示为根据本发明的一个实施例的方法，具体包括：
[0079] 步骤S21 :预处理，即对采集的主、辅SAR图像的数据进行图像配准；
[0080] 步骤S22:进行差异图像提取，即利用根据本发明实施例的改进的多视极化SAR乘 性模型拟合极化主、辅SAR图像，并基于根据本发明的相似性试题进行差异图提取；
[0081] 步骤S23 :阈值分割，即对提取的差异图进行自适应的KI阀值分割，并得到二值化 掩膜；
[0082] 步骤S24 :图像整合，即将得到的二值化掩膜与主和/或辅SAR图像进行图像整 合，然后处理输出。
[0083] 以下对极化SAR乘性模型进行定义，以说明本发明中使用的基于极不均匀纹理的 乘性模型的概率密度函数的原理：
[0084] 单视（即，L = 1)情况下，SIRV模型（球不变随机矢量模型）认为全极化数据的 散射矢量k是由复圆高斯矢量z与一个正随机参数τ的平方根相乘得到的，因而可以称其 为乘性模型。此处τ表示图像的纹理参数，ζ表征图像的斑噪。根据这一定义，可以把L视 下未标准化的极化相干矩阵W表示为
[0088] 式中，矩阵R是对角线元素为图像纹理的是实对角矩阵，H代表共轭转秩，矩阵Zh 是由斑噪矢量组成的均值为O多元复高斯矩阵，记为
其roF(概率密度 函数）为
[0090] 为了推导矩阵W的概率密度函数，首先假设矩阵R是给定的常数矩阵，即假设图像 纹理是已知的。定义
那么矩阵A到矩阵Zh的雅可比行列式为
[0092] 这样，可以得到矩阵R给定情况下，矩阵A的PDF为
[0094] 从式⑴可以得到W = AHA，因此矩阵A到矩阵W的雅可比行列式为
[0096] 此处Γ p (L)为复数形式的多元Gamma函数。
[0097] 这样，联合式（6)和式（5)，就可以得到在矩阵R给定情况下，未标准化的极化相干 矩阵W的PDF为
[0099] 上式pw (WI R ;L，Σ)对矩阵R在从0到〇〇进行积分，就可以得到W的PDF为
[0101] 式⑶就是多视极化SAR乘性模型的表达式，或者叫做SIRV模型到多视情况的扩 展，记为W e s (L，R，Σ)。当图像纹理为Dirac delta函数时，即pR (R) = δ (R-I1)，其中矩 阵込为LXL的单位矩阵。这时图像中每个像素的纹理相等，多视极化SAR乘性模型就退 化为复Wishart分布模型。
[0103] 记为W e W(L，Σ)。这说明复Wishart分布模型可以看出多视极化SAR乘性模型 的一个特例。
[0104] 在式（8)中，由于纹理τ的分布是不确定的，因此矩阵R的分布很难求解。为了 简化计算，需要对R进行近似处理。当视数不大时，窗口内像素的纹理有很强的相关性。这 时，可以认为在窗口内，纹理是固定值，即τ L。这样，矩阵W变成
[0105] W= τ * ZZh (10)
[0106] 上式认为未标准化的相干矩阵是由纹理τ和未标准化的噪声协方差矩阵相乘得 到的。在这种情况下，式（8)变为
[0108] 以下基于上述基于极不均匀纹理的乘性模型的概率密度函数对极化SAR乘性模 型的参数进行估计，包括：
[0109] 在式（11)中，需要估计的参数由三类：纹理τ、斑噪协方差矩阵Σ和^(〇中 的纹理分布参数。其中斑噪协方差矩阵的估计是参数估计的关键，首先把W的PDF重写为
[0111] 其中
[0113] 对于窗口内的N个像素极化相干矩阵，假设他们之间是独立且同分布 (Independent and Identically Distributed, IID)的，这样N个相干矩阵的似然函数可以 表示为
[0115] 求解式中Σ的似然估计值，就是要找到式中似然函数Lw( ·)最大时Σ的值。这 个过程等价于令In Lw( ·)对Σ的一阶导数为0,即
[0117] 此处
[0119] 这样，可以得到Σ的似然估计值
[0121] 上式为Σ的实际似然估计值。但是，该式仍然包含有hM(x)，而hM(x)中依赖于 ρτ(τ)。也就是说，式（17)中其实有2个未知数，Σ的估计需要事先知道纹理，这显然是 不可能的。
[0122] 为了解决上面这个问题，需要采用逐个最大化的估计方法，分两步来求解Σ的似 然估计值。首先，假设图像纹理是已知的，即认为式（14)中仅有Σ -个参数，这时求解Σ 变成了解算下面的方程
[0124] 在这个假设下，Σ的估计值为
[0126] 然后，再假设斑噪协方差矩阵Σ式已知的，即认为式（14)中仅有τ 一个参数，这 时求解τ等价于解算下面的方程
[0128] 在这个假设下，τ的估计值为
[0130] 这样，联立式（19)和式（21)，就得到了 Σ的似然估计值
[0132] 可以通过下面的迭代处理解