CN103595981B - 基于非局部低秩的色彩滤波阵列图像去马赛克方法 - Google Patents

Abstract

本发明公开了一种基于非局部低秩的色彩滤波阵列图像去马赛克方法，主要解决现有技术易出现边缘模糊、拉链效应和块效应的问题。其实现步骤为：1.输入一幅色彩滤波阵列图像；2.对绿色通道作方向插值；3.取绿色插值图的待修正块；4.取待修正块的图像块集合并从集合中找出相似块，组成相似块矩阵；5.定义相似块矩阵的拉格朗日优化函数并求解，得到修正的相似块矩阵；6.将修正的相似块矩阵的中间列向量重组，得到修正的块；7.对红色和蓝色通道作方向插值；8.取红色和蓝色插值图的待修正块，重复执行步骤4-6，修正红色和蓝色插值图的所有块；9.输出全彩色图像。本发明能避免边缘模糊，抑制拉链效应和块效应，用于对色彩滤波阵列图像的恢复。

Description

基于非局部低秩的色彩滤波阵列图像去马赛克方法

技术领域

本发明属于图像处理技术领域，具体的说是一种基于非局部低秩的色彩滤波阵列图像去马赛克方法，可用于恢复单传感器芯片相机中色彩滤波阵列图像的完整彩色信息，从而弥补由于减少相机硬件成本所带来的图像彩色信息的丢失。

背景技术

随着数码相机的广泛使用，其成本和体积是一个不可忽视的问题，为此现有的大部分数码相机都采用单块电荷耦合器件CCD或互补金属氧化物半导体CMOS作为图像的采样传感器，其上面覆盖有一层按Bayer方式排列的色彩滤波阵列，这种阵列使同一个像素点只能采样三种基色，即红R，绿G，蓝B中的一种，而另外两个颜色值则需要估计出来，从而得到所需要的全彩色图像，这种处理称做去马赛克。去马赛克是数码相机产品中的核心技术。

现有的色彩滤波阵列图像去马赛克方法包括最近邻复制、双线性插值、三次样条插值，目前较好的方法包括方向线性最小均方误差算法，局部方向插值和非局部均值滤波去马赛克的方法和基于局部方向插值和自适应阈值的去马赛克方法。

方向线性最小均方误差去马赛克算法，是Lei Zhang在文献“Color Demosaicking ViaDirectional Linear Minimum Mean Square-Error Estimation.IEEE Trans.on Image Processing,vol.14,no.12,pp.2167-2178,Dec.2005.”中提出来的。该方法首先通过最小均方误差LMMSE将水平和竖直方向的颜色差值信号估计出来，然后通过混合方向加权得到最终的颜色差值信号，最终将各个通道的缺失分量估计出来。这种方法由于只考虑水平和竖直方向的差值估计，所以存在拉链效应及虚假色彩的现象。

局部方向插值和非局部均值滤波去马赛克的方法，是Lei Zhang在文献“ColorDemosaicking by Local Directional Interpolation and Nonlocal Adaptive Thresholding,Journal ofElectronic Imaging20(2),023016(Apr-Jun),2011.”中提出的。该方法利用图像的非局部冗余来提高图像的局部色彩恢复效果。首先，对像素缺失的分量使用多个方向混合估计出来，然后分别对各个通道使用非局部均值来修复局部估计的结果。此方法会模糊一些边缘。

基于局部方向插值和自适应阈值的去马赛克方法，是Lei Zhang在文献“ColorDemosaicking by Local Directional Interpolation and Nonlocal Adaptive Thresholding,Journal ofElectronic Imaging20(2),023016(Apr-Jun),2011.”中提出的。此方法利用图像的非局部结构相似性来提高图像的局部色彩恢复效果，首先，对像素缺失的分量使用方向插值估计出来，然后分别对各个通道使用非局部自适应阈值来修复局插值结果。此方法可以有效的提高滤波阵列图像的去马赛克效果，并消除高饱和度图像边缘处产生的虚假结构，但仍会出现块效应。

发明内容

本发明的目的在于针对上述已有技术的不足，提出一种基于非局部低秩的色彩滤波阵列图像去马赛克方法，通过利用图像的局部和非局部冗余信息，并利用低秩来描述图像的几何结构相似性，以有效抑制拉链效应和块效应的出现，并避免了虚假颜色，提高色彩滤波阵列图像的去马赛克效果，得到高质量的全彩色图像。

为实现上述目的，本发明的技术方案包括如下步骤：

（1）输入一幅色彩滤波阵列图像I；

（2）通过方向插值法对色彩滤波阵列图像I中绿色通道缺失的像素进行估计，得到绿色通道图像的插值图像

（3）在绿色通道图像的插值图像中逐像素取一个5×5大小的图像块，作为当前待修正图像块X；

（4）在当前待修正图像块X的中心像素的32×32大小的邻域中取所有5×5的块，组成当前待修正图像块X的图像块集合Ω；

（5）计算当前待修正块X与其图像块集合Ω中每个块之间的像素欧式距离d，将欧式距离从小到大排序，取前100个块作为当前待修正块的相似块；

（6）将当前待修正块X和其相似块拉成列向量，组成相似块矩阵I_X，并使待修正块X拉成的待修正列向量Y作为I_X的中间列向量；

（7）定义相似块矩阵I_X的拉格朗日优化函数F(I′_X)：

F(I′_X)=||I_X-I′_X||_F+λ₁||I_X-I′_X||₁+λ₂||I′_X||_*，

其中，I′_X为函数F(I′_X)的自变量，c为常数，为相似块矩阵I_X中插值误差的方差，为相似块矩阵I_X的奇异值；

（8）最小化约束求解相似块矩阵I_X的拉格朗日优化函数F(I′_X)，得到修正后的相似块矩阵

${\hat{I}}_X=\underset{{I\text{'}}_X}{\min F\left({I\text{'}}_X\right),}$

其中，I′_X为函数F(I′_X)的自变量；

（9）根据修正后的相似块矩阵求修正后的列向量

$\hat{Y}=\mathrm{mid}\left({\hat{I}}_X\right),$

其中，为取修正后的相似块矩阵的中间列向量；

（10）将修正后的列向量重组为5×5的矩阵，得到修正后的图像块

（11）对绿色通道图像的插值图像的所有图像块重复执行步骤（3）-（10），完成对绿色通道图像的最终估计；

