CN1786735A - 一种雷达图像相干斑抑制技术 - Google Patents

Abstract

本发明提供了一种雷达图像相干斑抑制技术，它从雷达图像与相干斑噪声间相互独立这一重要统计特征出发，利用基于高阶统计量的分析方法，将相干斑噪声从雷达图像中分离出来，以实现雷达图像的相干斑抑制。采用本发明的方法，可以使雷达图像的相干斑噪声得到明显的降低、图像质量得到明显改善，在边缘信息、目标特征上具有更高的保持效果，很好地解决了抑制相干斑噪声与图像边缘特征保持这一对矛盾。而且，本发明具有结构简单、收敛速度快的特点。

Description

一种雷达图像相干斑抑制技术

技术领域

本发明属于图像处理领域，它特别涉及雷达图像相干斑抑制技术。

背景技术

当雷达波束照射在一个相对于雷达波长而言比较粗糙的表面时，接收的回波信号是由在一个分辨单元的许多散射体所反射的电磁波所合成的，由于表面粗糙的原因，在分辨单元内的散射体与接收机之间的距离是变化的，所以所接收到的电磁波在频率上是相关的，而在相位上是不相关的，因此回波信号时强时弱。雷达图像正是对连续回波信号的相关处理而形成的，因此在雷达图像中常出现很大的像素点与像素点之间的灰度级波动，这种波动在视觉上呈现出一种类似颗粒的噪声，称为相干斑噪声。

相干斑噪声相当于在图像上叠加了一种非高斯乘性噪声，它直接影响雷达图像的后处理。相干斑抑制是雷达图像后处理多年来一直研究的一个问题，相干斑抑制能力的强弱直接影响后处理的效果，它是雷达图像后处理的首要问题。

到目前为止，学者们提出了许多抑制相干斑噪声的方法，但是没有一种滤波器能够显著地去掉相干斑噪声而保留边缘特征。早期的去斑方法是在SAR成像处理中，采用多视平均处理技术，即平均几幅由同一合成孔径的不同分段形成的、不相干的SAR图像，以得到相干斑减少的SAR图像。这种技术是以牺牲空间分辨率为代价的，且相干斑噪声的抑制也并不理想。从20世纪80年代开始，以空域滤波为基础的斑点抑制技术逐渐成为主流。它是利用图像像素之间的空间相关性来对相干斑进行滤波，以减少相干斑。空域的处理方法常见的有均值滤波器、中值滤波器、J.S.Lee.的局部统计滤波器、最大后验(MAP)滤波器、最小均方误差滤波器(MMSE)、sigma滤波器、形态滤波器等。空域方法一般是利用一个滑动窗口，对窗口的像素进行加权以得到窗口中心点的像素值的滤波方法。该类方法虽可以较好地进行斑点噪声的去除，但效果不够理想，图像细节信息丢失较多。近年来，小波变换方法被用于信号的噪声抑制，它是将图像分成代表不同尺度信息的一系列图像，对代表低频成分的低分辨率图像进行滤噪处理，对代表高频成分的高分辨率图像进行适当的阈值处理以保留主要的边缘信息，然后再重建图像，其关键问题是小波域门限的选取。较为典型的例子是Donoho提出的小波变换软阈值理论和方法，该方法是在加性高斯白噪声条件下推导的，可以较好地保持图像的细节信息，但在实际应用中需要人工干预。总之，从多视处理、基本的滤波法到引入小波变换等，虽然在一定程度上抑制了相干斑，但降噪效果还不够理想，且在边缘特征、目标保持及自适应性等方面还存有不足。

发明内容

为了更有效地抑制相干斑，更好地保留雷达图像的细节成分及边缘特征，本发明方法提供了一种雷达图像相干斑抑制技术，采用本发明的方法能显著地降低雷达图像的相干斑噪声。

为了方便地描述本发明地内容，首先作术语定义：

