【计算机视觉】图像增强----图像的傅立叶变换

**个人简介： **

📦个人主页：赵四司机
🏆学习方向：JAVA后端开发
📣种一棵树最好的时间是十年前，其次是现在！
⏰往期文章：SpringBoot项目整合微信支付
💖喜欢的话麻烦点点关注喔，你们的支持是我的最大动力。

一：相关概念

1.什么是傅里叶变换

    傅里叶变换是线性系统分析的一个有力工具，它能够定量地分析诸如数字化系统、采样点、电子放大器、卷积滤波器、噪音和显示点等的作用。通过实验培养这项技能，将有助于解决大多数图像处理问题。对任何想在工作中有效应用数字图像处理技术的人来说，把时间用在学习和掌握博里叶变换上是很有必要的。 

2.傅里叶变换的定义

对于二维信号，二维Fourier变换定义为：

二维离散傅立叶变换为：

图像的傅立叶变换与一维信号的傅立叶变换变换一样，有快速算法，具体参见参考书目，有关傅立叶变换的快速算法的程序不难找到。实际上，现在有实现傅立叶变换的芯片，可以实时实现傅立叶变换。

二：傅里叶变换

a) 绘制一个二值图像矩阵,并将其傅立叶函数可视化。

f=zeros(30,30);

f(5:24,13:17)=1;

imshow(f,'InitialMagnification','fit')

图2.5.1-1 矩阵可视化结果

F=fft2(f);

F2=log(abs(F));

figure,imshow(F2,[-1 5],'InitialMagnification','fit');colormap(jet);

图2.5.1-2 矩阵频谱图

    分析：fft2表示将二维矩阵转换为频域，可以看到得到的像谱图不是中心对称的， 因此需要进行下一步处理。可以看到，由于矩阵的零频分量分布在矩阵的四周，而1 频分量分布在中心，所以在矩阵中间的梯度较大，对应的频谱图越暗。

F=fft2(f,256,256); %零填充为256×256矩阵

figure,imshow(log(abs(F)),[-1 5],'notruesize');colormap(jet);

图2.5.1-3 256*256矩阵频谱图

F2=fftshift(F); %将图像频谱中心由矩阵原点移至矩阵中心

figure,imshow(log(abs(F2)),[-1 5],'notruesize');colormap(jet);

图2.5.1-4 频谱中心化处理结果

    分析：根据快速傅里叶变换的计算要求，需要图像的行数、列数均满足2的n次方，所以进行转换前要进行补零操作。一个M行N列的二维图像f(x,y)，先按行队列变量y做一次长度为N的一维离散傅里叶变换，再将计算结果按列向对变量x做一次长度为M傅里叶变换就可以得到该图像的傅里叶变换结果。fft2将时域信号转换为频域信号，低频分布在四周，fftshift用于将fft2的 DC分量移动到频谱中心，低频区域移动到中央，这样便于处理，可以看到得到的频谱图是中心对称的。

b)利用傅立叶变换分析两幅图像的相关性，定位图像特征。读入图像‘text.png'，抽取其中的字母‘a’。

bw=imread('text.png');

a=bw(32:46,85:102);

imshow(bw);

图2.5.1-5 bw可视化结果

figure,imshow(a);

图2.5.1-6 截取’a’

C=real(ifft2(fft2(bw).*fft2(rot90(a,2),256,256)));%求相关性

figure,imshow(C,[]);

图2.5.1-7 计算相关性

thresh=max(C(:));

figure,imshow(C>thresh-10)

图2.5.1-8 相关值大于thresh-10

figure,imshow(C>thresh-15)

图2.5.1-9 相关值大于thresh-15

    在空域中f(x,y)与 与 g(x,y) 的相关等价于频域中 F(u,v) 的共轭与 G(u,v)  相乘,利用傅里叶变换，可以将积分求和过程转换成频域相乘（程序里求的是相乘后的实部），简化了计算过程，相关定理与卷积定理类似。相关定理可用于匹配当中，用于寻找感兴趣的物体，例如题中寻找字母a，匹配之后函数相关值在匹配点达到最大值（见图2.5.1-7，匹配点亮度最亮）。此外在寻找其他区域匹配点时，先找出了相关矩阵的最大值，然后显示出相关值大于最大值减10的匹配点，可以看到在第一第二行基本匹配成功（见2.5.1-8），显示出相关值大于最大值减15的匹配点时误差比较大。其实求的也即实部误差。

三：离散余弦变换(DCT)

a) 使用dct2对图像‘autumn.tif’进行DCT变换

RGB=imread('autumn.tif');

imshow(RGB)

图2.5.2-1 autumn.tif

I=rgb2gray(RGB); %转换为灰度图像

figure,imshow(I)

图2.5.2-2 autumn.tif灰度处理