（12）通过方向插值法对色彩滤波阵列图像I中红色和蓝色通道缺失的像素进行估计，得到红色通道图像的插值图像和得到蓝色通道图像的插值图像

（13）分别在红色通道图像的插值图像和蓝色通道图像的插值图像中逐像素取一个5×5大小的图像块，作为当前待修正的图像块X，重复执行步骤（4）-（10），完成对红色和蓝色通道图像的最终估计；

（14）输出含有绿色、红色和蓝色的全彩色图像。

本发明与现有的技术相比具有以下优点：

第一，本发明通过使用精确的方向插值方法，使得图像在纹理处的插值效果明显提高，有效的抑制了拉链效应。

第二，本发明利用拉格朗日优化函数来修复相似块矩阵，进一步提高了像素恢复的准确性，有效抑制块效应的出现，并避免了虚假颜色和边缘模糊，提高色彩滤波阵列图像的去马赛克效果。

附图说明

图1是本发明的实现流程图；

图2是本发明仿真使用McMaster数据库中图像的放大图；

图3是用现有技术与本发明方法对色彩阵列图像去马赛克的仿真结果放大图。

具体实施方式

以下参照附图对本发明的具体实现及效果做进一步的详细说明。

参照图1，本发明的实施步骤如下：

步骤1，输入一幅色彩滤波阵列图像I。

本步骤输入的色彩滤波阵列图像I为Bayer模式的色彩滤波阵列图像，该图像中每个像素点仅存在红、绿、蓝三基色中的一种颜色，另外两个颜色缺失，需要采用去马赛克方法进行估计。

步骤2，对色彩滤波阵列图像I中的绿色通道图像进行方向插值。

（2.1)计算以色彩滤波阵列图像I的像素点R(i,j)为中心的北、南、东、西、水平、竖直六个方向上的绿色分量与红色分量的颜色差值：

${\mathrm{\Δ}}_{\mathrm{gr}}^n=G(i-1,j)-(R(i,j)+R(i-2,j))/2,$

${\mathrm{\Δ}}_{\mathrm{gr}}^s=G(i+1,j)-(R(i,j)+R(i+2,j))/2$

${\mathrm{\Δ}}_{\mathrm{gr}}^e=G(i,j+1)-(R(i,j)+R(i,j+2))/2$

${\mathrm{\Δ}}_{\mathrm{gr}}^w=G(i,j-1)-(R(i,j)+R(i,j-2))/2,$

${\mathrm{\Δ}}_{\mathrm{gr}}^h=(G(i,j-1)+G(i,j+1))/2-R(i,j),$

${\mathrm{\Δ}}_{\mathrm{gr}}^v=(G(i-1,j)+G(i+1,j))/2-R(i,j),$

其中，(i,j)表示像素点的位置，即该像素位于色彩滤波阵列图像I的第i行第j列，R为红色通道图像，G为绿色通道图像，分别为像素点R(i,j)在北、南、东、西、水平、竖直方向上的绿色分量与红色分量的颜色差值；

（2.2)计算像素点R(i,j)沿北、南、东、西、水平、竖直方向的梯度： ${\▿}_n=|G(i-1,j)-G(i+1,j)|+|R(i,j)-R(i-2,j)|+\frac{1}{2}|G(i,j-1)-G(i-2,j-1)|+\frac{1}{2}|G(i,j+1)-G(i-2,j+1)|+\mathrm{\ϵ},$ ${\▿}_s=|G(i-1,j)-G(i+1,j)|+|R(i,j)-R(i+2,j)|+\frac{1}{2}|G(i,j+1)-G(i+2,j+1)|+\frac{1}{2}|G(i,j-1)-G(i+2,j-1)|+\mathrm{\ϵ},$ ${\▿}_e=|G(i,j-1)-G(i,j+1)|+|R(i,j)-R(i,j+2)|+\frac{1}{2}|G(i-1,j)-G(i-1,j+2)|+\frac{1}{2}|G(i+1,j)-G(i+1,j+2)|+\mathrm{\ϵ},$ ${\▿}_w=|G(i,j-1)-G(i,j+1)|+|R(i,j)-R(i,j-2)|+\frac{1}{2}|G(i-1,j)-G(i-1,j-2)|+\frac{1}{2}|G(i+1,j)-G(i+1,j-2)|+\mathrm{\ϵ},$ ${\▿}_h=\frac{1}{4}|G(i-1,j-2)-G(i-1,j)|+\frac{1}{4}|G(i-1,j)-G(i-1,j+2)|+\frac{1}{4}|G(i+1,j-2)-G(i+1,j)|+\frac{1}{4}|G(i+1,j)-G(i+1,j+2)|$ $+|G(i,j-1)-G(i,j+1)|+\frac{1}{2}|R(i,j)-R(i,j-2)|+\frac{1}{2}|R(i,j)-R(i,j+2)|+\mathrm{\ϵ},$ ${\▿}_v=\frac{1}{4}|G(i-2,j-1)-G(i,j-1)|+\frac{1}{4}|G(i,j-1)-G(i+2,j-1)|+\frac{1}{4}|G(i-2,j+1)-G(i,j+1)|+\frac{1}{4}|G(i,j+1)-G(i+2,j+1)|$ $+|G(i-1,j)-G(i+1,j)|+\frac{1}{2}|R(i-2,j)-R(i,j)|+\frac{1}{2}|R(i,j)-R(i+2,j)|+\mathrm{\ϵ},$

其中，下标n、s、e、w、h、v分别为北、南、东、西、水平、竖直方向，▽_n、▽_s、▽_e、▽_w、▽_h、▽_v分别为像素点R(i,j)沿北、南、东、西、水平、竖直方向的梯度，ε为一个很小的常数，以避免梯度为0，其值为：ε＝0.1；

（2.3)根据步骤（2.2）所述的像素点R(i,j)沿北、南、东、西、水平、竖直方向的梯度▽_n、▽_s、▽_e、▽_w、▽_h、▽_v，计算各个方向的颜色差值权重值，即：

${\tilde{w}}_n=\frac{1}{{\▿}_n},{\tilde{w}}_s=\frac{1}{{\▿}_s},{\tilde{w}}_e=\frac{1}{{\▿}_w},{\tilde{w}}_w=\frac{1}{{\▿}_w},{\tilde{w}}_h=\frac{1}{{\▿}_h},{\tilde{w}}_v=\frac{1}{{\▿}_v};$

