数字图像处理大作业

1图像变换

1.1实验背景

在数字图像应用领域，图像需要进行分析、变换、压缩或者增强来提高图像的可处理性和视觉效果。其中，图像变换是将图像从空间域变换到频率域，变换的目的是根据图像在变换域的某些性质对其进行处理，而这些性质在空间域难以获取，通常在频率域才能获取，在变换域处理完后再反变换到空间域，恢复图像。图像变换可以减少图像的计算量，提高处理性能。

图像变换常用的有三种变换方法：傅里叶变换，沃尔什-哈达码变换，离散余弦变换。此外，还有近年来兴起的小波变换。其中应用最广泛最重要的是傅里叶变换。它的变换核是复指数函数，转换域图像是原空间域图像的二维频谱，其直流项与原图像亮度的平均值成比例，高频项表征图像中边缘变化的强度和方向，快速傅里叶变换[3]是为了提高原傅里叶变换的运算速度孕育而生的。本文重点分析傅里叶变换的原理以及它在数字图像中的应用效果。

1.2实验目的

熟悉并掌握图像变换中的傅里叶变换原理，完成傅里叶变换的简单实例。

1.3实验原理

1807年，傅里叶首先提出傅里叶级数的概念，即任一周期信号可以分解为复正弦信号的叠加。在此基础上傅里叶在1822年提出傅里叶变换，它是数字图像中应用最重要最广泛的正交变换。函数)(x f 的一维连续傅里叶变换和反变换由下式定义：

dx e

x f u F ux

j ?

∞

∞

--=

π2)()( （1）

du e

u F x f ux

j ?

∞

∞

-=

π2)()( （2）

傅里叶变换的实部，虚部，振幅，能量和相位分别用下式表示： 实部 ?

∞

∞-=

dx ux x f u R )2cos()()(π （3） 虚部 ?

∞

∞

-=dx ux x f u I )2sin()()(π （4）

振幅 []21

2

2

)()()(u I u R u F += （5） 能量 2

)()(u F u E = （6）

相位 )

()(arctan )(u R u I u =φ （7）

傅里叶变换很容易扩展到二维情况，它的定义式如（8）所示：

??

∞∞

-∞

∞

-+-=

dxdy e

y x f v u F vy ux j )

