時域與頻域

非周期性函數的頻率成份是連續的，我們可使用傅立葉轉換來計算其頻率密度函數，假設時域函數為 $s(t)$，其傅立葉轉換公式如下：

$$S(f) = \int_{-\infty}^{\infty}f(x)\ e^{-j2\pi ft}\ dt$$

一般在信號處理時，我們多使用快速傅立葉轉換 (FFT) 來計算其頻譜，使用快速傅立葉轉換時，應牢記以下幾個重點：

1. 假設用來分析的信號段有 $N$ 個點，且其取樣間隔為 $\Delta t$，則此信號段總長為 $T=N\Delta t$，且取樣速度為 $R=1/\Delta t$。

2. 根據取樣定理，可以觀察的頻率範圍為 $[-f_0,f_0]$，$f_0=R/2$。

3. $N$ 個時域點做 FFT 之後變成 $N$ 個頻域點，故頻率解析度，亦即相鄰頻率點的間隔為 $2f_0/N = R/N = 1/(N\Delta t) = 1/T$。
   簡單說，取樣時間總長的倒數就是頻率取樣點之間的距離，也就是頻率解析度。

4. 要記得取樣速度的一半，是可觀察的最大頻率；取樣總長的倒數是頻率解析度。

5. 如果是實函數的話，做完 FFT 之後，其正頻率和負頻率相對的值會是共軛複數。

我們現在來觀察帽子函數的頻譜，定義帽子函數 rect(t) 如下：

使用 MATLAB 觀察其頻譜步驟如下：

1. 假設取樣範圍從-1.5到1.5，每隔0.1取一點，先計算時域信號：

t1 = -1.5:0.1:-0.6;
t2 = -0.4:0.1:0.4;
t3 =  0.6:0.1:1.5;
t = [t1 -0.5 t2 0.5 t3];
s = [zeros(size(t1)) 0.5 ones(size(t2)) 0.5 zeros(size(t3))];
plot(t, s);

這個程式片段寫得有點繁瑣，不過很容易理解，看一下結果是否正確。

2. 計算頻譜，直接使用 fft 函數將上面的時域信號轉到頻域，應注意轉出來的結果為複數，我們先觀察其振幅。

sf = fft(s);
plot(abs(sf));

注意在上面的程式中，畫圖時只給了y軸，橫軸會是點數，得到的結果如下：

3. 上面提過實函數的 FFT，其正頻率和負頻率相對的值是共軛複數。實際上我們現在看到的圖形中，圖的右半邊應該移到負的那邊去才對。MATLAB 里面有一個 fftshift 就做這事。另外，最大頻率為取樣率一半=5，間隔為總長倒數=1/3，我們要把橫軸 (頻率) 加上去。

sf = fft(s);
sf = abs(fftshift(sf));
f = -5:1/3:5;
plot(f, sf);
grid on;

最後得到修正的圖如下，這個圖形就是 sinc 函數的振幅。

完整程式碼如下：

t1 = -1.5:0.1:-0.6;
t2 = -0.4:0.1:0.4;
t3 =  0.6:0.1:1.5;
t = [t1 -0.5 t2 0.5 t3];
s = [zeros(size(t1)) 0.5 ones(size(t2)) 0.5 zeros(size(t3))];
sf = fft(s);
sf = abs(fftshift(sf));
f = -5:1/3:5;
plot(f, sf);
grid on;

練習 4

請分析時域函數 rect$((t-1)/2)\cdot\cos(20\pi t)$ 的頻譜。