CN103873743A - 一种基于结构张量和卡尔曼滤波的视频去噪方法 - Google Patents

Abstract

本发明提供了一种基于结构张量和卡尔曼滤波的视频实时去噪方法，包括以下步骤：获取当前时刻待处理图像帧，以及保存的已完成去噪处理的当前帧之前的n帧图像；对当前时刻待处理图像帧采用均值滤波器进行预滤波处理；基于结构张量，充分利用当前时刻待处理图像帧与前面相邻图像帧之间紧密的时空联系，对当前图像帧进行运动估计；基于运动估计结果，采用卡尔曼滤波方法在时域进行降噪处理；采用维纳滤波在空域进行降噪处理；综合两个去噪图像，加权获得最终的去噪图像。本发明通过上述方法能够实现大噪声视频的去噪处理，具有较好的去噪效果，而且由于没有复杂的迭代计算，易于FPGA等硬件实现，进而能够实现大噪声视频的实时去噪。

Description

一种基于结构张量和卡尔曼滤波的视频去噪方法

技术领域：

本发明属于视频处理领域，主要涉及视频去噪，特指一种基于结构张量和卡尔曼滤波的视频实时去噪方法，可用于自然大噪声视频的实时去噪。

背景技术：

随着数字光电成像技术的快速发展，数字光电成像设备已广泛应用于计算摄影、安防监控、机器人导航以及军事侦察等领域。通常，数字光电成像设备的传感器都是由CCD或者CMOS构成，在成像过程中，受光学噪声、电阻和电容等元器件噪声、传感器噪声、电路噪声等影响，输出的图像会不可避免的包含有许多噪声，这些噪声不仅破坏了图像的真实信息，还严重影响了图像的视觉效果。随着实际应用要求的提高，数字光电成像设备也大量应用于低照度环境，典型的应用如安防监控领域，在夜间等低照度环境下进行监控。然而，受低照度环境的影响，获取的视频受到严重的噪声污染，视频图像如图1所示，这极大地影响了实际的应用。因此，对于这些包含噪声的视频图像，需要进行降噪处理，还原得到清晰的视频图像，以便进行更高层次的处理操作。

目前，视频的去噪方法主要按空域、时域以及变换域来进行划分。空域去噪方法只针对视频的各帧图像进行去噪处理，较好的空域滤波方法有双边滤波方法[1]、非局部均值滤波方法[2]、稀疏表示下的去噪方法[3]等，对各帧图像均能得到较好的去噪效果。其中，双边滤波是一种非线性的滤波方法，是结合图像的空间邻近度和像素值相似度的一种折中处理，同时考虑空域信息和灰度相似性，达到保边去噪的目的。双边滤波方法具有简单、迭代、局部的特点，能很好地保存图像边缘，但是由于保存了过多的高频信息，对于彩色图像里的高频噪声，双边滤波器不能够干净的滤掉，只能够对于低频信息进行较好的滤波。非局部均值滤波方法是一种利用自然图像自身存在大量的重复冗余信息特点，用非局部自相似性来抑制图像噪声的方法。计算图像块之间相似性的匹配过程是非局部均值方法的关键技术。虽然非局部均值滤波有优异的去噪性能，但是过高的计算复杂度极大的限制了它的发展和应用。由于信号稀疏表示的优良特性，基于稀疏分解的图像去噪也越来越受到人们的重视。由于稀疏分解可以保留最匹配原始图像的有用信息而去除与信息无关的噪声，因而可以用此方法来去除图像中的随机噪声。虽然稀疏分解的理论研究很成功，但在实际应用中很难推广，影响其发展的最重要因素为稀疏分解的计算量十分巨大。

上述空域滤波方法除了自身存在的缺点外，共有的缺点是没有充分利用时域信息，所以不能得到理想的滤波效果。但是传统的时域滤波虽然考虑了视频帧间的相关性，但只适用于静止的环境，对于运动物体会产生伪影、拖影等现象，新的时域滤波在原有算法的基础上加入了运动估计，基于视频各帧图像的运动相关性，可以有效解决这些现象。

在变换域去噪方法中，较好的滤波方法有小波去噪方法[4]、三维块匹配方法BM3D[5]等。从信号学的角度看，小波去噪是一个信号滤波的问题，而且尽管在很大程度上小波去噪可以看成是低通滤波，但是由于在去噪后还能成功地保留信号特征，所以又优于传统的低通滤波器。小波去噪实际上是特征提取和低通滤波功能的综合。小波去噪方法包括三个基本的步骤：对含噪信号进行小波变换；对变换得到的小波系数进行某种处理，以去除其中包含的噪声；对处理后的小波系数进行小波逆变换，得到去噪后的信号。BM3D算法是基于块匹配的三维变换域滤波方法。将图像分成若干块，然后对每个参考块进行搜索，在整幅图像中搜索其相似块。并根据块的相似程度将匹配块进行组合形成一个三维矩阵。该矩阵的数据具有很高的相关性，通过三维酉变换可以有效地降低其相关性。同时，通过对变换域系数的滤波，可以大幅降低噪声。

