深入理解傅里叶变换
这原本是我在知乎上对傅立叶变换、拉普拉斯变换、Z变换的联系?为什么要进行这些变换。研究的都是什么?问题的回答,实际上是我在本科学习数学和信号处理期间的思考,知乎上的答案因为写得仓促,只写了一些大致思想,没有具体展开,也没有图,比较难以理解,这里重新整理了一下,汇成此文。
本文要求读者需要在对傅里叶变换有一定的了解的基础之上阅读,至少要知道怎么算傅里叶变换。 此外部分地方要求读者有一定的微分方程基础,至少会求简谐振子的二阶常微分方程吧。
什么是傅里叶变换
高等数学中一般是从周期函数的傅里叶级数开始介绍的,这里也不例外。简单的说,从高中我们就学过一个理想的波可以用三角函数来描述,但是实际上的波可以是各种奇形怪状的。首先我们来看具有固定周期的波,下图中展示了4种常见的周期波。傅里叶级数告诉我们,这些周期信号都可以分解为有限或无限个正弦波或余弦波的叠加,且这些波的频率都是原始信号频率\(f_0\)的整数倍。
\[s(x) = \frac{A_0}{2} + \sum_{n=1}^{\infty} A_n\cdot \sin(2\pi n f_0 x+\phi_n).\]这里\(f_0\)被称为这些波的基频,\(A_0/2\)代表直流系数,系数\(A_n\)被称为幅度,\(\phi_n\)被称作相位。根据幅度和相位可以利用反变换恢复信号的波形,因此幅度和相位包含了信号的全部信息。这里的幅度关于频率的函数,我们称之为频谱,相位关于频率的函数,称之为相位谱。
下图是矩形波分解为多个正弦波的示意图,随着正弦波数目的增加,可以无限地逼近矩形波。 对于非周期信号,我们不能简单地将它展开为可数个正弦波的叠加,但是可以利用傅里叶变换展开为不可数的正弦波的叠加,其表达式可以通过\(f_0 \rightarrow \infty\)简单得到。
\[\displaystyle {\hat {f}}(\xi )=\int_{-\infty }^{\infty }f(x)\ e^{-2\pi ix\xi }\,dx,\]我们日常遇到的琴音、震动等都可以分解为正弦波的叠加,电路中的周期电压信号等信号都可以分解为正弦波的叠加。 那么问题来了,为什么我们要将信号分解为正弦波的叠加呢?这里面包含两个问题,为什么要分解?为什么是正弦波(或余弦波),可不可以是其他的波?另一个问题是对通信的同学的,我们学过多个变换那么这些变换之间有哪些关系? 在下面的篇章中,我将回答这三个问题。
为什么要分解为正弦波的叠加
这个问题可以追溯到傅里叶变换的创始人傅里叶解热传导方程的时候,因为热传导方程要求读者对热力学有一定了解,这里我以简谐振子系统为例来说明这个问题。没有阻尼的简谐振子系统可以用下面这个微分方程来描述
\[{\frac {\mathrm {d} ^{2}x}{\mathrm {d} t^{2}}}+2 \omega_{0}{\frac {\mathrm {d} x}{\mathrm {d} t}}+\omega_{0}^{2}x=F(t).\]\(x, t, \omega_0,F\)分别代表位移、时间、系统固有频率和外界驱动力。当没有外界驱动力\(F\)时,这个系统有通解
\[x(t) = A \sin\left(\omega_0 t + \phi \right)\]现在我们考虑存在外界驱动力\(F\)的场景,熟悉常微分方程理论的可以知道此时的通解是上述其次方程的通解(\(F\)恒为0)加上一个特解,所谓特解就是某个满足上述非齐次方程(\(F\)不恒为0)的任意一个接!那为什么能做这种分解呢?原因在于这是一个线性系统,或者说这个方程是一个线性方程,因此遵循叠加原理,可以简单的证明这个一般性结论。假设线性系统可以由线性微分方程来描述
\[\hat{L} x(t) = F(t)\]\(\hat{L}\)是线性算子,你可以简单地理解为谐振子方程中的左边操作。如果\(C(t), x_0(t)\)分别是其次方程通解和非齐次方程特解,即他们满足
\[\hat{L} C(t) = 0, \\ \hat{L} x_0(t) = F(t).\]那么将这两个式子相加,就可以得到
\[\hat{L} \left( x_0(t) + C(t) \right) = F(t)\]因此,只剩下一个问题,对于给定的驱动力\(F(t)\),怎么找特解的问题了。 也许你还记得在高数的书上,对\(F(t)\)为三角函数和指数函数时,可以有和\(F(t)\)形式相同的特解。 例如\(F(t)=f_0 \sin(wt)\)时,可以假定非齐次方程也有这种形式的特解\(B \sin(wt)\),代入原方程,求出待定常数可得特解\(\frac{A}{(w - w_0)^2} \sin(wt)\)。指数形式的驱动力也类似,那么对于其他形式的驱动力,怎么求特解呢?很简单,利用线性叠加原理,我如果求出很多个\(F\)为正弦驱动力\(\sin(w_n t)\)下的特解\(x_n(t)\),并且如果\(F\)可以表达为这些正弦波的叠加,那么特解不就可以用这些特解的叠加得到了么?用数学语言表述就是
\[\hat{L} x_n(t) = \sin(w_n t), n=0,1... \\ \hat{L} \sum_n A_n x_n(t) = \sum_n A_n \sin(w_n t).\]上面第二个式子右边如果等于\(F(t)\),那么左边的\(\sum_n A_n x_n(t)\)就是原齐次方程的特解。 简单地说,就是将驱动力做傅里叶变换(如果是周期驱动力则展开为傅里叶级数),求出每个基驱动力的特解,然后叠加得到特解。当然实际求解不用那么绕,以简谐振动方程为例,直接对方程左右两边做傅里叶变换即得
\[w^2 \hat{X}(w) - 2w\omega_0 \hat{X}(w) + \omega_0^2 \hat{X}(w) = \hat{F}(w)\]上式带尖头的函数代表时域函数的傅里叶变换,这是一个代数方程,容易求得
\[\hat{X}(w) = \frac{\hat{F}(w)}{(w - \omega_0)^2}\]通过上述描述,我们可以看到,将一个函数做傅里叶变换或者展开为傅里叶级数,可以帮助我们求解线性微分方程,或者从实际意义来说,可以帮助我们分析一个线性系统对外界做出如何响应!之所以能这样展开,是因为我们分析的是线性系统,如果是非线性系统就不能这样操作了。至于为什么是三角函数,我在下面将会回答,接下来我们先来看看更多的例子。
傅里叶变换与信号系统
这里,我们对通信相关的领域再举一个例子,来说明展开为三角函数(或者复指数函数)的重要性。这种分析,我们称之为傅里叶分析,或者叫频谱分析。
一个信号,通常用一个时间的函数\(x(t)\)来表示,这样简单直观,因为它的函数图像可以看做信号的波形,比如声波和水波等等。很多时候,对信号的处理是很特殊的,比如说线性电路会将输入的正弦信号处理后,输出仍然是正弦信号,只是幅度和相位有一个变化。这是因为线性电路都可以用常系数线性微分方程来描述,输入信号可以看做外界驱动力,输出可以看做系统地响应,这和上面的谐振子方程类似。因此,如果我们将信号全部分解成正弦信号的线性组合(傅里叶变换)\(x(t)=\Sigma_\omega X(\omega ) e^{i \omega t}\),那么就可以用一个传递函数\(H(w)=Y(w)/X(w)\)来描述这个线性系统。倘若这个信号很特殊,例如\(e^{2t} sin(t)\),傅里叶变换在数学上不存在,这个时候就引入拉普拉斯变换来解决这个问题\(x(t)=\Sigma_s X(s) e^{st}\)。这样一个线性系统都可以用一个传递函数\(H(s)=Y(s)/X(s)\)来表示。所以,从这里可以看到将信号分解为正弦函数(傅里叶变换)或者 复指数函数(拉普拉斯变换)对分析线性系统也是至关重要的。
傅里叶变换与量子力学
量子力学的波函数可以用多种不同的表象来描述,例如坐标表象、动量表象、能量表象等,不同表象之间的变换实际上是希尔伯特空间的一个幺正变换,其中坐标表象和动量表象之间的变换就是傅里叶变换。
\[\Phi (p)=\frac{1}{\sqrt {2\pi \hbar }}\int \Psi (x)e^{-{\frac{i}{\hbar }}px}dx, \\ \Psi (x)=\frac{1}{\sqrt {2\pi \hbar }}\int \Phi (p)e^{\frac{i}{\hbar }px}dp.\]傅里叶变换、拉普拉斯、Z变换、离散傅里叶变换的关系
信号处理中经常要对信号做各种变换,其中傅里叶变换、拉普拉斯、Z变换、离散傅里叶变换是最基础的几个变换。 他们都是为了对信号做频谱分析而采用的变换,只不过被变换的信号会有一些差异。
如果只关心信号本身,不关心系统,这几个变换的关系可以通过下面这样一个过程联系起来。
从模拟信号\(x(t)\)开始,如果模型信号能量是有限的,那么我们可以对它做傅里叶变换,把它用频域表达为\(X(w)\)。如果信号的能量是无限的,那么傅里叶变换将不会收敛,这种时候可以对它做拉普拉斯变换\(X(s)\)。 如果我们将拉普拉斯的\(s=\sigma+j w\)域画出来,他是一个复平面,拉普拉斯变换\(X(s)\)是这个复平面上的一个复变函数。而这个函数沿虚轴\(j w\)的值\(X(jw)\)就是傅里叶变换。
拉普拉斯变换和傅里叶变换广泛应用在模拟电路分析当中,下图就是对模拟电路中基本元件的\(s\)域建模示意图,当\(s=jw\)时,就是傅里叶变换了。
需要明确一个观点,不管使用时域还是频域(或s域)来表示一个信号,他们表示的都是同一个信号!也就是说,上面的时域表达、频域表达和\(s\)域表达都表示的是同一个模拟信号。关于这一点,你可以从线性空间的角度理解。同一个信号,如果采用不同的坐标框架(或者说基向量),那么他们的坐标就不同。例如,采用\(\{\delta(t-\tau )|\tau \in R\}\)作为坐标,那么信号就可以表示为\(x(t)\),而采用\(\{e^{i w t}|w\in R\}\)则表示为傅里叶变换的形式\(X(w)\)。 两个不同坐标框架下,同一个向量的坐标可以通过一个线性变换联系起来,如果是有限维的空间,则可以表示为一个矩阵,在这里是无限维,这个线性变换就是傅里叶变换。
到现在,对信号的形式还没有多少假定,如果信号是带宽受限信号,也就是说\(X(jw)\)只在一个小范围内(如\(-B<w<B\))不为0。之所以要做这个假定以及这个假定的合理性是根据实际需要而定的。在一个通信系统或者信号处理系统中,无限带宽的信号是无法处理的,而且一般接受信号的期间都会有一定的带宽,所以这是对实际中的信号的一种理想假设。现代的信号处理系统多是数字信号处理系统,即使是模拟系统,现在也多将复杂的处理放到数字信号处理子系统端进行处理,这两个系统之间通过 AD、DA 连接起来。根据采样定理,只要采样的频率足够高(大于两倍带宽),就可以无失真地将信号还原出来。那么采样对信号的影响是什么呢?从s平面来看,时域的采样将\(X(s)\)沿虚轴方向作周期延拓!这个性质从数学上可以很容易验证。下图显示的是就是采样对信号频谱的影响,只画出虚轴上的图像。这个性质也很好的解释了为什么要两倍的采样频率,这样才能使得周期延拓后频谱不会重叠到一起。设\(f_s = w_s/2\pi\)是采样频率,则采样后信号在\(s\)域可以表达为
\[X_{\text{sampling}}(s) = X(s) \sum_{n=-\infty}^{\infty} e^{n s/f_s}\]对于采样后的信号,可以利用指数变换将\(s\)域的带状区域变换到单位圆内。这就是z变换,它可以看做拉普拉斯变换的一种特殊形式,即做了一个代换\(z=e^ {-s/ f_s}\),\(f_s\)是采样频率。这个变换将信号从s域变换到z域。请注意,s域和z域表示的是同一个信号,即采样完了之后的信号,只有采样才会改变信号本身!从复平面上来看,这个变换将与\(\sigma\)轴平行的条带变换到z平面的一个单叶分支\(2k\pi\le\theta \le 2(k+1)\pi\),并且将虚轴映射到单位圆。\(z = e^{-jw/f_s}\)时也称作离散时间傅里叶变换(DTFT)。你会看到前面采样导致的周期延拓产生的条带重叠在一起了,因为具有周期性,所以z域不同的分支的函数值\(X(z)\)是相同的。换句话说,如果没有采样,直接进行z变换,将会得到一个多值的复变函数!所以一般只对采样完了后的信号做z变换!
\[X(z) = X_{sampling}(z = e^ {-s/f_s})\\ = X(s = f_s \ln z)\sum_{n=-\infty}^{\infty} z^n\]这里讲了时域的采样,时域采样后,信号只有\(-f_s/2\rightarrow f_s/2\)间的频谱,即最高频率只有采样频率一半,但是要记录这样一个信号,仍然需要无限大的存储空间,可以进一步对频域进行采样。如果时间有限(实际上这与频率受限互相矛盾,但大多数信号近似成立)的信号,那么通过频域采样(时域做周期扩展)可以不失真地从采样的信号中恢复原始信号。并且信号长度是有限的,这就是离散傅里叶变换(DFT),它有著名的快速算法快速傅里叶变换(FFT)。为什么DFT这么重要呢,因为计算机要有效地对一般的信号做傅里叶变换,都是用DFT来实现的,除非信号具有简单的解析表达式!利用上述关系,可以推导出DFT在第k个频点的值为
\[X(k) = X(z = e^{-j\frac{2\pi}{N}k}) \\ = X(s= -j\frac{2\pi}{N}kf_s) \sum_{n=-\infty}^{\infty} e^{-j\frac{2\pi}{N}nk} \\ = X(s= -j\frac{2\pi}{N}kf_s) \\ = \int_{-\infty}^{\infty} x(t) e^{-j\frac{2\pi}{N}kf_s t} dt \\ = \sum_n x_n e^{-j\frac{2\pi}{N}nk}\]上述推导利用到两个基本公式
\[\sum_{n=-\infty}^{\infty} e^{-j\frac{2\pi}{N}nk} = 1 \\ \int_{-\infty}^{\infty} x(t) e^{-j\frac{2\pi}{N}kf_s t} dt = \sum_n x_n e^{-j\frac{2\pi}{N}nk}\]总结起来说,就是对于一个线性系统,输入输出是线性关系的,不论是线性电路还是光路,只要可以用一个线性方程或线性微分方程(如拉普拉斯方程、泊松方程等)来描述的系统,都可以通过傅里叶分析从频域来分析这个系统的特性,比单纯从时域分析要强大得多!两个著名的应用例子就是线性电路和傅里叶光学(信息光学)。甚至非线性系统,也在很多情况里面使用线性系统的东西!所以傅里叶变换才这么重要!你看最早傅里叶最早也是为了求解热传导方程(那里其实也可以看做一个线性系统)!
傅里叶变换的思想还在不同领域有很多演变,比如在信号处理中的小波变换,它也是采用一组基函数来表达信号,只不过克服了傅里叶变换不能同时做时频分析的问题。
傅里叶变换特殊的原因解释
最后,我从纯数学的角度说一下傅里叶变化到底是什么。 如果我们把函数\(f(t)\)看做向量,那么这些函数在加法和数乘两种运算下构成一个线性空间。 如果我们定义内积
\[<f, g> = \int_{-\infty}^{\infty} f(t)g(t) dt\]并且限定该集合是有界函数的子集,所谓有界是指内积\(<f, f>\)有界。 那么上述线性空间就是一个希尔伯特空间。这里我们忽略这些严格的泛函分析中的定义,就简单地与欧式空间中的向量和内积进行类比即可。 在这种类比下,一个函数就是一个向量。
在这种类比下(严格的证明需要用泛函分析那一套,这里我们只关注直观的图像理解),傅里叶变换就是这个向量空间中的一个幺正变换! 我们知道,欧式空间中的线性变换都可以用一个矩阵A来表示,即变换
\[b = A x\]表示把向量x通过变换A变换为b!傅里叶变换就把时域函数f(t)变换到频域函数F(w)! 利用傅里叶变换的基本性质,容易验证这个变换是一个幺正变换。
我们知道,线性变换的本质就是选取的基向量不同。向量的每一个坐标就是对应基向量前面的系数!
\[[x_1, x_2, ..., x_n] = \sum_{i=1}^n x_i \mathbf{e}_i\]那么在函数空间中,基向量是什么呢?在时域基向量可以看成delta函数\(\{\delta(t - s); s \in R\}\)
\[f(t) = \int_{-\infty}^{\infty} f(s) \delta(t -s) ds\]这里的积分可以类比于前面的求和,\(\delta(t - s)\)可以类比于基向量\(\mathbf{e}_i\)成为基函数,f(s)可以类比于\(x_i\)是基函数前面的系数!
同样的类比,傅里叶变换到频域选取的基函数是\(\{e^{-i wt}; w \in R\}\)
\[f(t) = \int_{-\infty}^{\infty} F(w) e^{-i wt} dw\]F(w)就是基函数\(e^{-i wt}\)前面的系数。傅里叶变换就是在这两组基函数间的线性变换!
那么,问题来了,线性变换这么多,为什么傅里叶变换这么特殊?
还记得线性代数中的线性方程Ax=b吗?解这个方程的方法很多,高斯消元法是最常用的方法之一。 但是如果A是一个对角方阵,那么这个向量版的线性方程可以变为多个独立的代数方程!
\[a_{ii} x_i = b_i, i=1,...,n\]这种情况下,很容易求得\(x_i = b_i/a_{ii}\)!
上述情况过于特殊,我们考虑更一般的情况,如果A是对称方阵,那么根据线性空间的特征值理论,可以找到矩阵A的所有互相正交的特征向量\(\{\mathbf{v}_i,i=1..n\}\)和特征值\(\{\lambda_i,i=1..n\}\),然后将向量x和b表示成特征向量的组合\(x=\Sigma_i x_i' \mathbf{v}_i, b=\Sigma_i b_i' \mathbf{v}_i\)。由于特征向量的正交关系,矩阵的代数方程可以化为n个标量代数方程
\[\lambda_i x_i' = b_i', i=1,...,n\]是不是很神奇!!一个向量版的线性方程通过重新选取了一组基向量变成多个独立的代数一次方程!
你会问这跟傅里叶变换有毛关系?别急,我们再来看非齐次线性常微分方程
\[(\frac{d}{dt} + a)y(t)=z(t)\]如果把左边的线性算子部分看做线性变换,那么这个方程完全可以和上述向量版的线性方程进行类比! 把算子\(\Lambda = \frac{d}{dt} + a\)看做线性变换,那么我们可以采用上述类似的思路,把这个方程变成多个独立的代数方程吗? 答案是肯定的,利用该算子的特征函数作为基函数重新选取基函数即可!可以验证指数函数\(y=e^{iwt}\)是该的特征函数,对应的特征值是\(iw+a\)
\[(\frac{d}{dt} + a) e^{st} = (iw + a) e^{st}\]利用相似的思路,我门把函数\(y(t), z(t)\)都表示为基函数的线性组合
\[y(t) = \int Y(w) e^{iwt} dw \\ z(t) = \int Z(w) e^{iwt} dw\]那么这样一来,前述微分方程变成了多个标量线性代数方程!
\[(iw + a) Y(w) = Z(w), w\in R\]其实这个过程也可以看做对原始方程左右两边同时做傅里叶变换!这也是傅里叶变换求解常系数微分方程的理论基础!
在常系数线性偏微分方程中也有类似结论!例如,考虑有源的拉普拉斯方程
\[\nabla^2 \phi(\mathbf{x}) = \rho(\mathbf{x})\]容易验证基函数(其实就是格林函数)
\[G(\mathbf{k},\mathbf{x}) = e^{i \mathbf{k} \cdot \mathbf{x}}\]是拉普拉斯算子的特征函数!将场\(\phi\)和源\(\rho\)按照基函数展开,就可以将原来的拉普拉斯方程变为多个标量代数线性方程
\[- (k_x^2 + k_y^2 + k_z^2) \phi(\mathbf{k}) = \rho(\mathbf{k})\]上述傅里叶变换也可以用拉普拉斯变换替换,结论一样! 以上是我在上数理方程课程的时候体会到的。归纳起来,就是说傅里叶变换就是线性空间中的一个特殊的正交变换!他之所以特殊是因为指数函数是常系数微分算子的特征函数!而自然界常见的规律大多是用常系数微分方程描述,信号系统中更是常见,线性时不变系统都可以利用常系数微分方程描述,这使得傅里叶变换应用十分广泛!
其他微分算子的特征函数举例
对于常系数线性微分算子,可以用指数函数作为基函数展开,而对于变系数线性微分算子,其基函数就不再是简单的指数函数了。 但是上述思想仍然可以利用,只不过基函数是一些特殊函数,如贝塞尔函数、勒让德多项书函数等等!
所谓常系数微分算子就是具有这种形式的微分算子
\[\hat{L} = \sum_{k=0}^n a_k \frac {\mathrm {d} ^{k}x}{\mathrm {d} t^{k}}, a_k\in R\]对于变系数的微分算子,\(a_k\)是自变量\(t\)的函数,这种算子的特征函数并没有一般性的结论。 基本上每一类算子都会有自己特殊的特征函数, 这里列举几个我遇到过多次的特征函数及变系数算子。
柱坐标下的贝塞尔函数是下述微分算子的特征函数
\[\hat{L} = x^{2}{\frac {d^{2}y}{dx^{2}}}+x{\frac {dy}{dx}}+(x^{2}-\alpha ^{2})\]球坐标下的勒让德多项式
\[{\displaystyle P_{n}(x)={1 \over 2^{n}n!}{d^{n} \over dx^{n}}\left[(x^{2}-1)^{n}\right].}\]它是下述微分算子的特征函数,这是一个变系数的微分算子
\[\hat{L} = {d \over dx}\left[(1-x^{2}){d \over dx}\right]+n(n+1)\] \[L_{n}(x)={\frac {e^{x}}{n!}}{\frac {d^{n}}{dx^{n}}}\left(e^{-x}x^{n}\right)={\frac {1}{n!}}\left({\frac {d}{dx}}-1\right)^{n}x^{n} \\ \hat{L} = x{d^2 \over dx^2} + (1 - x){d \over dx} + n\]这样的例子还有很多,这些函数实际上都是一个函数族,这些函数互相正交,这和实对称阵的本征向量互相正交的性质一样,这里的线性算子也是其泛函空间上的对称轭米算子。这些函数族构成一组完备正交基,可以表达对应泛函空间中的任意函数。这和傅里叶变换的基函数——复指数函数一样。