数字图像处理(六)——Matlab实现频域图像分析、FFT实现4:1的图像压缩

实验内容

实验一:

自选一幅图像,并对其分别添加一定强度的周期噪声和高斯噪声,然后分别采用高斯模板、中值滤波的时域方法以及傅里叶变换和小波变换的频率滤波方法对该含噪图像进行去噪处理,并基于PSNR值和视觉效果这两个指标来比较这四种滤波的方法对两种不同噪声的去噪能力。

实验二:

1、编写一个程序,要求实现下列算法:首先将图像分割为8×8的子图像,对每个子图像进行FFT,对每个子图像中的64个系数,按照每个系数的方差排序后,舍去小的变换系数,只保留16个系数,实现4:1的图像压缩。
2、小波变换是如何定义的?小波分析的主要优点是什么?
3、给定一幅行和列都为2的整数次幂图像,用Haar小波基函数对其进行二维小波变换,试着将最低尺度近似分量置零再反变换,结果是什么?如果把垂直方向的细节分量置零,反变换后结果又是什么呢?试解释一下原因。
4、基于小波变换对图像进行不同压缩比的压缩。在同压缩比情况下,对于基于小波变换和基于傅里叶变换的压缩结果,比较二者保留原图像能量百分比情况。

部分Matlab函数解释

fft2函数:图像的二维快速傅里叶变换,能够将图像从空间域变换到频率域。Ifft函数:图像的二维快速傅里叶变换逆变换,能够将图像从频率域变换到空间域。Blkproc函数:能够将每个显示块从图像中提取出来,然后将其作为参数传递给任何用户函数。
Im2col函数:将图像块排列成向量,转换成一列,重新组合成矩阵。
col2im函数:将向量排列成图像块。
sort(A)函数:如果A是向量,对A中元素按照升序排列,如果A是矩阵,对A按每一列元素按照升序排列。

1.1实验一:

自选一幅图像,并对其分别添加一定强度的周期噪声和高斯噪声,然后分别采用高斯模板、中值滤波的时域方法以及傅里叶变换和小波变换的频率滤波方法对该含噪图像进行去噪处理,并基于PSNR值和视觉效果这两个指标来比较这四种滤波的方法对两种不同噪声的去噪能力。

1.1添加噪声

通过imnoise函数对图像添加高斯噪声,通过for循环对图像添加sin和cos的信号40sin(40i)+40sin(40j),可以实现对图像添加周期噪声。实验结果如下图所示:
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

1.2通过高斯加权滤波模板对图像中的每个像素点进行处理,得出去掉噪声后的图像,通过调用sort函数对每个模块进行排序,取中间的像素点作为该像素点的灰度值,实现中值滤波去噪。去噪后的图像如图所示:

在这里插入图片描述
在这里插入图片描述

1.3通过调用fft2函数对该图像进行快速傅里叶变换,得出频谱图后进行低通滤波去噪,将绝对值小于50的频域值置零,在进行快速傅里叶反变换ifft2,得出傅里叶变换去噪后的图像。用Harr小波函数对图像进行分解,通过wrcoef2()函数对变换后的低频部分进行重构,得到小波变换去噪后的图像。实验结果如下图所示:

在这里插入图片描述
在这里插入图片描述
通过观察可以发现,使用不同方法对图像进行去噪处理后,图像都变得模糊了,但同时相应的噪声也被减弱,图像看起来更加平滑。

附录代码如下:

clear;
t=imread('lenna.jpg');
t=rgb2gray(t);%灰度化
[m,n]=size(t);
%%
%高斯噪声
I1=im2double(t);
V=0.01;
I2=imnoise(I1,'gaussian',0,V);

imshow(I2),title('高斯噪声');
%%
%周期噪声
tt=t;
for i=1:m
for j=1:n
tt(i,j)=t(i,j)+40*sin(40*i)+40*sin(40*j);
end
end
figure;
imshow(t),title('原图');
tt=double(tt)/255;%归一化 便于计算
t1=tt;
figure,imshow(tt),title('加入周期噪声后');
%%
%傅里叶变换去噪
I=fft2(tt);
I=fftshift(I);
% figure;
% imshow(log(abs(I)+1),[0 10]);
I(abs(I)<50)=0;
% figure;
% imshow(log(abs(I)+eps),[ ]);
I=ifftshift(I);
I=ifft2(I);
figure;
I=I*255;
subplot(121),
imshow(I,[0 255]),
title('傅里叶去噪周期噪声'); 
I1=fft2(I2);
I1=fftshift(I1);