（2.4)根据步骤（2.3）得到的各个方向的颜色差值权重值，计算各个方向的权重值之和C：

$C={\tilde{w}}_n+{\tilde{w}}_s+{\tilde{w}}_w+{\tilde{w}}_e+{\tilde{w}}_h+{\tilde{w}}_v;$

（2.5)根据步骤（2.4）得到的权重之和C，计算归一化后的各个方向的权重值：

$w_n=\frac{{\tilde{w}}_n}{C},w_s=\frac{{\tilde{w}}_s}{C},w_e=\frac{{\tilde{w}}_e}{C},w_w=\frac{{\tilde{w}}_w}{C},w_h=\frac{{\tilde{w}}_h}{C},w_v=\frac{{\tilde{w}}_v}{C},$

（2.6)根据步骤（2.1）得到的颜色差值和步骤（2.5）得到的权重值，计算像素点R(i,j)处绿色分量值与红色分量值的颜色差值，即：

${\widehat{\mathrm{\Δ}}}_{\mathrm{gr}}=w_n{\mathrm{\Δ}}_{\mathrm{gr}}^n+w_s{\mathrm{\Δ}}_{\mathrm{gr}}^s+w_w{\mathrm{\Δ}}_{\mathrm{gr}}^w+w_e{\mathrm{\Δ}}_{\mathrm{gr}}^e+w_h{\mathrm{\Δ}}_{\mathrm{gr}}^h+w_v{\mathrm{\Δ}}_{\mathrm{gr}}^v;$

（2.7)根据步骤（2.6）得到的颜色差值计算像素点R(i,j)处缺失的绿色分量值即：

$\hat{G}(i,j)=R(i,j)+{\widehat{\mathrm{\Δ}}}_{\mathrm{gr}};$

（2.8)对绿色通道图像上的所有缺失像素执行步骤（2.1)-（2.7)，得到绿色通道图像的插值图像

步骤3，从绿色通道图像的插值图像中逐像素取块，作为当前待修正块X。

从绿色通道图像的插值图像的第17行第17列像素起，至倒数第17行第17列像素结束，以每个像素点为中心取5×5大小的图像块作为当前待修正图像块X。

步骤4，获取当前待修正图像块X的图像块集合Ω。

对当前待修正图像块X的中心像素的32×32大小的邻域，从第3行第3列像素起，到倒数第3行第3列像素结束，以每个像素点为中心取5×5大小的图像块组成当前待修正图像块X的相似块集合Ω。

步骤5，在图像块集合Ω中取当前待修正块X的相似块。

（5.1)计算当前待修正块X与其图像块集合Ω中每个块Yi之间的像素欧式距离d：

$d={\left|\right|X-Y_i\left|\right|}_2=\sqrt{\frac{1}{t^2}\underset{m=1}{\overset{t}{\mathrm{\Σ}}}\underset{n=1}{\overset{t}{\mathrm{\Σ}}}{|X(m,n)-Y_i(m,n)|}^2},$

其中，Y_i为图像块集合Ω中的第i个图像块，t为图像块的行像素数，其值为t＝5，(m,n)为像素点在图像块中的位置；

（5.2）将欧式距离d由小到大排序，取前100个块作为当前待修正块X的相似块。

步骤6，取当前待修正块X的相似块矩阵I_X。

将当前待修正块X和其相似块拉成列向量，组成相似块矩阵I_X，并使当前待修正块X拉成的待修正列向量Y作为相似块矩阵I_X的中间列向量。

步骤7，定义相似块矩阵I_X的拉格朗日优化函数F(I′_X)。

利用相似块矩阵I_X的核范数||I′_X||_*来描述其结构相似性，并定义相似块矩阵I_X的拉格朗日优化函数F(I′_X)：

$F\left({I\text{'}}_x\right)={\left|\right|I_X-{I\text{'}}_X\left|\right|}_F+{\mathrm{\λ}}_1{\left|\right|I}_X-{I\text{'}}_X{\left|\right|}_1+{\mathrm{\λ}}_2{\left|\right|{I\text{'}}_X\left|\right|}_{*},$

其中，I′_X为函数F(I′_X)的自变量，c为常数，为相似块矩阵I_X中插值误差的方差，为相似块矩阵I_X的奇异值。

步骤8，最小化约束求解相似块矩阵I_X的拉格朗日优化函数F(I′_X)，得到修正后的相似块矩阵