此外，根据视频的时空特性，一些时空结合的滤波方法被提了出来，如时空双边滤波方法、视频三维块匹配方法VBM3D等，在一定程度上取得了较好的效果。

然而，受算法复杂度的影响，上述的这些滤波方法大部分都无法进行实时应用，如非局部均值滤波方法、稀疏表示下的去噪方法、BM3D以及VBM3D等。其他能满足实时应用的滤波方法，如时空双边滤波等，受去噪效果的限制，无法满足大噪声视频的去噪要求。因此，针对安防监控等领域的应用需求，对低照度环境下的大噪声进行实时的降噪处理，已成为急需解决的问题。

发明内容：

针对上述现有技术存在的不能既实时又有效地对大噪声视频进行降噪处理的问题，本发明提出一种基于结构张量和卡尔曼滤波的视频实时去噪方法，能够对视频进行实时地降噪处理，并且具有较好的去噪效果。

为实现上述发明目的，本发明采用的技术方案是：

一种基于结构张量和卡尔曼滤波的视频实时去噪方法，包括如下步骤：

步骤一：获取当前时刻待处理图像帧，以及保存的已完成去噪处理的当前帧之前的n帧图像，n为预先设置的整数值；

步骤二：采用均值滤波器对当前时刻待处理图像帧进行预滤波处理；

步骤三：针对步骤二预滤波处理后的图像和保存的已完成去噪处理的当前帧之前的n帧图像，计算每帧图像的结构张量，然后分别计算当前帧图像与之前n帧图像中每一帧图像之间的结构张量距离，再将这n个结构张量的距离求平均值，这个平均值可作为当前时刻待处理图像帧的运动估计值；

步骤四：根据步骤三得到的当前时刻待处理图像帧的运动估计值建立卡尔曼滤波的运动方差矩阵，在时间域上采用卡尔曼滤波方法对当前时刻待处理图像帧进行降噪处理，得到卡尔曼滤波后的去噪图像；

步骤五：在空间域上采用维纳滤波器对当前时刻待处理图像帧进行降噪处理，得到维纳滤波后的去噪图像；

步骤六：根据步骤三得到的当前时刻待处理图像帧中各像素点的运动估计值计算其高斯权重值，结合步骤四得到的卡尔曼滤波后的去噪图像以及步骤五得到的维纳滤波后的去噪图像，通过加权获得最终的去噪图像，并保存为当前帧已完成去噪处理的图像。

进一步地，所述步骤一中保存的已完成去噪处理的当前帧之前的n帧图像的初始值，即对于视频图像的第1帧至第n帧，将每一帧图像的原始含噪图像保存为其对应的已完成去噪处理的图像，其对应的已完成去噪处理的图像采用如下方式确定：

（1）对于视频图像的第1帧，将该第1帧图像的原始含噪图像保存为其对应的已完成去噪处理的图像；

（2）对于视频图像的第k帧，k为整数且1<k≤n，将第1帧至第k-1帧已完成去噪处理的图像作为所述步骤三中保存的已完成去噪处理的当前帧之前的n帧图像，依次执行步骤二至步骤六，得到第k帧对应的已完成去噪处理的图像。

进一步地，所述步骤二中采用均值滤波器进行滤波处理的具体方法为：对于当前时刻待处理图像帧中的每一个像素点(x,y)，选择以该像素点为中心的大小为N×N的模板，求模板中所有像素的均值，再把该均值赋予当前像素点(x,y)，作为处理后图像在该点上的像素值V₁(x,y)，即：

$V_1(x,y)=\frac{\underset{i-0}{\overset{N-1}{\mathrm{\&Sigma;}}}\underset{j=0}{\overset{N-1}{\mathrm{\&Sigma;}}}{V}_0(i,j)}{N^2}---\left(1\right)$

其中V₀(x,y)为待处理图像在该模板中处的像素值。

进一步地，所述步骤三中的每帧图像的结构张量采用的是黎曼空间（Riemannian space），它包括了欧式空间、球面、双曲空间为其特例。如果黎曼流行M上任何点处的任何二维切平面，其相应的截面曲率均为常数K，则称此黎曼流行为常曲率黎曼空间，由称常曲率空间。黎曼空间中线性结构张量的表达式为

$J_{{\mathrm{\&rho;}}_2}\left(p\right)=K_{\mathrm{\&rho;}2}*(\&dtri;p_{\mathrm{\&sigma;}\&prime;}\&CircleTimes;\&dtri;p_{\mathrm{\&sigma;}\&prime;}^T)=K_{\mathrm{\&rho;}2}*\left(\begin{array}{cc}{\left(I_x\left(p_{\mathrm{\&sigma;}\&prime;}\right)\right)}^2& I_x\left(p_{\mathrm{\&sigma;}\&prime;}\right)I_y\left(p_{\mathrm{\&sigma;}\&prime;}\right)\\ I_x\left(p_{\mathrm{\&sigma;}\&prime;}\right)I_y\left(p_{\mathrm{\&sigma;}\&prime;}\right)& {\left(I_y\left(p_{\mathrm{\&sigma;}\&prime;}\right)\right)}^2\end{array}\right)---\left(2\right)$