I1(abs(I1)<50)=0;
I1=ifftshift(I1);
I1=ifft2(I1);
I1=I1*255;
subplot(122),
imshow(I1,[0 255]),
title('傅里叶去噪高斯噪声'); 
%%
%小波变换的频域滤波
%用小波函数haar对x进行2层小波分解  
[c,s]=wavedec2(tt,2,'sym4');  
%提取小波分解中第一层的低频图像,即实现了低通滤波去噪  
a1=wrcoef2('a',c,s,'sym4');  % a1为 double 型数据;
[c,s]=wavedec2(I2,2,'sym4');  
%提取小波分解中第一层的低频图像,即实现了低通滤波去噪  
a2=wrcoef2('a',c,s,'sym4');  % a1为 double 型数据;
%画出去噪后的图像  
figure;
subplot(121),
imshow(a1),
title('小波变换去噪周期噪声图像');   
subplot(122),
imshow(a2),
title('小波变换去噪高斯噪声图像');  
%%
I4=I2;
P4=tt;
for i=2:m-1
    for j=2:n-1
        I4(i,j)=(I2(i,j)*4+I2(i-1,j)*2+I2(i+1,j)*2+I2(i,j-1)*2+I2(i,j+1)*2+I2(i-1,j-1)+I2(i-1,j+1)+I2(i+1,j-1)+I2(i+1,j+1))/16;
        P4(i,j)=(tt(i,j)*4+tt(i-1,j)*2+tt(i+1,j)*2+tt(i,j-1)*2+tt(i,j+1)*2+tt(i-1,j-1)+tt(i-1,j+1)+tt(i+1,j-1)+tt(i+1,j+1))/16;
    end
end    

a=zeros(1,9);
c=zeros(1,9);
I5=I2;
P5=tt;
for i=2:m-1
    for j=2:n-1
        a=[I2(i,j),I2(i-1,j),I2(i+1,j),I2(i,j-1),I2(i,j+1),I2(i-1,j-1),I2(i-1,j+1),I2(i+1,j-1),I2(i+1,j+1)];
        c=[tt(i,j),tt(i-1,j),tt(i+1,j),tt(i,j-1),tt(i,j+1),tt(i-1,j-1),tt(i-1,j+1),tt(i+1,j-1),tt(i+1,j+1)];
        b=sort(a);
        d=sort(c);
        I5(i,j)=b(5);
        P5(i,j)=d(5);
    end
end
figure;
subplot(121),imshow(I4),title('高斯滤波处理高斯噪声');    
subplot(122),imshow(I5),title('中值滤波处理高斯噪声');  
figure;
subplot(121),imshow(P4),title('高斯滤波处理周期噪声');    
subplot(122),imshow(P5),title('中值滤波处理周期噪声');   

实验二

1、编写一个程序,要求实现下列算法:首先将图像分割为8×8的子图像,对每个子图像进行FFT,对每个子图像中的64个系数,按照每个系数的方差排序后,舍去小的变换系数,只保留16个系数,实现4:1的图像压缩。

以256×256的图像为例,将该图像分割成8×8的子图像,经傅里叶变换后,每个子图像有64个傅里叶系数。按照每个系数的方差来排次序。由于图像是实值得,原来得64个复系数中只有32个是有差别的,可以保留高位次得系数,舍弃一些低位次得系数,从而实现数据压缩。使用系数的多少实际上决定了压缩比的大小,如可以保留最高位次的16个系数,来实现4:1的图像压缩。实验结果如下图所示:
在这里插入图片描述
首先将图像裁剪为256×256个像素点,然后调用blkproc函数将图像分割成8×8的子图像并进行快速傅里叶变换,然后调用im2col函数将图像每个8×8块系数排成一列,然后调用sort函数对其进行排序,然后将前48个系数置零,即实现图像的4:1压缩。然后调用col2im函数对其进行重新排列系数块,再调用blkproc函数对子图像块进行FFT反变换获得各子图像的恢复图像。4:1压缩效果不是很明显,分成16×16各像素块进行4:1压缩后可以发现图像变得模糊且粗糙,实现了图像压缩。

