加伯轉換是窗函數為高斯函數的短時距傅立葉變換。
數學定義
將短時距傅立葉轉換中的窗函數代入高斯函數,即可得下面的標準定義:
:G_x(t,f)=\int_{-\infty}^{\infty}e^{-\pi (\tau-t)^2} e^{-j2 \pi f \tau}x(\tau)\, d\tau
以下是幾種常見的替代定義:
Gx_1(t,f)=\int_{-\infty}^{\infty}e^{-\pi (\tau-t)^2} e^{-j2 \pi f (\tau-\frac{t}{2})}x(\tau)\, d\tau
#Gx_2(t,f)=\sqrt[4]{2}\int_{-\infty}^{\infty}e^{-\pi (\tau-t)^2} e^{-j2 \pi f \tau}x(\tau)\, d\tau
#Gx_3(t,\omega)=\int_{-\infty}^{\infty}e^{-\frac{(\tau-t)^2}{2}} e^{-j \omega\tau}x(\tau)\, d\tau
#Gx_4(t,\omega)=\sqrt{\frac{1}{2\pi}}\int_{-\infty}^{\infty}e^{-\frac{(\tau-t)^2}{2}} e^{-j \omega(\tau-\frac{t}{2})}x(\tau)\, d\tau
*註:在文獻上可能會看到不同形式的加伯轉換,但本質上都是一樣的。
由於實作時,不能計算無限大的積分式子,所以根據高斯函數會從兩側遞減的性質,我們可以將上式進一步化簡:
:\because \begin{cases}
\ e^{-\pi a^2} 1.9143 \\
\ e^{-a^2 /2} 4.7985
\end{cases}
:\therefore \begin{cases}
\ G_x(t,f) \simeq \int_{t-1.9143}^{t+1.9143}e^{-\pi (\tau-t)^2} e^{-j2 \pi f \tau}x(\tau)\, d\tau \\
\ Gx_4(t,\omega) \simeq \sqrt{\frac{1}{2\pi}}\int_{t-4.7985}^{t+4.7985}e^{-\frac{(\tau-t)^2}{2}} e^{-j\omega (\tau-t/2)}x(\tau)\, d\tau
\end{cases}
為何選擇高斯函數作為窗函數
其他窗函數的短時距傅立葉變換,如利用方型窗函數的短時距傅立葉變換,無法同時兼顧時間軸和頻率軸的解析度;一者解析度提升,另一者解析度必定下降。但高斯函數由海森堡測不準原理可得知,是最能同時讓兩軸兼顧解析度的窗函數(將於下面章節詳述)。
高斯函數為傅立葉轉換的特徵函數:
- \int_{-\infty}^{\infty}e^{-\pi t^2}e^{-j2 \pi f \tau}dt=e^{-\pi f^2}
- \int_{-\infty}^{\infty}e^{-\frac{t^2}{2}}e^{-j \omega t}dt=e^{-\frac{f^2}{2}}
因此經過轉換後其性質不變。因此可讓加伯轉換後在時間軸和頻率軸的性質相互對稱。
由測不準原理了解高斯函數的性質
上述提到,高斯函數是最能兼顧時間與頻率解析度的窗函數。我們利用這個章節來詳細討論。
*海森堡測不準原理:
:對於一個信號 x(t),當 |t|\to \infty,若 \sqrt{t} x(t)=0 ,則
: \sigma_t \sigma_f \ge 1/4 \pi \,
:其中 \sigma_t^2 = \int (t-\mu_t)^2 P_x(t)dt, \sigma_f^2 = \int (f-\mu_f)^2 P_X(f)df \,
:\qquad \mu_t = \int tP_x(t)dt\quad \quad \quad\ ,\mu_f = \int fP_X(f)df
:\qquad P_x(t) = \frac
從圖中可以發現方形窗函數的短時傅立葉轉換會有能量擴散的情形,而加伯轉換則是清晰的時頻圖。
加伯轉換的縮放
:由於高斯窗函數的寬度可以由一常數做調整,因此我們將這個參數加入加伯轉換的數學式子中,讓轉換更加彈性,如下式:
:G_x(t,f)=\sqrt[4]{\sigma} \int_{-\infty}^{\infty} e^{-\sigma \pi (\tau-t)^2} e^{-j2 \pi f\tau}x(\tau)d\tau
:而根據前面章節所述。實作時,不能計算無限大的積分式子,所以根據高斯函數會從兩側遞減的性質,我們可以將上式進一步化簡:
:G_x(t,f) \simeq \sqrt[4]{\sigma} \int_{t-\frac{1.9143}{\sqrt{\sigma}}}^{t+\frac{1.9143}{\sqrt{\sigma}}} e^{-\sigma \pi (\tau-t)^2} e^{-j2 \pi f\tau}x(\tau)d\tau
*根據傅立葉轉換的縮放公式,假設w(t)=e^{-\pi \sigma t^2},則傅立葉轉換後為W(f)=\frac{1}{\sqrt{\sigma}}e^{-\frac{\pi f^2}{\sigma}},使其能根據需求而調整時域解析度或頻域解析度
*改變高斯函數的寬度,和改變方形窗函數短時距傅立葉變換的效果類似。若選取較大的\sigma,時域的高斯窗函數較窄,則時域有較高的解析度,而頻域的高斯窗函數較寬,所以頻域的解析度會下降(通常用於需要時域解析度較高的應用,例如:音樂訊號);反之,若選取較小的\sigma,時域的高斯窗函數較寬,則時域的解析度下降,而頻域的高斯窗函數較窄,所以頻域的解析度會上升(通常運用在需要頻域解析度較高的應用,例如:氣候)。雖然還是有兩軸之間的解析度的犧牲,但比起其他無法滿足測不準原理下限的窗函數,加伯轉換的兩軸還是能相對維持較高的解析度。
*若應用於瞬時頻率改變較劇烈的應用,則可考慮使用窗寬度隨時間而變動的加伯轉換數學式子,如下
:G_x(t,f)=\sqrt[4]{\sigma (t)} \int_{-\infty}^{\infty} e^{-\sigma (t)\pi (\tau-t)^2} e^{-j2 \pi f\tau}x(\tau)d\tau
:當瞬時頻率變動非常快時,使用較大的\sigma 值,使其時域解析度能較高;當瞬時頻率變動很慢時,使用較小的\sigma 值,使其頻域解析度能較高。
實現方法及注意事項
Direct Implementation
X(t,f)=\int_{-\infty}^{\infty}w(t-\tau)x(\tau)e^{-j2 \pi f\tau}d\tau
w(t)=e^{- \pi \sigma t^2}
*Discrete Form:
令t = n\Delta_t , f = m\Delta_f ,\tau= p\Delta_t
可將式子改寫為離散形式:
{X}\left( {n{\Delta _t},m{\Delta _f}} \right) = \sum\limits_{p = - \infty }^\infty {w\left( {(n - p){\Delta _t}} \right){x}\left( {p{\Delta _t}} \right)}{e^{ - j2\pi \,mp{\Delta _t}{\Delta _f}}}{\Delta _t}
w(t) \cong 0 \qquad for\left| t\right| >B , \frac{B}{\Delta _t} = Q
w((n-p)\Delta_t)\cong 0 \qquad
for \left| n-p\right| > \frac{B}{\Delta _t}
,\left| p-n\right| > Q
therefore,only when
-Q
w((n-p)\Delta_t) is nonzero
可改寫為:
{X}\left( {n{\Delta _t},m{\Delta _f}} \right) = \sum\limits_{p = n-Q }^{ n+Q} {w\left( {(n - p){\Delta _t}} \right){x}\left( {p{\Delta _t}} \right)}{e^{ - j2\pi \,mp{\Delta _t}{\Delta _f}}}{\Delta _t} 按照此式即可實現
e^{- \pi \sigma a^2} when\left| a\right| > 1.9143
Q=\frac{1.9143}{\sqrt{\sigma}\Delta t}
B=\frac{1.9143}{\sqrt{\sigma}}
限制
*避免贋頻效應(aliasing effect)
(1){\Delta_t}
時間複雜度
O(TFQ) T:時間取樣點數 F:頻率取樣點數 Q:Q=\frac{1.9143}{\sqrt{\sigma}\Delta t}
優缺點
:優點:簡單實現,限制條件少
:缺點:時間複雜度高
FFT-Based Method(快速傅立葉轉換)
由Direct Implementation可得下式
{X}\left( {n{\Delta _t},m{\Delta _f}} \right) = \sum\limits_{p = n-Q }^{ n+Q} {w\left( {(n - p){\Delta _t}} \right){x}\left( {p{\Delta _t}} \right)}{e^{ - j2\pi \,mp{\Delta _t}{\Delta _f}}}{\Delta _t}
令q=p-(n-Q) \to p=(n-Q)+q且離散傅立葉轉換標準式Y[m]=\sum\limits_{n = 0 }^{ N-1}y[n]e^{-j\frac{2\pi mn}{N}}
可將式子整理為:
{X}\left( {n{\Delta _t},m{\Delta _f}} \right) = {\Delta _t}{e^{ j{\textstyle{{2\pi \,(Q-n)m} \over N}}}}\sum\limits_{q = 0}^{N-1} {x_1\left( {q} \right){e^{ - j{\textstyle{{2\pi \,qm} \over N}}}}}按照此式將 {x_1}以fft()算出帶入即可實現
其中 {x_1}\left( q \right) = w\left( {(Q - q ){\Delta _t}} \right) x\left( {(n - Q + q){\Delta _t}} \right) ,0 \le q \le 2Q,w(t)=e^{- \pi \sigma t^2}
:{x_1}\left( q \right) = 0,2Q
:Q=\frac{1.9143}{\sqrt{\sigma}\Delta t}
:B=\frac{1.9143}{\sqrt{\sigma}}
*Matlab及python 皆可呼叫fft函式完成Y[m]=\sum\limits_{n = 0 }^{ N-1}y[n]e^{-j\frac{2\pi mn}{N}}
*演算法
假設t=n_0\Delta_t,(n_0+1)\Delta_t,\cdots \cdots ,(n_0+T-1)\Delta_t
: \,f=m_0\Delta_f,(m_0+1)\Delta_f,\cdots \cdots,(m_0+F-1)\Delta_f
:step 1:計算n_0,m_0,T,F,N,Q
:step 2:n=n_0
:step 3:決定x_1(q)
:step 4:X_1(m)=FFT[x_1(q)]
:step 5:轉換X_1(m)成X(n\Delta_t,m\Delta_f)
:step 6:設n=n+1 and return to Step 3 until n=n_0+T+1
限制
*避免贋頻效應(aliasing effect)
:(1){\Delta_t} (基本上任何實現方法都要避免贋頻效應)
:(2){\Delta _t}{\Delta _f} = {\textstyle{1 \over {N}}}
:(3)N=1/{\Delta _t}{\Delta _f} \ge 2Q+1
時間複雜度
O(TN{\log _2}N)
優缺點
:優點:時間複雜度低
:缺點:限制條件較直接實現法多
Chirp Z Transform
可改寫為:
由Direct Implementation可得下式
{X}\left( {n{\Delta _t},m{\Delta _f}} \right) = \sum\limits_{p = n-Q }^{ n+Q} {w\left( {(n - p){\Delta _t}} \right){x}\left( {p{\Delta _t}} \right)}{e^{ - j2\pi \,mp{\Delta _t}{\Delta _f}}}{\Delta _t}
e^{- \pi \sigma a^2} when\left| a\right| > 1.9143
Q=\frac{1.9143}{\sqrt{\sigma}\Delta t}
B=\frac{1.9143}{\sqrt{\sigma}}
令exp( - j2\pi \,mp{\Delta _t}{\Delta _f} ) = exp( -j\pi \, p^2{\Delta _t}{\Delta _f}) exp( j\pi \, {(p-m)}^2{\Delta _t}{\Delta _f}) exp( -j\pi \, m^2{\Delta _t}{\Delta _f})
可將式子改寫為:
{X}\left( {n{\Delta _t},m{\Delta _f}} \right) = {\Delta _t} \sum\limits_{p = n-Q }^{ n+Q} {w\left( {(n - p){\Delta _t}} \right){x}\left( {p{\Delta _t}} \right)}{e^{ - j2\pi \,mp{\Delta _t}{\Delta _f}}} \to {X}\left( {n{\Delta _t},m{\Delta _f}} \right) = {\Delta _t} {e^{ - j\pi \,m^2{\Delta _t}{\Delta _f}}}\sum\limits_{p = n-Q }^{ n+Q} {w\left( {(n - p){\Delta _t}} \right){x}\left( {p{\Delta _t}} \right)}{e^{ - j\pi \,p^2{\Delta _t}{\Delta _f}}}{e^{ j\pi \,{(p-m)}^2{\Delta _t}{\Delta _f}}}按此式即可實現
*演算法
:Step1:x_1[p] = w((n-p)\Delta_t)x(p\Delta_t)e^{-j\pi p^2 \Delta_t\Delta_f} \quad \quad n-Q \le p \le n+Q
:Step2:X_2[n,m] = \sum_{p=n-Q}^{n+Q}x_1[p]c[m-p] \quad \quad c[m]=e^{j\pi m^2 \Delta_t\Delta_f}
:Step3:X(n\Delta_t,m\Delta_f)=\Delta_t e^{-j\pi m^2 \Delta_t\Delta_f}X_2[m,n]
限制
*避免贋頻效應(aliasing effect)
(1){\Delta_t}
時間複雜度
O(TN{\log _2}N)
優缺點
:優點:限制條件與Direct Implementation法一樣基本上沒有限制
:缺點:時間複雜度與FFT-Based Method(快速傅立葉轉換)一樣
*但由於加伯轉換無法使用Recursive Method(遞迴法)所以此不能算是缺點
加伯轉換的性質
加伯轉換(Gabor Transform)是短時距傅立葉轉換(Short-Time Fourier Transform, STFT)的一種特殊形式,其定義為使用高斯函數(Gaussian function)作為窗函數(Window function)。定義如下:
:G_{x}(t,f)=\int_{-\infty}^{\infty}e^{-j2\pi f\tau}e^{-\pi(\tau-t)^{2}}x(\tau)d\tau
相比於其他窗函數,高斯函數在時頻分析中擁有獨特的地位。加伯轉換不僅繼承了 STFT 的特性,還因為高斯函數的數學完美性而具備了以下優良性質:
積分與還原性質
加伯轉換保留了訊號的所有資訊,這意味著我們可以從轉換後的結果中完全重建原始訊號。
- 還原性質(Recovery Property):當我們對加伯轉換結果在頻率軸上進行積分(並乘上適當係數)時,可以精確地還原出原始的時間域訊號 x(t)。這表示加伯轉換是一種可逆轉換(Invertible Transform),訊號在轉換過程中不會丟失任何成分。
:\int_{-\infty}^{\infty}G_{x}(t,f)e^{j2\pi tf}df=x(t)
- 廣義積分性質:更一般地,若引入參數 k,則對頻率的積分會得到經過高斯包絡調變後的訊號時移版本:
:\int_{-\infty}^{\infty}G_{x}(t,f)e^{j2\pi ktf}df=e^{-\pi(k-1)^{2}t^{2}}x(k t)
位移與調變性質
這兩個性質描述了加伯轉換對於時間位移和頻率位移的響應,說明了它在時頻平面上的「平移不變性」特徵。
- 位移性質(Shifting Property):若原始訊號在時間上延遲了 t_0,其加伯轉換的模(Magnitude)在時間軸上也會平移 t_0,但在相位上會產生一個與頻率相關的旋轉項 e^{-j2\pi ft_{0}}。
:若 y(t)=x(t-t_{0}),則 G_{y}(t,f)=G_{x}(t-t_{0},f)e^{-j2\pi ft_{0}}
- 調變性質(Modulation Property):若原始訊號在頻率上被調變(即乘上 e^{j2\pi f_{0}t}),其加伯轉換在頻率軸上會直接平移 f_0。這使得加伯轉換非常適合分析頻率隨時間變化的訊號(如調頻訊號)。
:若 y(t)=x(t)e^{j2\pi f_{0}t},則 G_{y}(t,f)=G_{x}(t,f-f_{0})
線性性質
加伯轉換滿足線性疊加原理。這意味著,如果一個訊號是由多個不同分量組成的(例如一個音樂和弦),那麼該訊號的加伯轉換等於各個分量加伯轉換的總和。這一性質使得我們能夠將複雜訊號分解為簡單的子訊號進行個別分析,而不會產生非線性的交叉項(Cross-terms)干擾。
:若 z(t)=\alpha x(t)+\beta y(t),則 G_{z}(t,f)=\alpha G_{x}(t,f)+\beta G_{y}(t,f)
能量衰減性質
這是加伯轉換最重要的特性之一,源於高斯窗函數的快速衰減特性。它保證了時頻分析的「局部化」能力。
- 時域局部化:若訊號在某個時間點 t_0 之後消失(即當 t>t_0 時,x(t)=0),那麼在該時間點之後,加伯轉換結果的能量平均值會隨著時間 t 遠離 t_0 而以高斯函數的形式快速衰減。
:|G_{x}(t,f)|^{2} 的平均值 的平均值
- 頻域局部化:同理,若訊號是頻帶受限的(Band-limited),其加伯轉換在頻帶外的能量也會快速衰減。這特性極大程度地減少了頻譜洩漏(Spectral Leakage),使得加伯轉換的時頻圖看起來非常清晰,干擾極小。
能量守恆性質
加伯轉換滿足帕塞瓦爾定理(Parseval's Theorem)的推廣形式。這表示訊號在時頻域中的總能量等於其在時域(或頻域)中的總能量。
:\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}G_{x}(t,f)G_{y}^{}(t,f)df dt=\int_{-\infty}^{\infty}x(\tau)y^{}(\tau)d\tau
當 y(t)=x(t) 時,上式說明了能量守恆:
:\iint |G_{x}(t,f)|^2 df dt = \int |x(\tau)|^2 d\tau
這意味著加伯轉換後的能量分佈真實反映了原始訊號的能量強度,適合用於能量分析。
最佳化測不準原理
這通常被認為是加伯轉換最強大的優點。根據海森堡測不準原理(Heisenberg Uncertainty Principle),一個信號不能在時間域和頻率域同時擁有無限小的寬度。時間寬度 \sigma_t 和頻率寬度 \sigma_f 的乘積有一個下限:
:\sigma_t \sigma_f \ge \frac{1}{4\pi}
在所有可能的窗函數中,只有高斯函數(Gaussian function)能夠達到這個下限的等號(\sigma_t \sigma_f = \frac{1}{4\pi})。這意味著使用高斯函數作為窗函數的加伯轉換,在理論上能夠提供最佳的時頻解析度,是在時間解析度和頻率解析度之間取得平衡的最優選擇。
特徵函數與對稱性
高斯函數是傅立葉轉換的特徵函數(Eigenfunction),即高斯函數的傅立葉轉換仍然是高斯函數。
:\mathcal{F}\{e^{-\pi t^2}\} = e^{-\pi f^2}
這使得加伯轉換在時間域和頻率域的數學形式具有高度的對稱性。這種對稱性在光學、雷達訊號處理等領域具有重要的物理意義,因為它確保了訊號在經過透鏡或傳播過程(類比於傅立葉轉換)後,其波包形狀保持不變。
優缺點
Gabor Transform 的優點
- 最佳時間-頻率局部化特性
** Gabor Transform 使用高斯窗函數,與其他常見窗函數(如Rectangle、Triangle、Hanning、Hamming)相比,滿足測不準定理的最小下限(Minimum Uncertainty Principle)。這意味著,高斯函數能夠在時間域和頻率域中同時提供最佳的解析度,避免信號特徵的模糊或失真。
*** 高時間分辨率:能捕捉信號的快速變化,對於瞬態信號(如語音中的短促音位或振動信號中的瞬時變化)尤為重要。
*** 高頻率分辨率:能精確分辨信號中的穩態頻率成分,特別適合於分析連續且平穩的周期信號。
- 算法穩健且實現簡單
** Gabor Transform 基於傅里葉變換的數學理論,其結構清晰、明了且實現相對簡單。現代數值計算技術(如快速傅里葉變換,FFT)的發展進一步提升了 Gabor Transform 的計算效率,使其能夠在高效實現的同時保持穩健性。
** 穩健性:由於其依賴於成熟的數學基礎,在實施中容易檢測和修正潛在錯誤。
** 實現便利性:現有的數學工具庫(如 MATLAB、Python 的 Scipy、Octave)提供了高度封裝的 Gabor Transform 函數,大幅降低了實現門檻,讓開發者能更專注於應用場景設計,而非底層算法調試。
- 廣泛的應用場景
*# 語音去噪:利用 Gabor Transform 可以有效提取語音信號的時頻特徵,通過將語音信號分解為多個頻帶,對噪聲進行有效抑制,從而提升語音的清晰度和識別準確度,特別是在低信噪比環境下
*# 圖像處理
# 紋理分析:有效捕捉圖像的方向與頻率特徵,用於紋理分類和圖像分割。
# 邊緣檢測:適用於醫學圖像和場景理解,改善邊緣檢測效果。
*# 機械振動信號分析
# 故障檢測:由於Gabor Transform能夠提供高時間和頻率解析度,它能有效捕捉非平穩信號中的瞬時頻率變化。這使得它特別適合用於檢測如軸承、齒輪等機械部件的故障。轉換後的信號圖像可以作為特徵輸入至卷積神經網絡(CNN),進行自動化分類和故障診斷。
Gabor Transform 的缺點
- 計算複雜度較高
** Gabor Transform 在處理高維數據(如圖像信號處理)時,計算複雜度可能大幅增加。每個窗函數的計算都需要執行一次傅立葉變換,這對於大數據集或實時應用場景來說,可能會成為系統性能的瓶頸。
** 在圖像處理中,Gabor 變換通常需要對圖像的不同尺度和方向應用一組 Gabor 濾波器,以提取豐富的特徵信息。這意味著每個尺度和方向都需要單獨進行濾波操作,隨著濾波器數量的增加,計算量會線性增長。此外,對於高分辨率圖像,每次濾波操作都需要處理大量像素,從而進一步增加了計算負擔。
**為了提高計算效率,基於離散傅立葉變換(DFT)的快速算法應運而生,快速算法用於二維離散 Gabor 變換。可以顯著降低了計算複雜度
- 解析度折衷的不可避免性
** 根據測不準定理,Gabor Transform 的時間和頻率分辨率達到了理論的最佳折衷,但這也意味著:
*** 受測不凖定理約束,當需要同時對信號的快速變化與細微頻率差異進行精確分析時,時間和頻率的分辨率會有可能不足以同時滿足所有需求。
*** 相較於 Gabor Transform, Wigner Distribution Function(WDF)等方法,因是對訊號的自相關函數做傅立葉轉換,可以超越測不準原理約束的下限,因此能提供更高的時頻解析度,尤其是對於結構複雜的信號。然而,WDF 的非線性特性容易引入交叉干擾項(cross-terms),而為了為了結合兩者的優點,Gabor Wigner Transform應運而生
参见
- 闵可夫斯基空间
- 柯西不等式
- 三角不等式
- 完备空间
參考書目、資料來源
评论 (0)