《数字图像处理》实验指导书(2011年用过)

实验一、图象的基本操作（2学时）

实验二、图像的傅立叶变换（2学时）

实验三、图像增强（2学时）

实验四、图像压缩（2学时）

实验五、图像融合（选做）

实验一 图像的基本操作

1．实验目的 

熟悉常用的图像文件格式与格式转换； 

熟悉图像矩阵的显示方法（灰度、索引、黑白、彩色）； 

熟悉图像矩阵的格式转换使学生熟悉图象文件格式、图象的载入、彩色图象的构成、图象的显示等。

2．实验内容 

练习图像读写命令 imread和imwrite并进行图像文件格式间的转换。特别是索引图像与1，4，8，16 比特图像的存储与转换。

Image file I/O. 

    imread     - Read image file. 

    imwrite    - Write image file. 

Image display. 

    colorbar   - Display colorbar. 

    getimage   - Get image data from axes. 

    image      - Create and display image 

    imagesc    - Scale data and display as image. 

    immovie    - Make movie from multiframe 

image. 

    imshow     - Display image. 

subimage- Display multiple images in single figure. 

truesize - Adjust display size of image. 

warp   - Display image as texture-mapped surface. 

zoom  - Zoom in and out of image or 2-D plot.

3． 实验步骤 

a.  从硬盘加载 cameraman.tif  图象 (用函数 imread). 

b.  在图象显示窗口显示图象 (using function image or imshow). 

c.  显示彩色图象的 3 基色图象。

参考程序： 

figure('Name', 'Football'); 

imshow(football); 

figure('Name', 'Football color planes');

subplot(2,2,1); imshow(football); 

redfootball = football; 

redfootball(:, :, 2:3) = 0; 

subplot(2,2,2); imshow(redfootball); 

greenfootball = football;

greenfootball(:, :, 1) = 0; 

greenfootball(:, :, 3) = 0; 

subplot(2,2,3); imshow(greenfootball); 

bluefootball = football; 

bluefootball(:, :, 1:2) = 0; 

subplot(2,2,4); imshow(bluefootball);

实验二   图象的基本操作与变换

一、实验目的 

通过实验了解并掌握图象的基本操作、图象的傅立叶变换及频谱、傅立叶变换可分离性

和哈达吗变换。 

二、实验内容 

（1） 进行图象的基本算术操作（+、-、*） 

  载入图象‘’,加入躁声,观察图象的变化. 

（2） 进行二值图象的逻辑运算（and、or、not） 

自编程序生成 2 个二值图象(256x256),进行逻辑操作. 

(3)  设简单的二值图象为: 

imag=zeros(256,256); 

imag(100:180, 100:150)=1; 

编制程序用两次一维傅立叶变换计算其二维傅立叶变换.并观察其幅值频谱和相位频谱. 

(fft,  abs, imshow, angle, fftshift) 

(4)  载入图象‘’,用 fft2 和 ifft2 观察图象的频谱,并进行逆傅立叶变换,观察还原

后的图象与原图象的差异. 

（5）观察图象的旋转和其傅立叶变换。

三  实验要求 

（1） 每位同学完成实验内容 

（2） 实验报告包括程序代码、运行结果和体会。

（3） 实验报告中变换前后的图像须打印

1,图像的基本操作

imshow(a) 

title ('original checker image'); 

figure(2) 

imshow(b) 

title('original castles image'); 

BW = im2bw(I, 0.0005);   % binary conversion 

BW2 = im2bw(I2, 0.0005); 

%compute A and B 

figure(3) 

A = and(BW, BW2);        % and operation 

imshow(A) 

title('output images (checker and castles) '); 

%compute A xor B 

figure(4) 

B = xor (BW, BW2); 

imshow(B) 

title('output images (checker xor castles) ');

%compute not_A xnor not_B 

figure(5) 

NI = not (I); 

NI2 = not (I2); 

c1 = and (I, I2); 

c2 = and (NI, NI2); 

C = or (c1 ,c2); 

imshow(C) 

title('output images (not checker xnor not castles) '); 

%compute E = max(C,D) 

figure(6) 

imshow(I3) 

title('original lenna image');

figure(7)

imshow(I4) 

figure(8) 

D = max(I3, I4); 

imshow(D) 

title('output image of max(lenna, lv40)'); 

%compute difference using imsubtract and imabsdiff 

figure(9) 

imshow(I5) 

title('seq1 frame 0001 original image');

figure(10) 

imshow(I6) 

title('seq1 frame 0002 original image'); 

figure(11) 

E = imsubtract(I5, I6); 

imshow(E) 

figure(12) 