()^T为转置操作；abs()为绝对值函数；diag()为求矩阵对角元素组成的向量的函数；E{}为均值函数；min()为最小值函数；max()为最大值函数；norm()为范数函数；x[i，j]，x₀[i，j]分别为矩阵X，X₀的第i行第j列的元素；W(k)为第k次循环时解混矩阵；W(k-1)为第k-1次循环时的解混矩阵；var()为求方差函数；β为相干斑指数。

假设s＝(s₁，s₂，…，s_n)^T是n个相互独立的源信号(SAR图像和噪声)组成的矩阵，也即各s_i中没有其它信号的污染，S即是待求的矩阵。

本发明方法提供的一种雷达图像相干斑抑制技术，其特征是采用下面的步骤实现：

步骤1图像数据预准备

对待处理的多幅不同极化方式的SAR图像，图像数据预准备是将每幅二维的图像数据利用行与行之间首尾相连的形式转化成一维向量形式，从而多幅图像数据就得到一个矩阵X₀，再利用公式 $x[i,j]=\frac{x_0[i,j]-\min \left(X_0\right)}{\max \left(X_0\right)-\min \left(X_0\right)}$ 对之进行归一化处理，这样预准备后的输出矩阵定义为矩阵X。

步骤2数据预处理

数据预处理是将经步骤1预准备后的矩阵X进行中心化和白化。

中心化是使矩阵X每一行向量均值为零，首先，求出X各行向量的均值，再将矩阵各行的元素分别减去对应的均值，得到的新矩阵X₁就是经中心化处理后的矩阵。

白化是将中心化后的矩阵X₁作为输入，通过线性变换矩阵V，输出白化后的矩阵Z＝VX₁。其中，白化矩阵V＝D^-1/2E^T，D为E{X₁X₁ ^T}的特征值组成的对角阵，E为E{X₁X₁ ^T}的特征向量组成的正交阵。

在白化的同时还要进行降维操作，将对角阵D的元素(特征值)从大到小排列成新的对角阵，选取左上角维数为原未知独立源个数(包含噪声)的子对角阵代替原对角矩阵D的值，由于矩阵E的每一个列向量对应矩阵D相应序数处的每一个对角阵元素，按维数减小后的矩阵D中保留的对角线元素及其位置，选取并排列矩阵E中的列向量，用这样新组成的矩阵作为矩阵E的值，从而计算出新的白化矩阵V。这样Z＝VX₁由具有单位方差且互不相关的各分量z_i组成，同时，Z的维数相比X₁得到了进一步降低。

步骤3独立分量分析(Independent Component Analysis，简称ICA)处理

将白化输出矩阵Z作为独立分量分析处理的输入，经过解混矩阵W线性变换后，输出Y＝W^TZ的各分量间相互独立。独立分量分析算法的核心任务就是确定解混矩阵W，其具体实现步骤如下：

(1)产生一个随机单位阵，其维数为矩阵Z的行向量个数，再将该产生的单位阵设为解混矩阵W的初值W(0)；

(2)利用白化输出矩阵Z及第(k-1)次循环时的解混矩阵W(k-1)(k从1开始取值)的值和迭代公式W(k)＝E{Z(W(k-1)^TZ)³}-3W(k-1)就可以计算解混矩阵W(k)；

(3)对步骤(2)得到的矩阵W(k)进行正交化处理，其实现公式为：W(k)＝(W(k)W(k)^T)^-1/2W(k)，这样经正交化后得到各向量间相互正交的矩阵W(k)；

(4)利用公式W(k)＝W(k)/norm(W(k))对正交化后的矩阵W(k)进行归一化处理；

(5)由步骤(4)得到的W(k)与步骤(2)中的W(k-1)判断步骤(4)得到的W(k)是否收敛，判断的原则是：若1-min(abs(diag(W(k)^T*W(k-1))))的值小于收敛条件(收敛条件根据需要人为设定，如10^-6，该值越小，收敛精度越高，但收敛速度变慢)，则认为解混合矩阵W(k)已收敛，转至步骤(6)，否则，返回步骤(2)；