(8.1）通过如下函数，初始化修正的相似块矩阵

${\widehat{I^0}}_X=\arg \underset{{I\text{'}}_X}{\min }{\left|\right|{I\text{'}}_X\left|\right|}_{*}+{\mathrm{\τ}}_0{\left|\right|{I\text{'}}_X-I_X\left|\right|}_F,$

${\mathrm{\τ}}_0=c_0\frac{S_{I_X}}{{\mathrm{\σ}}_{I_X}},$

其中，||·||_*为核范数，||·||_F为F范数，I′_X为自变量矩阵，c₀为常数，其值为：c₀＝65，为相似块矩阵I_X的奇异值，为相似块矩阵I_X中插值误差的方差；

（8.2）根据步骤(8.1)得到的初始化修正的相似块矩阵求稀疏条件修正的相似块矩阵H_X：

$H_X=\underset{{H\text{'}}_X}{\arg \min }{\left|\right|{H\text{'}}_X-I_X\left|\right|}_l+{\mathrm{\τ}}_1{\left|\right|{H\text{'}}_X-I_X\left|\right|}_F+{\mathrm{\τ}}_2{\left|\right|{H\text{'}}_X-{\widehat{I^0}}_X\left|\right|}_F,$

${\mathrm{\τ}}_1=c_1/{\mathrm{\σ}}_{I_X},$

${\mathrm{\τ}}_2=c_2/{\mathrm{\σ}}_{{\widehat{I^0}}_X},$

其中，||·||₁为l¹范数，H′_X为自变量矩阵，c₁和c₂为不同的常数，其值分别为：c₁＝0.1，c₂＝0.2，为初始化修正的相似块矩阵中插值误差的方差；

(8.3）根据步骤（8.2）得到的稀疏条件修正的相似块矩阵H_X，求低秩条件修正的相似块矩阵

${\stackrel{\_}{I}}_X=\arg \underset{{I\text{'}\text{'}}_X}{\min }{\left|\right|I_X^{\text{'}\text{'}}\left|\right|}_{*}+\mathrm{\τ}{\left|\right|I_X^{\text{'}\text{'}}-H_X\left|\right|}_F,$

$\mathrm{\τ}=c\frac{S_{H_X}}{{\mathrm{\σ}}_{I_X}},$

其中，c为常数，其值为：c＝65，I″_X为求自变量矩阵，为稀疏条件修正的相似块矩阵H_X的奇异值，为稀疏条件修正的相似块矩阵H_X中插值误差的方差；

(8.4）根据步骤（8.3）得到的低秩条件修正的相似块矩阵从而求最终修正后的相似块矩阵

${\hat{I}}_X=\underset{{I\text{'}}_X}{\min }F\left({I\text{'}}_X\right)$

的估计值

步骤9，根据修正后的相似块矩阵求修正后的列向量

修正后的列向量通过如下公式计算：

$\hat{Y}=\mathrm{mid}\left({\hat{I}}_X\right),$

其中，为取修正后的相似块矩阵的中间列向量。

步骤10，将修正后的列向量重组为5×5的矩阵，得到修正后的图像块

步骤11，对绿色通道图像的插值图像的所有图像块，重复执行步骤3到步骤10，完成对绿色通道图像的最终估计。

步骤12，对色彩滤波阵列图像I中的蓝色通道图像和红色通道图像进行方向插值。

（12.1)计算以像素点G(i,j)为中心的北、南、东、西、水平、竖直六个方向上的蓝色分量与绿色分量的颜色差值：

${\mathrm{\Δ}}_{\mathrm{bg}}^n=B(i-1,j)-G(i-1,j),$

${\mathrm{\Δ}}_{\mathrm{bg}}^s=B(i+1,j)-G(i+1,j),$

${\mathrm{\Δ}}_{\mathrm{bg}}^w=B(i,j-1)-G(i,j-1),$

${\mathrm{\Δ}}_{\mathrm{bg}}^e=B(i,j+1)-G(i,j+1),$

${\mathrm{\Δ}}_{\mathrm{bg}}^h=(B(i,j-1)+B(i,j+1))/2-G(i,j),$

${\mathrm{\Δ}}_{\mathrm{bg}}^v=(B(i-1.j)+B(i+1,j))/2-G(i,j),$

其中，(i,j)表示像素点的位置，B为蓝色通道图像，G为绿色通道图像， 分别为像素点G(i,j)在北、南、东、西、水平、竖直方向上的蓝色分量与绿色分量的颜色差值；

（12.2)计算像素点G(i,j)沿北、南、东、西、水平、竖直方向的梯度：

${\▿}_n=|G(i-2,j)-G(i,j)|+|B(i-1,j)-B(i+1,j)|+\frac{1}{2}|R(i-2,j-1)-R(i,j-1)|+\frac{1}{2}|R(i-2,j+1)-R(i,j+1)|+\mathrm{\ϵ},$

${\▿}_s=|G(i+2,j)-G(i,j)|+|B(i-1,j)-B(i+1,j)|+\frac{1}{2}|R(i+2,j-1)-R(i,j-1)|+\frac{1}{2}|R(i+2,j+1)-R(i,j+1)|+\mathrm{\ϵ},$

${\▿}_e=|R(i,j-1)-R(i,j+1)|+|G(i,j)-G(i,j+2)|+\frac{1}{2}|B(i-1,j)-B(i-1,j+2)|+\frac{1}{2}|B(i+1,j)-B(i+1,j+2)|+\mathrm{\ϵ},$

${\▿}_w=|R(i,j-1)-R(i,j+1)|+|G(i,j)-G(i,j-2)|+\frac{1}{2}|B(i-1,j-2)-B(i-1,j)|+\frac{1}{2}|B(i+1,j-2)-B(i+1,j)|+\mathrm{\ϵ},$

${\▿}_h=\frac{1}{4}|G(i-1,j-1)-G(i-1,j+1)|+\frac{1}{4}|G(i+1,j-1)-G(i+1,j+1)|+\frac{1}{4}|B(i-1,j-2)-B(i-1,j)|+\frac{1}{4}|B(i-1,j)-B(i-1,j+2)|+$

$\frac{1}{4}|B(i+1,j-2)-B(i+1,j)|+\frac{1}{4}|B(i+1,j)-B(i+1,j+2)|+|R(i,j-1)-R(i,j+1)|+\frac{1}{2}G(i,j)-G(i,j-2)|+\frac{1}{2}G(i,j)-G(i,j+2)|+\mathrm{\ϵ},$

${\▿}_v=\frac{1}{4}|G(i-1,j-1)-G(i+1,j-1)|+\frac{1}{4}|G(i-1,j+1)-G(i+1,j+1)|+\frac{1}{4}|R(i,j-1)-R(i-2,j-1)|+\frac{1}{4}|R(i,j-1)-R(i+2,j-1)|+$

$\frac{1}{4}|R(i-2,j+1)-R(i,j+1)|+\frac{1}{4}|R(i,j+1)-R(i+2,j+1)|+|B(i-1,j)-B(i+1,j)|+\frac{1}{2}|G(i,j)-G(i-2,j)|+\frac{1}{2}G(i,j)-G(i+2,j)|+\mathrm{\ϵ},$

其中，下标n、s、e、w、h、v分别为北、南、东、西、水平、竖直方向，▽_n、▽_s、▽_e、▽_w、▽_h、▽_v分别为像素点G(i,j)沿北、南、东、西、水平、竖直方向的梯度，ε为一个很小的常数，以避免梯度为0，其值为：ε＝0.1；

（12.3)根据步骤（12.2）所述的像素点G(i,j)沿北、南、东、西、水平、竖直方向的梯度▽_n、▽_s、▽_e、▽_w、▽_h、▽_v，计算各个方向颜色差值的权重值，即：

${\tilde{w}}_n=\frac{1}{{\▿}_n},{\tilde{w}}_s=\frac{1}{{\▿}_s},{\tilde{w}}_e=\frac{1}{{\▿}_e},{\tilde{w}}_w=\frac{1}{{\▿}_w},{\tilde{w}}_h=\frac{1}{{\▿}_h},{\tilde{w}}_v=\frac{1}{{\▿}_v};$

（12.4)根据步骤（12.3）所述的各个方向颜色差值的权重值，计算各个方向的权重值之和C：

$C={\tilde{w}}_n+{\tilde{w}}_s+{\tilde{w}}_w+{\tilde{w}}_e+{\tilde{w}}_h+{\tilde{w}}_v;$

（12.5)根据步骤（12.4）得到的权重之和C，计算归一化后的各个方向的权重值：

$w_n=\frac{{\tilde{w}}_n}{C},w_s=\frac{{\tilde{w}}_s}{C},w_e=\frac{{\tilde{w}}_e}{C},w_w=\frac{{\tilde{w}}_w}{C},w_h=\frac{{\tilde{w}}_h}{C},w_v=\frac{{\tilde{w}}_v}{C};$

(12.6)根据步骤（12.1）得到的颜色差值和步骤（12.5）得到的权重值，计算像素点G(i,j)处蓝色分量与绿色分量的颜色差值，即：

${\widehat{\mathrm{\Δ}}}_{\mathrm{bg}}=w_n{\mathrm{\Δ}}_{\mathrm{bg}}^n+w_s{\mathrm{\Δ}}_{\mathrm{bg}}^s+w_e{\mathrm{\Δ}}_{\mathrm{bg}}^e+w_w{\mathrm{\Δ}}_{\mathrm{bg}}^w+w_h{\mathrm{\Δ}}_{\mathrm{bg}}^h+w_v{\mathrm{\Δ}}_{\mathrm{bg}}^v;$

（12.7)根据步骤（12.6）得到的像素点G(i,j)处的颜色差值计算像素点G(i,j)处缺失的蓝色分量

