数字图像处理上机图像平滑锐化

数字图像处理上机实习

（第二专题）

学生姓名：                         

班    级：                      

学    号：                           

指导老师：                           

实验题目

一、图象灰度变换

0,  显示灰度图象p02-01~p02-06及直方图；

1，对灰度图象p02-04的灰度范围进行适当展宽；

二、图象平滑

4，对p02-04jy、p02-04gs进行5*5方形窗口的最大均匀性平滑滤波，并比较其效果；

三、图象锐化

4，利用3*3的Krisch算子对p02-04实施图象锐化，二维梯度模板为：

H1=， H2=

H3=， H4=

实验内容

一 图象灰度变换

要求（1-0）  显示灰度图象p02-01~p02-06及直方图；

1.程序代码

 subplot(621)

 imshow(I)

 title('图像')

 subplot(622),imhist(I,256)

 title('图像直方图')

 I=imread('D:\\matlab2011\\work\\p02-02.tif');

 subplot(623)

 imshow(I)

 subplot(624),imhist(I,256)

 I=imread('D:\\matlab2011\\work\\p02-03.tif');

 subplot(625)

 imshow(I)

 subplot(626),imhist(I,256)

 I=imread('D:\\matlab2011\\work\\p02-04.tif');

 subplot(6,2,7)

 imshow(I)

 subplot(6,2,8),imhist(I,256)

 subplot(6,2,9)

 imshow(I)

 subplot(6,2,10),imhist(I,256)

 I=imread('D:\\matlab2011\\work\\p02-06.tif');

 subplot(6,2,11)

 imshow(I)

 subplot(6,2,12),imhist(I,256)

所用函数：

I=imread(‘path’)：函数imread用于读取图片文件中的数据，path为图像文件存放路径。

figure('Name', '图像')：打开一个命名为“图像显示类”图表显示图片。

subplot(m,n,p)：subplot是将多个图画到一个平面上的工具。其中，m表示是图排成m行，n表示图排成n列，也就是整个figure中有n个图是排成一行的，一共m行，如果m=2就是表示2行图。p表示图所在的位置，p=1表示从左到右从上到下的第一个位置。本程序中为6行2列。

imshow(I)：imshow是matlab中显示图像的函数，I为图像中的信息。其调用方式有很多，比如imshow(BW);imshow(I,[low high])等。

imhist(I,256):n为灰度图像灰度级，I为灰度图像，缺省值为256。

2.运行结果

要求（1-4） 对灰度图象p02-04的灰度范围进行适当展宽

1.程序代码

J=imread('D:\\matlab2011\\work\\p02-04.tif');

subplot(131)

imshow(J,[100,150]);

title('灰度范围100-150')

J=imread('D:\\matlab2011\\work\\p02-04.tif');

subplot(132)

imshow(J,[50,200]);

title('灰度范围50-200')

J=imread('D:\\matlab2011\\work\\p02-04.tif');

subplot(133)

imshow(J,[0,255]);

title('灰度范围0-255')

相关函数：

imshow(J,[50,200]):J代表所显示图像的灰度矩阵 ，[50,200]为图像数据的值域。

2.运行结果

二  图象平滑

要求（2-2） 对p02-02jy、p02-02gs进行3*3方形窗口的灰度最相近的K个邻点平均法滤波

1.程序代码

A=imread('D:\\matlab2011\\work\\p02-04.tif');

figure('Name', '图像显示');

subplot(3,2,1);

imshow(A);  %显示灰度图像

title('原图象P01-01');

J1=imnoise(A,'salt & pepper',0.02); %加均值为0，方差为0.02的椒盐噪声

subplot(3,2,3);

imshow(J1);  %显示灰度图像

title('椒盐噪声图象');

K1=filter2(fspecial('average',3),J1)/255;

subplot(3,2,4);

imshow(K1);  %显示灰度图像

title('椒盐噪声平滑图象');

J2=imnoise(A,'gaussian',0.02); %加均值为0，方差为0.02的高斯噪声。

subplot(3,2,5);

imshow(J1);  %显示灰度图像

title('高斯噪声图象');

K2=filter2(fspecial('average',3),J1)/255;

subplot(3,2,6);

imshow(K2);  %显示灰度图像

title('高斯噪声平滑图象');

2.运行结果

对p02-02jy、p02-02gs进行3*3方形窗口的灰度最相近的K个邻点平均法滤波，滤波后椒盐噪声和高斯噪声均减少，但图像变得模糊。

要求（2-4） 对p02-04jy、p02-04gs进行方形窗口的最大均匀性平滑滤波

1.算法设计

