基于MATLAB的腦電信號帶通濾波器的仿真與比較

張井想++梁雨++許朋

摘 要 腦電信號的頻率范圍一般在0.5～35Hz。針對腦電信號的頻率特性，用MATLAB的信號處理工具箱的函數設計巴特沃斯、切比雪夫和橢圓函數帶通濾波器。通過仿真及結果的分析得出更加適合腦電信號軟件處理的帶通濾波器。

【關鍵詞】腦電信號MATLAB 帶通濾波器 仿真

1 引言

腦部疾病長期以來一直威脅著人類的健康，因此對其預防和及時發現在減少腦部疾病危害中極為重要。腦電信號受到工頻50Hz信號等噪聲及生理信號等干擾，故能夠設計一種可以分離出腦電信號的帶通濾波器就顯得尤為重要。本文從巴特沃斯帶通濾波器、切比雪夫帶通濾波器和橢圓帶通濾波器的MATLAB仿真及結果的分析得出更加適合腦電信號軟件編程處理的帶通濾波器，從而獲取更加純凈的腦電信號。

2 MATLAB簡介

MATLAB語言是一種面向工程與科學的計算語言。MATLAB信號處理工具箱提供了設計巴特沃斯、切比雪夫和橢圓函數濾波器等函數，本文利用這些函數，進行了巴特沃斯、切比雪夫和橢圓函數濾波器的程序設計，并通過仿真結果并分析這三種濾波器的優缺點及適用場合。

2.1 巴特沃斯帶通濾波器的設計與仿真

2.1.1 butter函數

butter函數是用于設計巴特沃斯的濾波器

[k， l] = butter（n， Wn）；

當Wn = [W1W2]時，它可以設計2n 階的巴特沃斯帶通濾波器，其通帶為 W1 < W < W2。

2.1.2 巴特沃斯帶通濾波器設計

設計一個巴特沃斯帶通濾波器，繪制原始信號、濾波后信號FFT、原始信號FFT和歸一化的信噪比圖。20Hz和50Hz 正弦波組成原始信號，將噪聲 50Hz 的正弦波濾掉，通過函數 butter 設計一組帶通濾波器系數，其階數是2，0.5Hz 到 35Hz的通帶頻率，采樣率 1Kbps。程序代碼如下：

fc=1000； %設置1k的采樣頻率

N=1024； %采樣點數

n=0：N-1；

t=0：1/fc：1-1/fc； %時間序列

f=n*fc/N； %頻率序列

X1=sin（2*pi*50*t）； %50Hz的噪聲

X2=sin（2*pi*20*t）； %20Hz的信號

X=X1+X2； %信號混合

subplot（221）；

plot（t，X）； %繪制原始信號

xlabel（時間）；

ylabel（幅值）；

title（原始信號）；

grid on；

subplot（222）；

Y=fft（X，N）； %繪制原始信號的幅頻響應

plot（f，abs（Y））；

xlabel（頻率/Hz）；

ylabel（振幅）；

title（原始信號 FFT）；

grid on；

subplot（223）；

Wn=[0.5*2 35*2]/fs； %設置通帶 0.5Hz 到 35Hz

[k，l]=butter（1，Wn）； %注意第一個參數雖然是 1，但生成的卻是 2 階 IIR 濾波器系數

Y2=filtfilt（k，l，X）； %計算濾波后的波形 y2

Y3=fft（Y2，N）； %濾波后波形的幅頻響應

plot（f，abs（Y3））；

xlabel（頻率/Hz）；

ylabel（振幅）；

title（濾波后信號 FFT）；

grid on；

[H，F]=freqz（k，l，512）；

subplot（224）；

plot（F/pi，abs（H））；

xlabel（歸一化頻率）； %繪制絕對幅頻響應

ylabel（幅度）；

P1=sum（X2.^2）； %腦電信號的總功率

P2=sum（（Y2-X2）.^2）； %剩余噪聲的功率

SNR=10*log10（P1/P2）； %腦電信號的總功率和剩余噪聲的功率比值

title（[Order=，int2str（2）， SNR=，num2str（SNR）]）；

grid on；

Matlab 的仿真結果知，SNR=4.5311。

2.2 切比雪夫帶通濾波器設計與仿真

2.2.1 cheby1 函數

cheby1 函數是用來設計切比雪夫I 型的濾波器。

[k，l] = cheby1（n， Rp， Wn）；cheby1 函數可以設計帶通切比雪夫I 型數字濾波器，其阻帶內為單調，通帶內為等紋波。

[k，l] = cheby1（n， Rp， Wn）；Rp 用來確定通帶內的紋波，該濾波器的截止頻率是Wn，可以設計 n 階低通切比雪夫I 型數字濾波器。當 Wn=[W1， W2]時，cheby1 函數可以設計出其通帶為 W1

2.2.2 切比雪夫帶通濾波器設計

設計一個巴特沃斯帶通濾波器，繪制原始信號、濾波后信號FFT、原始信號FFT和歸一化的信噪比圖。20Hz和50Hz 正弦波組成原始信號，將噪聲 50Hz 的正弦波濾掉，通過函數 cheby1 設計一組帶通濾波器系數，其階數是2，0.5Hz 到 35Hz的通帶頻率，采樣率 1Kbps，通帶紋波 1db。 Matlab程序代碼同上，不同部分如下：

Wn=[0.5*2 35*2]/fs； %設置通帶 0.5Hz 到 35Hz

[k，l]=cheby1（1，1，Wn）； %注意第一個參數雖然是 1，但生成的卻是 2 階 IIR 濾波器系數

Matlab 的仿真結果知，SNR=8.4301。

2.3 橢圓帶通濾波器設計與仿真

2.3.1 ellip函數

ellip函數是用來設計橢圓型的濾波器

[k，l] = ellip（n， Rp， RS， Wn）；

2.3.2 橢圓帶通濾波器設計

設計一個巴特沃斯帶通濾波器，繪制原始信號、濾波后信號FFT、原始信號FFT和歸一化的信噪比圖。20Hz和50Hz 正弦波組成原始信號，將噪聲 50Hz 的正弦波濾掉，通過函數 ellip 設計一組帶通濾波器系數，其階數是2，0.5Hz 到 35Hz的通帶頻率，采樣率 1Kbps，阻帶40db，通帶紋波 1db。Matlab程序代碼同上，不同部分如下：

Wn=[0.5*2 35*2]/fs； %設置通帶 0.5Hz 到 35Hz

[k，l] = ellip（2， 1， 40， Wn）；

Matlab 的仿真結果知，SNR=8.4434。

3 結束語

通過三個帶通濾波器的仿真圖的信噪比比較可知，巴特沃斯函數設計的數字帶通濾波器具有最大的平坦幅度，但其截止頻率處的下降斜度會受到損失，使幅度響應衰減較慢。與切比雪夫和巴特沃斯濾波器相比，橢圓函數只需較低的階數可以設計出衰減更快、下降斜度更大的濾波器，但通帶和阻帶內均為等紋波。綜上，高階巴特沃斯帶通濾波器能滿足腦電信號對通帶內的幅度響應平坦及截止頻率處的下降斜率的要求。

參考文獻

[1]李鐘慎.基于MATLAB設計巴特沃斯低通濾波器[J].信息技術，2003（03）：49-50+52.

[2]劉凌云，趙鵬宇，弓美桃.基于MATLAB的低通巴特沃斯濾波器仿真[J].數字技術與應用，2013（02）：124.

作者簡介

張井想（1991-），男，工程碩士學位。主要研究方向為醫療電子。

作者單位

江蘇師范大學 江蘇省徐州市 221116