其中

表示图像梯度算子，p_σ'表示图像p经过方差为p_σ'的高斯滤波后的输出图像，

表示结构张量的积，结构张量乘积的计算可表示为图像在x和y两个方向的梯度I_x(p_σ')，I_y(p_σ')来计算，*表示方差为ρ2的高斯滤波器K_ρ2与结构张量的卷积操作，通常地ρ₂>σ'。因为高斯滤波器均匀地合成周围各个方向的结构张量信息，所以称为线性结构张量。

包含图像的结构信息，通过对

进行正交分解，可以得到特征值λ₁与λ₂，以及特征向量

与

特征值表示结构张量的大小，特征向量表示结构张量的方向。不同的特征值代表不同的图像结构，通常用来表示反应图像结构的大小，如图3所示，（b1）和（b2）表示从含噪图像帧（a1）和（a2）中所提取的图像结构信息。

当有运动发生时，结构张量会发生变化，计算不同图像帧之间的结构张量距离可以进行运动估计，结构张量的距离可以用下式表示

$d_{\mathrm{ST}}(p_{\mathrm{current}},p_{\mathrm{past},i})=\sqrt{\mathrm{Trace}\left({(\log \left(J_{{\mathrm{\&rho;}}_2}\left(p_{\mathrm{current}}\right)\right)-\log \left(J_{{\mathrm{\&rho;}}_2}\left(p_{\mathrm{past},i}\right)\right))}^2\right)}---\left(3\right)$

其中，Trace(·)表示矩阵的轨迹，log(·)是结构张量的对数操作因子，

代表当前含噪帧图像的结构张量，

代表当前图像帧之前的图像帧中第i-th帧图像的结构张量，分别计算出当前处理图像帧与之前n帧图像间的结构张量距离，对其求平均值

$d_{\mathrm{ST}}=\frac{\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{d}_{\mathrm{ST}}(p_{\mathrm{current}},p_{\mathrm{past},i})}{n}---\left(4\right)$

这个平均值代表当前时刻待处理图像帧的运动估计值。

进一步地，步骤四中采用的卡尔曼滤波方法包括以下五个步骤：

（1）假设当前待处理图像帧处于k时刻，基于k-1时刻采用卡尔曼滤波后的图像，初步估计当前待处理图像的去噪图像：

X(k|k-1)=AX(k-1|k-1)    （5）

其中，X(k|k-1)为对当前待处理图像进行去噪估计后的图像，X(k-1|k-1)为k-1时刻采用卡尔曼滤波后的图像，A为状态转移参数矩阵；

（2）基于待处理图像帧中各像素点的运动估计值，计算误差协方差矩阵：

P(k|k-1)=AP(k-1|k-1)A^T+Q(k)     （6）

其中，P(k|k-1)为进行降噪估计后的图像的协方差矩阵，P(k-1|k-1)为k-1时刻采用卡尔曼滤波后的图像的协方差矩阵，Q(k)为当前待处理图像中运动方差矩阵；

视频中的运动会产生噪声，所以对当前待处理图像帧中任意像素点(x,y)，定义

Q_k-1(x,y)=d_ST(x,y)

根据上述设定，使得运动区域的方差值相对静止区域而言更大；

（3）计算卡尔曼增益：

Kg(k)=P(k|k-1)/P(k|k-1)+R(k)      （7）

其中，Kg(k)即为当前k时刻的卡尔曼增益矩阵，R(k)为当前待处理图像的噪声协方差矩阵，视频中的噪声协方差为一恒定不变的值；

（4）结合当前待处理图像及其初步估计的降噪图像，计算当前k时刻的最优化降噪估计图像：

X(k|k)=X(k|k-1)+Kg(k)Z(k)-X(k|k-1)      （8）

其中，X(k|k)即为对当前待处理图像的最优化降噪估计图像，Z(k)为当前k时刻包含噪声的待处理图像；

（5）更新当前k时刻的协方差矩阵，为下一帧图像进行卡尔曼滤波做准备：

P(k|k)=(I-Kg(k))P(k|k-1)     （9）

进一步地，所述步骤六的具体方法为：

S61.根据步骤三得到的当前时刻待处理图像帧中各像素点的运动估计值计算其高斯权重值：

$w(i,j)=e^{\frac{{\mathrm{MEV}}_m^2}{{\mathrm{\&sigma;}}^2}}---\left(10\right)$

其中w(i,j)表示像素点(i,j)的高斯权重值，(i,j)表示该像素点的坐标位置，MEV_m表示该像素点所在的第m个块图像的运动估计值；σ为可调节参数；

S62.结合步骤四得到的卡尔曼滤波后的去噪图像以及步骤五得到的维纳滤波后的去噪图像，通过加权获得最终的去噪图像：

