傅里叶级数 傅里叶变换 及应用

傅里叶级数和傅立叶变换是傅里叶分析的两个主要工具,它们之间有密切的关系。

什么是傅里叶级数

傅里叶级数是将一个周期函数分解为一系列正弦和余弦函数的和。它适用于周期性信号,可以将周期函数表示为一组振幅和相位不同的谐波分量的和。傅里叶级数展示了一个周期函数在不同频率上的频谱内容。

我们平时理解的更多是一维中的加减成除,给你任意一个数A,你都可以找到一些数线性组合相加等于A, A = a n B + b n C + … A = a_n B+b_nC+dots A=anB+bnC+,同样的傅里叶就是多维的加法去拟合函数。

傅立叶级数可以用以下公式表示:

f T ( x ) = a 0 2 + ∑ n = 1 N ( a n cos ⁡ ( 2 π n T x ) + b n sin ⁡ ( 2 π n T x ) ) f_T(x) = frac{a_0}{2} + sum_{n=1}^{N} left(a_n cosleft(frac{2pi n}{T}xright) + b_n sinleft(frac{2pi n}{T}xright)right) fT(x)=2a0+n=1N(ancos(T2πnx)+bnsin(T2πnx))

展开写理解起来会容 易一些:


f ( x 0 ) = a 0 2 + a 1 cos ⁡ ( 2 π n 1 T x 0 ) + b 1 sin ⁡ ( 2 π n 1 T x 0 ) + a 2 cos ⁡ ( 2 π n 2 T x 0 ) + b 2 sin ⁡ ( 2 π n 2 T x 0 ) + … + a N cos ⁡ ( 2 π N T x 0 ) + N sin ⁡ ( 2 π N T x 0 ) f(x_0)= frac{a_0}{2}+a_1cosleft(frac{2pi n_1}{T}x_0right) + b_1 sinleft(frac{2pi n_1 }{T}x_0right) +a_2cosleft(frac{2pi n_2}{T}x_0right) + b_2 sinleft(frac{2pi n_2}{T}x_0right) \ +dots \ +a_Ncosleft(frac{2pi N}{T}x_0right) + N sinleft(frac{2pi N}{T}x_0right) f(x0)=2a0+a1cos(T2πn1x0)+b1sin(T2πn1x0)+a2cos(T2πn2x0)+b2sin(T2πn2x0)++aNcos(T2πNx0)+Nsin(T2πNx0)


f ( x 1 ) = a 0 2 + a 1 cos ⁡ ( 2 π n 1 T x 1 ) + b 1 sin ⁡ ( 2 π n 1 T x 1 ) + a 2 cos ⁡ ( 2 π n 2 T x 1 ) + b 2 sin ⁡ ( 2 π n 2 T x 1 ) + … + a N cos ⁡ ( 2 π N T x 1 ) + N sin ⁡ ( 2 π N T x 1 ) f(x_1)= frac{a_0}{2}+a_1cosleft(frac{2pi n_1}{T}x_1right) + b_1 sinleft(frac{2pi n_1 }{T}x_1right) +a_2cosleft(frac{2pi n_2}{T}x_1right) + b_2 sinleft(frac{2pi n_2}{T}x_1right) \ +dots \ +a_Ncosleft(frac{2pi N}{T}x_1right) + N sinleft(frac{2pi N}{T}x_1right) f(x1)=2a0+a1cos(T2πn1x1)+b1sin(T2πn1x1)+a2cos(T2πn2x1)+b2sin(T2πn2x1)++aNcos(T2πNx1)+Nsin(T2πNx1)


f ( x 2 ) = a 0 2 + a 1 cos ⁡ ( 2 π n 1 T x 2 ) + b 1 sin ⁡ ( 2 π n 1 T x 2 ) + a 2 cos ⁡ ( 2 π n 2 T x 2 ) + b 2 sin ⁡ ( 2 π n 2 T x 2 ) + … + a N cos ⁡ ( 2 π N T x 2 ) + N sin ⁡ ( 2 π N T x 2 ) f(x_2)= frac{a_0}{2}+a_1cosleft(frac{2pi n_1}{T}x_2right) + b_1 sinleft(frac{2pi n_1 }{T}x_2right) +a_2cosleft(frac{2pi n_2}{T}x_2right) + b_2 sinleft(frac{2pi n_2}{T}x_2right) \ +dots \ +a_Ncosleft(frac{2pi N}{T}x_2right) + N sinleft(frac{2pi N}{T}x_2right) f(x2)=2a0+a1cos(T2πn1x2)+b1sin(T2πn1x2)+a2cos(T2πn2x2)+b2sin(T2πn2x2)++aNcos(T2πNx2)+Nsin(T2πNx2)


