数字图像处理学习笔记4第四章 图像变换 附实验

第四章 图像变换 附实验

前言

图像变换：为达到某种目的将原始图像变换映射到另一个空间上，使得图像的某些特征得以突出，以便于后面的处理和识别。

4.1连续傅里叶变换

一维变换

用傅里叶变换表示的函数特征完全可以通过傅里叶反变换来重建，而不会丢失任何信息。若把一个一维输入信号做一维傅里叶变换，该信号就被变换到频域上，即得到了构成该输入信号的频谱，频谱反映了该输入信号由哪些频率构成。

函数f(x)的一维连续傅里叶变换是，及其傅里叶变换F(x)的反变换是：

1. 其中f(x)与F(x)是一个傅里叶变换对，并且对于任意一个f(x)它的F(x)是唯一的，反之亦然
2. f(x)是实函数，而F(x)是复函数，因为有 j2=-1。
3. x称为时域变量，u称为频域变量。

F(u)的实部、虚部、振幅、能量和相位分别表示如下：

注意之后就会用R(u)代替实部，I(u)带图虚部。

欧拉公式

上述过程中R(u)与I(u)的得到用到欧拉公式

空间频率

空间频率是指在一定方向上的单位空间（距离）波动的周期数。它不仅具有大小而且具有方向，是一个矢量。一幅图像有明暗和色彩的差别，是一种光的强度和颜色按空间的分布。这种空间分布的特征可以用空间频率来表示。

引进空间频率这一概念的目的：

1. 研究输入图像由哪些空间频率成分构成。
2. 在空间频率域中进行各种处理。

相对于空间频率域，有时也把图像本身叫做空间域。

空间频率可以表示单位长度上的正弦波形变化的重复次数。用空间频率为u, v的二维平面（空间频率域）来表示（分别对应于x轴方向和y轴方向）。

推广到二维

设函数f(x,y)是连续可积，且F(u,v)可积，则存在如下傅里叶变换对：

存在限制

傅里叶变换的f(x)为连续（模拟）信号，而计算机处理的是数字信号（图像数据）。

数学上采用无穷大概念，而计算机只能进行有限次计算。

4.2离散傅里叶变换DFT

一维离散傅里叶变换

离散傅里叶变换是傅里叶变换在时间域和空间域上都呈离散形式，实际上是对变换两端（在时间域和空间域）的原本无限长的序列进行采样，是其主值序列。

一维离散傅里叶变换对：

其中x，u=0，1，2，…, N－1

推广到二维

2-D DFT 与 2-D IDFT。

其中u, x=0, 1, 2, …, M-1；v, y=0, 1, 2, …, N-1；

4.3快速傅里叶变换FFT

变换过程

快速傅里叶变换是离散傅里叶变换的一种。

快速傅里叶变换

在分析离散傅里叶变换中的多余运算的基础上

，进而在消除这些重复工作的思想指导下得到的，所以在运算中大大节省了工作量， 达到了快速的目的。

对于一个有限长序列，它的（离散）傅里叶变换由下式表示:

令

其傅里叶变换转换为

实际的快速傅里叶变换过程比较麻烦，这里就不细写了。

对于原本的傅里叶变换可以看出要得到每一个频率分量，需进行N次乘法和 (N-1)次加法运算。要完成整个变换需要N2次乘法和N(N-1)次加法运算，计算复杂度很高。

库利和图基于1965年发表了一篇论文，提出了快速傅里叶变换算法。该算法把原始的N点序列依次分解成一系列短序列，求出这些短序列的离散傅里叶变换，以此来减少乘法运算量。

傅里叶变换的计算量为N2次复数乘法运算,快速傅里叶变换的计算量为N/2 * log2N次复数乘法运算。

快速傅里叶变换对一维离散傅里叶变换的直接实现的计算优势：

C(N)=N2 / (N/2 * log2N) = N / log2N

通常假设N=2n，所以可以用n来表示：C(N)=2n / n。

4.4傅里叶变换的性质

可分离性