F = imabsdiff(I5, I6); 

imshow(F) 

%Connected-Components Analysis  

board = checkerboard(1)>0.5;

four_conn = bwlabel(board, 4) 

eight_conn = bwlabel(board, 8)

clear; 

clc; 

close all; 

pack; 

% read images into memory 

%no need to convert these images to binary -- they are already binary 

%but if they were not they would be thresholded to binary here 

%do the logical operations 

%1.a.:  AND operation on two binary images 

and_img = and(img_checkers, img_castles); 

figure; 

imshow(and_img); 

%1.b.: XOR operation on two binary images 

xor_img = xor(img_checkers, img_castles); 

figure;

imshow(xor_img); 

%1.c.: NOT...XNOR...NOT operation on two binary images: 

% A' XNOR B' = (A' XOR B')'  

not_img = not(xor(not(img_checkers), not(img_castles))); 

figure; 

imshow(not_img); 

%1.d.: Find max value for two images

img_max = max(img_lenna, img_lv); 

figure; 

imshow(img_max); 

2 傅立叶变换

I=imread(‘原图像名读入原图像文件

imshow(I);        %显示原图像

fftI=fft2(I);        %二维离散傅立叶变换

sfftI=fftshift(fftI);   %直流分量移到频谱中心

RR=real(sfftI);     %取傅立叶变换的实部

II=imag(sfftI);     %取傅立叶变换的虚部

A=sqrt(RR.^2+II.^2);%计算频谱幅值

A=(A-min(min(A)))/(max(max(A))-min(min(A)))*225;

                 %归一化

figure;           %设定窗口

imshow(A);       %显示原图像的频谱

//选做  

%1.e.: Find difference image between two animated motion sequence frames 

figure; 

imshow(img_seq1); 

title('Frame One'); 

figure; 

imshow(img_seq2); 

title('Frame Two');

img_diff= imabsdiff(img_seq2, img_seq1); 

figure; 

imshow(img_diff); 

title('Difference between Frame One and Frame Two--Shows Motion'); 

%Connected Components Analysis 

%generate a test image  

img_con_test = checkerboard(1)>0.5

figure; 

imshow(img_con_test); 

title('Checkerboard to Show Connectivity');  

%compute matrices of 4-connectivity and 8-connectivity in the test image 

[four_con_matrix, val_4] = bwlabeln(img_con_test, 4)  

figure; 

imshow(four_con_matrix); 

title('Image of 4-Connectivity Matrix');

[eight_con_matrix, val_8] = bwlabeln(img_con_test, 8) 

figure; 

imshow(eight_con_matrix); 

title('Image of 8-Connectivity Matrix');

实验三 图像增强（2 学时）

1．   实验目的 

1掌握灰度直方图的概念及其计算方法；

2熟练掌握直力图均衡化和直方图规定化的计算过程；

3熟练掌握空域滤波中常用的平滑和锐化滤波器；

4掌握色彩直方图的概念和计算方法

5利用MATLAB程序进行图像增强。

2． 实验内容 

练习图像增强的 Matlab命令。 

熟悉下列模块函数 

Image enhancement. 

    histeq      - Read image file. 

    imadust    - Adust image intensity values or colormap. 

Image noising. 

imnoise     - Add noise to an image. 

Image filtering 

    medfilt2     - Perform 2-D median filtering. 

    ordfilt2      - Perform 2-D order-statistic filtering. 

    wiener2      - Perform 2-D adaptive noise-removal filtering. 

3． 实验步骤

实验一

a.  从硬盘加载 cameraman.tif  图象 (using function imread). 

b.  显示图象. 

c.  显示图象的直方图 (using function imhist). 

d.  用直方图均衡方法进行图象增强. 

e.  对处理后的图象显示其直方图. 

f.  比较图象的质量并且进行讨论. 

g.  对图象加入躁声，比较躁声均值和方差的影响。. 

h.  去除图象躁声（用函数 medfilt2, ordfilt2  和 wiener2.比较结果。

imag=im2double(imag); 

subplot(221); imshow(imag); 

subplot(222); imhist(imag); 

imag1=histeq(imag); 

subplot(223); imshow(imag1); 

subplot(224); imhist(imag1); 

figure(2) 

a=min(min(imag)); 

b=max(max(imag)); 

J = imadjust(imag,[a b] ,[0.3 0.7],1.7); 

subplot(121); imshow(imag), 

subplot(122); imshow(J) 

imag3=imnoise(imag,'salt', 0.02); 

figure; imshow(imag3) 

实验二 下面给出直方图均衡化增强图像对比度的MATLAB程序：

I=imread(‘原图像名读入原图像

J=histeq(I); %对原图像进行直方图均衡化处理

Imshow(I); %显示原图像

Title(‘原图像’);  %给原图像加标题名

Figure；imshow(J);  %对原图像进行屏幕控制；显示直方图均衡化后的图像

Title(‘直方图均衡化后的图像’) ; %给直方图均衡化后的图像加标题名

Figure; subplot(1,2,1) ;%对直方图均衡化后的图像进行屏幕控制；作一幅子图作为并排两幅图的第1幅图

Imhist(I,);  %将原图像直方图显示为级灰度

Title(‘原图像直方图’) ; %给原图像直方图加标题名

Subplot(1,2,2);  %作第2幅子图

Imhist(J,) ; %将均衡化后图像的直方图显示为级灰度

Title(‘均衡变换后的直方图’) ; %给均衡化后图像直方图加标题名

处理后的图像直方图分布更均匀了，图像在每个灰度级上都有像素点。从处理前后的图像可以看出，许多在原始图像中看不清楚的细节在直方图均衡化处理后所得到的图像中都变得十分清晰。

一、思考题

1． 直方图是什么概念？它反映了图像的什么信息？

　　2． 直方图均衡化是什么意思？它的主要用途是什么？

实验四、图像压缩

一、实验目的

1．理解有损压缩和无损压缩的概念；

2．理解图像压缩的主要原则和目的；　

3．了解几种常用的图像压缩编码方式。

4．利用MATLAB程序进行图像压缩。

二、实验原理

1.图像压缩原理

图像压缩主要目的是为了节省存储空间，增加传输速度。图像压缩的理想标准是信息丢失最少，压缩比例最大。不损失图像质量的压缩称为无损压缩，无损压缩不可能达到很高的压缩比；损失图像质量的压缩称为有损压缩，高的压缩比是以牺牲图像质量为代价的。压缩的实现方法是对图像重新进行编码，希望用更少的数据表示图像。

信息的冗余量有许多种，如空间冗余，时间冗余，结构冗余，知识冗余，视觉冗余等，数据压缩实质上是减少这些冗余量。高效编码的主要方法是尽可能去除图像中的冗余成分，从而以最小的码元包含最大的图像信息。

编码压缩方法有许多种，从不同的角度出发有不同的分类方法，从信息论角度出发可分为两大类。

（1）.冗余度压缩方法，也称无损压缩、信息保持编码或嫡编码。具体说就是解码图像和压缩编码前的图像严格相同，没有失真，从数学上讲是一种可逆运算。

（2）信息量压缩方法，也称有损压缩、失真度编码或烟压缩编码。也就是说解码图像和原始图像是有差别的，允许有一定的失真。

应用在多媒体中的图像压缩编码方法，从压缩编码算法原理上可以分为以下3类：

（1）无损压缩编码种类

哈夫曼（Huffman）编码，算术编码，行程（RLE）编码，Lempel zev编码。

（2）有损压缩编码种类

预测编码，DPCM，运动补偿；

频率域方法：正交变换编码(如DCT)，子带编码；

空间域方法：统计分块编码；

模型方法：分形编码，模型基编码；

基于重要性：滤波，子采样，比特分配，向量量化；

（3）混合编码。

有JBIG，H261，JPEG，MPEG等技术标准。

本实验主要利用MATLAB程序进行离散余弦变换（DCT）压缩和行程编码（Run Length Encoding， RLE）。

1)离散余弦变换(DCT)图像压缩原理

离散余弦变换DCT在图像压缩中具有广泛的应用，它是JPEG、MPEG等数据压缩标准的重要数学基础。

    和相同图像质量的其他常用文件格式(如GIF(可交换的图像文件格式)，TIFF(标签图像文件格式)，PCX(图形文件格式))相比，JPEG是目前静态图像中压缩比最高的。JPEG比其他几种压缩比要高得多，而图像质量都差不多(JPEG处理的图像只有真彩图和灰度图)。正是由于其高压缩比，使得JPEG被广泛地应用于多媒体和网络程序中。JPEG有几种模式，其中最常用的是基于DCT变换的顺序型模式，又称为基本系统(Baseline)。

    用DCT压缩图像的过程为：

    (1)首先将输入图像分解为8×8或16×16的块，然后对每个子块进行二维DCT

变换。

    (2)将变换后得到的量化的DCT系数进行编码和传送，形成压缩后的图像格

式。

用DCT解压的过程为：

    (1)对每个8×8或16×16块进行二维DCT反变换。

    (2)将反变换的矩阵的块合成一个单一的图像。

    余弦变换具有把高度相关数据能量集中的趋势，DCT变换后矩阵的能量集中在矩阵的左上角，右下的大多数的DCT系数值非常接近于0。对于通常的图像来说，舍弃这些接近于0的DCT的系数值，并不会对重构图像的画面质量带来显著的下降。所以，利用DCT变换进行图像压缩可以节约大量的存储空间。压缩应该在最合理地近似原图像的情况下使用最少的系数。使用系数的多少也决定了压缩比的大小。

在压缩过程的第2步中，可以合理地舍弃一些系数，从而得到压缩的目的。在压缩过程的第2步，还可以采用RLE和Huffman编码来进一步压缩。

2)行程编码（RLE）原理：