附录代码如下:

clear;
f1=imread('lenna.jpg');
f1=double(f1)/255;%uint8转换为【01】内的双精度数
for m=1:256
    for n=1:256
        f(m,n)=f1(m,n);   %将图像转换为256×256个像素点
    end
end
F1=blkproc(f,[8 8],'fft2');%将图像分割成8×8的子图像进行fft快速傅里叶变换
F2=blkproc(f,[16 16],'fft2');%将图像分割成16×16的子图像进行fft快速傅里叶变换
C=im2col(F1,[8 8],'distinct');%将每个8×8块系数排成一列
C1=im2col(F2,[16 16],'distinct');
coe1=C;
coe2=C1;
[Y,Ind]=sort(C);%对系数C按列排序号,(Ind为对应元素在该列由小到大的序号)
[Y1,Ind1]=sort(C1);
rat=4;%设置压缩比
Snum1=64-64/rat;%8×8个系数中保留1/416个系数
Snum2=256-256/rat;
for i=1:n
    coe1(Ind(1:Snum1),i)=0;%令每个8×8块系数中排序靠前的3/4个系数为0,其余保留
    coe2(Ind1(1:Snum2),i)=0;
end
F1=col2im(coe1,[8 8],[256 256],'distinct');%重新排列系数块
F2=col2im(coe2,[16 16],[256 256],'distinct');%重新排列系数块
S1=blkproc(F1,[8 8],'ifft2');%对子图像块进行FFT反变换获得各子图像的恢复图像    
S2=blkproc(F2,[16 16],'ifft2');
figure;
subplot(221);imshow(f);title('原始图像','fontsize',16);
subplot(222);imshow(abs(F1));title('原始图像8×8分块FFT频谱','fontsize',16);
subplot(223);imshow(S1);title('8×8的4:1FFT压缩图像','fontsize',16);
subplot(224);imshow(S2);title('16×16的4:1FFT压缩图像','fontsize',16);

2.小波变换是如何定义的?小波分析的主要优点是什么?
定义:信号分析是为了获得时间和频率之间的相互关系。傅里叶变换提供了有关频率的信息,但有关时间的局部化信息却基本丢失。与傅里叶变换不同,小波变换是通过缩放母小波的宽度来获得信号的频率特征,通过平移母小波来获得信号的时间信息。对母小波的缩放和平移操作是为了计算小波系数,这些小波系数反映了小波和局部信息之间的想过系数。小波变换可以理解为用经过缩放和平移的一系列小波函数代替傅里叶变换的正弦波和预先波进行傅里叶变换的结果。
优点:小波变换的主要优点之一就是提供局部分析与细化的能力。小波分析在时域和频域都具有良好的局部化特性,而且,由于对高频采取逐渐精细的时域或空域步长,从而可以聚焦到分析对象的任意细节,与传统的信号分析技术相比,小波分析能在无明显损失的情况下,对信号进行压缩和去噪。小波变换可以把信号分解为多个具有不同的时间和频率分辨率的信号,从而可在一个变换中同时研究信号的低频和高频信息,因此,用小波变换对图像这种不平稳的复杂信号源进行处理,能有效克服傅里叶分析和其他分析方法的不足。因此,小波变换在图像编码,图像去噪、图像检测和图像复原等方面得到了广泛应用。

3.给定一幅行和列都为2的整数次幂图像,用Haar小波基函数对其进行二维小波变换,试着将最低尺度近似分量置零再反变换,结果是什么?如果把垂直方向的细节分量置零,反变换后结果又是什么呢?试解释一下原因。
首先调用dwt2函数对图像进行一层haar小波变换,分别得到分解后的低频近似信息和水平、竖直和对角方向的高频信息,然后分别将低频信息部分和垂直方向细节分量质量,调用idwt2进行反变换得到结果图如下图所示:

在这里插入图片描述
通过观察实验结果可以发现,将最低尺度近似分量置零后,图像原本的低频部分消失不见,只留下了高频的细节部分。而将垂直方向细节分量置零后,图像没有太大的变化,原因则是因为人眼对高频部分的信息识别很弱,而对低频信息识别性很强,这也正是进行图像压缩的原理,即去除掉图像的高频部分,只留下低频,图像不会有明显的失真。

附录代码如下:

