[发明专利]基于treelet融合和水平集分割的遥感图像变化检测

权利要求书

1.一种基于treelet融合和水平集分割的遥感图像变化检测方法，包括如下步骤：

（1）对输入的两幅已配准的大小均为mxn的多时相遥感图像I¹和I²分别进行均 值漂移滤波，得到滤波后图像X¹和X²；

（2）对滤波后图像X¹和X²分别进行3次二维平稳小波分解，且每一幅滤波后 图像3次分解的最高分解层数分别为s=1，...3，一次二维平稳小波分解得到四个小波 系数矩阵，即一个低频系数矩阵和三个分别表示水平、垂直、对角方向的高频小波系 数矩阵；

（3）对滤波后图像X¹和X²相同分解层数下对应方向子带的小波系数矩阵做差， 得到每一个分解层s的低频小波系数差矩阵和三个分别表示水平、垂直、对角方向的 高频小波系数差矩阵；

（4）对步骤（3）中的水平方向小波系数差矩阵和垂直方向小波系数差矩阵利用 sobel算子进行增强，保持低频小波系数差矩阵和对角方向小波系数差矩阵不变；

（5）使用（3）中低频小波系数差矩阵、对角方向小波系数差矩阵和（4）中增 强后的水平、垂直方向小波系数差矩阵，对每一个分解层s进行二维逆平稳小波变换， 得到每个分解层s的重构图像RIs；

（6）对每个分解层s的重构图像RIs使用treelet变换进行融合，得到融合后的 差异图D；

（7）对融合后的差异图D进行水平集分割，得到变化检测结果图Z；

其中步骤（6）所述的采用treelet算法对重构后的图像RIs进行融合，按如下步 骤进行：

（6a）将重构后的图像RIs均变换为mxn大小的列向量RI1′，RI2′，RI3′，组 成初始样本X=[RI1’,RI2’,RI3′]；

（6b）定义treelet变换的逐层聚类层数l＝0,1，...L-1，L为初始样本X中列向量的 个数，L=3，在第0层，每个变量采用初始样本X中的列向量表示，初始化基矩阵 B⁰为LxL的单位矩阵及和变量下标集合Ω={1,2，...L}，计算样本X的初始协方差矩阵 C⁰和初始相关系数矩阵M⁰，计算公式如下：

C_ij＝E[(s_p-Es_p)(s_q-Es_q)^T]

Mij=CijCiiCjj]]>

式中C_ij表示初始协方差矩阵C⁰第i行第j列的计算值，s_p，s_q表示样本X中两 个不同的列向量，E表示求数学期望，T表示转置操作，M_ij表示初始相关系数矩阵 M⁰第i行第j列的计算值；

（6c）当逐层聚类层数l≠0时，寻找相关系数矩阵M^l中最大的两个值，将最大 值和次大值的对应位置序号分别记为α和β：

(α,β)=argmaxi,j∈ΩMijl-1]]>

这里i＜j，分别表示相关系数矩阵M^l中任意值的行和列,且只在和变量下标集合Ω内 进行；

（6d）计算Jacobi旋转矩阵J：

其中c_n=cos(θ₁)，s_n=sin(θ₁)；旋转角θ_l由C^l＝J^TC^l-1J及计算得到， J^T表示旋转矩阵J的转置：

θl=tan-1[Cααl-Cββl±(Cαβl-Cββl)2+4(Cαβl)22Cαβl]]]>且|θl|≤π4;]]>

式中，表示协方差矩阵C^l中行、列位置序号均为α的元素值，表示协方 差矩阵C^l中行、列位置序号均为β的元素值，表示协方差矩阵C^l中行位置序号为 α、列位置序号均为β的元素值；

（6e）由矩阵J更新去相关后的正交基B^l＝B^l-1J，协方差矩阵C^l＝J^TC^l-1J，其 中B^l-1为更新前的正交基，B^l为更新后的正交基，C^l-1为更新前的协方差矩阵，C^l为 更新后的协方差矩阵，并利用C^l更新相关系数矩阵M^l-1为M^l，M^l-1为更新前的相关 系数矩阵，M^l为更新后的相关系数矩阵：

Mijl=CijlCiilCjjl]]>

其中为更新后相关系数矩阵M^l中第i行j列的更新值，为更新后协方差矩 阵C^l中第i行j列的值，为更新后协方差矩阵C^l中第i行i列的值，为更新后协 方差矩阵C^l中第j行j列的值；

（6f）定义更新后正交基矩阵B^l中位置序号为α和位置序号为β的两个列向量分 别为尺度函数φ_l和细节函数ψ_l，定义当前l层的尺度向量集合{φ_l}是尺度函数φ_l和上一 层的尺度向量集合{φ_l-1}的合集，将差变量的位置序号β从和变量的下标集合Ω中去 除，即Ω=Ω/{β}；

（6g）重复步骤（6c）至（6f）直至聚类的最高层l＝L-1，可得最终的正交基矩 阵B=[φ_L-1,ψ₁，...ψ_L-1]，其中φ_L-1∈{ψ_L-1}，ψ₁为第一层treelet变换得到的细节函数， ψ_L-1为最高层treelet变换得到的细节函数；

（6h）将初始样本X沿正交基矩阵B转置的方向进行投影，即P=XxB^T，并将投 影向量P变回mxn大小的图像，得到融合后的差异图D。