　　例如如下这幅 的二值图像，

　　如果采用行程编码可以按如下格式保存

　　其中10和8表示图像的宽和高。在这个小例子中行程编码并没有起到压缩图像的作用。这是由于这个图的尺寸过小，当图像尺寸较大时行程编码还是不错的无损压缩方法。对于灰度图像和二值图像，用行程编码—般都有很高的压缩率。行程编码方法实现起来很容易，对于具有长重复值的串的压缩编码很有效，例如：对于有大面积的阴影或颜色相同的图像，使用这种方法压缩效果很好。很多位图文件格式都采用行程编码，如TIFF，PCX，GEM，BMP等。

3. 图像压缩编码的MATLAB程序语句

1）利用DCT变换进行图像压缩的MATLAB程序

RGB = imread('原图像名.tif');

I = rgb2gray(RGB);

J = dct2(I);

imshow(log(abs(J)),[]), colormap(jet()), colorbar

J(abs(J) < 10) = 0;

K = idct2(J);

figure,imshow(I)

figure,imshow(K,[0 255])

2）利用离散余弦变换进行JPEG图像压缩

I=imread(‘原图像名.tif’);  %读入原图像；

I=im2double(I);  %将原图像转为双精度数据类型；

T=dctmtx(8);  %产生二维DCT变换矩阵