J=dct2(I);

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

图2.5.2-3 autumn图像dct系数

    对原始图像进行离散余弦变换，首先要把彩色图片进行灰度化处理，可以看到变换后的DCT系数能量主要集中在左上角，其余大部分DCT系数接近于0，说明离散余弦变换适用于对图像进行压缩。低频系数体现的是图像中目标的轮廓和灰度分布特性，高频系数体现的是目标形状的细节信息。DCT变换之后，能量主要集中在低频分量处，这也是DCT变换去相关性的一个体现。

b) 将上述DCT变换结果中绝对值小于10的系数舍弃，使用idct2重构图像并与原图像比较。

K=idct2(J);

figure,imshow(K,[0 255])

图2.5.2-4 DCT逆变换

J(abs(J)<10)=0; %舍弃系数

K2=idct2(J);

figure,imshow(K2,[0 255]);

图2.5.2-5 DCT压缩处理

    DCT变换和FFT变换都属于变换压缩方法，变换压缩的一个特点是将从前密度均匀的信息分布变换为密度不同的信息分布。在图像中，低频部分的信息量要大于高频部分的信息量，尽管低频部分的数据量比高频部分的数据量要小的多。对DCT系数设置门限（将DCT系数低于10的灰度值置0），这就是图像压缩的量化过程，再经过逆变换将频域转换成空域便可达到压缩效果，可以看到，虽然图像经过了压缩处理，但是图像质量肉眼看不出变化（见图2.5.2-5）。

c) 利用DCT变换进行图像压缩。

I=imread('cameraman.tif');

I=im2double(I);

T=dctmtx(8);

% B=blkproc(I,[8,8],'P1xP2',T,T');

fun = @(block_struct)Tblock_struct.dataT';

B=blockproc(I,[8 8],fun);

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];

% B2=blkproc(B,[8 8],'P1.*x',mask);

fun = @(block_struct)Mask.*block_struct.data;

B2=blockproc(B,[8 8],fun);

% I2=blkproc(B2,[8 8],'P1xP2',T',T);

fun = @(block_struct)T'block_struct.dataT;

I2=blockproc(B2,[8 8],fun);

imshow(I)

图2.5.2-6 cameramen.tif

figure,imshow(I2);

图2.5.2-5 图片压缩

    在运行参考代码时候，出现警告说建议使用blockproc来替换掉blkproc，经过查 阅得知两者的区别是blkproc必须一次性把要处理的矩阵全部导入内存中，对一些超大 型图像就无法使用了，随后matlab推出了可以处理任意大的图像blockproc函数。把 blkproc改成blockproc后，出现了错误，查阅帮助文档得知blockproc函数的原型为B = blockproc(A,[M N],FUN)，A是要处理的图像矩阵，[M,N]是每次分块处理矩阵的大小， FUN是函数句柄，即对每块矩阵的处理函数。不同于blkproc，FUN须事先定义好，且 输入必须是结构体，其data域是我们的矩阵数据，这里由blockproc分块后的机制决定。

    fun = @(block_struct)T*block_struct.data*T'的作用为进行离散余弦变换，以便于后 面进行压缩处理，其中T表示8*8的DCT变换矩阵。由于DCT系数能量值主要分布 于左上角，因此mask将左上角的系数设置为1，其余位置设置为0，这样就能保留每 个块中的低频系数，而舍弃高频系数，以达到图像压缩的目的。最后再通过fun = @(block_struct)T'*block_struct.data*T;进行反余弦变换，将频域转换为空域，得到 压缩后的图像。

四：反变换

①对(1)中的矩阵进行反变换，代码如下：

F3 = iff2(F);

figure,imshow(log(abs(F3)),[-1 5],'InitialMagnification','fit');colormap(jet);

图2.5.3-1 对矩阵进行傅里叶反变换

    由于进行傅里叶变换的是矩阵，将空域转换成了频域，而傅里叶反变换则将空域 转换回频域，重新显现出原来的图形。

②对(2)中的图像进行反变换

K=idct2(J);

figure,imshow(K,[0 255])

图2.5.3-2 DCT逆变换

    DCT可以用于对图片进行压缩，DCT过程是将空域转换为频域的过程，而经过DCT逆变换则可以将频域转换为空域，将图像还原。

五：不同的图像内容与FFT、DCT频谱之间的对应关系

    进行傅里叶变换（FFT2）后，频谱的低频部分主要集中在四周，中间则为高频部分； 而DCT则从左上角到右下角从低频到高频，系数逐渐减小。相同点为图像能量均集中在 低频部分，低频系数体现出图像中目标的轮廓和灰度分布特性，高频系数体现了目标形 状细节。需要注意的是，频谱图上各点和原图像不存在一一对应的关系，即使是在不移 频的情况下也没有。