⋮ vdots


f ( x N ) = a 0 2 + a 1 cos ⁡ ( 2 π n 1 T x N ) + b 1 sin ⁡ ( 2 π n 1 T x N ) + a 2 cos ⁡ ( 2 π n 2 T x N ) + b 2 sin ⁡ ( 2 π n 2 T x N ) + … + a N cos ⁡ ( 2 π N T x N ) + N sin ⁡ ( 2 π N T x N ) f(x_N)= frac{a_0}{2}+a_1cosleft(frac{2pi n_1}{T}x_Nright) + b_1 sinleft(frac{2pi n_1 }{T}x_Nright) +a_2cosleft(frac{2pi n_2}{T}x_Nright) + b_2 sinleft(frac{2pi n_2}{T}x_Nright) \ +dots \ +a_Ncosleft(frac{2pi N}{T}x_Nright) + N sinleft(frac{2pi N}{T}x_Nright) f(xN)=2a0+a1cos(T2πn1xN)+b1sin(T2πn1xN)+a2cos(T2πn2xN)+b2sin(T2πn2xN)++aNcos(T2πNxN)+Nsin(T2πNxN)


这里的合成的概念是时域上的叠加的概念,图片来源wikipedia

动态合成图

其中, f ( x ) f(x) f(x)是周期为 T T T的函数, a 0 a_0 a0是直流分量(平均值), a n a_n an b n b_n bn是傅立叶系数,表示在频率为 2 π n T frac{2pi n}{T} T2πn的正弦和余弦函数的振幅。这个公式表示了周期函数 f ( x ) f(x) f(x)可以由一系列不同频率的正弦和余弦函数的和来逼近。

3dft

请添加图片描述

为什么可以这样展开

傅里叶级数的证明是一个相对复杂的过程,涉及到数学分析和傅里叶变换的理论。这里我将提供一个简要的概述,介绍傅里叶级数的基本思想和证明思路。
傅里叶级数的证明基于以下两个主要思想: 正交性:傅里叶级数利用了三角函数的正交性质。正弦和余弦函数在一个周期内是正交的,即它们的内积在不同频率下为零,而在相同频率下非零。这意味着不同频率的正弦和余弦函数是线性无关的基函数,可以用来表示不同频率的分量。 逼近性:傅里叶级数的另一个关键思想是逼近函数的概念。根据逼近定理,任何一个具有有限能量的周期函数,都可以用无限多个正弦和余弦函数的线性组合来逼近。通过增加正弦和余弦函数的振幅和相位,可以越来越接近原始函数。基于这两个思想,傅里叶级数的证明过程可以大致概括为以下步骤: 首先,我们考虑一个周期为T的函数f(x)。我们要证明,可以将它表示为一系列正弦和余弦函数的和。 假设f(x)可以表示为以下形式的级数:
f ( x ) = a 0 2 + ∑ n = 1 ∞ ( a n cos ⁡ ( ω n x ) + b n sin ⁡ ( ω n x ) ) f(x) = frac{a_0}{2} + sum_{n=1}^{infty} (a_n cos(omega_n x) + b_n sin(omega_n x)) f(x)=2a0+n=1(ancos(ωnx)+bnsin(ωnx))
其中, a 0 a_0 a0 a n a_n an b n b_n bn是待定的系数, ω n omega_n ωn是频率。

接下来,我们将使用正交性质来计算各个系数。通过将f(x)与正弦和余弦函数进行内积运算,并利用正交性质,我们可以得到各个系数的表达式。

通过计算内积,我们可以得到以下公式

a 0 = 1 π ∫ − π π f ( x ) cos ⁡ ( ω n x ) d x , ω n = 0 a_0 = frac{1}{pi} int_{-pi}^{pi} f(x) cos(omega_n x)dx, omega_n = 0 a0=π1ππf(x)cos(ωnx)dx,ωn=0

a n = 1 π ∫ − π π f ( x ) cos ⁡ ( ω n x ) d x a_n = frac{1}{pi} int_{-pi}^{pi} f(x) cos(omega_n x) dx an=π1ππf(x)cos(ωnx)dx