B=blkproc(I,[8 8],’P1*x*P2’,T,T’);  %计算二维DCT，矩阵T及其转置T’是DCT函数P1*x*P2的参数

Mask=[ 1 1 1 1 0 0 0 0

             1 1 1 0 0 0 0 0

             1 1 0 0 0 0 0 0

             1 0 0 0 0 0 0 0

             0 0 0 0 0 0 0 0

             0 0 0 0 0 0 0 0

             0 0 0 0 0 0 0 0

                        0 0 0 0 0 0 0 0];      %二值掩膜，用来压缩DCT系数，只留下DCT系数中左上角的10个

      B2=blkproc(B,[8 8],’ P1.*x’,Mask);   %只保留DCT变换的10个系数

      I2= blkproc(B2,[8,8],’P1*x*P2’,T’,T); %逆DCT，重构图像

      Subplot(1,2,1);

        Imshow(I);title(‘原图像’)； %显示原图像

        Subplot(1,2,2);  

        Imshow(I2);title(‘压缩图像’)；%显示压缩后的图像。对比原始图像和压缩后的图像，虽然舍弃了85%的DCT系数，但图像仍然清晰（当然有一些质量损失）

3）利用行程编码（RLE）进行图像压缩

I()=checkerboard(10,2);  %调入原图像

[m n]=size(I);

J=[];

for i=1:m

    value=I(i,1);

    num=1;

    for j=2:n

        if I(i,j)==value

            num=num+1;

        else

            J=[J num value];

            num=1;

            value=I(i,j);

        end

end

I=[J num value 0 0];    %添加的行判断位 0 0

end

disp(‘原图像大小：’)

whos(‘I’);

disp(‘压缩图像大小：’)

whos(‘J’);

disp(‘图像的压缩比：’)

disp(m*n/length(J))

三、实验步骤

1打开计算机，启动MATLAB程序；

2调入“实验一”中获取的数字图像，并进行数据的行程（RLE）编码压缩处理；

3 将原图像在Photoshop软件中打开，分别以不同的位图文件格式进行“另保存”，比较它们的数据量。

4记录和整理实验报告

四、思考题

1．图像中哪些信息是主要的，哪些信息是次要的？

　　2．简述离散余弦变换（DCT）和行程编码（RLE）的原理。

实验六 图像融合（选作）

一、实验目的

1熟悉图像融合的意义和用途，理解图像融合的原理；

2掌握图像融合的一般方法；

3掌握运用MATLAB软件进行图像融合的操作。

二、实验原理

图像融合的目的把来自多传感器的数据互补信息合并成一幅新的图像,以改善图像的质量。图像融合最简单的理解就是两个（或多个）图像间的相加运算。这一技术广泛应用于多频谱图像理解和医学图像处理等领域。主要分为空域和频域相加。

