數字世界的“棱鏡”:離散傅里葉變換(DFT)完全解析
在上一章,我們探索了Z變換和離散時間傅里葉變換(DTFT)。我們知道,DTFT是一個無比強大的理論工具,它能將一個時域離散序列的“基因圖譜”——也就是它的頻譜——完整地揭示出來。理論上,只要我們知道了信號的DTFT,就知道了它包含的所有頻率成分。
但這里有一個巨大的“但是”。DTFT是理論上的完美,卻在實踐中遇到了一個無法逾越的鴻溝:
- DTFT的定義是針對無限長序列的。而在現實世界中,我們用計算機處理的信號,無論是音頻、圖像還是一段傳感器數據,都必然是有限長度的。
- DTFT計算出的頻譜 X ( e j ω ) X(e^{j\omega}) X(ejω) 是一個關于頻率 ω \omega ω 的連續函數。計算機天生只能處理離散的數據點,它無法存儲和計算一個連續函數。
這就好比我們擁有一張描繪了整個地球的、精度無限的完美地圖(DTFT),但我們的工具只是一臺內存有限的普通電腦。我們不可能把整張地圖裝進電腦。我們該怎么辦?
答案是:對地圖進行采樣!我們只在有限的區域內(信號的有限長度),選取有限個有代表性的地點(頻譜的有限個點)來觀察。這個為計算機量身定做的、從理論走向實踐的“采樣版”傅里葉變換,就是我們今天的主角——離散傅里葉變換(Discrete Fourier Transform, DFT)。
DFT是數字信號處理領域應用最廣泛的工具,沒有之一。從你的手機音樂播放器里的頻譜顯示,到通信系統中的調制解調,再到各種高效的濾波算法,背后都有DFT的身影。理解它,是你從理論學習邁向工程實踐的決定性一步。
從理想到現實:DFT的誕生
讓我們再次明確一下從完美的DTFT到實用的DFT,我們必須解決的兩個核心問題:
- 信號長度問題:我們無法處理無限長的信號 x ( n ) x(n) x(n)。解決方案:我們只截取信號中一段有限的長度,比如 N N N 個點,來進行分析。
- 頻譜連續問題:我們無法計算連續的頻譜 X ( e j ω ) X(e^{j\omega}) X(ejω)。解決方案:我們只在頻譜的一個周期(通常是 [ 0 , 2 π ) [0, 2\pi) [0,2π))內,等間隔地計算 N N N 個點的頻譜值。
這兩個解決方案合在一起,就構成了DFT的本質。
DFT的核心定義
對于一個長度為 N N N 的有限長序列 x ( n ) x(n) x(n)(其中 n = 0 , 1 , … , N ? 1 n=0, 1, \dots, N-1 n=0,1,…,N?1),它的N點DFT定義為:
X ( k ) = ∑ n = 0 N ? 1 x ( n ) W N k n , k = 0 , 1 , … , N ? 1 X(k) = \sum_{n=0}^{N-1} x(n) W_N^{kn}, \quad k=0, 1, \dots, N-1 X(k)=n=0∑N?1?x(n)WNkn?,k=0,1,…,N?1
相應地,從頻域序列 X ( k ) X(k) X(k) 恢復時域序列 x ( n ) x(n) x(n) 的逆變換(IDFT)為:
x ( n ) = 1 N ∑ k = 0 N ? 1 X ( k ) W N ? k n , n = 0 , 1 , … , N ? 1 x(n) = \frac{1}{N} \sum_{k=0}^{N-1} X(k) W_N^{-kn}, \quad n=0, 1, \dots, N-1 x(n)=N1?k=0∑N?1?X(k)WN?kn?,n=0,1,…,N?1
這里的 W N W_N WN? 是一個非常重要的符號,它被稱為旋轉因子(Twiddle Factor),其定義為:
W N = e ? j 2 π N W_N = e^{-j\frac{2\pi}{N}} WN?=e?jN2π?
W N k W_N^k WNk? 代表了在復平面單位圓上,從(1, 0)點開始,按順時針方向轉過 k k k 個 2 π / N 2\pi/N 2π/N 角度后到達的點。整個DFT的計算,本質上就是將時域信號 x ( n ) x(n) x(n) 與這些在單位圓上旋轉的點進行內積運算,從而衡量出信號在各個特定頻率上的分量大小。
DFT與Z變換/DTFT的血緣關系
DFT并不是一個憑空冒出來的全新概念,它與我們在上一章學習的Z變換和DTFT有著深刻的血緣關系。
一個有限長序列 x ( n ) x(n) x(n) 的Z變換為:
X ( z ) = ∑ n = 0 N ? 1 x ( n ) z ? n X(z) = \sum_{n=0}^{N-1} x(n) z^{-n} X(z)=n=0∑N?1?x(n)z?n
它的DTFT則是Z變換在單位圓上的取值:
X ( e j ω ) = X ( z ) ∣ z = e j ω = ∑ n = 0 N ? 1 x ( n ) e ? j ω n X(e^{j\omega}) = X(z)|_{z=e^{j\omega}} = \sum_{n=0}^{N-1} x(n) e^{-j\omega n} X(ejω)=X(z)∣z=ejω?=n=0∑N?1?x(n)e?jωn
現在,讓我們在DTFT的頻譜上,等間隔地采樣N個點。采樣點的頻率為 ω k = 2 π k N \omega_k = \frac{2\pi k}{N} ωk?=N2πk?,其中 k = 0 , 1 , … , N ? 1 k=0, 1, \dots, N-1 k=0,1,…,N?1。將這個 ω k \omega_k ωk? 代入DTFT的公式:
X ( e j ω k ) = ∑ n = 0 N ? 1 x ( n ) e ? j ( 2 π k N ) n = ∑ n = 0 N ? 1 x ( n ) ( e ? j 2 π N ) k n = ∑ n = 0 N ? 1 x ( n ) W N k n X(e^{j\omega_k}) = \sum_{n=0}^{N-1} x(n) e^{-j (\frac{2\pi k}{N}) n} = \sum_{n=0}^{N-1} x(n) (e^{-j\frac{2\pi}{N}})^{kn} = \sum_{n=0}^{N-1} x(n) W_N^{kn} X(ejωk?)=n=0∑N?1?x(n)e?j(N2πk?)n=n=0∑N?1?x(n)(e?jN2π?)kn=n=0∑N?1?x(n)WNkn?
看!這個結果不就是我們定義的DFT X ( k ) X(k) X(k) 嗎?
所以,我們得出了一個至關重要的結論:
一個N點有限長序列的N點DFT,就是其DTFT(或Z變換)在單位圓上的N個等間隔采樣點。
(上圖展示了Z平面上的單位圓,DFT的計算就相當于在這個單位圓上等間隔地取N個點,并計算Z變換在這些點上的值。)
這個結論完美地將理論(DTFT/Z變換)與實踐(DFT)聯系在了一起。DFT通過在頻域進行采樣,成功地將一個連續的頻譜問題轉化為了一個計算機可以處理的離散問題。
DFT的“怪癖”:循環的魔咒
將無限和連續的世界強行壓縮到有限和離散的世界,必然會帶來一些“副作用”。DFT最核心、最獨特,也最容易讓初學者困惑的特性,就是它的隱含周期性(Implied Periodicity)。
當你對一個長度為N的序列 x ( n ) x(n) x(n) 進行N點DFT時,DFT算法“看待”這個序列的方式,并不是一個孤立的、長度為N的片段。它會固執地認為,你給它的這個 x ( n ) x(n) x(n),其實是一個周期為N的無限長周期序列中的一個主值周期。
(上隱含周期性:一個有限長序列被DFT看作是首尾相連的圓形序列,并進行周期延拓。)
這個“自作多情”的假設,給DFT的性質帶來了深刻的影響,其中最關鍵的就是循環移位和循環卷積。
線性卷積 vs. 循環卷積:一個不小心就踩的“坑”
在上一章我們知道,LTI系統中信號的傳遞是通過線性卷積完成的,并且在頻域上有一個優美的性質:時域的線性卷積對應頻域的乘積。
x ( n ) ? h ( n ) ? X ( e j ω ) H ( e j ω ) x(n) \ast h(n) \longleftrightarrow X(e^{j\omega}) H(e^{j\omega}) x(n)?h(n)?X(ejω)H(ejω)
我們自然希望DFT也能保持這個美好的性質。但由于隱含周期性的存在,事情發生了變化。在DFT的世界里,時域的運算不再是線性的,而是循環的。
線性卷積的公式是:
y L ( n ) = ∑ k = ? ∞ ∞ x ( k ) h ( n ? k ) y_L(n) = \sum_{k=-\infty}^{\infty} x(k) h(n-k) yL?(n)=k=?∞∑∞?x(k)h(n?k)
而N點循環卷積的公式是:
y C ( n ) = x ( n ) ? h ( n ) = ∑ m = 0 N ? 1 x ( m ) h ( ( n ? m ) N ) y_C(n) = x(n) \otimes h(n) = \sum_{m=0}^{N-1} x(m) h((n-m)_N) yC?(n)=x(n)?h(n)=m=0∑N?1?x(m)h((n?m)N?)
注意那個 ( ( n ? m ) N ) ((n-m)_N) ((n?m)N?),它表示對 n ? m n-m n?m 的結果進行模N運算(即取余數)。這就像在一個有N個刻度的時鐘上做減法。例如,在8點循環卷積中,從第1個位置后退3個位置,不是到了-2,而是到了第6個位置( ( 1 ? 3 ) 8 = ? 2 ( m o d 8 ) = 6 (1-3)_8 = -2 \pmod 8 = 6 (1?3)8?=?2(mod8)=6 )。
DFT的卷積定理說的是:
N點循環卷積對應N點DFT的乘積。
x ( n ) ? h ( n ) ? X ( k ) H ( k ) x(n) \otimes h(n) \longleftrightarrow X(k) H(k) x(n)?h(n)?X(k)H(k)
這帶來了一個巨大的問題。我們通常想計算的是線性卷積,但如果我們直接計算兩個序列的DFT,在頻域相乘,再做IDFT,我們得到的將是循環卷積的結果。如果兩個序列的長度不合適,循環卷積的結果會因為“周期性”的重疊(稱為時間混疊)而變得面目全非,根本不是我們想要的線性卷積。
那么,我們還能用DFT來幫助我們計算線性卷積嗎?
答案是肯定的,只要我們使用一個絕妙的技巧——補零(Zero-Padding)。
假設我們有兩個序列, x ( n ) x(n) x(n) 長度為 M M M, h ( n ) h(n) h(n) 長度為 P P P。我們知道,它們的線性卷積結果 y L ( n ) y_L(n) yL?(n) 的長度將是 L = M + P ? 1 L = M+P-1 L=M+P?1。
快速線性卷積的正確步驟:
- 確定長度: 計算出線性卷積結果的長度 L = M + P ? 1 L = M+P-1 L=M+P?1。
- 補零: 將序列 x ( n ) x(n) x(n) 和 h ( n ) h(n) h(n) 的尾部都補上若干個零,使它們的長度都達到 至少為L 的新長度 N N N(為了后續FFT算法的效率,通常會選擇一個大于等于L的2的冪次方作為N)。
- DFT: 對補零后的兩個序列做N點DFT,得到 X ( k ) X(k) X(k) 和 H ( k ) H(k) H(k)。
- 頻域相乘: 計算 Y ( k ) = X ( k ) ? H ( k ) Y(k) = X(k) \cdot H(k) Y(k)=X(k)?H(k)。
- IDFT: 對 Y ( k ) Y(k) Y(k) 做N點IDFT,得到最終的結果 y ( n ) y(n) y(n)。
通過補零,我們相當于把循環卷積的“時鐘”刻度 N 拉得足夠長,使得在一個周期內,翻轉后的 h ( n ) h(n) h(n) 與 x ( n ) x(n) x(n) 進行卷積時,其“尾部”不會“從時鐘的另一頭繞回來”污染到“頭部”。這樣,在一個周期內計算出的循環卷積結果,就和線性卷積的結果完全一樣了。
這個利用DFT(及其快速算法FFT)來計算線性卷積的方法,被稱為快速卷積,是數字信號處理中最高效的卷積計算方法。
DFT的威力:從頻譜分析到快速計算
掌握了DFT的定義和核心特性后,我們來看看它在實際工程中的兩大應用。
應用一:數字頻譜分析
DFT最直接的應用,就是作為一種頻譜分析儀。只要我們將一段信號(比如一段錄音)送入計算機,對其進行DFT,得到的 X ( k ) X(k) X(k) 的幅度 ∣ X ( k ) ∣ |X(k)| ∣X(k)∣ 就告訴了我們,這段信號在各個離散頻率點 f k = k ? f s / N f_k = k \cdot f_s / N fk?=k?fs?/N 上的能量強度。這就像一個棱鏡,將混合在一起的白光(時域信號)分解成赤橙黃綠青藍紫的彩色光譜(頻域信號)。
然而,我們必須清醒地認識到,這個“數字棱鏡”并非完美,它會帶來一些由“有限”觀測所導致的必然誤差:
- 頻譜泄漏(Spectral Leakage): 我們對信號的觀測是有限的,這相當于給無限長的原始信號乘上了一個矩形“窗戶”。這個操作在頻域上會引起卷積,使得原本單一頻率的尖銳譜線,“泄漏”成一個主瓣和多個旁瓣。這會導致我們探測到一些本不存在的頻率分量,并可能掩蓋掉一些微弱的真實信號。選擇不同的窗函數(如漢寧窗、哈明窗)可以在主瓣寬度和旁瓣抑制之間進行權衡,以減小泄漏。
(上圖展示了頻譜泄漏效應。一個理想單頻正弦波的頻譜(左)應該是一個尖峰,但由于加窗截斷,其實際DFT頻譜(右)變成了一個有主瓣和多個旁-瓣的形狀。)
- 柵欄效應(Picket-Fence Effect): DFT只在N個離散的頻率點上對真實頻譜進行采樣。如果一個信號的真實頻率峰值恰好落在了兩個采樣點之間,那么我們在兩個采樣點上觀測到的幅度都將低于真實峰值。這就像我們透過一個柵欄看風景,我們只能看到柵欄縫隙里的景象,可能會錯過最精彩的部分。增加DFT的點數N(通過補零)可以加密“柵欄”的縫隙,更精確地找到峰值位置和大小。
應用二:快速卷積的基石
如前文所述,利用DFT進行快速卷積是其另一大核心應用。直接計算M點和P點序列的線性卷積,其計算復雜度大約是 O ( M ? P ) O(M \cdot P) O(M?P)。而利用DFT(特別是下一章要講的FFT算法)進行快速卷積,其計算復雜度大約是 O ( N log ? N ) O(N \log N) O(NlogN),其中 N N N 是補零后的長度。當序列很長時,FFT方法的效率優勢是壓倒性的。
想象一下,在圖像處理中,一個512x512的模糊核要對一張1920x1080的圖像進行卷積濾波,如果用直接法計算,那將是天文數字。而利用2D-FFT進行快速卷積,則可以在幾秒鐘內完成。這背后全都是DFT的功勞。
DFT的出現,為我們打開了在頻域對信號進行處理的大門,也為下一章的快速傅里葉變換(FFT)——這個被譽為20世紀最重要的算法之一——鋪平了道路。FFT并不是一種新的變換,它只是實現DFT的一種極其高效的算法。沒有DFT的理論基礎,FFT就無從談起。
總結:掌握DFT,掌握實踐
今天,我們踏上了一段從理想到現實的旅程。我們理解到:
- DFT是為計算機而生的傅里葉變換。它通過對信號時域和頻域的同時“離散化”,解決了DTFT無法在計算機上實現的問題。
- DFT是DTFT的采樣。N點DFT的本質是在單位圓上對序列的Z變換(或DTFT)進行N點等間隔采樣。
- DFT的核心“怪癖”是隱含周期性。這導致了時域運算的“循環”特性,最典型的就是循環卷積。
- 補零是關鍵技巧。通過對序列進行恰當的補零,我們可以利用DFT的循環卷積定理來高效地計算線性卷積,這就是快速卷積的原理。
- DFT是強大的分析工具。作為數字頻譜分析儀,它能揭示信號的頻率構成,但我們也需警惕頻譜泄漏和柵欄效應等實踐中的問題。
DFT是連接理論與工程的堅固橋梁。掌握了它,你就不再僅僅是看著數學公式,而是真正擁有了用計算機分析和處理信號的利器。在下一章,我們將為這把利器裝上一個強力馬達——FFT算法,讓它的威力發揮到極致。
習題
第1題 (概念題)
一位工程師說:“我用計算機直接計算了一個信號的DTFT頻譜,發現它是一個連續的平滑曲線。” 這句話在嚴格意義上正確嗎?為什么?
第2題 (計算題)
有兩個長度為4的序列, x ( n ) = { 1 , 2 , 2 , 1 } x(n)=\{1, 2, 2, 1\} x(n)={1,2,2,1} 和 h ( n ) = { 1 , ? 1 , 0 , 0 } h(n)=\{1, -1, 0, 0\} h(n)={1,?1,0,0}。請計算它們的4點循環卷積 y C ( n ) = x ( n ) ? h ( n ) y_C(n) = x(n) \otimes h(n) yC?(n)=x(n)?h(n)。
第3題 (應用題)
你希望使用DFT(以及未來的FFT)來計算兩個有限長序列的線性卷積。第一個序列 x 1 ( n ) x_1(n) x1?(n) 的長度為30個點,第二個序列 x 2 ( n ) x_2(n) x2?(n) 的長度為20個點。為了得到正確的結果,你需要對這兩個序列進行補零,然后進行DFT。請問,你選擇的DFT點數N最小應該是多少?
答案
第1題答案:
這句話在嚴格意義上不正確。
原因: 計算機無法直接計算或存儲一個真正的連續函數。DTFT的頻譜 X ( e j ω ) X(e^{j\omega}) X(ejω) 是一個關于連續變量 ω \omega ω 的函數。這位工程師看到的“連續平滑曲線”實際上是計算了大量點的DFT(例如,對原信號大量補零后做DFT),然后在繪圖時將這些密集的離散點用直線連接起來,從而在視覺上產生了連續的效果。其本質仍然是離散點,而不是真正的連續函數。
第2題答案:
我們可以使用定義法計算。循環卷積 y C ( n ) = ∑ m = 0 3 x ( m ) h ( ( n ? m ) 4 ) y_C(n) = \sum_{m=0}^{3} x(m) h((n-m)_4) yC?(n)=∑m=03?x(m)h((n?m)4?)。
-
y C ( 0 ) = x ( 0 ) h ( 0 ) + x ( 1 ) h ( ? 1 ) 4 + x ( 2 ) h ( ? 2 ) 4 + x ( 3 ) h ( ? 3 ) 4 y_C(0) = x(0)h(0) + x(1)h(-1)_4 + x(2)h(-2)_4 + x(3)h(-3)_4 yC?(0)=x(0)h(0)+x(1)h(?1)4?+x(2)h(?2)4?+x(3)h(?3)4?
= x ( 0 ) h ( 0 ) + x ( 1 ) h ( 3 ) + x ( 2 ) h ( 2 ) + x ( 3 ) h ( 1 ) = x(0)h(0) + x(1)h(3) + x(2)h(2) + x(3)h(1) =x(0)h(0)+x(1)h(3)+x(2)h(2)+x(3)h(1)
= ( 1 ) ( 1 ) + ( 2 ) ( 0 ) + ( 2 ) ( 0 ) + ( 1 ) ( ? 1 ) = 1 ? 1 = 0 = (1)(1) + (2)(0) + (2)(0) + (1)(-1) = 1 - 1 = \mathbf{0} =(1)(1)+(2)(0)+(2)(0)+(1)(?1)=1?1=0 -
y C ( 1 ) = x ( 0 ) h ( 1 ) + x ( 1 ) h ( 0 ) + x ( 2 ) h ( ? 1 ) 4 + x ( 3 ) h ( ? 2 ) 4 y_C(1) = x(0)h(1) + x(1)h(0) + x(2)h(-1)_4 + x(3)h(-2)_4 yC?(1)=x(0)h(1)+x(1)h(0)+x(2)h(?1)4?+x(3)h(?2)4?
= x ( 0 ) h ( 1 ) + x ( 1 ) h ( 0 ) + x ( 2 ) h ( 3 ) + x ( 3 ) h ( 2 ) = x(0)h(1) + x(1)h(0) + x(2)h(3) + x(3)h(2) =x(0)h(1)+x(1)h(0)+x(2)h(3)+x(3)h(2)
= ( 1 ) ( ? 1 ) + ( 2 ) ( 1 ) + ( 2 ) ( 0 ) + ( 1 ) ( 0 ) = ? 1 + 2 = 1 = (1)(-1) + (2)(1) + (2)(0) + (1)(0) = -1 + 2 = \mathbf{1} =(1)(?1)+(2)(1)+(2)(0)+(1)(0)=?1+2=1 -
y C ( 2 ) = x ( 0 ) h ( 2 ) + x ( 1 ) h ( 1 ) + x ( 2 ) h ( 0 ) + x ( 3 ) h ( ? 1 ) 4 y_C(2) = x(0)h(2) + x(1)h(1) + x(2)h(0) + x(3)h(-1)_4 yC?(2)=x(0)h(2)+x(1)h(1)+x(2)h(0)+x(3)h(?1)4?
= x ( 0 ) h ( 2 ) + x ( 1 ) h ( 1 ) + x ( 2 ) h ( 0 ) + x ( 3 ) h ( 3 ) = x(0)h(2) + x(1)h(1) + x(2)h(0) + x(3)h(3) =x(0)h(2)+x(1)h(1)+x(2)h(0)+x(3)h(3)
= ( 1 ) ( 0 ) + ( 2 ) ( ? 1 ) + ( 2 ) ( 1 ) + ( 1 ) ( 0 ) = ? 2 + 2 = 0 = (1)(0) + (2)(-1) + (2)(1) + (1)(0) = -2 + 2 = \mathbf{0} =(1)(0)+(2)(?1)+(2)(1)+(1)(0)=?2+2=0 -
y C ( 3 ) = x ( 0 ) h ( 3 ) + x ( 1 ) h ( 2 ) + x ( 2 ) h ( 1 ) + x ( 3 ) h ( 0 ) y_C(3) = x(0)h(3) + x(1)h(2) + x(2)h(1) + x(3)h(0) yC?(3)=x(0)h(3)+x(1)h(2)+x(2)h(1)+x(3)h(0)
= ( 1 ) ( 0 ) + ( 2 ) ( 0 ) + ( 2 ) ( ? 1 ) + ( 1 ) ( 1 ) = ? 2 + 1 = ? 1 = (1)(0) + (2)(0) + (2)(-1) + (1)(1) = -2 + 1 = \mathbf{-1} =(1)(0)+(2)(0)+(2)(?1)+(1)(1)=?2+1=?1
所以,4點循環卷積的結果是 y C ( n ) = { 0 , 1 , 0 , ? 1 } y_C(n) = \{0, 1, 0, -1\} yC?(n)={0,1,0,?1}。
第3題答案:
線性卷積結果的長度為 L = M + P ? 1 L = M + P - 1 L=M+P?1。
在這里, M = 30 M = 30 M=30, P = 20 P = 20 P=20。
所以,線性卷積結果的長度 L = 30 + 20 ? 1 = 49 L = 30 + 20 - 1 = 49 L=30+20?1=49。
為了使用DFT/FFT計算線性卷積而不產生時間混疊,補零后的長度(即DFT的點數N)必須大于或等于線性卷積結果的長度L。
所以,DFT點數N的最小值是 49。
(在實際應用中,為了利用FFT算法的高效性,通常會選擇大于等于49的第一個2的冪次方,即 N = 64 N=64 N=64。)