数字信号处理中通常是取其有限的时间片段进行分析,而不是对无限长的信号进行测量和运算。具体做法是从信号中截取一个时间片段,然后对信号进行傅里叶变换、相关分析等数学处理。信号的截断产生了能量泄漏,而用FFT算法计算频谱又产生了栅栏效应,从原理上讲这两种误差都是不能消除的。在FFT分析中为了减少或消除频谱能量泄漏及栅栏效应,可采用不同的截取函数对信号进行截短,截短函数称为窗函数,简称为窗。
泄漏与窗函数频谱的两侧旁瓣有关,对于窗函数的选用总的原则是,要从保持最大信息和消除旁瓣的综合效果出发来考虑问题,尽可能使窗函数频谱中的主瓣宽度应尽量窄,以获得较陡的过渡带;旁瓣衰减应尽量大,以提高阻带的衰减,但通常都不能同时满足这两个要求。频谱中的如果两侧瓣的高度趋于零,而使能量相对集中在主瓣,就可以较为接近于真实的频谱。不同的窗函数对信号频谱的影响是不一样的,这主要是因为不同的窗函数,产生泄漏的大小不一样,频率分辨能力也不一样。信号的加窗处理,重要的问题是在于根据信号的性质和研究目的来选用窗函数。图1是几种常用的窗函数的时域和频域波形,其中矩形窗主瓣窄,旁瓣大,频率识别精度最高,幅值识别精度最低,如果仅要求精确读出主瓣频率,而不考虑幅值精度,则可选用矩形窗,例如测量物体的自振频率等;布莱克曼窗主瓣宽,旁瓣小,频率识别精度最低,但幅值识别精度最高;如果分析窄带信号,且有较强的干扰噪声,则应选用旁瓣幅度小的窗函数,如汉宁窗、三角窗等;对于随时间按指数衰减的函数,可采用指数窗来提高信噪比。表1 是几种常用的窗函数的比较。 如果被测信号是随机或者未知的,或者是一般使用者对窗函数不大了解,要求也不是特别高时,可以选择汉宁窗,因为它的泄漏、波动都较小,并且选择性也较高。但在用于校准时选用平顶窗较好,因为它的通带波动非常小,幅度误差也较小。
表1 几种常用的窗函数的比较
下面是几种窗函数归一化DTFT幅度的MATLAB程序: 附上DTFT函数(dtft.m): function [ X ] = dtft( x,n,w ) % Computes Discrete-time Fourier Transform % [X] = dtft(x,n,w) % X = DTFT values computed at w.frequencies % x = finite duration sequence over n % n = sample position vector % w = frequency location vector X = x*exp(-j*n'*w); % end
矩形窗: %DTFT of a Rectangular Window, M=10,25,50,101 clc; close all; Hf_1=figure; set(Hf_1,'NumberTitle','off','Name','P0304a'); w=linspace(-pi,pi,501); wtick=[-1:0.5:1]; magtick=[0:0.5:1.1]; % M=10 M=10; n=0:M; x=ones(1,length(n)); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); ylabel('|X|'); title(['M=10']); set(gca,'XTick',wtick,'YTick',magtick); % M=25 M=25; n=0:M; x=ones(1,length(n)); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,2); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); title(['M=25']); set(gca,'XTick',wtick,'YTick',magtick); % M=50 M=50; n=0:M; x=ones(1,length(n)); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); xlabel('\omega/\pi'); ylabel('|X|'); title('M=50'); set(gca,'XTick',wtick,'YTick',magtick); % M=101 M=101; n=0:M; x=ones(1,length(n)); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,4); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); xlabel('\omega/\pi'); ylabel('|X|'); title(['M=101']); ![]()
三角窗: % Triangular Window: % DTFT of a Triangular Window,M = 10,25,50,101 clc; close all; Hf_1=figure; set(Hf_1,'NumberTitle','off','Name','P0304b'); w=linspace(-pi,pi,501); wtick=[-1:0.5:1]; magtick=[0:0.5:1.1]; % M = 10 M=10; n=0:M; x=(1-(abs(M-1-(2*n))/(M+1))); x=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); ylabel('|X|'); title(['M = 10']); set(gca,'XTick',wtick,'YTick',magtick); % M = 25 M=25; n=0:M; x=(1-(abs(M-1-(2*n))/(M+1))); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,2); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); title(['M = 25']); set(gca,'XTick',wtick,'YTick',magtick); % M = 50 M=50; n=0:M; x=(1-(abs(M-1-(2*n))/(M+1))); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); xlabel('\omega/\pi'); ylabel('|X|'); title(['M = 50']); set(gca,'XTick',wtick,'YTick',magtick); % M = 100 M=101;n=0:M; x=(1-(abs(M-1-(2*n))/(M+1))); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,4); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); xlabel('\omega/\pi'); title(['M = 101']); set(gca,'XTick',wtick,'YTick',magtick); ![]()
海宁窗: % Hann Window % DTFT of a Hann Window, M = 10,25,50,101 clc;close all; Hf_1 = figure; set(Hf_1,'NumberTitle','off','Name','P0304c'); w=linspace(-pi,pi,501); wtick=[-1:0.5:1]; magtick=[0:0.5:1.1]; % M = 10 M=10;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1))); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); ylabel('|X|'); title(['M = 10']); set(gca,'XTick',wtick,'YTick',magtick); % M = 25 M=25;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1))); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,2); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); title(['M = 25']); set(gca,'XTick',wtick,'YTick',magtick); % M = 50 M=50;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1))); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); xlabel('\omega/\pi'); ylabel('|X|'); title(['M = 50']); set(gca,'XTick',wtick,'YTick',magtick); % M = 101 M=101;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1))); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,4); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); xlabel('\omega/\pi'); title(['M = 101']); set(gca,'XTick',wtick,'YTick',magtick); ![]()
哈明窗: % Hamming Window: clc; close all; Hf_1=figure; set(Hf_1,'NumberTitle','off','Name','P0304d'); w=linspace(-pi,pi,501); wtick=[-1:0.5:1]; magtick=[0:0.5:1.1]; % M = 10 M=10; n=0:M; x=(0.54-0.46*cos((2*pi*n)/(M-1))); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); ylabel('|X|'); title(['M = 10']); set(gca,'XTick',wtick,'YTick',magtick); % M = 25 M=25; n=0:M; x=(0.54-0.46*cos((2*pi*n)/(M-1))); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,2); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); title(['M = 25']); set(gca,'XTick',wtick,'YTick',magtick); % M = 50 M=50; n=0:M; x=(0.54-0.46*cos((2*pi*n)/(M-1))); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); xlabel('\omega/\pi'); ylabel('|X|'); title(['M=50']); set(gca,'XTick',wtick,'YTick',magtick); % M = 101 M=101; n=0:M; x=(0.54-0.46*cos((2*pi*n)/(M-1))); X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX); subplot(2,2,4); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]); xlabel('\omega/\pi');title(['M=101']); set(gca,'XTick',wtick,'YTick',magtick); ![]() |
|