b n = 1 π ∫ − π π f ( x ) sin ⁡ ( ω n x ) d x b_n = frac{1}{pi} int_{-pi}^{pi} f(x) sin(omega_n x) dx bn=π1ππf(x)sin(ωnx)dx

ω n = 2 π n T omega_n = frac{2pi n}{T} ωn=T2πn

其中, a 0 a_0 a0是直流分量, a n a_n an b n b_n bn是频率为 ω n omega_n ωn的余弦和正弦分量的振幅。

由欧拉公式:

cos ⁡ θ = e i θ + e − i θ 2 costheta = frac{e^{itheta}+e^{-itheta}}{2} cosθ=2eiθ+eiθ

sin ⁡ θ = − i e i θ − e − i θ 2 sintheta =-i frac{e^{itheta}-e^{-itheta}}{2} sinθ=i2eiθeiθ

下面的证明摘抄自复变函数与积分变换的教材:

请添加图片描述
请添加图片描述
在这里插入图片描述
在这里插入图片描述

F ( ω ) = ∫ − ∞ ∞ f ( t ) e − j ω t d t F(omega)=int_{-infty}^{infty}f(t)e^{-jomega t}dt F(ω)=f(t)etdt

就是傅里叶变换了

什么是傅里叶变换

傅里叶变换则是将非周期性函数或信号分解为一组连续的正弦和余弦函数(复指数函数)的积分。傅立叶变换将一个时域函数转换为频域函数,显示了信号在不同频率上的频谱特征。

F ( ω ) = ∫ − ∞ ∞ f ( t ) e − j ω t d t F(omega)=int_{-infty}^{infty}f(t)e^{-jomega t}dt F(ω)=f(t)etdt
1 2 π ∫ − ∞ ∞ F ( ω ) e j ω t d ω frac{1}{2pi}int_{-infty}^{infty}F(omega)e^{jomega t}domega 2π1F(ω)etdω


f T ( x ) = a 0 2 + ∑ n = 1 N ( a n cos ⁡ ( 2 π n T x ) + b n sin ⁡ ( 2 π n T x ) ) f_T(x) = frac{a_0}{2} + sum_{n=1}^{N} left(a_n cosleft(frac{2pi n}{T}xright) + b_n sinleft(frac{2pi n}{T}xright)right) fT(x)=2a0+n=1N(ancos(T2πnx)+bnsin(T2πnx))

两者之间的关系可以用以下方式表达:

傅里叶级数是傅立叶变换的一种特殊情况。当一个周期为T的函数被表示为傅里叶级数时,其傅立叶变换是一个离散的频谱,只包含一系列离散的频率分量。

对于一个连续非周期函数,傅立叶变换将其表示为一个连续的频谱,包含了所有可能的频率分量。这相当于将傅里叶级数中的频率分量离散化为连续的频谱。

傅立叶变换可以视为傅里叶级数的极限情况,当周期趋向于无穷大时,傅里叶级数的频率间隔趋向于零,从而得到了连续的频谱。

总的来说,傅里叶级数和傅立叶变换是傅里叶分析的两种表示方法,它们可以相互转换,但适用于不同的函数类别和分析需求。傅里叶级数适用于周期函数的频谱分析,而傅立叶变换适用于非周期函数的频谱分析。

FFT快速傅里叶变换

DFT定义:

X k = ∑ n = 0 N − 1 x n e − i 2 π n N k , k = 0 , 1 , … , N − 1 X_k = sum_{n=0}^{N-1}x_ne^{-ifrac{2pi n}{N}k},k = 0,1,dots,N-1 Xk=n=0N1xneiN2πnkk=0,1,,N1

展开各项:

X 0 = x 0 e − i 2 π 0 N k + x 1 e − i 2 π 1 N k + ⋯ + x N − 1 e − i 2 π ( N − 1 ) N k , k = 0 X_0=x_0e^{-ifrac{2pi 0}{N}k}+x_1e^{-ifrac{2pi 1}{N}k}+dots+x_{N-1}e^{-ifrac{2pi (N-1)}{N}k},k=0 X0=x0eiN2π0k+x1eiN2π1k++xN1eiN2π(N1)k,k=0
X 1 = x 0 e − i 2 π 0 N k + x 1 e − i 2 π 1 N k + ⋯ + x N − 1 e − i 2 π ( N − 1 ) N k , k = 1 X_1=x_0e^{-ifrac{2pi 0}{N}k}+x_1e^{-ifrac{2pi 1}{N}k}+dots+x_{N-1}e^{-ifrac{2pi (N-1)}{N}k},k=1 X1=x0eiN2π0k+x1eiN2π1k++xN1eiN2π(N1)k,k=1
… dots
X N − 1 = x 0 e − i 2 π 0 N k + x 1 e − i 2 π 1 N k + ⋯ + x N − 1 e − i 2 π ( N − 1 ) N k , k = N − 1 X_{N-1}=x_0e^{-ifrac{2pi 0}{N}k}+x_1e^{-ifrac{2pi 1}{N}k}+dots+x_{N-1}e^{-ifrac{2pi (N-1)}{N}k},k=N-1 XN1=x0eiN2π0k+x1eiN2π1k++xN1eiN2π(N1)k,k=N1

还是不够直观?举个离散数据的例子来看:

x 0 = 1 x_0 = 1 x0=1
x 1 = 2 x_1 = 2 x1=2
x 3 = 3 x_3 = 3 x3=3
x 4 = 4 x_4 = 4 x4=4
X 0 = x 0 e − i 2 π 0 N k + x 1 e − i 2 π 1 N k + ⋯ + x N − 1 e − i 2 π ( N − 1 ) N k , k = 0 , N = 4 X_0=x_0e^{-ifrac{2pi 0}{N}k}+x_1e^{-ifrac{2pi 1}{N}k}+dots+x_{N-1}e^{-ifrac{2pi (N-1)}{N}k},k=0,N=4 X0=x0eiN2π0k+x1eiN2π1k++xN1eiN2π(N1)k,k=0,N=4
X 1 = x 0 e − i 2 π 0 N k + x 1 e − i 2 π 1 N k + ⋯ + x N − 1 e − i 2 π ( N − 1 ) N k , k = 1 , N = 4 X_1=x_0e^{-ifrac{2pi 0}{N}k}+x_1e^{-ifrac{2pi 1}{N}k}+dots+x_{N-1}e^{-ifrac{2pi (N-1)}{N}k},k=1,N=4 X1=x0eiN2π0k+x1eiN2π1k++xN1eiN2π(N1)k,k=1,N=4
X 2 = x 0 e − i 2 π 0 N k + x 1 e − i 2 π 1 N k + ⋯ + x N − 1 e − i 2 π ( N − 1 ) N k , k = 2 , N = 4 X_2=x_0e^{-ifrac{2pi 0}{N}k}+x_1e^{-ifrac{2pi 1}{N}k}+dots+x_{N-1}e^{-ifrac{2pi (N-1)}{N}k},k=2,N=4 X2=x0eiN2π0k+x1eiN2π1k++xN1eiN2π(N1)k,k=2,N=4
X 3 = x 0 e − i 2 π 0 N k + x 1 e − i 2 π 1 N k + ⋯ + x N − 1 e − i 2 π ( N − 1 ) N k , k = 3 , N = 4 X_3=x_0e^{-ifrac{2pi 0}{N}k}+x_1e^{-ifrac{2pi 1}{N}k}+dots+x_{N-1}e^{-ifrac{2pi (N-1)}{N}k},k=3,N=4 X3=x0eiN2π0k+x1eiN2π1k++xN1eiN2π(N1)k,k=3,N=4