V(i,j)=w(i,j)V_k(i,j)+[1-w(i,j)]V_b(i,j)     （11）

其中，V_k(i,j)为步骤四中经过卡尔曼滤波得到的去噪图像在(i,j)处的像素值，V_b(i,j)为步骤五中经过维纳滤波得到的去噪图像在(i,j)处的像素值，V(i,j)即为加权求得的最终的去噪图像在处的像素值；

S63.将加权求得的最终的去噪图像保存为当前帧已完成去噪处理的图像。

进一步地，所述步骤一中n的值为3-6。

本发明通过上述步骤，结合了卡尔曼滤波和维纳滤波对大视频进行降噪，并在卡尔曼滤波中引入结构张量得出运动估计值，最终能够实现大噪声视频的降噪处理，去噪后的图像兼具卡尔曼滤波和维纳滤波的去噪效果，即在静止区域具有卡尔曼滤波的去噪效果，在运动区域具有维纳滤波的去噪效果；同时，采用结构张量进行运动估计，能更好地体现相邻图像帧之间的时空联系；由于整个方法没有任何复杂的迭代计算，易于FPGA等硬件实现，因此能够避免伪影、拖影等现象的产生，对于视频监控领域以及其他相关应用领域具有非常重要的意义。

参考文献：

[1]余博,郭雷,赵天云,钱晓亮.红外图像的自适应混合双边滤波算法.红外与激光工程,Vol.5,pp.517-523,2012.

[2]张丽果.快速非局部均值滤波图像去噪.信号处理,Vol.8,pp.1043-1049,2013.

[3]张晓阳,柴毅,李华锋.基于K-SVD和残差比的低信噪比图像稀疏表示去噪算法.光学技术,Vol.1,pp.23-29,2012.

[4]王香云.基于FPGA的小波去噪的硬件实现.激光技术,Vol.6,pp.786-790,2013.

[5]李政,刘文江,戎蒙恬,刘太智.BM3D视频去噪算法实现与评估.信息技术,Vol.4,pp.30-32,2012.

附图说明：

图1为低照度环境下拍摄的大噪声视频图像；

图2为本发明所述方法的流程图；

图3为本发明所述预滤波方法原理图；

图4为含噪图像及其提取的结构张量图；

图5为图4所示的两个结构张量的距离图；

图6为本发明所述高斯权重曲线图；

图7是本发明所述大噪声视频图像经过去噪后的图像。

具体实施方式：

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

如图2所示，本发明所述一种基于结构张量和卡尔曼滤波的视频实时去噪方法，包括如下步骤：

步骤一：获取当前时刻待处理图像帧，以及保存的已完成去噪处理的当前帧之前的n帧图像，n的值优选为3-6，本实施例中取值为4。

本发明中对于所述保存的已完成去噪处理的当前帧之前的n帧图像的初始值，即对于视频图像的第1帧至第n帧，将每一帧图像的原始含噪图像保存为其对应的已完成去噪处理的图像，其对应的已完成去噪处理的图像采用如下方式确定：

（1）对于视频图像的第1帧，将该第1帧图像的原始含噪图像保存为其对应的已完成去噪处理的图像；

（2）对于视频图像的第k帧，k为整数且1<k≤n，将第1帧至第k-1帧已完成去噪处理的图像作为所述步骤三中保存的已完成去噪处理的当前帧之前的n帧图像，依次执行步骤二至步骤六，得到第k帧对应的已完成去噪处理的图像。

步骤二：采用均值滤波器对当前时刻待处理图像帧进行滤波处理；

由于在大噪声条件下，视频图像受到噪声的严重污染，现有的传统运动检测方法很难提取出视频图像帧中运动物体的完整轮廓，因此，只能对运动进行估计。但是，由于大噪声的影响，可能会出现严重的错检情况，将噪声检测成运动物体。因此，在进行运动估计前，首先需要对当前时刻待处理图像帧进行预滤波处理，消弱噪声的影响，提高运动估计的准确性。本发明中对视频当前时刻待处理图像帧进行预滤波处理，可以消弱噪声的影响，提高运动估计的准确性。

本发明中进行预滤波处理的基本原理是用均值代替图像的各个像素值，如图3所示，对于待处理的当前像素点(x,y)，选择以该像素点为中心的大小为N×N的模板，求模板中所有像素的均值，再把该均值赋予当前像素点(x,y)，作为处理后图像在该点上的像素值V₁(x,y)，即：

$V_1(x,y)=\frac{\underset{i=0}{\overset{N-1}{\mathrm{\&Sigma;}}}\underset{j=0}{\overset{N-1}{\mathrm{\&Sigma;}}}{V}_0(i,j)}{N^2}---\left(1\right)$

其中V₀(x,y)为待处理图像在该模板中处的像素值。

