中南大学数学院数字图像处理报告1.docx

中南大学数学院数字图像处理报告1.docx

《中南大学数学院数字图像处理报告1.docx》由会员分享，可在线阅读，更多相关《中南大学数学院数字图像处理报告1.docx（20页珍藏版）》请在冰豆网上搜索。

中南大学数学院数字图像处理报告1

数据图像处理报告

学生姓名

学号

学院数学科学与计算技术学院

专业班级

本科生院制

图像处理实验

（一）图像的幅度谱图及频域滤波

快速傅氏变换（FFT），是离散傅氏变换的快速算法，它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的，能克服时间域与频率域之间相互转换的计算障碍，在光谱、大气波谱分析、数字信号处理等方面有广泛应用。

产生下图所示图象f1（x,y）（128×128大小，暗处=0，亮处=255）

（1）同屏显示原图f1和FFT（f1）的幅度谱图

f1=ones（128）*0.;

fori=25:

100

forj=60:

70

f1（i,j）=255;

end

end

figure

（1）;

subplot（2,3,1）;

imshow（f1）;

title（'f1原始图'）;

f3=fft2（f1）;

subplot（2,3,2）;

imshow（abs（f3））;

title（'FFT（f1）幅谱图'）;

（2）令f2（x,y）=

，重复以上过程,并简述二者幅度谱异同的理由

fori=1:

128

forj=1:

128

f2（i,j）=（（-1）^（i+j））*f1（i,j）;

end

end

subplot（2,3,3）;

imshow（f2）;

title（'f2的原始图'）;

f=fft2（f2）;

subplot（2,3,4）;

imshow（abs（f））;

title（'FFT（f2）幅谱图'）;

理由：

因为f2是f1通过平移得到，而且正好平移到f1中部，若再做快速傅里叶

变换，则f2相比于f1,变换后主要能量（低频分量）集中在频率平面的中心。

（3）对该图进行动态范围压缩以增强图象

利用DCT变化进行动态压缩

figure;

I=im2double（f1）;

T=dctmtx（8）;//求8*8模板的dct变换矩阵

B=blkproc（I,[88],'P1*x*P2',T,T'）;//施DCT变换B=TIT’

mask=[11110000

11100000

11000000

10000000

00000000

00000000

00000000

00000000];

B2=blkproc（B,[88],'P1.*x',mask）;//只取10个有效点

I2=blkproc（B2,[88],'P1*x*P2',T',T）;//IDCT变换

subplot（1,2,1）;

imshow（f1）;title（'原图像'）;

subplot（1,2,2）;

imshow（I2）;title（'变化后原图像'）;//显示原图与idct变换后的图

运行结果

（4）将f2（x,y）顺时针旋转45度得到f3（x,y），试显示FFT（f3）的幅度谱，并与FFT（f2）的幅度谱进行比较

f4=imrotate（f1,315,'bilinear','crop'）;

//双线性插值法旋转图像45度,并剪切图像,使其和原图像大小一致

subplot（2,3,5）

imshow（f4）

title（'顺时针旋转图像45^{0}'）

f5=fft2（f4）;

subplot（2,3,6）;

imshow（abs（f5））;

title（'旋转图像后的幅谱图'）;

以上程序运行的结果显示如下：

图像处理实验

（二）数字图像增强

灰度变换是图像增强的一种重要手段，使图像对比度扩展，图像更加清晰，特征更加明显。

灰度级的直方图给出了一幅图像概貌的描述，通过修改灰度直方图来得到图像增强。

（1）用直接灰度变换方法进行对比度增强，设计分段线性变换函数，做增强处理

由于本文采用的是张彩色图片，首先需通过算法将其变成黑白

a=imread（'trees.tif'）;//读入图片

ifndims（a）==3//如果图片是3维（彩图）

a=rgb2gray（a）;//转成灰图

end//结束

imshow（a）;//展示图片

接着再对该图片进行灰度变换原图像f（m,n）的灰度范围[a,b]线形变换为图像g（m,n），灰度范围[c,d]

公式：

g（m,n）=c+（d-c）*f（m,n）/（b-a）

figure;

b=imadjust（a,[0.30.7],[0,1],1）;

subplot（1,2,1）,imshow（a）;

title（'原图像'）;

subplot（1,2,2）,imshow（b）;

title（'灰度变换后图像'）;

运行结果如下：

（2）统计原图的灰度直方图，并利用直方图均衡方法进行图象增强

figure

subplot（1,2,1）,imhist（a,64）;

title（'原图像直方图'）;

c=histeq（a）;

subplot（1,2,2）,imshow（c）;

title（'直方图均衡化所得图像'）;

运行结果

（3）利用Matlab函数IMNOISE（），在原图上分别叠加高斯噪声和椒盐噪声，对比线性平滑滤波器和非线性平滑滤波器（中值滤波）的性能；

figure

（1）

d=imnoise（a,'salt&pepper',0.02）;//添加椒盐噪声

subplot（1,2,1）,imshow（d）;

title（'添加椒盐噪声图像'）;

e=imnoise（a,'gaussian',0,0.01）;//添加高斯噪声

subplot（1,2,2）,imshow（e）;

