-
实用标准文案
从头到尾彻底理解傅里叶变换算法、上
前言
第一部分、
DFT
第一章、傅立叶变换的由来
第二章、实数形式离散傅立叶变换（
Real
DFT
）
从头到尾彻底理解傅里叶变换算法、下
第三章、复数
第四章、复数形式离散傅立叶变换
/**************************************
***********************************************
**************/
这一片的傅里叶变换算法，讲解透
彻，希望对大家会有所帮助。感谢原作者们
（
July
、
dznlong
）的精心编写。
/********************************
**************************************************
***
*************/
前言
：
“<
/p>
关于傅立叶变换，
无论是书本还是在网上可以很容易找到关于傅立
叶变换的描述，
但是大
文档
实用标准文案
都是些故弄玄虚的文章
，
太过抽象，
尽是一些让人看了就望而生畏的公式的罗列，
让人很难
能够从感性上得到理解
”
---dznlong
，
那么，到底什么是傅里叶变换算法列
?
傅里叶变换所涉及到
的公式具体有多复杂列
?
傅里叶变换
（
Fourier tra
nsform
）是一种线性的积分变换。因其基本思想首先由法国学者
< br>傅里叶系统地提出，所以以其名字来命名以示纪念。
哦，
傅里叶变换原来就是一种变换而已，
只是这种变换是从时间转换为频率的变化。
这下，
你就知道了
，傅里叶就是一种变换，一种什么变换列
?
就是一种从时间到频
率的变化或其相
互转化。
ok
，咱们再来总体了解下傅里叶变换，让各位对其有个总体大概的印象，也顺便看看傅里
叶变换所涉及到的公式，究竟有多复杂：
以下
就是傅里叶变换的
4
种变体（摘自，维基百科）
连续傅里叶变换
<
/p>
一般情况下，若“傅里叶变换”一词不加任何限定语，则指的是“连续傅里叶变换”。连<
/p>
续傅里叶变换将平方可积的函数
f
（
t
）表示成复指数函数的积分或级数形式。
文档
实用标准文案
这是将频率域的函数
F(
ω
)
表示为时间域的函数
f
（
< br>t
）的积分形式。
连续傅里叶变换的逆变换
(inverse Fourier
transform)
为：
即将时间域的函数
f
（
t
）表示为频率域的函数
F(
ω
)
的积分。
一般
可称函数
f
（
t
）
为原函数，而称函数
F(
ω
)
为傅里叶变换的像函数，原函数和像函数构
成一个傅里叶变换对（
transform
pair
）。
文档
实用标准文案
除此之外，还有其它型
式的变换对，以下两种型式亦常被使用。在通信或是信号处理方面，
常以
或者是因系数重分配而得到新的变换对：
文档
来代换，而形成新的变换对
:
实用标准文案
一种对连续傅里叶变换的推广称为分数傅里叶变换（
Fractional Fourier Transform
）。分<
/p>
数傅里叶变换
(fractional Fourier tra
nsform,FRFT)
指的就是傅里叶变换
(Fourie
r
transform,FT)
的广义化。
< br>
分数傅里叶变换的物理意义即做傅里叶变换
a
次，其中
a
不一定要为整数；而做了分数
傅里
叶变换之后，信号或输入函数便会出现在介于时域
(time
domain)
与频域
(frequency
domain)
之间的分数域
(fractional
domain)
。
当
f
（
t
）为偶函数（或奇函数
）时，其正弦（或余弦）分量将消亡，而可以称这时的变换为
余弦变换（
cosine
transform
）或正弦变换（
sine
transform
）
.
另一个值得
注意的性质是，当
f
（
t
）为纯实函数时，
F(
?ω
) = F*(
ω
)
成立
.
傅里叶级数
连续形式的傅里叶变换其实是傅里叶级数
(Fourier series)
的推广，因为积分其实是一
种
极限形式的求和算子而已。对于周期函数，其傅里叶级数是存在的：
< br>
文档
实用标准文案
其中
Fn
为复幅度。对于实值函数，函数的傅里叶级数可以
写成：
其中
an
和
bn
是实频率分量的幅度。<
/p>
离散时域傅里叶变换
离散傅里叶变换是离散时间傅里叶变换（
DTFT
）的特例（有时作为后者的近似）。
DTFT
在时域上离散，在频域上则是周期的。
DTFT
可以被看作是傅里叶级数的逆变换。
文档
实用标准文案
离散傅里叶变换
离散傅里叶变换（
DFT
），是连续傅里叶变换在时
域和频域上都离散的形式，将时域信号
的采样变换为在离散时间傅里叶变换（
DTFT
）频域的采样。在形式上，变换两端（时域和
频域上）
的序列是有限长的，
而实际上这两组序列都应当被
认为是离散周期信号的主值序列。
即使对有限长的离散信号作
D
FT
，也应当将其看作经过周期延拓成为周期信号再作变换。
在
实际应用中通常采用快速傅里叶变换以高效计算
DFT
。
为了在科学计算和数字信号处理等领
域使用计算机进行傅里叶变换，
必须将函数
xn
定义
在离散点而非连续域内，
且须满足有限性或周期性
条件。
这种情况下，
使用
离散傅里叶变
换
（
DFT
），
将函数
xn
表示为下面的求和形式：
其中
Xk
是
傅里叶幅度。直接使用这个公式计算的计算复杂度为
O
（
n*n
），而
快速傅里叶
变换（
FFT
）可以将复杂度改进为
O
（
n*lgn
）。（
后面会具体阐述
FFT
是如何将复杂度降
为
O
（
n*lgn
）
的。
）
计算复杂度
的降低以及数字电路计算能力的发展使得
DFT
成为在信号
处理领域十分实用且重要的方法。
文档
实用标准文案
下面，比较下上述傅立叶变换的
4
种变体，
如上，容易发
现：函数在时（频）域的离散对应于其像函数在频（时）域的周期性。反之
连续则意味着
在对应域的信号的非周期性。
也就是说，
时间上的离散性对应着
频率上的周期
性。同时，注意，离散时间傅里叶变换，时间离散，频率不离散，它在频域
依然是连续的。
如果，读到此
，你不甚明白，大没关系，不必纠结于以上
4
种变体，继续往下
看，你自
会豁然开朗。（有什么问题，也恳请提出，或者批评指正）
ok
，
本文，
接下来，由傅里叶变换入手，
后重点阐述离散傅
里叶变换、
快速傅里叶算法，
到最后彻底实现
< br>FFT
算法，全篇力求通俗易懂、阅读顺畅，教你从头到尾彻底理解傅里叶
变换算法。
由于傅里叶变换，也称傅立叶变换，下文所称为
傅立叶变换
，同一个变换，不同
叫法，读者不必感
到奇怪。
第一部分、
DFT
第一章、傅立叶变换的由来
要理解傅立叶变换，先得知道傅立叶变换是怎么变换的，当然，也需要一定的高等数学<
/p>
文档
实用标准文案
基础，最基本的是级数
变换，其中傅立叶级数变换是傅立叶变换的基础公式。
一、傅立叶变换的提出
傅立叶是一位法国数学家和物理学家，原名是
Jean
Baptiste Joseph
Fourier(1768-1830), Fou
rier
于
1807
年在法国科学学会
上发表了一篇论文，论文里描述
运用正弦曲线来描述温度分布，
论文里有个在当时具有争议性的决断：
任何连续周期信号都
可以
由一组适当的正弦曲线组合而成。
< br>当时审查这个论文拉格朗日坚决反对此论文的发表，而后在近
50
年的时间里，拉格朗
日坚持认为傅立叶的方法无法表示带有棱角的信号，
如在方波中出现非连续变化斜率。
直到
拉格朗
日死后
15
年这个论文才被发表出来。
谁是对的呢？拉格朗日是对的：正弦曲线无法组合成一
个带有棱角的信号。但是，我们
可以用正弦曲线来非常逼近地表示它，
< br>逼近到两种表示方法不存在能量差别，
基于此，
傅立
叶是对的。
为什么我们要用正弦曲线来代替原来的曲线呢？如我们也还可以用方波或三角波来代替
呀，分解信号的方法是无穷多的，但分解信号的目的是为了更加简单地处理原来的信号。
用正余弦来表示原信号会更加简单，因为正余弦拥
有原信号所不具有的性质：正弦曲线
保真度。
一个正余弦曲线信
号输入后，
输出的仍是正余弦曲线，
只有幅度和相位可能发生变
化，
但是频率和波的形状仍是一样的。
且只有正余弦曲线才拥有这样的性质，
正因如此我们
才不用方波
或三角波来表示。
文档
实用标准文案
二、傅立叶变换分类
根据原信号的不同类型，我们可以把傅立叶变换分为四种类别：
1
、非周期性连续信号
傅立叶变换（
Fourier
Transform
）
2
、周期性连续信号
傅立叶级数
(Fourier Series)
3
、非周期性离散信号
离散
时域
傅立叶变换（
Discrete Time Fourier
Transform
）
4
、周期性离散信号
离散傅立叶变换
(Discrete Fourier
Transform)
下图是四种原信号图
例（从上到下，依次是
FT
，
FS
，
DTFT
，
DF
T
）：
这四种傅立叶变换都是针对正无穷大和负无穷大的信号
，即信号的的长度是无穷大的，
我们知道这对于计算机处理来说是不可能的，
那么有没有针对长度有限的傅立叶变换呢？
没
有。
因为正余弦波被定义成从负无穷小到正无穷大，
我们无法把一个
长度无限的信号组合成
长度有限的信号。
面对这种困难，方法是：把长度有限的信号表示成长度
无限的信号。如，可以把信号无
限地从左右进行延伸，
延伸的部
分用零来表示，
这样，
这个信号就可以被看成是
非周期性
离
散
信号
，
我们可以用到
离散时域傅立叶变换（
DTFT
）
的方法。也可以把信号用复制的方法<
/p>
文档
实用标准文案
进行延伸，这样信号就
变成了
周期性
离散
信号，这时我们就可
以用
离散傅立叶变换方法
（
DFT
）
进行变换。本章我们要讲的是
离散信号<
/p>
，对于连续信号我们不作讨论，因为计算机
只能处理离散的数值信
号，我们的最终目的是运用计算机来处理信号的。
但是对于非周期性的信号，我们需要用无穷多不同频率
的正弦曲线来表示，这对于计算
机来说是不可能实现的。
所以对
于离散信号的变换只有
离散傅立叶变换
（
DFT
）
才能被适用，
对于计算机来
说只有离散的和有限长度的数据才能被处理，
对于其它的变换类型只有在数学
演算中才能用到，
在计算机面前我们只能用
DFT
方法，
后面我们要理解的也正是
DFT
方法。
这里要理解的是我们使用周期性的信号目的是为了能够用数学方法来解决问题，至于考
虑周期性信号是从哪里得到或怎样得到是无意义的。
每种傅立叶变换都分成实数和复数两种方法，对于实数
方法是最好理解的，但是复数方
法就相对复杂许多了，
需要懂得
有关复数的理论知识，
不过，
如果理解了实数离散傅立叶变
换
(real DFT)
，再去理解复数傅
立叶变换就更容易了，所以我们先把复数的傅立叶变换放到
一边去，
先来理解实数傅立叶变换，
在后面我们会先讲讲关于复数的基本理论，
然后在理解
了实数傅立叶变换的基础上再来理解复数傅立叶变换。
还有，这
里我们所要说的变换
(transform)
虽然是数学意义上
的变换，但跟函数变换是不
同的，函数变换是符合一一映射准则的，对于离散数字信号处
理（
DSP
），有许多的变换：
傅立叶
变换、拉普拉斯变换、
Z
变换、希尔伯特变换、离散余弦变换等
，这些都扩展了函数
变换的定义，
允许输入和输出有多种的值，
简单地说变换就是把一堆的数据变成另一堆的数
据的方法。
文档
实用标准文案
三、一个关于实数离散傅立叶变换
(Real
DFT)
的例子
先来看一个变换实例，下图是一个原始信号图像：
这个信号的长度是<
/p>
16
，
于是可以把这个信号分解
9
个余弦波和
9
个正弦
波
（一个长度
为
N
的信号可以分解成
N/2+1
个正余弦信号，
这是为什么呢？结合下面的
18
个正余弦图
,
我想从计算机处理精度上就不难理解，
一个长
度为
N
的信号，
最多只能有
N/2+1
个不同频
率，再多的频率就超过了计算
机所能所处理的精度范围），如下图：
9
个余弦信号：
文档
实用标准文案
9
个正弦信号：
把以上所有信号相加
即可得到原始信号，
至于是怎么分别变换出
9
< br>种不同频率信号的，
我们先不急，
先看看对于以上的变换
结果，
在程序中又是该怎么表示的，
我们可以看看下面
这个示例图：
文档
实用标准文案
上图中左边表示时域中的信号，右边是频域信号表示方法，
<
/p>
从左向右，
-->
，
表示
正向转换
(Forward DFT)
，
从右向左，
<--
，
表示
逆向转换
(Inverse
DFT)
，
用小写
< br>x[]
表示信号在每个时间点上的幅度值数组
,
用大写
X[]
表示每种频率的副度值数组
（即时间
x-->
频率
X
）
,
因为有
N/2+1
种频率，所以该数组长度为
N/2+1
，
X[
]
数组又分两种，一种是表示余弦波的不同频率幅度值：
Re
X[]
，
另一种是表示正弦波的不同频率幅度值：
Im
X[]
，
< br>Re
是实数
(Real)
的意思
，
Im
是虚数
(Imagine)
的意思，采用复数的表示方法把正余弦波
组合起来进行表示，
但这里我们不考虑复数的其它作用，
只记住是一种组合方法而已，
目的
是为了便于表达（在后面我们会知道，复数形式的傅立叶变换长
度是
N
，而不是
N/2+1
）。
如此，再回过头去，看上面的正余弦各
9
种频率的变化，相信，问题不大了。
文档
实用标准文案
第二章、实数形式离散傅立叶变换（
Real
DFT
）
上一章，我们看到了一个实数形式离散傅立叶变换的例子，通过这个例子能够让我们
先对傅立叶变换有一个较为形象的感性认识，
现在就让我们来看看实数形式离散
傅立叶变换
的正向和逆向是怎么进行变换的。在此，我们先来看一下频率的多种表示方法
。
一、
频域中关于频率的四种表示方法
1
、序号表示方法，根据时域中信号的样本数取
0 ~ N/2
，用这种方法在程序中使用起来可
以更
直接地取得每种频率的幅度值，因为频率值跟数组的序号是一一对应的
:
X[k]
，取值范
围是
0 ~
N/2
；
2
、分数表示方法，根据时域中信号的样本数的比例值取
0 ~ 0.5:
X[
?
]
，?
=
k/N
，取值范
围是
0 ~
1/2
；
3
、
用弧度值来表示，
把?乘以一个
2<
/p>
π得到一个弧度值，
这种表示方法叫做自然频率
< br>(natural
frequency)
：
X[
ω
]
，ω
= 2
π?
= 2
π
k/N
，取值范围是
0 ~
π；
4
、以赫兹
(Hz)
为单位来表示，这个一般是应用于一些特殊应用，如取样率为
10 kHz
表示
每秒有
10,000
个
样本数：取值范围是
0
到取样率的一半。
二、
DFT
基本函数
ck[i] =
cos(2
π
ki/N)
sk[i]
= sin(2
π
ki/N)
文档
实用标准文案
其中
k
表示每个正余弦波的频率，
< br>如为
2
表示在
0
到
N
长度中存在两个完整的周期，
10
即有
10
个周期，如下图：<
/p>
上图中至于每个波的振幅
(amplitude)
值
(Re X[k],Im X[k])
是怎么算出来的
,
这个是
DFT
的核心，也是
最难理解的部分，我们先来看看如何把分解出来的正余弦波合成原始信
号
(Inverse DFT)
。
三、
合成运算方法
(Real Inverse DFT)
DFT
合成等式（合成原始
时间
信号，频率
-->
时间，逆向变换
）：
文档
实用标准文案
如果有学过傅立叶级数，
对这个等式就会有似曾相识的感觉，
不错！
这个等式跟傅立叶级数
是非常相似的：
当然，差别是肯定是存在的，因为这两个等式是在两个不同条件下运用的，至于怎
么证明
DFT
合成公式，这个我想需要非常强的高等数学理论知识了，
这是研究数学的人的
工作，
对于普通应用者就不需要如此的追根
究底了，
但是傅立叶级数是好理解的，
我们起码
可以从傅立叶级数公式中看出
DFT
合成公式的合理性
。
_
_
文档
实用标准文案
DFT
合成等式中的
Im
X[k]
和
Re
X[k]
跟之前提到的
Im
X[k]
和
Re X[k]
是不一样的
，
下
面是转换方法（
关于此公式的解释
，见下文
）
:
但<
/p>
k
等于
0
和
N/2
时
,
实数部
分的计算要用下面的等式
:
上面四个式中的
N
是时域中点的总数，
k
是从
0
到
N/2
的序号。
文档
实用标准文案
为什么要这样进行转换呢？这个可以从频谱密度
(spectral
density)
得到理解，如下图
就是个频谱图：
这是一个频谱图，
横坐标表示频率大小，
纵坐标表示振幅大小，
原始信号长度为
N
（这
里是
32
），经
DFT
转换后得到的
17
个
频率的频谱，频谱密度表示每单位带宽中为多大的
振幅，
那么带
宽是怎么计算出来的呢？看上图，
除了头尾两个，
其余点的所占
的宽度是
2/N
，
这个宽度便是每个点
的带宽，头尾两个点的带宽是
1/N,
而
Im X[k]
和
Re X[k]
表
示的是频谱
密度，即每一个单位带宽的振幅大小，但
文档
表
实用标准文案
示
2/N
（或
1/N
）带宽的振幅大小，
所以
别应当
是
Im X[k]
和
Re X[k]<
/p>
的
2/N
（或
1
/N
）。
分
频谱密度就象物理中物质密度，
原始信号中的每一个点就象
是一个混合物，
这个混合物是由
不同密度的物质组成的，
混合物中含有的每种物质的质量是一样的，
除了最大和最小两个密
度的物质外，
这样我们只要把每种物质的密度加起来就可以得到该混
合物的密度了，
又该混
合物的质量是单位质量，所以得到的密度
值跟该混合物的质量值是一样的。
至于为什么虚数部分是负数，这是为了跟复数<
/p>
DFT
保持一致，这个我们将在后面会知
道这是数学计算上的需要（
Im X[k]
在计算时就已经加上
了一个负号（
稍后，由下文，便可
文档
实用标准文案
知
），
有变化）。
再加上负号，结果便是正的，等于没
如果已经得到了
DFT
结果，这时要进行
逆转换，即合成原始信号
，则可按如下步骤进
行
转换：
1
、先根据上面四个式子计算
得出
2
、再根据
DFT
合成等式得到原始信号数据。
下面是用
BASIC
语言来实现的转换源代码：
100
‘
DFT
逆转换方法
的值；
文档
实用标准文案
110
‘
/XX[]
数组存储计算结果（时域中的原始信号
）
120
‘
/REX[]
数组存储频域中的实数分量，
IMX[]
为虚分量
130
‘
140 DIM XX[511]
150 DIM REX[256]
160 DIM
IMX[256]
170
‘
180 PI = 3.14159265
190 N% =
512
200
‘
210 GOSUB XXXX
‘转到子函数去获取
REX[]
和
IMX[]
数据
220
‘
230
‘
240
‘
250 FOR K% = 0
TO 256
260
REX[K%] =
REX[K%] / (N%/2)
270
IMX[K%] = -IMX[K%] / (N%/2)
280 NEXT k%
290
‘
300 REX[0] =
REX[0] / N
310 REX[256] = REX[256] / N
320
‘
文档
实用标准文案
330
‘
初始化
X
X[]
数组
340 FOR I%
= 0 TO 511
350
XX[I%] = 0
360 NEXT I%
370
‘
380
‘
390
‘
400
‘
410
‘
420 FOR K% =0
TO 256
430
FOR I%=0 TO 511
440
‘
450
XX[I%] = XX[I%] +
REX[K%] * COS(2 * PI * K% * I% / N%)
460
XX[I%] = XX[I%] +
IMX[K%] * SIN(2 * PI * K% * I% / N%)
470
‘
480
NEXT I%
490
NEXT K%
500
‘
510 END
上面代码中
420
至
490
换
成如下形式也许更好理解，但结果都是一样的：
420
FOR I% =0 TO 511
文档
-
-
-
-
-
-
-
-
-
上一篇：SCI收录期刊——土壤科学学科
下一篇：SCI收录的园艺类期刊