为避免消除噪声时引起边缘模糊，最大均匀性平滑算法先找出环绕每像素的灰度最均匀 窗口，然后用该窗口的灰度均值代替该像素原来的灰度值。

具体来说，对图像中任一像素（x,y）的5个有重叠的3*3邻域，用梯度衡量它们的灰度变换大小。把其中灰度变换最小的邻域作为最均匀的窗口，用其平均灰度代替像素（x,y）的灰度值。

2.程序代码

A=imread('D:\\matlab2011\\work\\p02-04.tif');

J1=imnoise(A,'salt & pepper',0.01);%加椒盐噪声

G1=imnoise(A,'gaussian',0.01);%加高斯噪声

subplot(2,2,1);

imshow(J1); 

title('p02-04jy');

subplot(2,2,3);

imshow(G1); 

title('p02-04gs');

%%%%%%%%%%%3*3最大均匀性平滑滤波%%%%%%%%%%%%

[W,L]=size(J1);

J2=double(J1);

G2=double(G1);

G3=double(G1);

J3=double(J1);

Com=zeros(1,5);

COM=double(Com);

COM1=double(Com);

%%%%%%%%%%%%%%%%计算梯度%%%%%%%%%%%%%%%

for i=3:L-2  

    for j=3:W-2

        COM(1)=abs(J2(i-2,j-2)-J2(i,j))+abs(J2(i-2,j-1)-J2(i,j))+abs(J2(i-2,j)-J2(i,j))+abs(J2(i-1,j-2)-J2(i,j))+abs(J2(i-1,j-1)-J2(i,j))+abs(J2(i-1,j)-J2(i,j))+abs(J2(i,j-2)-J2(i,j))+abs(J2(i,j-1)-J2(i,j));

        COM(2)=abs(J2(i+2,j-2)-J2(i,j))+abs(J2(i+2,j-1)-J2(i,j))+abs(J2(i+2,j)-J2(i,j))+abs(J2(i+1,j-2)-J2(i,j))+abs(J2(i+1,j-1)-J2(i,j))+abs(J2(i+1,j)-J2(i,j))+abs(J2(i,j-2)-J2(i,j))+abs(J2(i,j-1)-J2(i,j));

        COM(3)=abs(J2(i-2,j+2)-J2(i,j))+abs(J2(i-2,j+1)-J2(i,j))+abs(J2(i-2,j)-J2(i,j))+abs(J2(i-1,j+2)-J2(i,j))+abs(J2(i-1,j+1)-J2(i,j))+abs(J2(i-1,j)-J2(i,j))+abs(J2(i,j+2)-J2(i,j))+abs(J2(i,j+1)-J2(i,j));

        COM(4)=abs(J2(i+2,j+2)-J2(i,j))+abs(J2(i+2,j+1)-J2(i,j))+abs(J2(i+2,j)-J2(i,j))+abs(J2(i+1,j+2)-J2(i,j))+abs(J2(i+1,j+1)-J2(i,j))+abs(J2(i+1,j)-J2(i,j))+abs(J2(i,j+2)-J2(i,j))+abs(J2(i,j+1)-J2(i,j));

        COM(5)=abs(J2(i-1,j-1)-J2(i,j))+abs(J2(i-1,j)-J2(i,j))+abs(J2(i-1,j+1)-J2(i,j))+abs(J2(i,j+1)-J2(i,j))+abs(J2(i,j-1)-J2(i,j))+abs(J2(i+1,j-1)-J2(i,j))+abs(J2(i+1,j)-J2(i,j))+abs(J2(i+1,j+1)-J2(i,j));

        COM1(1)=abs(G2(i-2,j-2)-G2(i,j))+abs(G2(i-2,j-1)-G2(i,j))+abs(G2(i-2,j)-G2(i,j))+abs(G2(i-1,j-2)-G2(i,j))+abs(G2(i-1,j-1)-G2(i,j))+abs(G2(i-1,j)-G2(i,j))+abs(G2(i,j-2)-G2(i,j))+abs(G2(i,j-1)-G2(i,j));

        COM1(2)=abs(G2(i+2,j-2)-G2(i,j))+abs(G2(i+2,j-1)-G2(i,j))+abs(G2(i+2,j)-G2(i,j))+abs(G2(i+1,j-2)-G2(i,j))+abs(G2(i+1,j-1)-G2(i,j))+abs(G2(i+1,j)-G2(i,j))+abs(G2(i,j-2)-G2(i,j))+abs(G2(i,j-1)-G2(i,j));

        COM1(3)=abs(G2(i-2,j+2)-G2(i,j))+abs(G2(i-2,j+1)-G2(i,j))+abs(G2(i-2,j)-G2(i,j))+abs(G2(i-1,j+2)-G2(i,j))+abs(G2(i-1,j+1)-G2(i,j))+abs(G2(i-1,j)-G2(i,j))+abs(G2(i,j+2)-G2(i,j))+abs(G2(i,j+1)-G2(i,j));

        COM1(4)=abs(G2(i+2,j+2)-G2(i,j))+abs(G2(i+2,j+1)-G2(i,j))+abs(G2(i+2,j)-G2(i,j))+abs(G2(i+1,j+2)-G2(i,j))+abs(G2(i+1,j+1)-G2(i,j))+abs(G2(i+1,j)-G2(i,j))+abs(G2(i,j+2)-G2(i,j))+abs(G2(i,j+1)-G2(i,j));

        COM1(5)=abs(G2(i-1,j-1)-G2(i,j))+abs(G2(i-1,j)-G2(i,j))+abs(G2(i-1,j+1)-G2(i,j))+abs(G2(i,j+1)-G2(i,j))+abs(G2(i,j-1)-G2(i,j))+abs(G2(i+1,j-1)-G2(i,j))+abs(G2(i+1,j)-G2(i,j))+abs(G2(i+1,j+1)-G2(i,j));