clc;
clear;
% 读取原始图像
X = rgb2gray(imread('lenna.jpg'));
% 进行二维小波分解
[C, S] = wavedec2(X, 2, 'haar');
% 获得分解以后的低频近似信息

L1 = appcoef2(C, S, 'haar', 2);
L2 = appcoef2(C, S, 'haar', 1);
% 分别获得各层级的高频细节信息
[H2, V2, D2] = detcoef2('all', C, S, 2);
[H1, V1, D1] = detcoef2('all', C, S, 1);
% 去掉第一层的高频信息(替换成0),然后进行小波重建
% 注意这里乘以3是有HVD三种高频信息
D = [C(1: end - 3*size(H1, 1)*size(H1, 2)), zeros(1, 3*size(H1, 1)*size(H1, 2))];
CD1 = waverec2(D, S, 'haar');
% 去掉第一和第二层的高频信息,然后进行小波重建
D = [C(1: end - 3*size(H1, 1)*size(H1, 2) - 3*size(H2, 1)*size(H2, 2)), ...
    zeros(1, 3*size(H1, 1)*size(H1, 2) + 3*size(H2, 1)*size(H2, 2))];
CD2 = waverec2(D, S, 'haar');
%按照分解层级将分解系数排列拼接为一副图像
DD3 = [L2, H1; V1, D1];
% 结果显示
subplot(2, 2, 1), imshow(X, []), title('原始图像');
subplot(2, 2, 2), imshow(DD3, []), title('小波分解系数');
subplot(2, 2, 3), imshow(CD1, []), title('去掉第一层高频信息后');
subplot(2, 2, 4), imshow(CD2, []), title('去掉第二层高频信息后');
%%
A = X;
A = A(1:256,1:256);
[m, n]=size(A);
[LL, LH, HL, HH]=dwt2(A,'haar'); 
A=[LL LH;HL HH];      %一层分解
temp = zeros(128);
result1 = idwt2(temp,LH,HL,HH,'Haar');
figure;
subplot(1, 2, 1),
imshow(result1,[]),
title('将最低尺度近似分量置零');
result2= idwt2(LL,LH,temp,HH,'Haar');
subplot(1, 2, 2),
imshow(result2,[]),
title('将垂直方向细节分量置零');

4.基于小波变换对图像进行不同压缩比的压缩。在同压缩比情况下,对于基于小波变换和基于傅里叶变换的压缩结果,比较二者保留原图像能量百分比情况。
通过分解后的图像的高频信息,将其替换成零后进行小波重建,即可完成小波变换对图像的压缩,实验结果如下图所示:
在这里插入图片描述

而通过傅里叶变换后得到的压缩图像如下图所示:

在这里插入图片描述

附录代码如下:

clc;
clear;
%小波变换的频域滤波
I=imread('lenna.jpg');%读取图像信息
I=rgb2gray(I);
%I=double(I)/255;
subplot(2,2,1); imshow(I); 
title('原图');     
%进行小波分解
[c,s]=wavedec2(I,2,'haar',1);  
% 获得分解以后的低频近似信息
ca1=appcoef2(c,s,'haar',1);
%分别获得各层级的高频细节信息
ch1=detcoef2('h',c,s,1);      %水平方向
cv1=detcoef2('v',c,s,1);      %垂直方向
cd1=detcoef2('d',c,s,1);      %对角方向
%各频率成份重构
a1=wrcoef2('a',c,s,'haar',1);
h1=wrcoef2('h',c,s,'haar',1);
v1=wrcoef2('v',c,s,'haar',1);
d1=wrcoef2('d',c,s,'haar',1);
ca1=[a1,h1;v1,d1];
subplot(2,2,2); imshow(ca1,[]); 
title('第一次小波变换后的低频和高频信息');     
%%
ca1=appcoef2(c,s,'haar',1);
cal=wcodemat(ca1,440,'mat',0);
%改变图像高度并显示
ca1=0.5*ca1;
subplot(223);imshow(cal,[]);
title('第一次压缩后的图像');
%保留小波分解第二层低频信息进行压缩
ca2=appcoef2(c,s,'haar',2);
%首先对第二层信息进行量化编码
ca2=wcodemat(ca2,440,'mat',0);
%改变图像高度并显示
ca2=0.25*ca2;
subplot(224);imshow(ca2,[]);
title('第二次压缩后的图像');  

评论 2
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值