傅立叶变换，在图形学上貌似应用不怎么多，海水波浪的运动算法方面需要吧，然后不怎么见了。而在图像处理上，频域处理这个大领域，傅立叶变换可是算法基础口牙~——ZwqXin.com

其实啊,这个傅立叶的实现同是为了作业,但就在形态学处理之前。之所以迟迟不记录在博，主要还是因为自己对其原理还是云里雾里的。嘛~有时候写多了日志就有种“什么东西要说就得说明白啊”的错觉。其实嘛，本博客也不是什么传道授业解惑的地方，只是写着写着就有点忘记初衷了。

废话太多了。其实我只想说，找傅立叶原理算法的同学，接下来这篇文章没啥太大价值的。姑且当作一个被傅立叶实现搞晕了一两天的傻子的自言自语，恩，这里是我的记录簿。

本文来源于 ZwqXin (http://www.zwqxin.cn/), 转载请注明
      原文地址：http://www.zwqxin.cn/archives/image-processing/dif-fast-fourier.html

最初被傅立叶缠上应该是上自动控制原理的时候吧，印象大概就是时域下的东西转去频率域吧，然后扔出一个复杂的公式。都知道，连续波都可以由正弦搅成。那么，代表任意连续波的函数由正弦函数搅成如何呢？正弦波状态特征是频率，振幅，相位吧，那么，一段时间内的波的函数形式分成几种不同周期（不同频率）的正弦表示的“有机组合”，就是所谓的时域 - 频域转换了。这很有用吧，在很多领域，包括自动控制.....不过其实我不太清楚具体起来是什么用，包括自动控制领域。

而在图像上的应用，其实也就是空间域 - 频率域。一张500*400的灰度图像，不妨假想为一张有500*400个方格的地形图（啥叫3D中毒...），某坐标处的灰度自然就是地形在该点的高度了（当然，没可能就像一张地形灰度图那样能形成比较连续的地形）。图像是按像素记数的，每个相邻像素的灰度或多或少都有个差值，形成一个梯度——关键就在于这个“或多或少”上。现在，忽略相邻像素灰度间的插值等问题，就有一个地形出来了，地形当然有陡坡有缓坡（正是由梯度造成）——这就是我们要处理的空间波了。傅立叶变换功效在于，变换结果是频率（梯度）——经过变换所得的是一堆按灰度大小聚集的新像素（总的像素数量不变，但意义改变了，反映的是梯度而不是原来的“灰度”，而且与原来的空间位置无关）。

该看什么呢？就看“聚集的情况”（如，图中的白色区域和黑色区域，注意，它们没有明显的分界，因此只能直观去感受~）对于中部高亮区域，反映了与相邻像素梯度比较大的像素量；周围的黑色区域，反映了与相邻像素梯度比较小的像素量。这么说，就是图像中的高频部分和低频部分了——图像中的边缘处就是最直观的“高频部分”了。至于要分离出高频和低频有什么用，那就是频率域处理的事情了。傅立叶变换作用就只是生成一张这样的图以把高低频部分分开。

图像中离散的像素可以认为是二维波的采样点。因此对每行进行一维离散傅立叶变换，再在此基础上对每列进行一维离散傅立叶变换，则能得出两个纬度结合的傅立叶频率图（如上）。离散傅立叶变换的计算形式对计算机而言太残酷了（前面说了，式子复杂得很），因此就有人发明出快速傅立叶变换这种东西。

最容易理解的是“基-2”的快速傅立叶，因为它的采样点尺寸是2的N次方（N> 0），能直接执行蝶形算法。一般图像长宽都不一定满足——因此比较受限。补或截都是迫不得已的方法。N为合数的 FFT 算法 （混合基）可以不受限，不过复杂死了，我不想去深究。对“基-2”的快速傅立叶，又分为“时域抽取基-2”和“频域抽取基-2”，其实差不多，不过中间的步骤顺序有点差异而已（前者还要倒码之类的）。

我用的是“频域抽取基-2快速傅立叶”（DIF-FFT）。算法推导，诸君在网络上找找吧，反正结论就是一个符合蝶形算法的形式。

//	DIT – FFT [频域抽取基-2FFT] 蝶形算法：double angle;	for(i  =  0;  i  <  Num / 2;  i++)    	{    	  angle  =  -i  *  PI  *  2  /  Num;    	  if(!InverseTranform)		  W[i]  =  complex <double>  (cos(angle),  sin(angle)); 	  else		  W[i]  =  complex <double>  (cos(angle),  -sin(angle)); 	 } 	int Apartsize = 0, halfPartsize = 0, Partoffset = 0;	for(k = 0; k < r; k++)	{		for(j = 0; j < (1<<k); j++)		{			Apartsize = 1 << (r - k);			halfPartsize = Apartsize / 2;			for(i = 0; i < halfPartsize; i++)			{				Partoffset = j * Apartsize;				X2[i + Partoffset] 					= X1[i + Partoffset] + X1[i + Partoffset + halfPartsize];				X2[i + Partoffset + halfPartsize] 				    = (X1[i + Partoffset] - X1[i + Partoffset + halfPartsize]) * W[i * (1<<k)];			}		}        X = X1;//		X1 =X2;        X2 = X;//	}//到二维，横竖方向各转一次：	for(j = 0; j <Exheight; ++j)		DifFastFourierTransform1D(&input[j * Exwidth], &output[j * Exwidth], Exwtimes);	for(j = 0; j <Exheight; ++j)	  for(i = 0; i <Exwidth; ++i)          input[i + j * Exwidth] = output[j + i * Exheight];	 for(i = 0; i <Exwidth; ++i)		  DifFastFourierTransform1D(&input[i * Exheight], &output[i * Exheight], Exhtimes);

另外，对彩色图片，一般是针对每个颜色通道都做一个快速傅立叶变换哦。

再另外，其实傅立叶变换除了可以生成频率域图，还能生成相位图哦。只是输出的时候输出结果的振幅就可了，相位图能反映出梯度点在整张图片上的方位呢~

最后一句话：反变换跟正变换差不多的。不多得结合频率域图和相位图的信息哦，才能还原。

程序截图：
[若看不到本博客图片，可参照此帖，简单的复原法：[显示本站所有图片] ]

 
（原始图）

(Fourier)

(彩色图，原图)

     
[对非2幂大小的某图像，按需（见对话框）取其各个角域。图上为变换（上）和反变换]

[按需察看彩色图的不同通道下的频域图，图左为上图的蓝色通道；图右：上图的相位图]

本文来源于 ZwqXin (http://www.zwqxin.cn/), 转载请注明
      原文地址：http://www.zwqxin.cn/archives/image-processing/dif-fast-fourier.html