(6)根据步骤(5)求得的已收敛的矩阵W(k)及白化输出矩阵Z和公式Y＝W(k)^TZ，计算独立分量分析处理阶段的输出矩阵Y，矩阵Y即为待求矩阵S的估计。

步骤4图像恢复

矩阵Y的每个行向量就是一幅图像数据，利用公式 $\mathrm{\β}=\sqrt{\mathrm{var}\left(y\right)}/E\left(y\right)$ 计算各行向量的相干斑指数β，其中相干斑指数最小的就是利用ICA降噪后的SAR图像，其余则是噪声。

经过以上步骤后，就可以得到抑制相干斑后的雷达图像。

需要说明的是，对于极化雷达图像，不同的极化方式会包含一些不同的细节特征。极化雷达对同一场景可提供一组极化图像数据，目标特性分散在这些被相干斑污染的图像中，且与相干斑噪声相互独立，抑制相干斑的过程就是从观测图像中分离出期望的雷达图像信号和噪声，这样就初步具备了ICA的模型基础。ICA是基于高阶统计量的一种分析方法，当最多只有一个源服从高斯分布时，能够在接收信号瞬时混叠参数未知的情况下，仅根据源信号的统计独立性，由混合的观测信号恢复出源信号。

在对数量化的基础上，极化雷达图像(所得的图像数据已进行了对数量化处理)可以被看成是目标信号与一些相干斑噪声的线性混合，并且满足最多只有一个高斯源的ICA条件，因此利用ICA的分离性可从混合图像中分离出期望的雷达图像信号。传统的基于二阶统计量方法描述的是图像幅度谱，并不涉及相位信息，而ICA的基石——高阶统计量涵盖了重要的图像结构和相位特征，这就保证了分离出的独立分量具有很好的细节成分和边缘信息。

本发明的核心是利用ICA技术来对雷达图像进行相干斑抑制。采用本发明的方法，能显著地降低雷达图像的相干斑噪声。本发明提出利用ICA技术来分离极化雷达图像中的目标特性和相干斑噪声，它从雷达图像与相干斑噪声间相互独立这一重要统计特征出发，利用基于高阶统计量的分析方法，将相干斑噪声从雷达图像中分离出来，以达到从根本上降低雷达图像相干斑噪声和提高图像质量的作用。其有益效果是，可以使雷达图像的相干斑噪声得到明显的降低、图像质量得到明显改善，在边缘信息、目标特征上具有更高的保持效果，很好地解决了抑制相干斑噪声与图像边缘特征保持这一对矛盾。而且，本算法结构简单、收敛速度快。

附图说明

图1为ICA结构图

图中：X为观测信号矢量，V为白化矩阵，W为解混合矩阵，Y为输出矢量。

图2为本发明实现的流程图

图3为实施例中相干斑性指标比较表

具体实施方式：

下面结合附图和实施例对本发明进一步说明。

实验数据：采用航天飞机SIR-C雷达L波段香港地区的三幅极化图像数据，大约5.23视，是1994年10月4日所成的像。这里截取了一小块区域(像素250*250)，其中既有岛屿又有海洋。实验中，在对图像白化时进行了降维处理，即仅选用最大的两个特征值及其对应的特征向量，白化后输出矩阵Z的大小也因此变为2*62500。再将白化后的数据通过ICA处理，当设收敛条件为0.00001时，经过7次迭代后算法收敛，得到解混矩阵 $W=\left(\begin{array}{cc}-0.1242& -0.9923\\ 0.9923& -0.1242\end{array}\right).$ 由Y＝W^TZ可计算出经ICA处理后的输出矩阵Y，Y的两个行向量分别是独立分量1和独立分量2，它们的性能指标如图3所示。从图3可以看出，独立分量1具有最低的相干斑指数，也就是无相干斑噪声污染的雷达图像信号的估计，实验表明该算法在保证图像边缘、细节信息的同时，显著地降低了极化SAR图像的相干斑噪声，提高了图像质量。

Claims (1)