本实验主要应用MATLAB软件进行两幅图像的融合。方法有：

1图像直接融合；

2图像傅立叶变换融合；

3图像小波变换融合。

图像融合的MATLAB程序如下：

1）调入、显示两幅图像的程序语句

load tartan;

X1=X;map1=map;

Load sinsin;

X2=X;map2=map;  %打开图像

Subplot(1 2 1)

Image(X1),colormap(map1);

Title(‘图像map1’)

Subplot(1 2 2)

Image(X2),colormap(map2);

Title(‘图像map2’)  %显示两幅图像

2）两幅图像直接融合的程序语句

figure,subplot(1 3 1)

image((X1+X2)/2),colormap(map2); %在空域内直接融合

title(‘两图像直接相加融合’)  %显示融合后的图像，并命名为“两图像直接相加融合”

3）两幅图像傅立叶变换融合的程序语句

F1=fft2(X1);

F2=fft2(X2);   %分别计算两幅图像的快速傅立叶变换

X=abs(ifft2(F1+F2)/2); %两幅图像在频域内相加后的傅立叶逆变换

Subplot(1 3 2)

Image(X),colormap(map2);  %显示融合后的图像

Title(‘两幅图像傅立叶变换融合’)  %给融合后的图像命名并显示在图上

4）两幅图像小波变换融合的程序语句

[C1,L1]=wavedec2(X1,2, ‘sym4’);

[C2,L2]=wavedec2(X2,2, ‘sym4’); %分别对两幅原图像进行小波分解

C=C1+C2; %对分解系数进行融合

X=waverec2(C,L1, ‘sym4’);  %对融合后的信号进行图像重构

Subplot(1 3 3)

Image(X/2),colormap(map2);  %显示经过小波变换融合后的图像

Title(‘两图像小波变换融合’)  %给融合后的图像命名并显示在图上

三、实验步骤

1打开计算机，启动MATLAB程序；

2调入“实验一”中获取的两幅数字图像，并进行三种方法的图像融合；

3显示并保存图像融合的结果。

4记录和整理实验报告。

四、实验仪器

        1计算机；

2 MATLAB程序；

3移动式存储器（软盘、U盘等）；

4记录用的笔、纸。

五、实验报告内容

1叙述实验过程；

2提交实验的原始图像和结果图像。

六、思考题

1简述图像小波变换的实现过程。

2小波变换的分解层数对融合图像的光谱质量和空间质量有什么影响?

实验 七   频域滤波（选作）

一、实验目的 

通过实验了解频域高频和低频滤波器对图象处理的效果，了解离散余弦变换在图象变换编码

中的作用。 

二、实验内容 

a)  载入图象文件，用低通滤波器（一阶巴特沃斯滤波器）对图象进行变换。（imfilter）  

（lowpassfilter） 

b)  用高通滤波器对图象进行变换。(highpassfilter)  

c)  用离散余弦变换（dct2 和 idct2），对图象进行变换，然后对变换后的谱系数进

行阀值化；观察阀值化后的主系数谱的离散余弦逆变换。比较结果。 

d)  整理以上 4次实验的内容。 

3  实验要求 

） 完成实验内容，对出现的结果进行解释。复制他人结果者，零分处理。

） 整理完整的实验报告，打印，第 10周（下周）四下午交齐。

cutoff < 0 | cutoff > 0.5

 error('cutoff frequency must be between 0 and 0.5'); 

    end 

if rem(n,1) ~= 0 | n < 1

error('n must be an integer >= 1');

    end 

    rows = sze(1); cols = sze(2); 

    % X and Y matrices with ranges normalised to +/- 0.5 

    x =  (ones(rows,1) * [1:cols]  - (fix(cols/2)+1))/cols; 

    y =  ([1:rows]' * ones(1,cols) - (fix(rows/2)+1))/rows; 

    radius = sqrt(x.^2 + y.^2);        % A matrix with every pixel = radius relative to 

centre. 

f = fftshift( 1 ./ (1.0 + (radius ./ cutoff).^(2*n)) );   % The filter

function f = highpassfilter(sze, cutoff, n) 

if cutoff < 0 | cutoff > 0.5

 error('cutoff frequency must be between 0 and 0.5'); 

    end 

if rem(n,1) ~= 0 | n < 1

error('n must be an integer >= 1');

    end 

    f = 1.0 - lowpassfilter(sze, cutoff, n);

资源：

http://www.csse.uwa.edu.au/~pk/Research/MatlabFns/index.html