$\hat{B}(i,j)=G(i,j)+{\widehat{\mathrm{\Δ}}}_{\mathrm{bg}};$

（12.8)对蓝色通道图像上的所有缺失像素执行步骤（12.1)-（12.7)，得到蓝色通道图像的插值图像

（12.9)根据色彩滤波阵列图像I中红色像素与蓝色像素分布对称的特性，对红色通道图像上的所有缺失像素执行步骤（12.1)-（12.7)，得到红色通道图像的插值图像

步骤13，对红色通道的插值图像和蓝色通道的插值图像进行修正。

（13.1）分别在红色通道图像的插值图像和蓝色通道图像的插值图像中逐像素取一个5×5大小的图像块，作为当前待修正的图像块X；

（13.2）重复执行步骤4到步骤10，对红色通道的插值图像和蓝色通道的插值图像进行修正，完成对红色和蓝色通道图像的最终估计。

步骤14，输出含有绿色、红色和蓝色的全彩色图像，完成对色彩滤波阵列图像的去马赛克。

下面结合附图2和附图3对本发明的仿真效果做进一步的说明：

1．仿真条件：

在CPU为pentium(R)4处理器，主频1.86GHZ，内存2G，操作系统WINDOWS XP SP3，软件为Matlab7.10的平台上进行。

2．评价指标：

图像去马赛克的评价分为主观和客观两个方面。在主观上评价一幅色彩阵列图像去马赛克效果的优劣主要是通过人眼的视觉特性来衡量，图像质量好，图像没有拉链效应，没有模糊，感觉清晰则去马赛克效果好，反之效果则差，在客观上评价一幅色彩滤波阵列图像的去马赛克效果，本发明通过峰值信噪比PSNR来衡量。

3．仿真图像：

本发明采用McMaster数据库中的图像作为仿真图像，此数据库中的图像是由电影直接获取并数字化的图像，它含有8个高分辨的彩色图像，由于这些图像尺寸较大，在实验中，我们将它们剪裁为500×500的子图像来测量色彩滤波阵列图像去马赛克方法的性能。

4、仿真内容：

利用现有的方向线性最小均方误差去马赛克方法DLMMSE、局部方向插值和非局部均值滤波去马赛克的方法LDI-NLM、基于局部方向插值和自适应阈值的去马赛克方法LDI-NAT和本发明方法对McMaster数据库中的图像进行色彩滤波阵列图像去马赛克仿真，仿真结果如图3，其中，图3(a)是利用方向线性最小均方误差去马赛克方法DLMMSE得到的仿真结果放大图，图3(b)是利用局部方向插值和非局部均值滤波去马赛克的方法LDI-NLM得到的仿真结果放大图，图3(c)是利用基于局部方向插值和自适应阈值的去马赛克方法LDI-NAT得到的仿真结果放大图，图3(d)是用本发明方法得到的去马赛克仿真结果放大图，附图2为原图的放大图。

从图3和可见，现有方向线性最小均方误差去马赛克方法DLMMSE在图像边缘处出现了明显的拉链效应，并出现了虚假的颜色。现有局部方向插值和非局部均值滤波去马赛克的方法LDI-NLM虽然减弱了拉链效应，但在边缘处仍出现了模糊。现有基于局部方向插值和自适应阈值的去马赛克方法LDI-NAT出现了块效应。本发明的方法所处理的图像没有出现虚假颜色，避免了边缘模糊并抑制了拉链效应和块效应，更接近附图2中的原图。

为了客观上评价一幅色彩滤波阵列图像的去马赛克效果的好坏，本发明采用常用的峰值信噪比PSNR来衡量，令真实图像为A，色彩滤波阵列图像去马赛克后的图像为A?，则峰值信噪比按如下公式计算：

