一、實驗目的
1、掌握Matlah由各種窗數序列的生成方法;
2、掌握窗函數序列頻率特性的計算與畫圖方法;
3、掌握窗函數的相關參數對窗函數頻域性能的影響;
4、了解混合窗函數的定義、生成方法和頻域性能特點;
5、對各種常用的窗函數序列的時域及頻域性能特點有整體的認識;
6、了解窗函數序列實現頻率響應的線性相位特性需要滿足的基本條件。
二、實驗原理
1、常用窗函數序列
加窗處理是信號處理中常用的一種處理方式，在數字濾波器的設計、信號的離散譜分析中都有很廣泛的應用。常用的窗函數序列主要有:
(1）矩形窗函數序列
矩形窗序列是最基本形式的一種窗函數序列，其時域表示如
公式所示:
(2）漢寧窗函數序列
漢寧窗序列也稱為升余弦窗，在矩形窗序列的基礎上，增加一個余弦項，其時域表示如公式所示
(3）海明窗函數序列
海明窗序列和漢寧窗序列類似，只是窗函數的兩部分的系數有所不同,其時域表示如公式所示:
(4）布萊克曼窗序列
布萊克曼窗序列進一步增加一個余弦項,其時域表示如公式所示:
(5）凱瑟窗序列
凱瑟窗序列和前面4種類型的窗函數序列在原理上有所不同,
其建立在第一類修正的零階貝塞爾函數基礎上,其時域表示如公式所示:
公式中，I0(g)是第一類修正的零階貝塞爾函數，β為調整參數,用以調整窗函數序列的性能
在 Matlab中,各種窗函數序列的生成均有相對應的函數，對應上述的各窗函數序列的分別為;
(1）矩形窗序列:boxcar()
(2)漢寧窗序列:hanning()
(3）海明窗序列:hamming()(漢明窗)
(4）布萊克曼窗序列:blackman()
(5）凱瑟窗序列:kaiser()
這幾個窗函數的調用格式比較簡單,其中,前4個函數只要給出窗函數的長度N即可,如:
Window=hamming(N);
該函數產生一個長度為N的海明窗序列,結果保存在輸出向量Window中。
kaiser()函數在使用中不僅要給出窗序列的長度N,還需給出窗函數的參數β，即其調用格式為:
Window=kaiser(N,belta);
其中,輸入參數belta即為窗函數的參數β，一般在4~8之間取值,用以調整凱瑟窗序列的性能。
2、窗函數的頻域特性
將各種窗函數的時域表示式進行離散時間傅里葉變換,即可得到窗函數的頻域表示。各種窗函數的幅頻特性一般都具有下圖的形式:
衡量窗函數的性能,通常通過兩個指標:
（1）主瓣寬度。
即窗函數幅頻特性曲線第一至第二過零點之間的頻帶寬度。在濾波器設計中,主瓣越寬，則濾波器過渡帶越寬;在頻譜分析中，主瓣越寬，則頻率分辨力越差。所以在實際應用中,總是希望窗函數有盡量窄的主瓣寬度。
(2）副瓣衰減。
一般用第一副瓣相對于主瓣的衰減進行表示。在濾波器設計中,副瓣衰減越大，則濾波器幅頻響應在通帶之內的平穩程度和阻帶衰減特性就越好;在頻譜分析中，副瓣衰減越大,則頻譜泄漏效應就會減弱。所以在實際應用中，總是希望窗函數有盡量大的副瓣衰減。窗函數的主全寬度和副瓣衰減性能之間往往是相互矛盾的,提高了一個性能,另外一個性能就會有所下降。
在Matlab中,窗函數序列頻率特性的計算可由 freqz）函數來實現。以下以漢寧窗序列為例，給出窗函數序列的頻率特性計算及畫圖的仿真程序實現方法,如下:
N=10;%設置窗函數的寬度成漢寧窗
wn=hanning(N);%生成頻率特性
[H,m] = freqz(wn, [1],1024,'whole');%計算窗函數的算幅度特性
mag = abs(H);
db =20*log10((mag"+eps)/max(mag));%將幅度特性轉換為對數形式算相位特性
pha=angle(H);
subplot(121)%以下代碼為繪制幅度特性曲線plot(m/pi,db); axis([01 -100 0]);
xlabel('w/pi');ylabel('dB');title('幅度特性); grid on;
subplot(122)%以下代碼為繪制相位特性曲線plot(m/pi,pha); axis([0 1 -pi pi]);
xlabel ('w/pi');ylabel('rad');title('相位特性'); grid on
各種窗函數序列的相位特性一般都是線性的（或分段線性）,即為頻率的線性函數,也稱為線性相位特性。但是,要保證窗函數的線性相位特性,需要保證窗函數序列是中心對稱的,即需要滿足:
3、混合窗函數序列簡介
混合窗函數序列是目前在信號處理、濾波器設計等領域經常提及的一種窗函數序列,其幅度特性較常用的窗函數序列更為理想一些。濾合窗承數序列的定義如下所示:
公式中，窗函數的前半部分為矩形窗利疊加，而后半部分則是標準的海明窗的后一半。α、β和γ是可調整系數,會影響窗函數的主那寬度和副瓣衰減, 這三個參數在實際應用中可根據實際需求進行調整，但需滿足α+β+γ=1。由混稱的,因而其相合窗函數的定義也可看出混合窗函數是非中心對位特性是非線性的,從本質上講是一種以犧牲相位特性來提升幅度特性的窗函數序列
三、實驗步驟、數據記錄及處理
本實驗以Matlab信號處理工具箱所提供的窗函數為基礎，結合頻率特性分析,對常用窗函數的時域和頻域特性進行分析和對驗內容和實驗步比,并對相關的結論與參數進行驗證。具體的步驟如下:
1、設置窗函數的長度N=15,分別生成矩形窗、漢寧窗、海明窗和布萊克曼窗序列，繪制各窗函數序列的時域圖形,觀察各種窗函數的時域性能。
2、利用頻率響應計算的函數或傅立葉變換函數,分別計算實驗內容1中所生成的各窗函數序列的頻率特性,繪制幅頻和相頻時對各種窗函數響應曲線，驗證理論課程中的相關性能指標，同的頻域性能進行對比。
3、設置窗函數的長度N=33,重復上述過程, 分析窗序列的長度對窗函數頻域性能的影響。
4、生成N=33，β分別為B=4,5,6,7,8的凱瑟窗,分別計算并繪制幅頻特性曲線，分析調整參數β對窗函數頻域性能的影響。
5、設置窗序列長度N=31,生成漢寧窗。設置N=31,a=0.42,，生成混合窗序列。畫出這兩個窗序列的時域圖β=0.4，和y=0.18形，并進行對比。
6、計算并繪制5中生成的兩個窗函數的幅度特性和相位特性,對兩種窗函數的幅度性能進行對比。同時,觀察兩個窗函數的相位特性,驗證窗函數取得線性相位特性需滿足的條件。
實驗程序：
clc;clear all;
N=15;%設置窗函數序列長度；
%矩形窗序列
window1=boxcar(N);%生成窗函數序列；
[H1,m1]=freqz(window1,[1],1024,'whole');%計算函數的頻率特性
mag1=abs(H1);%計算函數的幅度特性；
pha1=angle(H1);%計算函數的相位特性；
db1=20*log10((mag1+eps)/max(mag1));%將幅度特性轉換為對數形式
%繪制矩形窗頻域波形
figure('name','矩形窗頻域波形');
subplot(311);stem(window1);grid on;xlabel('n');ylabel('x(n)');title('N=15矩形窗函數時域圖形');
subplot(312);plot(m1/(2*pi),db1);grid on;xlabel('w/(2*pi)');ylabel('ab');title('N=15矩形窗幅頻響應');
subplot(313);plot(m1/(2*pi),pha1);grid on;xlabel('w/(2*pi)');ylabel('rad');title('N=15矩形窗相頻響應');%漢寧窗序列
window2=hanning(N);%生成窗函數序列；
[H2,m2]=freqz(window2,[1],1024,'whole');%計算函數的頻率特性
mag2=abs(H2);%計算函數的幅度特性；
pha2=angle(H2);%計算函數的相位特性；
db2=20*log10((mag2+eps)/max(mag2));%將幅度特性轉換為對數形式
%繪制漢寧窗頻域波形
figure('name','漢寧窗頻域波形');
subplot(311);stem(window2);grid on;xlabel('n');ylabel('x(n)');title('N=15漢寧窗函數時域圖形');
subplot(312);plot(m2/(2*pi),db1);grid on;xlabel('w/(2*pi)');ylabel('ab');title('N=15漢寧窗幅頻響應');
subplot(313);plot(m2/(2*pi),pha1);grid on;xlabel('w/(2*pi)');ylabel('rad');title('N=15漢寧窗相頻響應');%漢明窗序列
window3=hamming(N);%生成窗函數序列；
[H3,m3]=freqz(window3,[1],1024,'whole');%計算函數的頻率特性
mag3=abs(H3);%計算函數的幅度特性；
pha3=angle(H3);%計算函數的相位特性；
db3=20*log10((mag3+eps)/max(mag3));%將幅度特性轉換為對數形式
%繪制漢明窗頻域波形
figure('name','漢明窗頻域波形');
subplot(311);stem(window3);grid on;xlabel('n');ylabel('x(n)');title('N=15漢明窗函數時域圖形');
subplot(312);plot(m3/(2*pi),db1);grid on;xlabel('w/(2*pi)');ylabel('ab');title('N=15漢明窗幅頻響應');
subplot(313);plot(m3/(2*pi),pha1);grid on;xlabel('w/(2*pi)');ylabel('rad');title('N=15漢明窗相頻響應');%布萊克曼窗序列
window4=blackman(N);%生成窗函數序列；
[H4,m4]=freqz(window4,[1],1024,'whole');%計算函數的頻率特性
mag4=abs(H4);%計算函數的幅度特性；
pha4=angle(H4);%計算函數的相位特性；
db4=20*log10((mag4+eps)/max(mag4));%將幅度特性轉換為對數形式
%繪制布萊克曼頻域波形
figure('name','布萊克曼頻域波形');
subplot(311);stem(window4);grid on;xlabel('n');ylabel('x(n)');title('N=15布萊克曼窗函數時域圖形');
subplot(312);plot(m4/(2*pi),db1);grid on;xlabel('w/(2*pi)');ylabel('ab');title('N=15布萊克曼窗幅頻響應');
subplot(313);plot(m4/(2*pi),pha1);grid on;xlabel('w/(2*pi)');ylabel('rad');title('N=15布萊克曼窗相頻響應');
clc;clear all;
N=33;%設置窗函數序列長度；
%矩形窗函數
window1=boxcar(N);%生成窗函數序列；
[H1,m1]=freqz(window1,[1],1024,'whole');%計算函數的頻率特性
mag1=abs(H1);%計算函數的幅度特性；
pha1=angle(H1);%計算函數的相位特性；
db1=20*log10((mag1+eps)/max(mag1));%將幅度特性轉換為對數形式
%繪制矩形窗頻域波形
figure('name','矩形窗頻域波形');
subplot(311);stem(window1);grid on;xlabel('n');ylabel('x(n)');title('N=33矩形窗時域圖形');
subplot(312);plot(m1/(2*pi),pha1);grid on;xlabel('w/2pi');ylabel('ψ(w)');title('N=33矩形窗的相頻特性');
subplot(313);plot(m1/(2*pi),db1);grid on;xlabel('w/2pi');ylabel('db');title('N=33矩形窗的幅度特性');%漢寧窗函數
window2=hanning(N);%生成窗函數序列；
[H2,m2]=freqz(window2,[1],1024,'whole');%計算函數的頻率特性
mag2=abs(H2);%計算函數的幅度特性；
pha2=angle(H2);%計算函數的相位特性；
db2=20*log10((mag2+eps)/max(mag2));%將幅度特性轉換為對數形式
%繪制漢寧窗頻域波形
figure('name','漢寧窗頻域波形');
subplot(311);stem(window2);grid on;xlabel('n');ylabel('x(n)');title('N=33的漢寧窗函數時域圖形');
subplot(312);plot(m2/(2*pi),pha2);grid on;xlabel('w/2pi');ylabel('ψ(w)');title('N=33的漢寧窗函數的相頻特性');
subplot(313);plot(m2/(2*pi),db2);grid on;xlabel('w/2pi');ylabel('db');title('N=33的漢寧窗函數的幅度特性');%漢明窗函數
window3=hamming(N);%生成窗函數序列；
[H3,m3]=freqz(window3,[1],1024,'whole');%計算函數的頻率特性
mag3=abs(H3);%計算函數的幅度特性；
pha3=angle(H3);%計算函數的相位特性；
db3=20*log10((mag3+eps)/max(mag3));%將幅度特性轉換為對數形式
%繪制漢明窗頻域波形
figure('name','漢明窗頻域波形');
subplot(311);stem(window3);grid on;xlabel('n');ylabel('x(n)');title('N=33的漢明窗函數時域圖形');
subplot(312);plot(m3/(2*pi),pha3);grid on;xlabel('w/2pi');ylabel('ψ(w)');title('N=33的漢明窗函數的相頻特性');
subplot(313);plot(m3/(2*pi),db3);grid on;xlabel('w/2pi');ylabel('db');title('N=33的漢明窗函數的幅度特性');%布萊克曼窗函數
window4=blackman(N);%生成窗函數序列；
[H4,m4]=freqz(window4,[1],1024,'whole');%計算函數的頻率特性
mag4=abs(H4);%計算函數的幅度特性；
pha4=angle(H4);%計算函數的相位特性；
db4=20*log10((mag4+eps)/max(mag4));%將幅度特性轉換為對數形式
%繪制布萊克曼頻域波形
figure('name','布萊克曼頻域波形');
subplot(311);stem(window4);grid on;xlabel('n');ylabel('x(n)');title('N=33的布萊克曼窗函數時域圖形');
subplot(312);plot(m4/(2*pi),pha4);grid on;xlabel('w/2pi');ylabel('ψ(w)');title('N=33的布萊克曼窗函數的相頻特性');
subplot(313);plot(m4/(2*pi),db4);grid on;xlabel('w/2pi');ylabel('db');title('N=33的布萊克曼窗函數的幅度特性');
clc;clear all;
N=33;%設置窗函數序列長度；
belta1=4;
window1=kaiser(N,belta1);%凱瑟窗函數
[H1,m1]=freqz(window1,[1],1024,'whole');%計算函數的頻率特性
mag1=abs(H1);%計算函數的幅度特性；
pha1=angle(H1);%計算函數的相位特性；
db1=20*log10((mag1+eps)/max(mag1));%將幅度特性轉換為對數形式
figure('name','=4時凱瑟窗函數頻域波形');
subplot(311);stem(window1);grid on;xlabel('n');ylabel('x(n)');title('β=4的凱瑟窗函數時域圖形');
subplot(312);plot(m1/(2*pi),pha1);grid on;xlabel('w/2pi');ylabel('ψ(w)');title('β=4的凱瑟窗函數的相頻特性');
subplot(313);plot(m1/(2*pi),db1);grid on;xlabel('w/2pi');ylabel('db');title('β=4的凱瑟窗函數的幅度特性');belta2=5;
window2=kaiser(N,belta2);%凱瑟窗函數
[H2,m2]=freqz(window2,[1],1024,'whole');%計算函數的頻率特性
mag2=abs(H2);%計算函數的幅度特性；
pha2=angle(H2);%計算函數的相位特性；
db2=20*log10((mag2+eps)/max(mag2));%將幅度特性轉換為對數形式
figure('name','=5時凱瑟窗函數頻域波形');
subplot(311);stem(window2);grid on;xlabel('n');ylabel('x(n)');title('β=5的凱瑟窗函數時域圖形');
subplot(312);plot(m2/(2*pi),pha2);grid on;xlabel('w/2pi');ylabel('ψ(w)');title('β=5的凱瑟窗函數的相頻特性');
subplot(313);plot(m2/(2*pi),db2);grid on;xlabel('w/2pi');ylabel('db');title('β=5的凱瑟窗函數的幅度特性');belta3=6;
window3=kaiser(N,belta3);%凱瑟窗函數
[H3,m3]=freqz(window3,[1],1024,'whole');%計算函數的頻率特性
mag3=abs(H3);%計算函數的幅度特性；
pha3=angle(H3);%計算函數的相位特性；
db3=20*log10((mag3+eps)/max(mag3));%將幅度特性轉換為對數形式
figure('name','=6時凱瑟窗函數頻域波形');
subplot(311);stem(window3);grid on;xlabel('n');ylabel('x(n)');title('β=6的凱瑟窗函數時域圖形');
subplot(312);plot(m3/(2*pi),pha3);grid on;xlabel('w/2pi');ylabel('ψ(w)');title('β=6的凱瑟窗函數的相頻特性');
subplot(313);plot(m3/(2*pi),db3);grid on;xlabel('w/2pi');ylabel('db');title('β=6的凱瑟窗函數的幅度特性');belta4=7;
window4=kaiser(N,belta4);%凱瑟窗函數
[H4,m4]=freqz(window4,[1],1024,'whole');%計算函數的頻率特性
mag4=abs(H4);%計算函數的幅度特性；
pha4=angle(H4);%計算函數的相位特性；
db4=20*log10((mag4+eps)/max(mag4));%將幅度特性轉換為對數形式
figure('name','=7時凱瑟窗函數頻域波形');
subplot(311);stem(window4);grid on;xlabel('n');ylabel('x(n)');title('β=7的凱瑟窗函數時域圖形');
subplot(312);plot(m4/(2*pi),pha4);grid on;xlabel('w/2pi');ylabel('ψ(w)');title('β=7的凱瑟窗函數的相頻特性');
subplot(313);plot(m4/(2*pi),db4);grid on;xlabel('w/2pi');ylabel('db');title('β=7的凱瑟窗函數的幅度特性');belta5=8;
window5=kaiser(N,belta5);%凱瑟窗函數
[H5,m5]=freqz(window5,[1],1024,'whole');%計算函數的頻率特性
mag5=abs(H5);%計算函數的幅度特性；
pha5=angle(H5);%計算函數的相位特性；
db5=20*log10((mag5+eps)/max(mag5));%將幅度特性轉換為對數形式
figure('name','=8時凱瑟窗函數頻域波形');
subplot(311);stem(window5);grid on;xlabel('n');ylabel('x(n)');title('β=8的凱瑟窗函數時域圖形');
subplot(312);plot(m5/(2*pi),pha5);grid on;xlabel('w/2pi');ylabel('ψ(w)');title('β=8的凱瑟窗函數的相頻特性');
subplot(313);plot(m5/(2*pi),db5);grid on;xlabel('w/2pi');ylabel('db');title('β=8的凱瑟窗函數的幅度特性');
clear all;clc;
N=31;n=0:30;
a=0.42;b=0.4;c=0.18;
x=hanning(N);%漢寧窗函數
for n=0:1:(N-1)/2 %混合窗函數w(n+1)=a-b*cos(2*pi*n/(N-1))+c*sin(pi*n/(N-1));
end
for n=(N+1)/2:1:N-1w(n+1)=0.54-0.46*cos(2*pi*n/(N-1));
end
%繪制漢寧窗頻域波形
figure('name','漢寧窗時域波形');
subplot(211);stem(x);grid on;xlabel('n');ylabel('x(n)');title('N=31的漢寧窗函數時域圖形');
subplot(212);stem(w);grid on;xlabel('n');ylabel('w(n)');title('N=31的混合窗函數時域圖形');
clear all;clc;
N=31;n=0:30;
a=0.42;b=0.4;c=0.18;
x=hanning(N);%漢寧窗函數
for n=0:1:(N-1)/2 %混合窗函數w(n+1)=a-b*cos(2*pi*n/(N-1))+c*sin(pi*n/(N-1));
end
for n=(N+1)/2:1:N-1w(n+1)=0.54-0.46*cos(2*pi*n/(N-1));
end
[H1,m1]=freqz(x,[1],1024,'whole');%計算函數的頻率特性
mag1=abs(H1);%計算函數的幅度特性；
pha1=angle(H1);%計算函數的相位特性；
db1=20*log10((mag1+eps)/max(mag1));%將幅度特性轉換為對數形式
%繪制漢寧窗頻域波形
figure('name','漢寧窗頻域波形');
subplot(221);plot(m1/(2*pi),pha1);grid on;xlabel('w/2pi');ylabel('ψ(w)');title('N=31的漢寧窗函數的相頻特性');
subplot(223);plot(m1/(2*pi),db1);grid on;xlabel('w/2pi');ylabel('db');title('N=31的漢寧窗函數的幅度特性');
[H2,m2]=freqz(w,[1],1024,'whole');%計算函數的頻率特性
mag2=abs(H2);%計算函數的幅度特性；
pha2=angle(H2);%計算函數的相位特性；
db2=20*log10((mag2+eps)/max(mag2));%將幅度特性轉換為對數形式
subplot(222);plot(m2/(2*pi),pha2);grid on;xlabel('w/2pi');ylabel('ψ(w)');title('N=31的混合窗函數的相頻特性');
subplot(224);plot(m2/(2*pi),db2);grid on;xlabel('w/2pi');ylabel('db');title('N=31的混合窗函數的幅度特性');
四、分析
分析略
五、思考題
(1)漢寧窗、海明窗和布萊克曼窗是如何做到比矩形窗更高的副瓣衰減的？
(2)窗函數序數序列的長度N是否會影響窗函數的副瓣衰減?試從理論上加以證明。
(3)通過實驗分析并說明凱瑟窗序列中參數β參數與窗函數頻率特性之間的關系。
(4）在對信號進行頻域分析前通常需要進行加窗處理。截取要求的點數,即相關于加矩窗。也可將被處理信號序列與其他各種類型的窗函數進行相乘后再進行頻域分析。通過理論分析或者仿真實例說明加窗處理對譜分析的頻譜泄漏效應和分辨力特性有何影響?
(5)常用的矩形窗、漢寧窗、海明窗、布萊克曼窗和凱瑟窗的時域圖形是否是中心對稱的,這幾個窗函數的相頻特性是否為線性相位的?混合窗函數的時域圖形是否中心對稱的，其相頻特性是否為線性相位的？
六、總結
通過“窗函數法設計FIR數字濾波器”的實驗，我對于濾波器的設計又了更加深入的理解。特別是對于濾波器類型的選擇，還有根據所給出的條件進行濾波器的確定還有濾波器窗口長度的確定。對于不同種類的窗口，如：矩形窗、三角窗、漢寧窗、海明窗等都更加熟悉，也知道了它們的表達式和性能特點，以及它們的函數圖形。我也通過實驗，知道了已知單位沖激響應，要在Matlab中求其幅頻特性，并進行作圖的方法。對于FIR濾波器的知識也有了更扎實的掌握。對我理解窗函數設計濾波器方法很有幫助，我相信，這次數字信號處理實驗的參與僅僅是打開了我利用matlab信號分析與處理的大門，在今后我一定會不斷努力加油讓matlab真正成為我的一項技能。
總結
以上是生活随笔為你收集整理的窗函数性能分析——MATLAB的全部內容，希望文章能夠幫你解決所遇到的問題。
如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。