步骤三：针对步骤二预滤波处理后的图像和保存的已完成去噪处理的当前帧之前的n帧图像，计算每帧图像的结构张量，然后分别计算当前帧图像与之前n帧图像中每一帧图像之间的结构张量距离，再将这n个结构张量的距离求平均值，这个平均值可作为当前帧图像的运动估计值。具体过程为：

在黎曼空间，即欧式空间的对数空间中计算结构张量更加简单快捷。黎曼空间中线性结构张量的表达式为

$J_{{\mathrm{\&rho;}}_2}\left(p\right)=K_{\mathrm{\&rho;}2}*(\&dtri;p_{\mathrm{\&sigma;}\&prime;}\&CircleTimes;\&dtri;p_{\mathrm{\&sigma;}\&prime;}^T)=K_{\mathrm{\&rho;}2}*\left(\begin{array}{cc}{\left(I_x\left(p_{\mathrm{\&sigma;}\&prime;}\right)\right)}^2& I_x\left(p_{\mathrm{\&sigma;}\&prime;}\right)I_y\left(p_{\mathrm{\&sigma;}\&prime;}\right)\\ I_x\left(p_{\mathrm{\&sigma;}\&prime;}\right)I_y\left(p_{\mathrm{\&sigma;}\&prime;}\right)& {\left(I_y\left(p_{\mathrm{\&sigma;}\&prime;}\right)\right)}^2\end{array}\right)---\left(2\right)$

其中

表示图像梯度算子，p_σ'表示图像p经过方差为p_σ'的高斯滤波后的输出图像，

表示结构张量的积，结构张量乘积的计算可表示为图像在x和y两个方向的梯度I_x(p_σ')，I_y(p_σ')来计算，*表示方差为ρ2的高斯滤波器K_ρ2与结构张量的卷积操作，通常地ρ₂>σ'。因为高斯滤波器均匀地合成周围各个方向的结构张量信息，所以称为线性结构张量。

包含图像的结构信息，通过对

进行正交分解，可以得到特征值λ₁与λ₂，以及特征向量

与特征值表示结构张量的大小，特征向量表示结构张量的方向。不同的特征值代表不同的图像结构，通常用来表示反应图像结构的大小，如图4所示，（b1）和（b2）表示从含噪图像帧（a1）和（a2）中所提取的图像结构信息。

当有运动发生时，结构张量会发生变化，计算不同图像帧之间的结构张量距离可以进行运动估计，结构张量的距离可以用下式表示

$d_{\mathrm{ST}}(p_{\mathrm{current}},p_{\mathrm{past},i})=\sqrt{\mathrm{Trace}\left({(\log \left(J_{{\mathrm{\&rho;}}_2}\left(p_{\mathrm{current}}\right)\right)-\log \left(J_{{\mathrm{\&rho;}}_2}\left(p_{\mathrm{past},i}\right)\right))}^2\right)}---\left(3\right)$

其中，Trace(·)表示矩阵的轨迹，log(·)是结构张量的对数操作因子，

代表当前含噪帧图像的结构张量，

代表当前图像帧之前的图像帧中第i-th帧图像的结构张量。图5展示了图4中（b1）和（b2）之间的结构张量距离。

分别计算出当前处理图像帧与之前n帧图像间的结构张量距离，对其求平均值

$d_{\mathrm{ST}}=\frac{\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{d}_{\mathrm{ST}}(p_{\mathrm{current}},p_{\mathrm{past},i})}{n}---\left(4\right)$

这个平均值代表当前图像帧中的运动情况，即运动估计值。运动估计值是对目标物体位置的估计，其值越大，运动的可能性就越高。

步骤四：根据步骤四得到的待处理图像帧中各像素点的运动估计值建立卡尔曼滤波的运动方差矩阵，在时间域上采用卡尔曼滤波方法对当前时刻待处理图像帧进行降噪处理，得到卡尔曼滤波后的去噪图像。

卡尔曼滤波是一种递归的估计，即只要获知上一时刻状态的估计值以及当前状态的观测值，即可以计算出当前状态的估计值。将卡尔曼滤波应用到视频降噪处理中，其过程主要有以下五个步骤：

（1）假设当前待处理图像帧处于k时刻，基于k-1时刻采用卡尔曼滤波后的图像，初步估计当前待处理图像的去噪图像。

X(k|k-1)=AX(k-1|k-1)    （5）

其中，X(k|k-1)为对当前待处理图像进行去噪估计后的图像，X(k-1|k-1)为k-1时刻采用卡尔曼滤波后的图像，A为状态转移参数矩阵，如果已知上一时刻状态的估计值以及当前状态的观测值，则状态转移参数矩阵可以计算出来。

（2）基于待处理图像帧中各像素点的运动估计值，计算误差协方差矩阵。

P(k|k-1)=AP(k-1|k-1)A^T+Q(k)      （6）

其中，P(k|k-1)为进行降噪估计后的图像的协方差矩阵，P(k-1|k-1)为k-1时刻采用卡尔曼滤波后的图像的协方差矩阵，Q(k)为当前待处理图像中运动方差矩阵。