一个二维傅里叶变换可分解为两步进行，其中每一步都是一个一维傅里叶变换。

先对f(x, y)按行进行傅里叶变换得到F(x, v)，再对F(x, v)按列进行傅里叶变换，便可得到f(x, y)的傅里叶变换结果。

平移性

f(x,y)与一个指数相乘 等于 将变换后的频率域中心移到新的位置，同理 F(x,y)与一个指数相乘 等于 将其反变换后的空间域中心移到新的位置。

f(x,y)的平移将不改变频谱的幅值。

移中性

当u0=M/2,v0=N/2时

如图，在平移前，坐标原点在窗口的左上角，窗口四角分布低频成分。经过移中之后，坐标原点移至频谱图窗口中央，因而围绕坐标原点的是低频，向外是高频。

由图可知，图像的能量主要集中在低频去，相对高频的位置的幅度值很小，可能接近于0，之后进行傅里叶变换之后都要进行移中。

周期性与共轭对称性

旋转性质

当变量x，y，u，v都用极坐标表示时：

上述表明,对f(x,y)旋转一个角度θ0对应于F(u,v) 将其傅里叶变换也旋转相同的角度θ0。

卷积

设:f(x),g(x)是R上的两个可积函数，作积分

可以证明，关于几乎所有的实数x，上述积分是存在的。这样，随着x的不同取值，这个积分就定义了一个新函数h(x)，称为函数f与g的卷积，记为h(x)=(f*g)(x)。

卷积定理表明两个二维连续函数在空间域中的卷积可用求其相应的两个傅里叶变换乘积的逆变换而得。反之，在频域中的卷积可用空间域中乘积的傅里叶变换而得。

目的：与其在一个域中作不直观的和难懂的卷积，不如在另外一个域中作乘法，可以达到相同的效果。

尺度变化

a * f(x,y) = a * F(u,v)

图片

如图，在乘以e-1之后，图像整体变暗。中央低频成分代表图像的平均亮度。

平均值

频谱直流成分的 1/(MN)0.5 倍等于图像平面的亮度平均值。

4.5图像傅里叶变换实例

例1

加噪图存在明显的颗粒噪声，变换后的高频成分增多。

例2

1. 从幅值谱图像中得到的信息比在相位谱图像中得到的信息多，但对幅值谱图像重构后，即忽略相位信息，将其设为常数，所得到的图像与原始图像相比，结果差别很大。
2. 对相位谱图像重构后，即忽略幅值信息，将其设置为常数，可以从中看出图像的基本轮廓。

4.6离散余弦变换DCT

离散余弦变换（Discrete Cosine Transform， DCT）是可分离的变换，其变换核为余弦函数。在对语音信号、图像信号的变换中，DCT变换被认为是一种准最佳变换。

定义

矩阵表示

离散余弦变换具有很强的“能量集中”特性，能量主要集中在左上角处。

4.7课程实验

实验内容

（1）对图像lena、cameraman和face进行傅里叶变换，观察图像能量在频谱图中的分布情况。

（2）利用Matlab生成下列图像，并对其进行旋转30度、90度和120度，然后对他们分别进行傅里叶变换。

（3）对图像lena、cameraman和face用DCT变换进行图像压缩，舍掉的变换系数分别小于0.01、0.03、0.05，求经压缩、解压后的图像。

实验1

I1 =imread('face.jpg');I1=im2double(I1); 
I2 =imread('cameraman.tif');I2=im2double(I2); 
I3 =imread('lena.jpg');I3=im2double(I3); 

%求傅里叶变换
I1f = fft2(I1);I2f = fft2(I2); I3f = fft2(I3);
I1f1 =abs(I1f);I2f1 = abs(I2f); I3f1 = abs(I3f);
I1f2 = fftshift(I1f1); I2f2 = fftshift (I2f1); I3f2 = fftshift(I3f1);
%显示图像
figure;
%原图
subplot(3,3,1);
imshow(I1);
xlabel('原图1 '); 

subplot(3,3,2);
imshow(I2);
xlabel('原图2 '); 