$\mathrm{PSNR}={10\log }_{10}\left(\frac{{255}^2}{\frac{1}{m\×n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\underset{j=1}{\overset{m}{\mathrm{\Σ}}}{(\hat{A}(i,j)-A(i,j))}^2}\right),$

其中，(i,j)为像素点在图像中的位置，m为图像列像素点的个数，n为图像行像素点的个数，

上述两种现有方法和本发明的方法对McMaster数据库中的图像进行色彩滤波阵列图像去马赛克的PSNR评价指标如表1。

表1不同的方法对McMaster数据库色彩滤波阵列图像去马赛克的峰值信噪比

从表1可以看出，本发明方法对色彩滤波阵列图像去马赛克的PSNR高于其他三种现有的方法，说明本发明方法对色彩滤波阵列图像的像素恢复准确性高，验证了本发明对色彩滤波阵列图像去马赛克的有效性。

Claims (4)

1.一种基于非局部低秩的色彩滤波阵列图像去马赛克方法，包括如下步骤：

（1）输入一幅色彩滤波阵列图像I；

（2）通过方向插值法对色彩滤波阵列图像I中绿色通道缺失的像素进行估计，得到绿色通道图像的插值图像；

（3）在绿色通道图像的插值图像中逐像素取一个5×5大小的图像块，作为当前待修正图像块X；

（4）在当前待修正图像块X的中心像素的32×32大小的邻域中取所有5×5的块，组成当前待修正图像块X的图像块集合Ω；

（5）计算当前待修正块X与其图像块集合Ω中每个块之间的像素欧式距离d，将欧式距离从小到大排序，取前100个块作为当前待修正块的相似块；

（6）将当前待修正块X和其相似块拉成列向量，组成相似块矩阵I_X，并使待修正块X拉成的待修正列向量Y作为I_X的中间列向量；

（7）定义相似块矩阵I_X的拉格朗日优化函数F(I′_X)：

F(I′_X)=||I_X-I′_X||_F+λ₁||I_X-I′_X||₁+λ₂||I′_X||_*，

其中，I′_X为函数F(I′_X)的自变量，c为常数，为相似块矩阵I_X中插值误差的方差，为相似块矩阵I_X的奇异值；

（8）最小化约束求解相似块矩阵I_X的拉格朗日优化函数F(I′_X)，得到修正后的相似块矩阵

${\hat{I}}_X=\underset{{I\text{'}}_X}{\min F\left({I\text{'}}_X\right),}$

其中，I′_X为函数F(I′_X)的自变量；

（9）根据修正后的相似块矩阵求修正后的列向量

$\hat{Y}=\mathrm{mid}\left({\hat{I}}_X\right),$

其中，为取修正后的相似块矩阵的中间列向量；

（10）将修正后的列向量重组为5×5的矩阵，得到修正后的图像块

（11）对绿色通道图像的插值图像的所有图像块重复执行步骤（3）-（10），完成对绿色通道图像的最终估计；

（12）通过方向插值法对色彩滤波阵列图像I中红色和蓝色通道缺失的像素进行估计，得到红色通道图像的插值图像和得到蓝色通道图像的插值图像；

（13）分别在红色通道图像的插值图像和蓝色通道图像的插值图像中逐像素取一个5×5大小的图像块，作为当前待修正的图像块X，重复执行步骤（4）-（10），完成对红色和蓝色通道图像的最终估计；

（14）输出含有绿色、红色和蓝色的全彩色图像。

2.根据权利要求1所述的基于非局部低秩的色彩滤波阵列图像去马赛克方法，其特征在于：步骤（2）中所述的通过方向插值法对色彩滤波阵列图像I中绿色通道缺失的像素进行估计，按如下步骤进行：

（2a)计算以像素点R(i,j)为中心的北、南、东、西、水平、竖直六个方向上的绿色分量值与红色分量值的颜色差值：

${\mathrm{\Δ}}_{\mathrm{gr}}^n=G(i-1,j)-(R(i,j)+R(i-2,j))/2,$

${\mathrm{\Δ}}_{\mathrm{gr}}^s=G(i+1,j)-(R(i,j)+R(i+2,j))/2$

${\mathrm{\Δ}}_{\mathrm{gr}}^e=G(i,j+1)-(R(i,j)+R(i,j+2))/2$

${\mathrm{\Δ}}_{\mathrm{gr}}^w=G(i,j-1)-(R(i,j)+R(i,j-2))/2,$

${\mathrm{\Δ}}_{\mathrm{gr}}^h=(G(i,j-1)+G(i,j+1))/2-R(i,j),$

${\mathrm{\Δ}}_{\mathrm{gr}}^v=(G(i-1,j)+G(i+1,j))/2-R(i,j),$

其中，(i,j)表示像素点的位置，R为红色通道图像，G为绿色通道图像， 分别为像素点R(i,j)在北、南、东、西、水平、竖直方向上的绿色分量值与红色分量值的颜色差值；

（2b)计算像素点R(i,j)沿北、南、东、西、水平、竖直方向的梯度：

${\▿}_n=|G(i-1,j)-G(i+1,j)|+|R(i,j)-R(i-2,j)|+\frac{1}{2}|G(i,j-1)-G(i-2,j-1)|+\frac{1}{2}|G(i,j+1)-G(i-2,j+1)|+\mathrm{\ϵ},$

${\▿}_s=|G(i-1,j)-G(i+1,j)|+|R(i,j)-R(i+2,j)|+\frac{1}{2}|G(i,j+1)-G(i+2,j+1)|+\frac{1}{2}|G(i,j-1)-G(i+2,j-1)|+\mathrm{\ϵ},$

${\▿}_e=|G(i,j-1)-G(i,j+1)|+|R(i,j)-R(i,j+2)|+\frac{1}{2}|G(i-1,j)-G(i-1,j+2)|+\frac{1}{2}|G(i+1,j)-G(i+1,j+2)|+\mathrm{\ϵ},$

${\▿}_w=|G(i,j-1)-G(i,j+1)|+|R(i,j)-R(i,j-2)|+\frac{1}{2}|G(i-1,j)-G(i-1,j-2)|+\frac{1}{2}|G(i+1,j)-G(i+1,j-2)|+\mathrm{\ϵ},$

${\▿}_h=\frac{1}{4}|G(i-1,j-2)-G(i-1,j)|+\frac{1}{4}|G(i+1,j-2)-G(i+1,j)|+\frac{1}{4}|G(i+1,j)-G(i+1,j+2)|$

$+|G(i,j-1)-G(i,j+1)|\frac{1}{2}|R(i,j)-R(i,j-2)|+\frac{1}{2}|R(i,j)-R(i,j+2)|+\mathrm{\ϵ},$

${\▿}_v=\frac{1}{4}|G(i-2,j-1)-G(i,j-1)|+\frac{1}{4}|G(i,j-1)-G(i+2,j-1)|+\frac{1}{4}|G(i-2,j+1)-G(i,j+1)|+\frac{1}{4}G(i,j+1)-G(i+2,j+1)|$

$+|G(i-1,j)-G(i+1,j)|+\frac{1}{2}|R(i-2,j)-R(i,j)|+\frac{1}{2}|R(i,j)-R(i+2,j)|+\mathrm{\ϵ},$

其中，下标n、s、e、w、h、v分别为北、南、东、西、水平、竖直方向，▽_n、▽_s、▽_e、▽_w、▽_h、▽_v分别为像素点R(i,j)沿北、南、东、西、水平、竖直方向的梯度，ε为一个很小的常数；

