前言
源于图像处理与模式识别课作业。
最初接触是在大一数学分析课上,小远姐的一次大作业要求是写一篇关于傅里叶变换的应用的文章,当时有去图书馆查阅一些信号处理的书籍,第一了解到傅里叶变换。之后暑假中参与学校实验室做图像处理项目的时候,看书看到傅里叶变换,一是当时觉得太复杂了,二是所做的项目对于实时性要求较高,即便是fft可能对于3000*2000以上的图片也难以达到实时性,主要用到的边缘检测方面Sobel算法也能满足要求,所以当时也是一笔带过,没有仔细看。
这次需要作业要求手写快速傅里叶变换,才真正开始认真看这一部分内容。因为前前后后包括思考和动手写大概花了两个周的时间,在看老师课件同时也看了一下前人写的博客,自己在实现过程中也遇到过一些阻碍,不过最后还是按时交上了作业,于是把一些思路、用到相关资料和最后的代码写成一个博客,供别人参考也供自己以后参考。由于是针对图像处理,所以还是以二维变换为主。
好下面是正文。
一、离散傅里叶DFT
二维离散傅里叶DFT的公式如下:
-----公式(1)
上课的时候没有好好听,再加上高数基本上已经忘光了,第一眼看到这个公式的是懵逼的,这怎么实现,还带着自然数的指数,而且还有个 -j 是什么鬼。好吧,滚去看书。
我用的是刚萨雷斯的《数字图像处理》,从头认真看会发现作者真的是从最基础的复数将到离散傅里叶变换是怎么得来的。
此处用了一个j是通过欧拉公式:
-----公式(2)
我想它的作用就是用来表示实轴分量和虚轴分量。我们编程的时候,其实还是需要使用欧拉公式将公式(1)重新展开成含三角函数的表达式。
对于离散傅里叶的得出属于信号处理中的详细内容,我在这里就不再赘述了,不过后面如果有时间我也会继续补充。
直接使用公式暴力求解离散傅里叶的算法复杂度是相当高的——O(n^2),对于一张512×512的图片,运算量是10^10量级的,我跑了一张100×100的图片大概估计了一下可能需要两个半小时。这也恰恰说明的fft的必要性,fft的算法复杂度是O(nlogn),同样是一张512×512的图片,大概就能降到10^6量级。对于我自己写的fft差不多在30s左右orz,而使用OpenCV自带的dct()函数几乎是秒出。当然我后面给出的fft还有很多问题,用了很多vector来存储在一定程度上拖慢了速度,我也会花时间不断改进。
离散傅里叶的代码就不给出了,能来看这篇文章的人应该都能写出来一个暴力求解离散傅里叶的程序吧。
二、快速离散傅里叶FFT
快速傅里叶使用的算法叫蝶形算法,如下图所示:
这样的一个运算单元被叫做蝶形单元,emmm不得不说程序员们的脑回路真的是清奇,我为啥觉得叫沙漏算法更形象呢。(手动划掉)。老师ppt的给出的蝶形运算的图真的困扰了我很久,然后我就转去看快速傅里叶变换的公式,想甩开图直接从公式入手进行实现,公式如下:
-----公式组(3)
不得不说,从看这一组公式开始