1、一种雷达图像相干斑抑制技术，其特征是采用下面的步骤实现：

步骤1  图像数据预准备

对待处理的多幅不同极化方式的SAR图像，图像数据预准备是将每幅二维的图像数据利用行与行之间首尾相连的形式转化成一维向量形式，从而多幅图像数据就得到一个矩阵X₀，再利用公式 $x[i,j]=\frac{x_0[i,j]-\min \left(X_0\right)}{\max \left(X_0\right)-\min \left(X_0\right)}$ 对之进行归一化处理，这样预准备后的输出矩阵定义为矩阵X；

步骤2  数据预处理

数据预处理是将经步骤1预准备后的矩阵X进行中心化和白化；

中心化是使矩阵X每一行向量均值为零，首先，求出X各行向量的均值，再将矩阵各行的元素分别减去对应的均值，得到的新矩阵X₁就是经中心化处理后的矩阵；

白化是将中心化后的矩阵X₁作为输入，通过线性变换矩阵V，输出白化后的矩阵Z＝VX₁；其中，白化矩阵V＝D^-1/2E^T，D为E{X₁X₁ ^T}的特征值组成的对角阵，E为E{X₁X₁ ^T}的特征向量组成的正交阵；

在白化的同时还要进行降维操作，将对角阵D的元素(特征值)从大到小排列成新的对角阵，选取左上角维数为原未知独立源个数(包含噪声)的子对角阵代替原对角矩阵D的值，由于矩阵E的每一个列向量对应矩阵D相应序数处的每一个对角阵元素，按维数减小后的矩阵D中保留的对角线元素及其位置，选取并排列矩阵E中的列向量，用这样新组成的矩阵作为矩阵E的值，从而计算出新的白化矩阵V；这样Z＝VX₁由具有单位方差且互不相关的各分量z_i组成，同时，Z的维数相比X₁得到了进一步降低；

步骤3  独立分量分析处理

将白化输出矩阵Z作为独立分量分析处理的输入，经过解混矩阵W线性变换后，输出Y＝W^TZ的各分量间相互独立；独立分量分析算法的核心任务就是确定解混矩阵W，其具体实现步骤如下：

(1)产生一个随机单位阵，其维数为矩阵Z的行向量个数，再将该产生的单位阵设为解混矩阵W的初值W(0)；

(2)利用白化输出矩阵Z及第(k-1)次循环时的解混矩阵W(k-1)(k从1开始取值)的值和迭代公式W(k)＝E{Z(W(k-1)^TZ)³}-3W(k-1)就可计算解混矩阵W(k)；

(3)对步骤(2)得到的矩阵W(k)进行正交化处理，其实现公式为：

W(k)＝(W(k)W(k)^T)^-1/2W(k)，这样经正交化后得到各向量间相互正交的矩阵W(k)；

(4)利用公式W(k)＝W(k)/norm(W(k))对正交化后的矩阵W(k)进行归一化处理；

(5)由步骤(4)得到的W(k)与步骤(2)中的W(k-1)判断步骤(4)得到的W(k)是否收敛，判断的原则是：若1-min(abs(diag(W(k)^T*W(k-1))))的值小于收敛条件(收敛条件根据需要人为设定)，则认为解混合矩阵W(k)已收敛，转至步骤(6)，否则，返回步骤(2)；

(6)根据步骤(5)求得的已收敛的矩阵W(k)及白化输出矩阵Z和公式

Y＝W(k)^TZ，

计算独立分量分析处理阶段的输出矩阵Y，矩阵Y即为待求矩阵S的估计；

步骤4  图像恢复

矩阵Y的每个行向量就是一幅图像数据，利用公式 $\mathrm{\β}=\sqrt{\mathrm{var}\left(y\right)}/E\left(y\right)$ 计算各行向量的相干斑指数β，其中相干斑指数最小的就是利用ICA降噪后的SAR图像，其余则是噪声；

经过以上步骤后，就可以得到经过抑制相干斑后的雷达图像。

Images