title（'添加高斯噪声图像'）;

运行结果

接着用领域平均域波和中值域波分别对添加了椒盐噪声的图像进行平滑，并给出图像结果进行对比

figure

（2）

K1=filter2（fspecial（'average',3）,d）/255;

%应用3×3邻域窗口法

subplot（1,2,1）,imshow（K1）;

title（'添加椒盐噪声后的3×3窗的邻域平均滤波图像'）;

K=medfilt2（d）;

%缺省3×3的邻域窗的中值滤波

subplot（1,2,2）,imshow（K）;

title（'添加椒盐噪声后的缺省3×3邻域窗的中值滤波图像'）;

运行结果

接着用领域平均域波和中值域波分别对添加了高斯噪声的图像进行平滑，并给出图像结果进行对比

figure（3）

K2=filter2（fspecial（'average',3）,e）/255;

%应用3×3邻域窗口法

subplot（1,2,1）,imshow（K2）;

title（'添加高斯噪声后的3×3窗的邻域平均滤波图像'）;

K3=medfilt2（e）;

%缺省3×3的邻域窗的中值滤波

subplot（1,2,2）,imshow（K3）;

title（'添加高斯噪声后的缺省3×3邻域窗的中值滤波图像'）;

运行结果

（4）利用线性锐化器和非线性锐化滤波器增强原图象的边缘信息

Laplacian算子

figure

b=fspecial（'laplacian'）;

c=filter2（b,a）;

subplot（1,2,1）;

subimage（[0,500],[0,500],c）;

title（'线性锐化器增强边缘'）;

Sobel算子

figure

b=fspecial（'sobel'）;

c=filter2（b,a）;

subplot（1,2,2）;

subimage（[0,500],[0,500],c）;

title（'非线性锐化器增强边缘'）;

运行结果

由实验所得的结论

从实验中可以知道，对于椒盐噪声，中值域波比领域平均要好，这是因为椒盐噪声是复制近似相等但随机分布在不同的位置上，图像中又干净点也有污染点。

中值域波是选择适当的点来代替污染点的值，所以处理效果好。

因为噪声的均值不为0，所以均值域波不能很好地去除噪声。

图像处理实验（三）平滑、锐化处理的模板运算

平滑是增强图像的方法之一，其主要是抑制或者消除图像的噪音，为其他处理做准备。

锐化是增强图像的方法，它主要是提取图像的边缘、轮廓、线条等信息，供进一步识别，加重图像轮廓克服降质，以达到更好的视觉效果。

（1）利用二个低通邻域平均模板（3×3和9×9）对一幅图象进行平滑，验证模板尺寸对图象的模糊效果的影响。

a=imread（'trees.tif'）;

K1=filter2（fspecial（'average',3）,a）/255;

%应用3×3邻域窗口法

K2=filter2（fspecial（'average',9）,a）/255;

%应用9×9邻域窗口法

figure

（1）;

subplot（1,3,1）,imshow（a）

title（'原图像'）

subplot（1,3,2）,imshow（K1）

title（'3×3窗的邻域平均滤波图像'）

subplot（1,3,3）,imshow（K2）

title（'9×9窗的邻域平均滤波图像'）

运行结果

由图所得的结论：

由图可见模板尺寸影响图像的模糊效果，模板的尺寸越大，图像模糊层度越高，这是因为，用领域平均域波在一定程度上抑制噪声，但是领域平均法的平均作用会引起模糊现象，模糊程度与领域半径成正比

（2）利用一个低通模板对一幅有噪图象（GAUSS白噪声）进行滤波，检验两种滤波模板（分别使用一个5×5的线性邻域平均模板和一个非线性模板：

3×5中值滤波器）对噪声的滤波效果。

figure

a=imread（'trees.tif'）;

b=imnoise（a,'gaussian',0,0.01）;

c=filter2（fspecial（'average',5）,b）/255;

subplot（1,2,1）,imshow（c）;

title（'5×5线性邻域平均模板'）;

d=medfilt2（b,[3,5]）;

subplot（1,2,2）,imshow（d）;

title（'3×5中值滤波器'）;

运行结果

（4）选择一个经过低通滤波器滤波的模糊图象，利用sobel和prewitt水平边缘增强高通滤波器（模板）对其进行高通滤波图象边缘增强，验证模板的滤波效果。

经过低通邻域平均模板（9×9）平滑后的图像，分别利用sobel和prewitt算子进行边缘增强。

领域平均模板（9×9）平滑图像

a=imread（'trees.tif'）;

b=filter2（fspecial（'average',9）,a）/255;%应用9×9邻域窗口法

subplot（1,1,1）,imshow（b）

title（'9×9窗的邻域平均滤波图像'）;

运行结果

Sobel算子进行边缘增强

figure

b=fspecial（'sobel'）;

c=filter2（b,a）;

subplot（1,1,1）;

subimage（[0,500],[0,500],c）;

title（'sobel增强边缘'）;

运行结果

prewitt算子进行边缘增强

figure

b=fspecial（'prewitt'）;

c=filter2（b,a）;

subplot（1,1,1）;

subimage（[0,500],[0,500],c）;