(2),(),(π （8）

同理得其傅里叶反变换、频率谱、相位谱和能量谱也可以从一维扩展到二维。

离散傅里叶变换（DFT ）是傅里叶变换在时间域和频率域上都呈离散的形式，将信号的时间域采样变换为其DFT 的频率域采样，在形式上，变换两端的序列是有限长的，都应当被认为是离散周期信号的主

值序列，在实际应用中通常用快速傅里叶变换（FFT ）来计算DFT 。 离散傅里叶变换和反变换的定义如下式所示： ∑

-=-=10

/2)(1

)(N x N

ux j e

x f N

u F π （9）

∑

-==

1

/2)()(N u N

ux j e

u F x f π （10）

其中N 为离散序列)(x f 的长度，)(u F 也为复函数形式，其频谱，相位谱和能量谱形式与连续傅里叶变换类似，也可以推广到二维情形。

对于一个有限长序列)(x f ，傅里叶变换由下式表示：

∑-==10

)(1)(N x ux

N

W x f N u F （11） 令N j N e W π2-=，从式子可以看出要得到每一个频率分量，需要进行N 次乘法和N-1次加法，完成整个变换需要N 2

乘法和N(N-1)次加法，花费时间太大。快速傅里叶变换就是利用包含在DFT 系数矩阵中的规律，将矩阵元素巧妙排列替换，以减少乘法计算次数。通过快速傅里叶变换，可以快速获取图像的频谱、相位谱、能量谱等信息，以便于在频域中对图像进行处理。

傅里叶变换的性质包括可分离性、平移性、周期性、共轭对称性、旋转性和尺度变换等。通过这些性质可以地对图像的各种谱进行相应的处理来获得预期的效果。

1.4实验过程与结果分析

图1显示了lena 图像进行傅里叶变换后的频谱图以及lena 图乘以1

-e 进行变暗处理和加入高斯噪声再傅里叶变换后的两幅频谱图的对比，它们的频谱图都通过傅里叶变换的平移特性将其中心移到图中央。a 为lena 图像的原图，b 为进行暗处理后的lena 图像，c 为加入了高丝噪声的lena 图像。

图1 lena 图像的转换

从a和d中可以看出，图像的能量主要集中在低频处，高频出的幅值很小；通过b和e看出变暗后的图像的频谱图中央低频分量变小，可知中央低频代表了图片的平均亮度，当图片平均亮度发生变化时，对应的频谱图中央的低频分量也改变，说明一幅图的平均亮度由低频区域的大小决定；通过c和f看出，加入高斯噪声后，得到有颗粒噪声的图，变换后的频谱图高频幅值数值增加，分布增多。由此得出，图像灰度变化缓慢的区域，对应它变换后的低频分量部分，图像灰度呈阶跃变换的区域，对应变换后的高频分量部分，除颗粒噪声外，图像细节的边缘，轮廓处都是灰度变化突变区域，他们都具有变换后的高频分量特征。

图2所示为lena图的幅值谱图与相位谱图以及分别利用幅值谱图与相位谱图重构的lena图。

图2 lena图的幅值谱图与相位谱图及重构图

通过对比看出，我们能从幅值谱图中获得比相位谱图中更多的信息，因为相比较而言，幅值谱层次更分明。但是利用幅值谱重构的图却比利用相位谱重构的图效果更差，对幅值谱图像重构，与原始图像相比，差别很大，而对相位谱图像重构，可以看出图像的基本轮廓。

1.图像分割

2.1实验背景

图像分割在图像工程中的位置它起着承上启下的作用，可以认为是介于低层次处理和高层次处理的中间层间。最近几年又出现了许多新思路、新方法、或改进算法。下面对一些经典传统方法作简要的概述。多年来人们对图像分割提出了不同的解释和表述，借助集合概念对图像分割可给出如下定义：令集合R

代表整个图像区域，对R的图像分割可以看做是将R分成N个满足以下条件的非空子集R

1，R

2

，R

3

，…，

R N ；

（1）在分割结果中，每个区域的像素有着相同的特性；

（2）在分割结果中，不同子区域具有不同的特性，并且它们没有公共特性； （3）分割的所有子区域的并集就是原来的图像； （4）各个子集是连通的区域；

图像分割是把图像分割成若干个特定的、具有独特性质的区域并提取出感兴趣目标的技术和过程，这些特性可以是像素的灰度、颜色、纹理等提取的目标可以是对应的单个区域，也可以是对应的多个区域。图像分割方法有许多种分类方式，在这里将分割方法概括为四类：(1)边缘检测方法(2)区域提取方法(3)阈值分割方法(4)结合特定理论工具的分割方法。

根据所给的实验图片flower.bmp ，要把花朵从背景中分离出来，最适合的方法是阈值分割方法，所以本实验重点分析用全局阈值、局部阈值方法处理图片的方法，实现分割的目的。

2.2实验目的

熟悉并掌握图像分割中的阈值分割方法，完成阈值分割方法的简单实例

2.3实验原理

阈值化图像分割是一种最基本的图像分割方法，经过半个多世纪的研究，现已提取了大量的算法。其基本原理就是选取一个或多个处于灰度图像范围之中的灰度阈值，然后将图像中各个像素的灰度值与阈值比较，并根据比较的结果将图像中的对应像素分成两类或多类，从而把图像划分成互不重叠的区域集合，达成图像分割的目的。

采用阈值化图像分割时通常需要对图像作一定的模型假设。利用图像模型尽可能了解图像有几个不同的区域组成。基于图像分割模型经常采用这样一种假设：目标或背景内相邻像素间的灰度值是相似的，但不同目标或背景的像素在灰度上存有差异。如果图像的是暗的背景上有较亮的对象组成的，那么目标和背景像素会有两种主要模式的灰度级。一种从背景上提取对象的明显方法是选取一个阈值T 来分离这两种模式。任何满足T y x f ≥),(的点),(y x 称为对象点，其他点则称为背景点。换而言之，阈值处理后的图像

),(y x g 定义为

?

?≤=T

y x f T y x f y x g ),(0

),(1

),(＞

标注为1的像素对应于对象，而标注为0的像素则对应于背景。T 为常数时，这种方法称为全局阈值处理。

阈值变化的方法，称为局部阈值处理。 2.3.1全局阈值处理

选取阈值的一种方法是目视检查直方图，双峰的直方图可以很容易的找到一个阈值来分开背景和目标。另一种选择阈值的方法是反复试验，挑选不同的阈值，知道观测者觉得产生了较好的结果时为止。这在交互式环境下特别有效。

为自动选择一个阈值，Gonzalez和Woods描述了如下迭代步骤：1．选择一个T的初始估计值。

2．用T分割图像。这样做会生成两组像素：G

1由所有灰度值大于T的像素组成，而G

2

由所有灰度值

小于或等于T的像素组成。

3．对区域G

1和G

2

中的所有像素计算平均灰度值μ

1

和μ

2

。

4．计算新的门限值：T=(μ

1+μ

2

)/2

5．重复步骤2到4，直到逐次迭代所得的T值之差小于事先定义的参数T

。

2.3.2局部阈值处理

仔细观察原图，发现花朵的左右光线亮度不一致，背景照明不均匀。在这种情况下，一种常用的处理方法是针对照明问题做预处理以补偿图像，然后再对预处理后的图像采用全局阈值处理。

2.4实验过程与结果分析

2.3.1全局阈值处理

对于实验图flower.bmp及其直方图如图3所示，此图的直方图并非是明显的双峰图，所以用全局阈值的方法不是很容易区分，但还是可以用迭代的算法进行实验。

(a) 原图

图3 flower原图及其直方图

按照上述算法对图像处理，得到阈值T=58，分割结果如图4所示，很明显花朵并没有完全分割出来，后半部分缺失。

图4 全局阈值分割结果

2.3.2局部阈值处理

050100150200250

图5改进的阈值处理结果是通过应用一个形态学顶帽算子，并对得到的结果使用函数graythresh 来计算的。这种处理等同于使用局部变化的阈值函数),(y x T 对),(y x f 进行阈值处理：

?

?≤=T

y x f T y x f y x g ),(0

),(1

),(＞，其中，o o T y x f y x T +=),(),(

图像),(y x f o 是f 的形态学开运算，常数o T 是对o f 应用函数graythresh 后的结果。实验计算得阈值

T =18，如图5(c)所示，视觉上感觉花朵的右半边还不是很明显，如果使阈值更小，T =10时，如图5(d)

所示，则把很多背景选了出来，虽然右半边明显了，但是特别是边缘、细小的背景内容都错误的变成了目标对象。所以图5(c)的效果是可以的，此阈值比较适合。

’

(a)经开运算后的图像 (b)顶帽变换 (c)经阈值处理后的图像 (d)T=10的阈值分割结果

图5 局部阈值处理结果

附录：

1.傅里叶变换

%读入显示原图

I=imread('lena.bmp'); subplot(2,3,1) imshow(I)

xlabel('(a)lena 图'); %快速傅里叶变换后读取其频谱图 K1=fft2(I); J1=fftshift(K1); subplot(2,3,4)

imshow(log(abs(J1)),[8,10]) xlabel('(d)lena 图的频谱图'); %将图像变暗后进行FFT 变换读取其频谱图 Dark=I*exp(-1); subplot(2,3,2) imshow(Dark)

xlabel('(b)变暗的lena 图');

K2=fft2(Dark);

J2=fftshift(K2);

subplot(2,3,5)

imshow(log(abs(J2)),[8,10])

xlabel('(e)变暗的lena图的频谱图');

%图像加入高斯噪声后进行FFT变换读取其频谱图

Noise=imnoise(I,'gaussian',0,0.01);

subplot(2,3,3)

imshow(Noise)

xlabel('(c)加入高斯噪声的lena图');

K3=fft2(Noise);

J3=fftshift(K3);

subplot(2,3,6)

imshow(log(abs(J3)),[8,10])

xlabel('(f)加入高斯噪声的lena图的频谱图');

%求出幅值谱和相位谱，并进行图像重构

fftI=fft2(I);

sfftI=fftshift(fftI);

RRfdp1=real(sfftI);

IIfdp1=imag(sfftI);

a=sqrt(RRfdp1.^2+IIfdp1.^2);

a=(a-min(min(a)))/(max(max(a))-min(min(a)))*255; figure(2);

subplot(2,2,1)

imshow(real(a))

xlabel('(a)lena图的幅值谱图');

b=angle(fftI);

subplot(2,2,3)

imshow(real(b))

xlabel('(c)lena图的相位谱图');

theta=pi/6;

RR1=a*cos(theta);

II1=a*sin(theta);

fftI1=RR1+1i.^II1;

C=ifft2(fftI1)*255;

subplot(2,2,2)

imshow(real(C))

xlabel('(b)利用幅值谱图重构lena图');

MM=150;

RR2=MM*cos(angle(fftI));

II2=MM*sin(angle(fftI));

fftI2=RR2+1i.^II2;

D=ifft2(fftI2);

subplot(2,2,4)

imshow(real(D))

xlabel('(d)利用相位谱图重构lena图');

2.全局阈值分割函数quanju.m

function X=quanju(f)

T = 0.5*(double(min(f(:))) + double(max(f(:))));

%初始估计值为最大亮度和最小亮度值的中间值

done = false;

while ~done

g1=f>=T;%亮度值大于等于阈值th的所有像素组成g1

g2=f

u1=mean(f(g1));%g1范围内的像素的平均亮度u1

u2=mean(f(g2));%g2范围内的像素的平均亮度u2

Tnext = 0.5*(u1+u2);%计算一个新的阈值

done = abs(T - Tnext) < 0.5;%判断连续迭代中阈值th的差比预先制定的参数t0=0.5为止T = Tnext;

end

T2=graythresh(f);%阈值

t=T2*255;%转换为0-255的阈值

[M N]=size(f);

X=zeros(M,N);

for i=1:M

for j=1:N

if f(i,j)>=t

X(i,j)=1;%大于等于阈值的置为白色

else

X(i,j)=0;%小于阈值的置为黑色

end

end

end

figure;imshow(X);title('全局阈值分割后的图像');%显示分割后的图像

3.图像分割fenge.m

f=imread('flower.bmp');%读入图像

figure;imshow(f);%显示原图

figure;imhist(f);%显示原图的直方图

%全局阈值分割

X1=quanju(f);%显示全局阈值分割结果

%局部阈值处理

se=strel('disk',10);%strel生成指定形状的形态学结构，disk碟形fo=imopen(f,se);%对图像进行开运算产生对整个图像背景的合理估计figure;imshow(fo);title('图像背景');

f2=imsubtract(f,fo);%顶帽变换也可用f2=imtophat(f,se);

%从原图像中减去背景，可以生成一幅具有合适背景且均匀的背景的图像figure;imshow(f2);title('均匀背景的图像');

X2=quanju(f2);%显示新的全局阈值分割结果