图像Demosaic算法实验报告.doc

图像Demosaic算法实验报告
图像Demosaic算法
一、 马赛克成因
一般的数码相机采用CDC或CMOS做传感器,约占数码相机总成本的10%一25%,是相机的最昂贵的部件之一"对于专业相机,如用在科学或医学中的,使用了三个传感器,每一个像素点都可以直接获得RGB三个值,每一个传感器需要精确的机械控制,以便对应着各自的色彩通道,这是最昂贵的解决方法"在这种情况下,每一个传感器得到一幅完整的单通道图像,最后再合成一幅全彩色图像"不出现马赛克现象"

为了减少成本,缩小数码相机的体积,一般市场上的相机大都是单传感器,在传
感器的上面覆盖了一层按照Bayer模式排列的色彩滤波阵列,每一个像素点只能够捕获一个色彩值,而缺失另外两个色彩值,这时候得到的是一幅马赛克图像,为了得到全彩色的图像,我们需要利用其周围的像素点的色彩值来估计出缺失的另外两种颜色,这种处理叫做去马赛克"
二、常见的demosaic方法
关于 demosaicing 方法在过去二十年中一直是研究热点，已有的文献中就报道了大量有价值的插值算法，一般来说算法可以分为两类。第一类是单个颜色通道独立插值算法，包括邻域插值法、双线性插值法、卷积插值法等，规律是未知的绿色分量值仅由已知的绿色像素的强度进行估算，对红和蓝色通道亦是如此。其中最典型、应用最多的是双线性法。这一类算法容易实现，在平滑区域内也可以得到比较满意的效果，但在高频区域，尤其是在边界区域却失真明显。第二类则是利用多通道的相关性进行插值，这一类自适应的算法常常结合了对图像细节的分析判断以及颜色通道之间的相关性，插值效果要明显优于前一类，大部分的
算法都属于这一类。例如基于边界的算法，加权系数法，交互式插值法，最优化恢复，也有部分算法利用了信号变换如小波变换，傅立叶变换等，这些算法通常能得到比较满意的效果，尤其是清晰的边界及细节特征，但往往算法复杂，无法应用于数码相机信号转换系统。
（1） 双线性插值法
双线性（bilinear）法属于单通道独立插值法，作为最传统、简单、基础的插
值算法之一，它仍然是理解并设计新型算法的基础，对参考评估其它算法的性能
仍具有广泛的意义，参考图 3-1，该图为标注了坐标位置的 Bayer 采样阵列，传感
器只获得相应位置所对应的单色分量的强度。

考虑坐标(3,3)，此处只有红色信息，为了得到全彩色信息，必须恢复出蓝色和
绿色分量，利用双线性法可以得到：

由于双线性法是单通道之间独立插值，并且始终取 3×3 滤波器的平均值，忽
略了细节信息，以及图片三个颜色通道之间的相关性，所以双线性法往往不能得
到满意的插值效果，容易在细线状态结构的边缘处产生锯齿形图案。
（2） 自适应插值法
在 Bayer 格式中，绿色像素占了总像素的一半，也反映了更多的细节，因此，
插值都会优先考虑对绿色通道的恢复，而红色和蓝色的比例较少，细节信息少，
需要估算的信息多，失真会明显，为改进算法的性能，多会从改善红色及蓝色通
道着手。
为合理解决边界的插值问题，研究人员提出了基于梯度的插值算法。为保证
插值是沿着边界进行，而没有跨越边界，在恢复绿色分量时，先通过计算梯度来
检测边界方向，从而选择合适的插值方向。例如，为恢复出采样点的绿色分量，
先计算出该点水平和竖直方向的梯度算子。参考图 3-1，以 B44为例，假设α 为该采样点水平方向的梯度，β 为竖直方向的梯度，
Laroche则通过计算二阶微分即色度分量来确定边界方向，其中：

如果 α < β，则认为水平方向存在边界的概率比较大，插值会沿着最可能的边
界进行，即在水平方向进行，而水平方向的邻域像素及被插值点的像素的相关性
也会大于竖直方向的像素。插值结果如下式所示：

如果 α > β，则认为在竖直方向存在边界的概率大，插值在竖直方向进行，插
值结果为：

如果两者相比较相等的话，则认为该像素点处在一个平滑区域内，插值在整
个邻域内进行，认为：

一旦亮度被确定，则色度分量通常是由亮度信号和色度信号的相关性来进行
插值。

第二步：在对红色和蓝色通道插值时，用到的基本思想是基于色差的插值算
法。色差思想及色比恒定的思想是一致的，它认为在图片的小平滑区域内，色差
是恒定的。假设 (i ,j)邻近的一个采样点是 ( m, n)，则有：

上面的等式成立即是利用了色差恒定的原则。

假设我们需要恢复红色通道信息，根据阵列的排列规律，有三种不同的情况，
分别如下：
(1) 假设要恢复 G34处的红色分量：

(2) 假设要恢复 G43处的红色分量：

(3)假设要恢复 B44处的红色分量：