subplot(3,3,3);
imshow(I3);
xlabel('原图3 ');  

%快速傅里叶变换的图
subplot(3,3,4);
imshow(I1f,[5,30]);
xlabel('频谱图1');  

subplot(3,3,5);
imshow(I2f,[5,30]);
xlabel('频谱图2');

subplot(3,3,6);
imshow(I3f,[5,30]);
xlabel('频谱图3');

%傅里叶变换频谱中间零频率的图
subplot(3,3,7);
imshow(I1f2,[5,30]);
xlabel('中心移到零点的频谱图1');

subplot(3,3,8);
imshow(I2f2,[5,30]);
xlabel('中心移到零点的频谱图2');

subplot(3,3,9);
imshow(I3f2,[5,30]);
xlabel('中心移到零点的频谱图3');
           

实验2

%构造原始图像
I = zeros(256,256); 
I(28:228,108:148) = 1;
%旋转图像
I1=imrotate(I,30,'bilinear');
I2=imrotate(I,90,'bilinear');
I3=imrotate(I,120,'bilinear');
%傅里叶变换
If=fft2(I);F=abs(If);If=fftshift(F);
I1f=fft2(I1);F=abs(I1f);I1f=fftshift(F);
I2f=fft2(I2);F=abs(I2f);I2f=fftshift(F);
I3f=fft2(I3);F=abs(I3f);I3f=fftshift(F);
%显示图像
subplot(2,4,1);imshow(I);xlabel('原图 ');
subplot(2,4,2);imshow(I1);xlabel('旋转30度 ');
subplot(2,4,3);imshow(I2);xlabel('旋转90度');
subplot(2,4,4);imshow(I3);xlabel('旋转120度');

subplot(2,4,5);imshow(If,[5,50]);=xlabel('原图的像傅里叶频谱 ');
subplot(2,4,6);imshow(I1f,[5,50]);=xlabel('旋转30度的傅里叶图谱 ');
subplot(2,4,7);imshow(I2f,[5,50]);=xlabel('旋转90度的傅里叶图谱 ');
subplot(2,4,8);imshow(I3f,[5,50]);=xlabel('旋转120度的傅里叶图谱 ');
           

实验3

I =imread('cameraman.tif');
[M,N]=size(I);%M=512,N=512
I=im2double(I); 

%生成标准DCT变化中的矩阵（8x8）
n=8;[cc,rr]=meshgrid(0:n-1);
C=sqrt(2/n)*cos(pi*(2*cc+1).*rr/(2*n));
C(1,:)=C(1,:)/sqrt(2);

%光亮度量化表
a=[16 11 10 16 24 40 51 61;
    12 12 14 19 26 58 60 55;
    14 13 16 24 40 57 69 56;
    14 17 22 29 51 87 80 62;
    18 22 37 56 68 109 103 77;
    24 35 55 64 81 104 113 92;
    49 64 78 87 103 121 120 101;
    72 92 95 98 112 100 103 99 ];

%分块做DCT变换（8x8） DCT变换公式 正变换Y=CIC'
for i=1:8:M
    for j=1:8:N
        P=I(i:i+7,j:j+7);
        K=C*P*C';
        I1(i:i+7,j:j+7)=K;
        K=K./a; %量化（按位除）
        K(abs(K)<0.05)=0;  %舍掉的变换系数分别小于0.01 0.03 0.05
        I2(i:i+7,j:j+7)=K;  
    end 
end

  %分块做DCT反变换（8x8），逆变换 P=C'YC
for i=1:8:M
    for j=1:8:N
        P=I2(i:i+7,j:j+7).*a;%反量化
        K=C'*P*C;
        I3(i:i+7,j:j+7)=K;     
    end
end

subplot(2,2,1);imshow(I);xlabel('原图 ');
subplot(2,2,2);imshow(I1);xlabel('DCT变换后的频域图像 ');
subplot(2,2,3);imshow(I2);xlabel('量化后的频域图像 ');
subplot(2,2,4);imshow(I3);xlabel('复原图像 ');