视频中的运动会产生噪声，所以对当前待处理图像帧中任意像素点(x,y)，定义

Q_k-1(x,y)=d_ST(x,y)

根据上述设定，使得运动区域的方差值相对静止区域而言更大。经过卡尔曼滤波器的递归运算，静止区域的方差值会越来越小，而运动区域由于有运动产生，其方差值将维持在一个较大的范围。

（3）计算卡尔曼增益。

Kg(k)=P(k|k-1)/P(k|k-1)+R(k)      （7）

其中，Kg(k)即为当前k时刻的卡尔曼增益矩阵，R(k)为当前待处理图像的噪声协方差矩阵，视频中的噪声协方差为一恒定不变的值。

（4）结合当前待处理图像及其初步估计的降噪图像，计算当前k时刻的最优化降噪估计图像。

X(k|k)=X(k|k-1)+Kg(k)Z(k)-X(k|k-1)     （8）

其中，X(k|k)即为对当前待处理图像的最优化降噪估计图像，Z(k)为当前k时刻包含噪声的待处理图像。

（5）更新当前k时刻的协方差矩阵，为下一帧图像进行卡尔曼滤波做准备。

P(k|k)=(I-Kg(k))P(k|k-1)     （9）

步骤五：在空间域上采用维纳滤波器对当前时刻待处理图像帧进行降噪处理，得到维纳滤波后的去噪图；

在时间域上采用卡尔曼滤波方法对视频图像进行降噪处理时，为了保留运动物体并防止其在去噪过程中产生虚影，去噪后视频图像中的运动区域往往保留了大量的噪声。因此，还需要针对每帧图像进行一次空间上的滤波，消弱噪声的影响。

本发明采用维纳滤波在空间域上对视频图像进行降噪处理。维纳滤波是一种线性的滤波方法，可以看成是一种线性估计问题。维纳滤波一般在小波域进行计算，由于小波变换具有正交性，因此零均值方差为

的高斯白噪声在小波域内仍然是零均值方差为

的高斯白噪声。含噪图像的小波系数表示为：

V₀(x,y)=V₁(x,y)+n(x,y)       （10）

其中，V₀(x,y)为未处理的含噪图像在点处的像素值，V₁(x,y)为处理后的图像在点处的像素值，n(x,y)为噪声成分。考虑到噪声和信号不相关信号小波方差值

可以由下式得到：

${\widehat{\mathrm{\&sigma;}}}_{(x,y)}^2=E{\left\{V_1(x,y)\right\}}^2=\max \left\{E\right\{V_0{(x,y)}^2\}-{\mathrm{\&sigma;}}_n^2,0\}---\left(11\right)$

其中，期望值E{V₀(x,y)²}可由局部窗内小波系数求平均值得到：

$E\left\{V_0{(x,y)}^2\right\}=\frac{1}{\mathrm{MN}}\underset{i,j\&Element;W}{\mathrm{\&Sigma;}}{V}_0{(x,y)}^2---\left(12\right)$

其中，W是图像中每个像素的M×N的邻域。M和N分别表示窗的行数和列数。

最后可得到信号的维纳滤波估计为：

${\hat{V}}_{1(x,y)}=\frac{{\widehat{\mathrm{\&sigma;}}}_{x,y}^2}{{\widehat{\mathrm{\&sigma;}}}_{x,y}^2+{\mathrm{\&sigma;}}_n^2}{V}_0(x,y)---\left(13\right)$

步骤六：根据步骤三得到的待处理图像帧中各像素点的运动估计值计算其高斯权重值，结合步骤四得到的卡尔曼滤波后的去噪图像以及步骤五得到的维纳滤波后的去噪图像，通过加权获得最终的去噪图像，并保存为当前帧已完成去噪处理的图像。

S61.根据步骤三得到的待处理图像帧中各像素点的运动估计值计算其高斯权重值：

$w(i,j)=e^{\frac{{\mathrm{MEV}}_m^2}{{\mathrm{\&sigma;}}^2}}---\left(14\right)$

其中w(i,j)表示像素点(i,j)的高斯权重值，(i,j)表示该像素点的坐标位置，MEV_m表示该像素点所在的第m个块图像的运动估计值；σ为可调节参数；

如图6所示，静止区域的运动估计值小，其高斯权重值就偏大；而运动区域的运动估计值大，其高斯权重值就偏小。

S62.结合步骤四得到的卡尔曼滤波后的去噪图像以及步骤五得到的维纳滤波后的去噪图像，通过加权获得最终的去噪图像：

V(i,j)=w(i,j)V_k(i,j)+[1-w(i,j)]V_b(i,j)      （15）

其中，V_k(i,j)为步骤四中经过卡尔曼滤波得到的去噪图像在(i,j)处的像素值，V_b(i,j)为步骤五中经过维纳滤波得到的去噪图像在(i,j)处的像素值，V(i,j)即为加权求得的最终的去噪图像在处的像素值；

S63.将加权求得的最终的去噪图像保存为当前帧已完成去噪处理的图像。