%%%%%%%%%%%%%%找出梯度值最小的3*3邻域（灰度变换最小的邻域）%%%%%%%%%%

        temp=min(COM);

            for k=1:5

                if temp==(COM(k)) t=k;

                end

            end

%%%%%%%%%%%%%%%%%计算该邻域平均灰度%%%%%%%%%%%%%%%

        if(k==1)

        J3(i,j)=J2(i-2,j-2)/9+J2(i-2,j-1)/9+J2(i-2,j)/9+J2(i-1,j-2)/9+J2(i-1,j-1)/9+J2(i-1,j)/9+J2(i,j-2)/9+J2(i,j-1)/9+J2(i,j-2)/9;

        end

        if(k==2)

        J3(i,j)=J2(i+2,j-2)/9+J2(i+2,j-1)/9+J2(i+2,j)/9+J2(i+1,j-2)/9+J2(i+1,j-1)/9+J2(i+1,j)/9+J2(i,j-2)/9+J2(i,j-1)/9+J2(i,j-2)/9;

        end

        if(k==3)

        J3(i,j)=J2(i-2,j+2)/9+J2(i-2,j+1)/9+J2(i-2,j)/9+J2(i-1,j+2)/9+J2(i-1,j+1)/9+J2(i-1,j)/9+J2(i,j+2)/9+J2(i,j+1)/9+J2(i,j+2)/9;

        end

        if(k==4)

        J3(i,j)=J2(i+2,j+2)/9+J2(i+2,j+1)/9+J2(i+2,j)/9+J2(i+1,j+2)/9+J2(i+1,j+1)/9+J2(i+1,j)/9+J2(i,j+2)/9+J2(i,j+1)/9+J2(i,j+2)/9;

        end

        if(k==5)

        J3(i,j)=J2(i-1,j-1)/9+J2(i,j-1)/9+J2(i+1,j-1)/9+J2(i-1,j)/9+J2(i,j)/9+J2(i+1,j)/9+J2(i-1,j+1)/9+J2(i,j+1)/9+J2(i+1,j+1)/9;

        end

        temp=min(COM1);

            for k=1:5

                if temp==(COM1(k)) t=k;

                end

            end

        if(k==1)

        G3(i,j)=G2(i-2,j-2)/9+G2(i-2,j-1)/9+G2(i-2,j)/9+G2(i-1,j-2)/9+G2(i-1,j-1)/9+G2(i-1,j)/9+G2(i,j-2)/9+G2(i,j-1)/9+G2(i,j-2)/9;

        end

        if(k==2)

        G3(i,j)=G2(i+2,j-2)/9+G2(i+2,j-1)/9+G2(i+2,j)/9+G2(i+1,j-2)/9+G2(i+1,j-1)/9+G2(i+1,j)/9+G2(i,j-2)/9+G2(i,j-1)/9+G2(i,j-2)/9;

        end

        if(k==3)

        G3(i,j)=G2(i-2,j+2)/9+G2(i-2,j+1)/9+G2(i-2,j)/9+G2(i-1,j+2)/9+G2(i-1,j+1)/9+G2(i-1,j)/9+G2(i,j+2)/9+G2(i,j+1)/9+G2(i,j+2)/9;

        end

        if(k==4)

        G3(i,j)=G2(i+2,j+2)/9+G2(i+2,j+1)/9+G2(i+2,j)/9+G2(i+1,j+2)/9+G2(i+1,j+1)/9+G2(i+1,j)/9+G2(i,j+2)/9+G2(i,j+1)/9+G2(i,j+2)/9;

        end

        if(k==5)

        G3(i,j)=G2(i-1,j-1)/9+G2(i,j-1)/9+G2(i+1,j-1)/9+G2(i-1,j)/9+G2(i,j)/9+G2(i+1,j)/9+G2(i-1,j+1)/9+G2(i,j+1)/9+G2(i+1,j+1)/9;

        end

    end