蓝色通道处理往往是及红色通道一致，必须注意的是在色度值确立之时，已
经获得了整个图片的亮度信号。经过这两个步骤的处理，可以获得 Bayer 图片的重构图像。
以上所述，即是自适应插值法的过程。

三、 实验结果
原图
Bayer阵列

非自适应插值（双线性）

自适应插值法

四、 总结
通过三类图像的比较，可以看出真彩色图像的三个颜色通道之间有很强的相
关性，插值最好结合其相关性，而不能依赖于单通道独立插值；其次，由于图片
的失真往往出现在高频区域，必须考虑基于边界的插值方案，对图像的区域细节
特征进行分析，使得插值尽可能沿着边界方向，或是使得边界方向插值核的加权
系数尽可能大。

附录：
Matlab代码
%双线性插值
clear;
clc;
I=imread('9.jpg');
figure
imshow(I);

[h,w,c]=size(I);
[m,n]=size(I);
m=h,
n=w;
for i=1:2:h
for j=1:w
if mod(j,2)==1
IMG(i,j,1)=I(i,j,1);
else
IMG(i,j,2)=I(i,j,2);
end
end
end
for i=2:2:h
for j=1:w
if mod(j,2)==1
IMG(i,j,2)=I(i,j,2);
else
IMG(i,j,3)=I(i,j,3);
end
end
end
figure
imshow(IMG);
%IMG11=rgb2gray(IMG);
%save('Bayer.png','IMG11')
IMG=double(IMG);

for i=2:2:m-1
for j=2:2:n-1
x11 = IMG(i-1,j-1,1);
x12 =IMG(i-1,j+1,1);
x21 = IMG(i+1,j-1,1);
x22 = IMG(i+1,j+1,1);
IMG(i,j,1) = uint8( ( x11+x21+x12+x22)/4 );
end
end

for i=2:m-1
for j=2:n-1
if IMG(i,j,1)==0
IMG(i,j,1)=0.25*(IMG(i,j+1,1)+IMG(i,j-1,1)+IMG(i+1,j,1)+IMG(i-1,j,1));
else
IMG(i,j,1)=IMG(i,j,1);
end
end
end

for i=2:2:m-1
for j=2:2:n-1
IMG(i,j,2)=0.25*(IMG(i,j+1,2)+IMG(i,j-1,2)+IMG(i+1,j,2)+IMG(i-1,j,2));
end
end

for i=2:m-1
for j=2:n-1
if IMG(i,j,2)==0
x11 = IMG(i-1,j-1,2);
x12 =IMG(i-1,j+1,2);
x21 = IMG(i+1,j-1,2);
x22 = IMG(i+1,j+1,2);
IMG(i,j,2) = uint8( 0.25*( x11+x21+x12+x22) );
else
IMG(i,j,2)=IMG(i,j,2);
end
end
end

for i=3:2:m-1
for j=3:2:n-1
x11 = IMG(i-1,j-1,3);
x12 = IMG(i-1,j+1,3);
x21 = IMG(i+1,j-1,3);
x22 = IMG(i+1,j+1,3);
IMG(i,j,3) = uint8( 0.25*( x11+x21+x12+x22) );
end
end

for i=2:m-1
for j=2:n-1
if IMG(i,j,3)==0
IMG(i,j,3)=0.25*(IMG(i,j+1,3)+IMG(i,j-1,3)+IMG(i+1,j,3)+IMG(i-1,j,3));
else
IMG(i,j,3)=IMG(i,j,3);
end
end
end

figure
imshow(uint8(IMG));

%梯度
clear;
clc;
I=imread('DD.jpg');
figure
imshow(I);

[h,w,c]=size(I);
m=h;
n=w;
R=zeros(m,n);
G=zeros(m,n);
B=zeros(m,n);
R1=zeros(m+4,n+4);
G1=zeros(m+4,n+4);
B1=zeros(m+4,n+4);
for i=1:2:h
for j=1:w
if mod(j,2)==1
IMG(i,j,1)=I(i,j,1);

else
IMG(i,j,2)=I(i,j,2);

end
end
end
for i=2:2:h
for j=1:w
if mod(j,2)==1
IMG(i,j,2)=I(i,j,2);

else
IMG(i,j,3)=I(i,j,3);