本实施例通过上述方法，最终能够实现对包含大噪声的视频图像的降噪处理，采用本发明对图1所示的低照度环境下拍摄的大噪声视频图像进行去噪处理，最终得到的去噪图像如图7所示，可以看出，本发明具有良好的去噪效果。同时，由于整个方法没有任何复杂的迭代计算，易于FPGA等硬件实现，能够满足视频图像的实时性要求，因此能够实现大噪声视频的实时去噪。

Claims (7)

1.一种基于结构张量和卡尔曼滤波的视频实时去噪方法，其特征在于，包括如下步骤：

步骤一：获取当前时刻待处理图像帧，以及保存的已完成去噪处理的当前帧之前的n帧图像，n为预先设置的整数值；

步骤二：采用均值滤波器对当前时刻待处理图像帧进行预滤波处理；

步骤三：针对步骤二预滤波处理后的图像和步骤一保存的已完成去噪处理的当前帧之前的n帧图像，计算每帧图像的结构张量，然后分别计算当前帧图像与之前n帧图像中每一帧图像之间的结构张量距离，再将这n个结构张量距离求平均值，该平均值作为当前时刻待处理图像帧的运动估计值；

步骤四：根据步骤三得到的当前时刻待处理图像帧的运动估计值建立卡尔曼滤波的运动方差矩阵，在时间域上采用卡尔曼滤波方法对当前时刻待处理图像帧进行降噪处理，得到卡尔曼滤波后的去噪图像；

步骤五：在空间域上采用维纳滤波器对当前时刻待处理图像帧进行降噪处理，得到维纳滤波后的去噪图像；

步骤六：根据步骤三得到的当前时刻待处理图像帧中各像素点的运动估计值计算其高斯权重值，结合步骤四得到的卡尔曼滤波后的去噪图像以及步骤五得到的维纳滤波后的去噪图像，通过加权获得最终的去噪图像，并保存为当前帧已完成去噪处理的图像。

2.根据权利要求1所述基于结构张量和卡尔曼滤波的视频实时去噪方法，其特征在于，所述步骤一中保存的已完成去噪处理的当前帧之前的n帧图像的初始值，即对于视频图像的第1帧至第n帧，将每一帧图像的原始含噪图像保存为其对应的已完成去噪处理的图像，其对应的已完成去噪处理的图像采用如下方式确定：

（1）对于视频图像的第1帧，将该第1帧图像的原始含噪图像保存为其对应的已完成去噪处理的图像；

（2）对于视频图像的第k帧，k为整数且1<k≤n，将第1帧至第k-1帧已完成去噪处理的图像作为所述步骤三中保存的已完成去噪处理的当前帧之前的n帧图像，依次执行步骤二至步骤六，得到第k帧对应的已完成去噪处理的图像。

3.根据权利要求1所述基于结构张量和卡尔曼滤波的视频实时去噪方法，其特征在于，所述步骤二中采用均值滤波器进行预滤波处理的具体方法为：对于当前时刻待处理图像帧中的每一个像素点(x,y)，选择以该像素点为中心的大小为N×N的模板，求模板中所有像素的均值，再把该均值赋予当前像素点(x,y)，作为处理后图像在该点上的像素值V₁(x,y)，即：

$V_1(x,y)=\frac{\underset{i-0}{\overset{N-1}{\mathrm{\&Sigma;}}}\underset{j=0}{\overset{N-1}{\mathrm{\&Sigma;}}}{V}_0(i,j)}{N^2}---\left(1\right)$

其中V₀(x,y)为待处理图像在该模板中处的像素值。

4.根据权利要求1所述基于结构张量和卡尔曼滤波的视频实时去噪方法，其特征在于，所述步骤三中的每帧图像的结构张量采用的是黎曼空间（Riemannian space），黎曼空间中线性结构张量的表达式为

$J_{{\mathrm{\&rho;}}_2}\left(p\right)=K_{\mathrm{\&rho;}2}*(\&dtri;p_{\mathrm{\&sigma;}\&prime;}\&CircleTimes;\&dtri;p_{\mathrm{\&sigma;}\&prime;}^T)=K_{\mathrm{\&rho;}2}*\left(\begin{array}{cc}{\left(I_x\left(p_{\mathrm{\&sigma;}\&prime;}\right)\right)}^2& I_x\left(p_{\mathrm{\&sigma;}\&prime;}\right)I_y\left(p_{\mathrm{\&sigma;}\&prime;}\right)\\ I_x\left(p_{\mathrm{\&sigma;}\&prime;}\right)I_y\left(p_{\mathrm{\&sigma;}\&prime;}\right)& {\left(I_y\left(p_{\mathrm{\&sigma;}\&prime;}\right)\right)}^2\end{array}\right)---\left(2\right)$

其中

表示图像梯度算子，p_σ'表示图像p经过方差为p_σ'的高斯滤波后的输出图像，