（2c)根据步骤（2b）所述的像素点R(i,j)沿北、南、东、西、水平、竖直方向的梯度▽_n、▽_s、▽_e、▽_w、▽_h、▽_v，计算各个方向颜色差值的权重值，即：

${\tilde{w}}_n=\frac{1}{{\▿}_n},{\tilde{w}}_s=\frac{1}{{\▿}_s},{\tilde{w}}_e=\frac{1}{{\▿}_e},={\tilde{w}}_w\frac{1}{{\▿}_w},{\tilde{w}}_h=\frac{1}{{\▿}_h},{\tilde{w}}_v=\frac{1}{{\▿}_v},$

（2d)根据步骤（2c）所述的各个方向颜色差值的权重值，计算各个方向的权重值之和C，并对各个方向颜色差值的权重值进行归一化，得到归一化后的各个方向的权重值：

$w_n=\frac{{\tilde{w}}_n}{C},w_s=\frac{{\tilde{w}}_s}{C},w_e=\frac{{\tilde{w}}_e}{C},w_w=\frac{{\tilde{w}}_w}{C},w_h=\frac{{\tilde{w}}_h}{C},w_v=\frac{{\tilde{w}}_v}{C},$

其中， $C={\tilde{w}}_n+{\tilde{w}}_s+{\tilde{w}}_w+{\tilde{w}}_e+{\tilde{w}}_h+{\tilde{w}}_v;$

（2e)根据步骤（2a）得到的颜色差值和步骤（2d）得到的权重值，计算像素点R(i,j)处绿色分量值与红色分量值的颜色差值，即：

${\widehat{\mathrm{\Δ}}}_{\mathrm{gr}}=w_n{\mathrm{\Δ}}_{\mathrm{gr}}^n+w_s{\mathrm{\Δ}}_{\mathrm{gr}}^s+w_w{\mathrm{\Δ}}_{\mathrm{gr}}^w+w_e{\mathrm{\Δ}}_{\mathrm{gr}}^e+w_h{\mathrm{\Δ}}_{\mathrm{gr}}^h+w_v{\mathrm{\Δ}}_{\mathrm{gr}}^v;$

（2f)根据步骤（2e）所述像素点R(i,j)处的颜色差值计算像素点R(i,j)处缺失的绿色分量值，即：

$\hat{G}(i,j)=R(i,j)+{\widehat{\mathrm{\Δ}}}_{\mathrm{gr}};$

（2g)对绿色通道图像上的所有缺失像素执行步骤（2a)-（2f)，得到绿色通道图像的插值图像

3.据权利要求1所述的基于非局部低秩的色彩滤波阵列图像去马赛克方法，其特征在于：步骤（8）中所述的最小化约束求解相似块矩阵I_X的拉格朗日优化函数F(I′_X)，得到修正后的相似块矩阵按如下步骤进行：

