摘要:本文要介绍的该方法完成了对彩色滤波阵列(color filter array,CFA)采样高分辨率RGB图像的过程进行了数学建模,通过对公式重组分离出了一个亮度通道和两个色差通道,通过频域分析发现两色差通道被搬移到了特定高频位置,其中一个色差通道的频域表达具有各向异性,即水平/垂直方向与亮度通道的串扰程度不同。因此可使用各向异性的滤波器进行该色差通道的恢复,对另一个色差通道使用各向同性滤波器进行恢复,再结合CFA图像获得亮度通道图像,这样就可以重构出高分辨率RGB图像了。
注:本文源于对论文“Least-squares luma-chroma demultiplexing algorithm for Bayer demosaicking”及论文“Frequency-Domain Methods for Demosaicking of Bayer-Sampled Color Images”的解读。论文及程序下载地址为:https://download.csdn.net/download/u014230360/12527273。“Least-squares luma-chroma demultiplexing algorithm for Bayer demosaicking”忘记包含在内,读者可通过Sci-Hub自行下载。
CFA图像中R、G、B通道的采样函数分别为:
由此,CFA信号可表示为:
即:
其中:
反过来有:
因此,求出、、即可重构出、、。
对做傅里叶变换,有:
对于,显然两个色度通道值为0。一般来说,色度通道具有更低的能量,占据更低的带宽。典型CFA图像频谱如下图所示:
设,,可以看出,垂直方向上受到更少的串扰,的水平频率受到更少的串扰,因此,可以使用各向异性的滤波器从中提取和,使用各向同性的滤波器来提取。有了、和,也就不难得到了,进而、、即可重构出来。图2是提取各分量滤波器的频域响应(无需使用滤波器提取)
fig2. 提取各分量的滤波器的频域响应 (a) (b) (c) (d)
各向异性的滤波器、及各向同性的可采用利用训练集获得的最小二乘估计解(根据测试集中高分辨率RGB图像生成CFA图像,其色度通道标准值已知,使滤波恢复结果与标准值的误差平方和最小即可)。
下面给出算法去马赛克效果与常规方法的对比:
下面贴上论文给出的'GRBG'格式bayer图像去马赛克源码(其他格式的可通过变量代换完成,含有滤波器包的完整源码可从我上传的资源中获得,下载地址:https://download.csdn.net/download/u014230360/12527273):
% Demosaicking of Noisy Bayer-Sampled Color Images with Least-Squares
% Luma-Chroma Demultiplexing and Noise Level Estimation
%
% IEEE Trans. Image Processing (submitted)
% This software is for provided for non-commercial and research purposes only
% Copyright Eric Dubois ([email protected]), University of Ottawa 2011
%
% [RGB] = demos_LSLCD_NE(CFA, fpkg, sigma, method)
% CFA(:,:) input Bayer CFA mosaicked image, double according to the
% pattern G R
% B G
% fpkg path of the filter package that contains the follwing
% filter: h1, h2a, h2b, hL, hG1, hG2
% sigma added noise (sigma_A) or estimated noise (sigma_E)
% method 1: h_LN, 2: BM3D
% RGB(:,:,3) output RGB image, double
% BM3D can be downloaded: http://www.cs.tut.fi/~foi/GCF-BM3D/
function [RGB] = demos_LSLCD_NE(CFA, fpkg, sigma, method)
% Load the filter package
load(fpkg);
% Use the filters to demosaic the CFA image
S = size(CFA); N1 = S(1); N2 = S(2);
yc = 0:N1-1; xc = 0:N2-1;
[XC,YC] = meshgrid(xc,yc);
% % Filter the input image with the two Gaussian filters to get energy terms
% eX = imfilter(CFA,hG1,'replicate','same').^2;
% eY = imfilter(CFA,hG2,'replicate','same').^2;
% % Average energy with moving average filter
% NMA=5; h_MA = ones(NMA,NMA)/(NMA^2);
% eX = imfilter(eX,h_MA,'replicate','same');
% eY = imfilter(eY,h_MA,'replicate','same');
% % Compute weighting coefficients
% w = eY./(eX+eY);
w=0.5*ones(N1,N2);
% Extract chrominance in corners using h1
C1mhat = imfilter(CFA,h1,'replicate','same');
% Extract chrominance on sides at f_y = 0
C2mahat = imfilter(CFA,h2a,'replicate','same');
% Extract chrominance on sides at f_x = 0
C2mbhat = imfilter(CFA,h2b,'replicate','same');
% Estimate the C2 component
C2hat = (w.*C2mahat.*(-1).^(XC) - (1-w).*C2mbhat.*(-1).^(YC));
% Estimate the C1 component
C1hat = C1mhat.*(-1).^(XC+YC);
% Estimate the Luma component
Lhatn = CFA - C1mhat - C2hat.*((-1).^XC - (-1).^YC);
if method == 1
if sigma<=0
Lhat = Lhatn;
else
Lhat = imfilter(Lhatn,hL,'replicate','same');
end
elseif method == 2
if sigma<=0
Lhat = Lhatn;
else
Lhat=BM3D(255*Lhatn,sigma)/255;
end
else
display('error');
end
% Reconstruct the RGB image
RGB(:,:,1) = Lhat - C1hat - 2*C2hat;
RGB(:,:,2) = Lhat + C1hat;
RGB(:,:,3) = Lhat - C1hat + 2*C2hat;
RGB(RGB>1)=1;
RGB(RGB<0)=0;