title（'prewitt增强边缘'）;

运行结果

（4）选择一幅灰度图象分别利用一阶Sobel算子和二阶Laplacian算子对其进行边缘检测，验证检测效果。

figure;

a=imread（'trees.tif'）;

S0=[-1-2-1;000;121];

b=imfilter（a,S0,'replicate'）;

subplot（1,2,1）;imshow（b）;

title（'sobel水平图像边缘检测'）;

S1=[0-10;-14-1;0-10];

c=imfilter（a,S1,'replicate'）;

subplot（1,2,2）;imshow（c）;

title（'laplacian边缘检测'）;

运行结果

图像处理实验（四）图像压缩与编码

实现基本JPEG的压缩和编码

（1）通过DCT变换去除数据冗余

（2）使用量化表对DCT系数进行量化

（3）对量化后的系数进行Huffman编码

I=imread（'trees.tif'）;

lightable=[1611101624405161;1212141926586055;

   1413162440576956;1417222951878062;

   182237566810910377;243555648110411392;

   49647887103121120101;7292959811210010399];//亮度标准化

colortable=[1718244799999999;1821266699999999;

    2426569999999999;4766999999999999;

    9999999999999999;9999999999999999;

    9999999999999999;9999999999999999];//色度标准化

I1=I（:

:

1）;I2=I（:

:

2）;//三维数组，源图像第一页的数据给I1,第二页数据给I2

[mn]=size（I1）;//将I1的行数赋给m,将I1的列数赋给n

t1=8;ti1=1;

while（t1

 t1=t1+8;ti1=ti1+1;

end

 t2=8;ti2=1;

while（t2

 t2=t2+8;ti2=ti2+1;

end

 times=0;

 fork=0:

ti1-2

forj=0:

ti2-2

  dct8x8（I1（k*8+1:

k*8+8,j*8+1:

j*8+8）,lightable,times*64+1）;

  dct8x8（I2（k*8+1:

k*8+8,j*8+1:

j*8+8）,colortable,times*64+1）;

  times=times+1;

end

 block（I2（k*8+1:

k*8+8,j*8+1:

t2）,[88],'dctmtx（8）'）;//二值掩模，用来压缩DCT系数，只留下DCT系数中左上

end

forj=0:

ti2-2

  dct8x8（I1（k*8+1:

t1,j*8+1:

j*8+8）,times*64+1）;

times=times+1;

end

dct8x8（I1（k*8+1:

t1,j*8+1:

t2）,times*64+1）;

2．functiondct8x8（I,m,s）%定义DCT量化子程序

T=inline（'dctmtx（8）'）;

y=blkproc（I,[88],T）;

y=round（y./m）;

p=1;te=1;

while（p<=64）

forq=1:

te

y1（s+p）=y（te-q+1,q）;p=p+1;

end

forq=te:

-1:

1

y1（s+p）=y（te-q+1,q）;p=p+1;

end

end

f=haffman（y1）;

c（s:

s+64,1）=f（:

1）;c（s:

s+64,2）=f（:

2）;c（s:

s+64,3）=f（:

3）

3．functionc=haffman（I） %定义Huffman编码子程序

[m,n]=size（I）;

p1=1;s=m*n;

fork=1:

m

forh=1:

n

 f=0;

forb=1:

p1-1

if（c（b,1）==I（k,h））

f=1;break;

End

end

if（f==0）

c（p1,1）=I（k,h）;p1=p1+1;

end

end

end

forg=1:

p1-1

 p（g）=0;c（g,2）=0;

fork=1:

m

forh=1:

n

if（c（g,1）==I（k,h））

p（g）=p（g）+1;

end

end

End

p（g）=p（g）/s;

end

pn=0;po=1;

while

（1）

if（pn>=1.0）break;

else

[pmp2]=min（p（1:

p1-1））;p（p2）=1.1;

[pm2,p3]=min（p（1:

p1-1））;p（p3）=1.1;

pn=pm+pm2;p（p1）=pn;

tree（po,1）=p2;tree（po,2）=p3;

po=po+1;p1=p1+1;

end

end

fork=1:

po-1

tt=k;m1=1;

if（or（tree（k,1）<9,tree（k,2）<9））

if（tree（k,1）<9）

c（tree（k,1）,2）=c（tree（k,1）,2）+m1;

m2=1;

while（tt

m1=m1*2;

forh=tt:

po-1

if（tree（h,1）==tt+g）

c（tree（k,1）,2）=c（tree（k,1）,2）+m1;

m2=m2+1;tt=h;break;

elseif（tree（h,2）==tt+g）

m2=m2+1;tt=h;break;

end

end

end

c（tree（k,1）,3）=m2;

end

tt=k;m1=1;

if（tree（k,2）<9）

m2=1;

while（tt

m1=m1*2;

forh=tt:

po-1

if（tree（h,1）==tt+g）

c（tree（k,2）,2）=c（tree（k,2）,2）+m1;

m2=m2+1;tt=h;break;

elseif（tree（l,2）==tt+g）

m2=m2+1,tt=h;break;

end

end

end

c（tree（k,2）,3）=m2;

end

end

end