end;

subplot(2,2,2);

imshow(uint8(J3));

title('jy经3*3最大均匀滤波');

subplot(2,2,4);

imshow(uint8(G3));

title('gs经3*3最大均匀滤波');

3.运行结果

三  图象锐化

要求（3-4） 利用3*3的Krisch算子对p02-04实施图象锐化

1.算法设计

1971年，R.Kirsch提出了一种能检测边缘方向的Kirsch算子新方法：它使用了8个模板来确定梯度幅度值和梯度的方向。  

算法流程图

读入图像，对边界点赋值，再利用krisch算子算出每个像素点梯度值，找出最大的梯度值并保存。算出梯度后，根据需要生成三种不同的增强图像。

2.程序代码

A=imread('D:\\matlab2011\\work\\p02-04.tif');

figure('Name', '图像显示');

subplot(2,2,1);

imshow(A);  %显示灰度图像

f = im2double(A);     %类型转换

title('原图象P02-04');

B=imkrisch( f , 0.7,0)

subplot(2,2,2);

imshow(B);  %显示灰度图像

title('锐化图象1 P02-04');

B=imkrisch( f ,0.7, 1)

subplot(2,2,3);

imshow(B);  %显示灰度图像

title('锐化图象2 P02-04');

B=imkrisch( f , 0.7,2)

subplot(2,2,4);

imshow(B);  %显示灰度图像

title('锐化图象3 P02-04');

相关函数：

function [g] = imkrisch( f , T, way)

%用imkrisch算子进行图像锐化

%f 为输入图像

%T 为比较的阈值

%way 为算子的输出方式

%way 0、1、2

in = nargin;

if in < 3

    error('Not enough input argument');

end

fs = size(f);

m = fs(1,1);

n = fs(1,2);

if (m <= 1) || (n <= 1)

    g = f;

    return;

end

%直接对整个图像赋值，简化算法

g = f; 

for i = 2:m-1

    for j = 2:n-1

        m_krisch = [f(i-1,j-1),f(i-1,j),f(i-1,j+1),f(i,j-1), ...

            f(i,j+1),f(i+1,j-1),f(i+1,j),f(i+1,j+1)];

        w_krisch =      [ 5, 5, 5,-3,-3,-3,-3,-3; ...

                         -3,-3, 5,-3, 5,-3,-3, 5; ...

                         -3,-3,-3,-3,-3, 5, 5, 5; ...

                          5,-3,-3, 5,-3, 5,-3,-3];

         con =abs( m_krisch*( w_krisch') );

         g(i,j) = max(con);

    end

end

%判定输出的形式

%使用4.3.9公式  大于T:g    小于T：0黑

if way == 2

    for i = 1:m-1

        for j = 1:n-1

            if g(i,j) > T

                g(i,j) = g(i,j);

            else

                g(i,j) = 0;

               % g(i,j) = g(i,j);

            end

        end

    end

end

%else    %使用4.3.8公式   大于T:g  小于T ：f

if way == 1

    for i = 1:m-1

        for j = 1:n-1

            if g(i,j) < T

                g(i,j) = f(i,j);

            end

        end

    end

end

if way == 0

    for i = 1:m-1

        for j = 1:n-1

            if g(i,j) > T

                g(i,j) = g(i,j);

            else

                g(i,j) = g(i,j);

            end

        end

    end

end

end

3.运行结果及分析

（1）采用Krisch算子，采用三种不同的增强图像方法

锐化图像1：

使原图像各点（x,y）的灰度g(x,y)等于梯度，即：g(x,y) = grad(x,y)

此法的缺点是增强的图像仅显示灰度变化比较陡的边缘轮廓，而灰度变化比较平缓或均匀的区域呈黑色。

锐化图像2：

T为阈值。适当选取T，可使边缘轮廓突出，又不会破坏原来灰度变化计较平缓的背景

锐化图像3：

T为阈值，适当选取T。背景用固定灰度级表示，便于观察边缘灰度变化。

（2）采用Krisch算子，采用增强图像方法一，阈值分别设置为0.4、0.7、0.95