end
end
end
R=IMG(:,:,1);
G=IMG(:,:,2);
B=IMG(:,:,3);
R1(3:m+2,3:n+2)=R;
R1(3:m+2,1)=R(:,2);
R1(3:m+2,2)=R(:,1);
R1(3:m+2,n+3)=R(:,n);
R1(3:m+2,n+4)=R(:,n-1);
R1(1,:)=R1(4,:);
R1(2,:)=R1(3,:);
R1(m+3,:)=R1(m+2,:);
R1(m+4,:)=R1(m+1,:);
G1(3:m+2,3:n+2)=G;
G1(3:m+2,1)=G(:,2);
G1(3:m+2,2)=G(:,1);
G1(3:m+2,n+3)=G(:,n);
G1(3:m+2,n+4)=G(:,n-1);
G1(1,:)=G1(4,:);
G1(2,:)=G1(3,:);
G1(m+3,:)=G1(m+2,:);
G1(m+4,:)=G1(m+1,:);
B1(3:m+2,3:n+2)=B;
B1(3:m+2,1)=B(:,2);
B1(3:m+2,2)=B(:,1);
B1(3:m+2,n+3)=B(:,n);
B1(3:m+2,n+4)=B(:,n-1);
B1(1,:)=B1(4,:);
B1(2,:)=B1(3,:);
B1(m+3,:)=B1(m+2,:);
B1(m+4,:)=B1(m+1,:);
X=cat(3,R1,G1,B1);
figure,imshow(uint8(X))
for i=4:2:m
for j=4:2:n
a=abs(2*B1(i,j)-B1(i,j-2)-B1(i,j+2));
b=abs(2*B1(i,j)-B1(i-2,j)-B1(i+2,j));
if a<b
G1(i,j)=(G1(i,j-1)+G1(i,j+1))/2;
elseif a>b
G1(i,j)=(G1(i-1,j)+G1(i+1,j))/2;
else
G1(i,j)=(G1(i-1,j)+G1(i,j-1)+G1(i,j+1)+G1(i+1,j))/4;

end
end
end
for i=3:2:m-1
for j=3:2:n-1
a=abs(G1(i,j-1)-G1(i,j+1));
b=abs(G1(i-1,j)-G1(i+1,j));
if a<b
G1(i,j)=(G1(i,j-1)+G1(i,j+1))/2;
else if a>b
G1(i,j)=(G1(i-1,j)+G1(i+1,j))/2;
else
G1(i,j)=(G1(i,j-1)+G1(i,j+1)+G1(i-1,j)+G1(i+1,j))/4;
end
end
end
end

for i=4:2:m
for j=4:2:n

a=abs(2* (G1(i,j)-G1(i-1,j+1)-G1(i+1,j-1)))+abs(G1(i-1,j+1)-G1(i+1,j-1));
b=abs(2*(G1(i,j)-G1(i-1,j-1)-G1(i+1,j+1)))+abs(G1(i-1,j-1)-G1(i+1,j+1));
if a<b
R1(i,j)=(R1(i-1,j+1)+R1(i+1,j-1))/2+(2*G1(i,j)-G1(i-1,j+1)-G1(i+1,j-1))/2;
else if a>b
R1(i,j)=(R1(i-1,j-1)+R1(i+1,j+1))/2+(2*(G1(i,j)-G1(i-1,j-1)-G1(i+1,j+1)))/2;
else
R1(i,j)=(R1(i-1,j-1)+R1(i-1,j+1)+R1(i+1,j-1)+R1(i+1,j+1))/4+(4*G1(i,j)-G1(i-1,j-1)-G1(i-1,j+1)-G1(i+1,j-1)-G1(i+1,j+1))/4;
end
end
end
end
for i=3:2:m-1;
for j=4:2:n;
R1(i,j)=(R1(i,j-1)+R1(i,j+1))/2-(G1(i,j-1)+G1(i,j+1))/2+G1(i,j);
end
end
for i=4:2:m;
for j=3:2:n-1;
R1(i,j)=(R1(i-1,j)+R1(i+1,j))/2-(G1(i-1,j)+G1(i+1,j))/2+G1(i,j);
end
end
%恢复蓝色分量
for i=3:2:m-1
for j=3:2:n-1
a=abs(2*(G1(i,j)-G1(i-1,j+1)-G1(i+1,j-1)))+abs(G1(i-1,j+1)-G1(i+1,j-1));
b=abs(2*(G1(i,j)-G1(i-1,j-1)-G1(i+1,j+1)))+abs(G1(i-1,j-1)-G1(i+1,j+1));
if a<b
B1(i,j)=(B1(i-1,j+1)+B1(i+1,j-1))/2+(2*G1(i,j)-G1(i-1,j+1)-G1(i+1,j-1))/2;
else if a>b
B1(i,j)=(B1(i-1,j-1)+B1(i+1,j+1))/2+(2*(G1(i,j)-G1(i-1,j-1)-G1(i+1,j+1)))/2;
else
B1(i,j)=(B1(i-1,j-1)+B1(i-1,j+1)+B1(i+1,j-1)+B1(i+1,j+1))/4+(4*G1(i,j)-G1(i-1,j-1)-G1(i-1,j+1)-G1(i+1,j-1)-G1(i+1,j+1))/4;
end
end
end
end
for i=3:2:m-1;
for j=4:2:n;
B1(i,j)=(B1(i,j-1)+B1(i,j+1))/2-(G1(i,j-1)+G1(i,j+1))/2+G1(i,j);
end
end
for i=4:2:m;
for j=3:2:n-1;
B1(i,j)=(B1(i-1,j)+B1(i+1,j))/2-(G1(i-1,j)+G1(i+1,j))/2+G1(i,j);
end
end
X=cat(3,R1,G1,B1);
X=uint8(X);
figure
imshow(X);

22 / 22