X 0 = 1 e − i 2 π 0 4 0 + 2 e − i 2 π 1 4 0 + 3 e − i 2 π ( 2 ) 4 0 + 4 e − i 2 π ( 3 ) 4 0 = 10 + 0 i X_0=1e^{-ifrac{2pi 0}{4}0}+2e^{-ifrac{2pi 1}{4}0}+3e^{-ifrac{2pi (2)}{4}0}+4e^{-ifrac{2pi (3)}{4}0}=10+0i X0=1ei42π00+2ei42π10+3ei42π(2)0+4ei42π(3)0=10+0i
X 2 = 1 e − i 2 π 0 4 1 + 2 e − i 2 π 1 4 1 + 3 e − i 2 π ( 2 ) 4 1 + 4 e − i 2 π ( 3 ) 4 1 = − 2 + 2 i X_2=1e^{-ifrac{2pi 0}{4}1}+2e^{-ifrac{2pi 1}{4}1}+3e^{-ifrac{2pi (2)}{4}1}+4e^{-ifrac{2pi (3)}{4}1}=-2+2i X2=1ei42π01+2ei42π11+3ei42π(2)1+4ei42π(3)1=2+2i
X 0 = 1 e − i 2 π 0 4 2 + 2 e − i 2 π 1 4 2 + 3 e − i 2 π ( 2 ) 4 2 + 4 e − i 2 π ( 3 ) 4 2 = − 2 + 0 i X_0=1e^{-ifrac{2pi 0}{4}2}+2e^{-ifrac{2pi 1}{4}2}+3e^{-ifrac{2pi (2)}{4}2}+4e^{-ifrac{2pi (3)}{4}2}=-2+0i X0=1ei42π02+2ei42π12+3ei42π(2)2+4ei42π(3)2=2+0i
X 0 = 1 e − i 2 π 0 4 3 + 2 e − i 2 π 1 4 3 + 3 e − i 2 π ( 2 ) 4 3 + 4 e − i 2 π ( 3 ) 4 3 = − 2 − 2 i X_0=1e^{-ifrac{2pi 0}{4}3}+2e^{-ifrac{2pi 1}{4}3}+3e^{-ifrac{2pi (2)}{4}3}+4e^{-ifrac{2pi (3)}{4}3}=-2-2i X0=1ei42π03+2ei42π13+3ei42π(2)3+4ei42π(3)3=22i

令: w = e − i 2 π N w = e^{frac{-i2pi}{N}} w=eNi2π
把这个案例写成矩阵的形式:
[ X 0 X 1 X 2 X 3 ] = [ 1 1 1 1 1 w w 2 w 3 1 w 2 w 4 w 6 1 w 3 w 6 w 9 ] [ x 0 x 1 x 2 x 3 ] begin{bmatrix} X_0 \ X_1 \ X_2 \ X_3 \ end{bmatrix} = begin{bmatrix} 1 & 1 &1 &1 \ 1 &w&w^2&w^3 \ 1 &w^2 & w^4 &w^6 \ 1 & w^3 & w^6 & w^9 end{bmatrix} begin{bmatrix} x_0 \ x_1 \ x_2 \ x_3 end{bmatrix} X0X1X2X3 = 11111ww2w31w2w4w61w3w6w9 x0x1x2x3
完成这4个点DFT,需计算 16次乘法,即矩阵F中每个元素都会参与一次乘法,共 4个元素;12次加法,即每完成一个点需要3次加法。直到现在DFT终于高懂了,但是计算起来太麻烦了,效率太低。数据越大,计算量越大。FFT出现解决了计算量大的问题。

w = e − i 2 π N w = e^{frac{-i2pi}{N}} w=eNi2π 具有对称性,这是FFT简化的计算的方法,即:
W N k + N 2 = − W N k W_N^{k+frac{N}{2}} = -W_N^k WNk+2N=WNk

X k = ∑ n = 0 N − 1 x n e − i 2 π n N k X_k = sum_{n=0}^{N-1}x_ne^{-ifrac{2pi n}{N}k} Xk=n=0N1xneiN2πnk
再看这个公式,变化即是信号值乘以对应的复数 e − i 2 π n N k e^{-ifrac{2pi n}{N}k} eiN2πnk再求累加。
根据对称性可以简化傅里叶的计算,这就是FFT的意义。

分治法实现

后续添加……

倍增法实现

后续添加……

位逆序置换

后续添加……

蝶形运算优化

后续添加……


X k = ∑ n = 0 N − 1 x n e − i 2 π n N k , k = 0 , 1 , … , N − 1 X_k = sum_{n=0}^{N-1}x_ne^{-ifrac{2pi n}{N}k},k = 0,1,dots,N-1 Xk=n=0N1xneiN2πnkk=0,1,,N1
X k = ∑ n = 0 N − 1 x ( n ) × cos ⁡ ( 2 π n N k ) − i ∑ n = 0 N − 1 x ( n ) × sin ⁡ ( 2 π n N k ) X_k = sum_{n=0}^{N-1}x(n)timescos(frac{2pi n}{N}k)-isum_{n=0}^{N-1}x(n)timessin(frac{2pi n}{N}k) Xk=n=0N1x(n)×cos(N2πnk)in=0N1x(n)×sin(N2πnk)

应用

等有时间再作整理……

C++代码实现

等有时间再作整理……