DSP筆記4_能量頻譜密度 Energy Spectral Density(ESD)_從時域訊號到幅度大小的平方

 

能量頻譜密度 Energy Spectral Density(ESD)

用於描述訊號在不同頻率下能量分布情形。

用 $S_E(f)$ 來表示，X(f)為原先時間域x(t)的訊號做傅立葉轉換的結果。

能量譜密度的總能量關係滿足如下式子，就稱作帕塞瓦爾定里Parseval's theorem)。

總能量會等同於對時域跟頻域函數從負無窮到正無窮取平方的積分

所以

練習程式ver1.矩形脈衝Energy Signal

import numpy as np
import matplotlib.pyplot as plt

fs = 1000 #Sampling Frequency
T = 1        # Duration (seconds)
t = np.linspace(0, T, int(T*fs), endpoint=False)

#建立一個全為0的訊號
x = np.zeros(len(t))
x[int(0.4*fs):int(0.6*fs)] = 1 # 0.4 ~ 0.6 秒設定為 1

# Plot
plt.figure(figsize=(8, 4))
plt.plot(t, x)

plt.title("Rectangular Pulse")
plt.xlabel("Time (s)")
plt.ylabel("Amplitude")

plt.grid(True)
plt.show()

運行效果

可以先得知目前此時間訊號的函數區間定義如下

$$x(t)= \begin{cases} 1, & 0.4\leq t<0.6\\ 0, & \text{其他時間} \end{cases}$$

也就是一個寬度

$$0.6-0.4=0.2\text{ 秒}$$

的矩形脈衝，是一個典型的 Energy Signal。。

它的能量為：

$$E=\int_{-\infty}^{\infty}|x(t)|^2dt$$

在這個例子裡：

$$E = \int_{0.4}^{0.6}1^2dt = 0.2$$

所以我們甚至還沒做 FFT，就已經知道：

練習程式ver2.做到 FFT Magnitude Spectrum，算出∣X(f)∣ 

將時域訊號 x(t)  ，透過 FFT 轉換到頻域 X(f)

做完FFT結果的 X 本身通常是複數，因此多做計算複數絕對值，也就是各個頻率成分的幅度大小（Magnitude）。

import numpy as np
import matplotlib.pyplot as plt
from numpy.fft import fft, fftfreq

fs = 1000 # Sampling Frequency
T = 1     # Duration (seconds)

t = np.linspace(0, T, int(T*fs), endpoint=False)
freq = np.fft.fftfreq(len(t), 1/fs)

# Rectangular pulse
x = np.zeros(len(t))
x[int(0.4*fs):int(0.6*fs)] = 1  # pulse from 0.4s to 0.6s (using sampling frequency)

# FFT
X = np.fft.fft(x) #頻域 X(f)
# Magnitude Spectrum
magnitude = np.abs(X)

# Plot
plt.figure(figsize=(8,4))
plt.plot(freq[:len(freq)//2], magnitude[:len(freq)//2])
plt.title("Magnitude Spectrum of Rectangular Pulse")
plt.xlabel("Frequency (Hz)")
plt.ylabel("Magnitude")
plt.grid(True)
plt.show()

運行效果

練習程式ver3.從 FFT 推導到 ESD，∣X(f)∣ 多做平方

import numpy as np
import matplotlib.pyplot as plt
from numpy.fft import fft, fftfreq

fs = 1000 #Sampling Frequency
T = 1        # Duration (seconds)
t = np.linspace(0, T, int(T*fs), endpoint=False)
freq = np.fft.fftfreq(len(t), 1/fs)

# Rectangular pulse
x = np.zeros(len(t))
x[int(0.4*fs):int(0.6*fs)] = 1  # pulse from 0.4s to 0.6s (using sampling frequency)

# FFT and ESD
X = np.fft.fft(x)
ESD = (np.abs(X)**2)

# Plot
plt.figure(figsize=(8,4))
plt.plot(freq[:len(freq)//2], ESD[:len(freq)//2])
plt.title("ESD of Rectangular Pulse")
plt.xlabel("Frequency (Hz)")
plt.ylabel("Energy")
plt.grid(True)
plt.show()

要觀察的不只是 Fourier Transform 的幅度，而是「能量在不同頻率上的分布」，根據 ESD 定義：

$$S_E(f)=|X(f)|^2$$

因此只需要將 Magnitude 再平方：

運行效果