表示结构张量的积，结构张量乘积的计算可表示为图像在x和y两个方向的梯度I_x(p_σ')，I_y(p_σ')来计算，*表示方差为ρ2的高斯滤波器K_ρ2与结构张量的卷积操作，通常地ρ₂>σ'；

上述

包含图像的结构信息，通过对

进行正交分解，得到表示结构张量的大小的特征值λ₁与λ₂，以及表示结构张量的方向的特征向量

与

结构张量的距离用下式表示

$d_{\mathrm{ST}}(p_{\mathrm{current}},p_{\mathrm{past},i})=\sqrt{\mathrm{Trace}\left({(\log \left(J_{{\mathrm{\&rho;}}_2}\left(p_{\mathrm{current}}\right)\right)-\log \left(J_{{\mathrm{\&rho;}}_2}\left(p_{\mathrm{past},i}\right)\right))}^2\right)}---\left(3\right)$

其中，Trace(·)表示矩阵的轨迹，log(·)是结构张量的对数操作因子，

代表当前含噪帧图像的结构张量，

代表当前图像帧之前的图像帧中第i-th帧图像的结构张量，分别计算出当前处理图像帧与之前n帧图像间的结构张量距离，对其求平均值

$d_{\mathrm{ST}}=\frac{\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{d}_{\mathrm{ST}}(p_{\mathrm{current}},p_{\mathrm{past},i})}{n}---\left(4\right)$

这个平均值代表当前时刻待处理图像帧的运动估计值。

5.根据权利要求1所述基于结构张量和卡尔曼滤波的视频实时去噪方法，其特征在于，所述步骤四中采用的卡尔曼滤波方法包括以下五个步骤：

（1）假设当前待处理图像帧处于k时刻，基于k-1时刻采用卡尔曼滤波后的图像，初步估计当前待处理图像的去噪图像：

X(k|k-1)=AX(k-1|k-1)       （5）

其中，X(k|k-1)为对当前待处理图像进行去噪估计后的图像，X(k-1|k-1)为k-1时刻采用卡尔曼滤波后的图像，A为状态转移参数矩阵；

（2）基于待处理图像帧中各像素点的运动估计值，计算误差协方差矩阵：

P(k|k-1)=AP(k-1|k-1)A^T+Q(k)        （6）

其中，P(k|k-1)为进行降噪估计后的图像的协方差矩阵，P(k-1|k-1)为k-1时刻采用卡尔曼滤波后的图像的协方差矩阵，Q(k)为当前待处理图像中运动方差矩阵；

视频中的运动会产生噪声，所以对当前待处理图像帧中任意像素点(x,y)，定义

Q_k-1(x,y)=d_ST(x,y)

根据上述设定，使得运动区域的方差值相对静止区域而言更大；

（3）计算卡尔曼增益：

Kg(k)=P(k|k-1)/P(k|k-1)+R(k)     （7）

其中，Kg(k)即为当前k时刻的卡尔曼增益矩阵，R(k)为当前待处理图像的噪声协方差矩阵，视频中的噪声协方差为一恒定不变的值；

（4）结合当前待处理图像及其初步估计的降噪图像，计算当前k时刻的最优化降噪估计图像：

X(k|k)=X(k|k-1)+Kg(k)Z(k)-X(k|k-1)     （8）

其中，X(k|k)即为对当前待处理图像的最优化降噪估计图像，Z(k)为当前k时刻包含噪声的待处理图像；

（5）更新当前k时刻的协方差矩阵，为下一帧图像进行卡尔曼滤波做准备：

P(k|k)=(I-Kg(k))P(k|k-1)     （9）

6.根据权利要求1所述基于结构张量和卡尔曼滤波的视频实时去噪方法，其特征在于，所述步骤六的具体实现步骤为：

S61.根据步骤三得到的当前时刻待处理图像帧中各像素点的运动估计值计算其高斯权重值：

$w(i,j)=e^{\frac{{\mathrm{MEV}}_m^2}{{\mathrm{\&sigma;}}^2}}---\left(10\right)$

其中w(i,j)表示像素点(i,j)的高斯权重值，(i,j)表示该像素点的坐标位置，MEV_m表示该像素点所在的第m个块图像的运动估计值；σ为可调节参数；

S62.结合步骤四得到的卡尔曼滤波后的去噪图像以及步骤五得到的维纳滤波后的去噪图像，通过加权获得最终的去噪图像：

V(i,j)=w(i,j)V_k(i,j)+[1-w(i,j)]V_b(i,j)       （11）

其中，V_k(i,j)为步骤四中经过卡尔曼滤波得到的去噪图像在(i,j)处的像素值，V_b(i,j)为步骤五中经过维纳滤波得到的去噪图像在(i,j)处的像素值，V(i,j)即为加权求得的最终的去噪图像在处的像素值；

S63.将加权求得的最终的去噪图像保存为当前帧已完成去噪处理的图像。

7.根据权利要求1-6之一所述基于结构张量和卡尔曼滤波的视频实时去噪方法，其特征在于，所述步骤一中n的值为3-6。

Images