(8a）通过如下函数，初始化修正的相似块矩阵

${{\hat{I}}^0}_X=\arg \underset{{I\text{'}}_X}{\min }{\left|\right|{I\text{'}}_X\left|\right|}_{*}+{\mathrm{\τ}}_0{\left|\right|{I\text{'}}_X-I_X\left|\right|}_F,$

${\mathrm{\τ}}_0=c_0\frac{S_{I_X}}{{\mathrm{\σ}}_{I_X}},$

其中，||·||_*为核范数，||·||_F为F范数，I′_X为自变量矩阵，c₀为常数，为相似块矩阵I_X的奇异值，为相似块矩阵I_X中插值误差的方差；

(8b）根据步骤(8a)得到的初始化修正的相似块矩阵求稀疏条件修正的相似块矩阵H_X：

$H_X=\underset{{H\text{'}}_X}{\arg \min }{\left|\right|{H\text{'}}_x-I_X\left|\right|}_l+{\mathrm{\τ}}_1{\left|\right|H\text{'}}_X-I_X{\left|\right|}_F+{\mathrm{\τ}}_2\left|\right|{H\text{'}}_X-{\widehat{I^0}}_X{\left|\right|}_F,$

${\mathrm{\τ}}_1=c_1/{\mathrm{\σ}}_{I_X},$

${\mathrm{\τ}}_2=c_2/{\mathrm{\σ}}_{{{\hat{I}}^0}_X},$

其中，||·||₁为l¹范数，H′_X为自变量矩阵，c₁和c₂为不同的常数，其值分别为：c₁＝0.1，c₂＝0.2，为初始化修正的相似块矩阵中插值误差的方差；

(8c）根据步骤（8b）得到的稀疏条件修正的相似块矩阵H_X，求低秩条件修正的相似块矩阵：

${\stackrel{\_}{I}}_X=\arg \underset{{I\text{'}}_X}{\min }{\left|\right|I_X^{\text{'}\text{'}}\left|\right|}_{*}+\mathrm{\τ}{\left|\right|I_X^{\text{'}\text{'}}-H_X\left|\right|}_F,$

$\mathrm{\τ}=c\frac{S_{H_X}}{{\mathrm{\σ}}_{H_X}},$

其中，c为常数，I″_X为自变量矩阵，为稀疏条件修正的相似块矩阵H_X的奇异值，为稀疏条件修正的相似块矩阵H_X中插值误差的方差；

(8d）根据步骤（8c）得到的低秩条件修正的相似块矩阵求最终修正后的相似块矩阵

${\hat{I}}_X=\underset{{I\text{'}}_X}{\min }F\left({I\text{'}}_X\right)$ 的估计值

4.根据权利要求1所述的基于非局部低秩的色彩滤波阵列图像去马赛克方法，其特征在于：步骤（12）所述的通过方向插值法对色彩滤波阵列图像I中蓝色通道缺失的像素和红色通道缺失的像素进行估计，按如下步骤进行：

（12a)计算以像素点G(i,j)为中心的北、南、东、西、水平、竖直六个方向上的蓝色分量值与绿色分量值的颜色差值：

${\mathrm{\Δ}}_{\mathrm{bg}}^n=B(i-1,j)-G(i-1,j),$

${\mathrm{\Δ}}_{\mathrm{bg}}^s=B(i+1,j)-G(i+1,j),$

${\mathrm{\Δ}}_{\mathrm{bg}}^w=B(i,j-1)-G(i,j-1),$

${\mathrm{\Δ}}_{\mathrm{bg}}^e=B(i,j+1)-G(i,j+1),$

${\mathrm{\Δ}}_{\mathrm{bg}}^h=(B(i,j-1)+B(i,j+1))/2-G(i,j),$

${\mathrm{\Δ}}_{\mathrm{bg}}^v=(b(i-1,j)+B(i+1,j))/2-G(i,j),$

其中，(i，j)表示像素点的位置，B为蓝色通道图像，G为绿色通道图像， 分别为像素点G(i,j)在北、南、东、西、水平、竖直方向上的蓝色分量值与绿色分量值的颜色差值；

（12b）计算像素点G(i,j)沿北、南、东、西、水平、竖直方向的梯度： ${\▿}_n=|G(i-2,j)-G(i,j)|+|B(i-1,j)-B(i+1,j)|+\frac{1}{2}|R(i-2,j-1)-R(i,j-1)|+\frac{1}{2}|R(i-2,j+1)-R(i,j+1)|+\mathrm{\ϵ},$ ${\▿}_s=|G(i+2,j)-G(i,j)|+|B(i-1,j)-B(i+1,j)|+\frac{1}{2}|R(i+2,j-1)-R(i,j-1)|+\frac{1}{2}|R(i+2,j+1)-R(i+2,j+1)|+\mathrm{\ϵ},$ ${\▿}_e=|R(i,j-1)-R(i,j+1)|+|G(i,j)-G(i,j+2)|+\frac{1}{2}|B(i-1,j)-B(i-1,j+2)|+\frac{1}{2}|B(i+1,j)-B(i+1,j+2)|+\mathrm{\ϵ},$ ${\▿}_w=|R(i,j-1)-R(i,j+1)|+|G(i,j)-G(i,j-2)|+\frac{1}{2}|B(i-1,j-2)-B(i-1,j)|+\frac{1}{2}|B(i+1,j-2)-B(i+1,j)|+\mathrm{\ϵ},$ ${\▿}_h=\frac{1}{4}|G(i-1,j-2)-G(i-1,j+1)|+\frac{1}{4}|G(i+1,j-1)-G(i+1,j+1)|+\frac{1}{4}|B(i-1,j-2)-B(i-1,j)|+\frac{1}{4}|B(i-1,j)-B(i-1,j+2)|+$ $\frac{1}{4}|B(i+1,j-2)-B(i+1,j)|+\frac{1}{4}|B(i+1,j)-B(i+1,j+2)|+|R(i,j-1)-R(i,j+1)|+\frac{1}{2}G(i,j)-G(i,j-2)|+\frac{1}{2}G(i,j)-G(i,j+2)|+\mathrm{\ϵ},$ ${\▿}_v=\frac{1}{4}|G(i-1,j-1)-G(i+1,j-1)|+\frac{1}{4}|G(i-1,j+1)-G(i+1,j+1)|+\frac{1}{4}|R(i,j-1)-R(i-2,j-1)+\frac{1}{4}|R(i.j-1)R(i+2,j-1)|+$ $\frac{1}{4}|R(i-2,j+1)-R(i,j+1)|+\frac{1}{4}|R(i,j+1)-R(i+2,j+1)|+|B(i-1,j)-B(i+1,j)|+\frac{1}{2}|G(i,j)-G(i-2,j)|+\frac{1}{2}|G(i,j)-G(i,j+2)|+\mathrm{\ϵ},$

其中，下标n、s、e、w、h、v分别为北、南、东、西、水平、竖直方向，▽_n、▽_s、▽_e、▽_w、▽_h、▽_v分别为像素点G(i,j)沿北、南、东、西、水平、竖直方向的梯度，ε为一个很小的常数；

（12c)根据步骤（12b）所述的像素点G(i,j)沿北、南、东、西、水平、竖直方向的梯度▽_n、▽_s、▽_e、▽_w、▽_h、▽_v，计算各个方向颜色差值的权重值，即：

${\tilde{w}}_n=\frac{1}{{\▿}_n},{\tilde{w}}_s=\frac{1}{{\▿}_s},{\tilde{w}}_e=\frac{1}{{\▿}_e},{\tilde{w}}_w=\frac{1}{{\▿}_w},{\tilde{w}}_h=\frac{1}{{\▿}_h},{\tilde{w}}_v=\frac{1}{{\▿}_v},$

（12d)根据步骤（12c）所述的各个方向颜色差值的权重值，计算各个方向的权重值之和C，并对各个方向颜色差值的权重值进行归一化，得到归一化后的各个方向的权重值：

$w_n=\frac{{\tilde{w}}_n}{C},w_s=\frac{{\tilde{w}}_s}{C},w_e=\frac{{\tilde{w}}_e}{C},w_w=\frac{{\tilde{w}}_w}{C},w_h=\frac{{\tilde{w}}_h}{C},w_v=\frac{{\tilde{w}}_v}{C},$

其中， $C={\tilde{w}}_n+{\tilde{w}}_s+{\tilde{w}}_e+{\tilde{w}}_w+{\tilde{w}}_h+{\tilde{w}}_v;$

(12e)根据步骤（12a）得到的颜色差值和步骤（12d）得到的权重值，计算像素点G(i,j)处蓝色分量值与绿色分量值的颜色差值，即：

${\widehat{\mathrm{\Δ}}}_{\mathrm{bg}}=w_n{\mathrm{\Δ}}_{\mathrm{bg}}^n+w_s{\mathrm{\Δ}}_{\mathrm{bg}}^s+w_e{\mathrm{\Δ}}_{\mathrm{bg}}^e+w_w{\mathrm{\Δ}}_{\mathrm{bg}}^w+w_h{\mathrm{\Δ}}_{\mathrm{bg}}^h+w_v{\mathrm{\Δ}}_{\mathrm{bg}}^v;$

（12f)根据步骤（12e）所述像素点G(i,j)处的颜色差值计算像素点G(i,j)处缺失的蓝色分量

$\hat{B}(i,j)=G(i,j)+{\widehat{\mathrm{\Δ}}}_{\mathrm{bg}};$

（12g)对蓝色通道图像上的所有缺失像素执行步骤（12a)-（12f)，得到蓝色通道图像的插值图像

（12h)根据色彩滤波阵列图像I中红色像素与蓝色像素分布对称的特性，对红色通图像上的所有缺失像素执行步骤（12a)-（12f)，得到红色通道图像的插值图像