怎样在MATLAB中实现这几种常见的窗函数

uperrua

怎样在MATLAB中实现这几种常见的窗函数


数字信号处理中通常是取其有限的时间片段进行分析，而不是对无限长的信号进行测量和运算。具体做法是从信号中截取一个时间片段，然后对信号进行傅里叶变换、相关分析等数学处理。信号的截断产生了能量泄漏，而用FFT算法计算频谱又产生了栅栏效应，从原理上讲这两种误差都是不能消除的。在FFT分析中为了减少或消除频谱能量泄漏及栅栏效应，可采用不同的截取函数对信号进行截短，截短函数称为窗函数，简称为窗。
, v0 h: p" y. Z7 O# `# T  h泄漏与窗函数频谱的两侧旁瓣有关，对于窗函数的选用总的原则是，要从保持最大信息和消除旁瓣的综合效果出发来考虑问题，尽可能使窗函数频谱中的主瓣宽度应尽量窄，以获得较陡的过渡带；旁瓣衰减应尽量大，以提高阻带的衰减，但通常都不能同时满足这两个要求。频谱中的如果两侧瓣的高度趋于零，而使能量相对集中在主瓣，就可以较为接近于真实的频谱。不同的窗函数对信号频谱的影响是不一样的，这主要是因为不同的窗函数，产生泄漏的大小不一样，频率分辨能力也不一样。信号的加窗处理，重要的问题是在于根据信号的性质和研究目的来选用窗函数。图1是几种常用的窗函数的时域和频域波形，其中矩形窗主瓣窄，旁瓣大，频率识别精度最高，幅值识别精度最低，如果仅要求精确读出主瓣频率，而不考虑幅值精度，则可选用矩形窗，例如测量物体的自振频率等；布莱克曼窗主瓣宽，旁瓣小，频率识别精度最低，但幅值识别精度最高；如果分析窄带信号，且有较强的干扰噪声，则应选用旁瓣幅度小的窗函数，如汉宁窗、三角窗等；对于随时间按指数衰减的函数，可采用指数窗来提高信噪比。表1 是几种常用的窗函数的比较。
0 P% @! E9 x7 `$ a& i) B
如果被测信号是随机或者未知的，或者是一般使用者对窗函数不大了解，要求也不是特别高时，可以选择汉宁窗，因为它的泄漏、波动都较小，并且选择性也较高。但在用于校准时选用平顶窗较好，因为它的通带波动非常小，幅度误差也较小。

1 Y9 c* V$ ?: P4 ~3 J3 C

表1 几种常用的窗函数的比较

名称	特点	应用
矩形窗 Rectangle	矩形窗使用最多，习惯上不加窗就是使信号通过了矩形窗。这种窗的优点是主瓣比较集中，缺点是旁瓣较高，并有负旁瓣，导致变换中带进了高频干扰和泄漏，甚至出现负谱现象。频率识别精度最高，幅值识别精度最低，所以矩形窗不是一个理想的窗。	如果仅要求精确读出主瓣频率，而不考虑幅值精度，则可选用矩形窗，例如测量物体的自振频率等，也可以用在阶次分析中。
汉宁窗 Hanning	又称升余弦窗。主瓣加宽并降低，旁瓣则显著减小，从减小泄漏观点出发，汉宁窗优于矩形窗．但汉宁窗主瓣加宽，相当于分析带宽加宽，频率分辨力下降。它与矩形窗相比，泄漏、波动都减小了,并且选择性也提高。	是很有用的窗函数。如果测试信号有多个频率分量，频谱表现的十分复杂，且测试的目的更多关注频率点而非能量的大小，需要选择汉宁窗。如果被测信号是随机或者未知的，选择汉宁窗。
海明窗 （汉明窗） Hamming	与汉宁窗都是余弦窗，又称改进的升余弦窗，只是加权系数不同，使旁瓣达到更小。但其旁瓣衰减速度比汉宁窗衰减速度慢。	与汉明窗类似，也是很有用的窗函数。
平顶窗 Flap Top	平顶窗在频域时的表现就象它的名称一样有非常小的通带波动。	由于在幅度上有较小的误差，所以这个窗可以用在校准上。
凯塞窗 Kaiser	定义了一组可调的由零阶贝塞尔Bessel 函数构成的窗函数，通过调整参数β可以在主瓣宽度和旁瓣衰减之间自由选择它们的比重。对于某一长度的Kaiser 窗，给定β，则旁瓣高度也就固定了。
布莱克曼窗 Blackman	二阶升余弦窗，主瓣宽，旁瓣比较低，但等效噪声带宽比汉宁窗要大一点，波动却小一点。频率识别精度最低，但幅值识别精度最高，有更好的选择性。	常用来检测两个频率相近幅度不同的信号。
高斯窗 Gaussian	是一种指数窗。主瓣较宽，故而频率分辨力低；无负的旁瓣，第一旁瓣衰减达一55dB。常被用来截短一些非周期信号，如指数衰减信号等。	对于随时间按指数衰减的函数，可采用指数窗来提高信噪比。
三角窗 （费杰窗） Fejer	是幂窗的一次方形式。与矩形窗比较，主瓣宽约等于矩形窗的两倍，但旁瓣小，而且无负旁瓣。	如果分析窄带信号，且有较强的干扰噪声，则应选用旁瓣幅度小的窗函数，如汉宁窗、三角窗等；
切比雪夫窗（Chebyshev）	在给定旁瓣高度下，Chebyshev窗的主瓣宽度最小，具有等波动性，也就是说，其所有的旁瓣都具有相等的高度。

. `- T, M' r2 b, f& K# x
$ P* ^' T+ T7 s* _" C8 h
2 }7 v" c# \$ f下面是几种窗函数归一化DTFT幅度的MATLAB程序:

附上DTFT函数（dtft.m）：
! v! d( V4 ]& T9 `+ d. r
  ^- ?- L& [% o5 C0 r7 v$ ~function [ X ] = dtft( x,n,w )
7 b  A4 X& h4 m$ ~- S
& ]& g2 l8 e( \" M% d%   Computes Discrete-time Fourier Transform

%   [X] = dtft(x,n,w)
3 d9 |9 ?3 `! g, K
0 u- \' A  m: L9 ]+ U$ I%    X = DTFT values computed at w.frequencies
. g& [! O; i& E' U
%    x = finite duration sequence over n
& {  l# c$ I0 z3 z# A' P! ^
2 ]) y' K3 c1 K. |( U5 j%    n = sample position vector

/ p0 t8 {* X2 Q* C%    w = frequency location vector
+ ~' [8 b4 x" z# `6 a1 M
X = x*exp(-j*n'*w);
2 r: a  i) P: n0 ^; `" S
%
* v) k4 G( g, H" {- |; y- E) k
/ I! O7 V& ^6 F) z7 ]end

* S% b/ t3 ^: u' T8 D* \5 _

矩形窗:
" W3 I5 `4 l+ ~9 S/ t8 T1 A8 Z6 m$ G
%DTFT of a Rectangular Window, M=10,25,50,101
3 f1 P, C% M/ r- S5 m- k2 ]* ?
clc; close all;
% U9 ~$ ]* k3 m% }1 T1 q, V
* J5 b! |* c/ J8 ^9 eHf_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));
; o& o; n9 W, i, N/ Y- Z6 n% O
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]);

! \+ `* a5 w2 I4 J$ {ylabel('|X|'); title(['M=10']);
( t, T. _+ H, D( h2 d
set(gca,'XTick',wtick,'YTick',magtick);
8 Y, G& C: ?$ f3 I9 J, E
1 {, X9 @6 o1 p$ A3 E7 ~* a, q/ p" W% M=25
! T; g2 I* [: D# B' c, y+ p
& ~( P& d3 Y8 L9 nM=25; n=0:M; x=ones(1,length(n));

X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
( e( Q3 `7 E9 q6 C2 G# b( P
subplot(2,2,2); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

" p+ G5 c; c6 V+ l5 E2 Jtitle(['M=25']); set(gca,'XTick',wtick,'YTick',magtick);
. r' k  \9 b4 u# O$ V7 P. d
% M=50

M=50; n=0:M; x=ones(1,length(n));

X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
6 e3 l3 p- y. K- g8 o
subplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
" @9 d5 j0 U4 p/ S  W
) b, L9 U: j1 Y. L$ g$ hxlabel('\omega/\pi'); ylabel('|X|');

title('M=50'); set(gca,'XTick',wtick,'YTick',magtick);
: O) Q# {/ U& H+ N3 J+ R
% M=101

0 O8 F* }  \- P7 z0 g* }# ?2 FM=101; n=0:M; x=ones(1,length(n));

X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
1 k  I6 I# V2 O( y/ u) H
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']);

/ d& S( q- x1 D
* j) b" E- f. f% ]: S6 n. S# E三角窗:
, q! j/ ^& F9 S5 O
% Triangular Window:

6 O( ]9 V8 U, q& ?4 d5 B0 c4 R" z% DTFT of a Triangular Window,M = 10,25,50,101
" F& T3 |& a, M: c! e
clc; close all;

1 D* w( m- _: e( |9 ~& `" BHf_1=figure; set(Hf_1,'NumberTitle','off','Name','P0304b');

: z5 ~9 {$ h9 Aw=linspace(-pi,pi,501); wtick=[-1:0.5:1]; magtick=[0:0.5:1.1];

1 c7 s3 _) i& v, ^7 g* B% M = 10

0 L8 Q- u4 U2 ^2 q* \8 n# ~; {- ?M=10; n=0:M; x=(1-(abs(M-1-(2*n))/(M+1)));

  L6 \& m" `5 K6 g% bx=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]);
: L& r2 r1 w# k8 m
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)));
; b7 Z0 e: m5 f8 `9 ^/ m) a
X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
9 F! L8 h0 d4 k  X9 O% |! T
" i  R* u) c  l6 e) usubplot(2,2,2); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
# X: g$ a  q( W$ q$ b
title(['M = 25']); set(gca,'XTick',wtick,'YTick',magtick);
% U; D+ Y4 B5 k9 p) L. I
' J# J0 u) b# {0 {% M = 50
) r! w, F5 b2 D: }! P' D1 ]: h! c5 K
; n4 C  `# D: R( E) r9 W; W5 q5 uM=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);

; [0 u2 ^+ I) i" psubplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
/ r1 q. I7 L8 Y- ?' e
- I: A9 p( E5 D$ a( `xlabel('\omega/\pi'); ylabel('|X|'); title(['M = 50']);
9 o; L! a( n8 f8 p( o
9 q7 s/ ]: F, qset(gca,'XTick',wtick,'YTick',magtick);
" K6 t- Y+ Q- e" X! q' G4 ]# A
6 F! d4 d2 @7 Q% M = 100

M=101;n=0:M; x=(1-(abs(M-1-(2*n))/(M+1)));
" O; S" e/ W& L, T& ^+ }( p
X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);

. H5 m9 r( u- G! a. Tsubplot(2,2,4); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
( V- \, `$ Z, U4 p
xlabel('\omega/\pi'); title(['M = 101']);

set(gca,'XTick',wtick,'YTick',magtick);

- S  A; P, f2 R( H! Y6 w6 ]

# k# o# R' n! R; O" g海宁窗:

$ R3 D' S- d5 `3 W& J0 N8 Y' a% Hann Window
* t" J5 ~+ _6 ^9 `; G  X0 o
5 p- Z/ q" l" ], j7 G# K4 A% DTFT of a Hann Window, M = 10,25,50,101

clc;close all;

% O( b* T5 l- sHf_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];
8 q* }: D9 ]6 I8 b2 i
% M = 10

M=10;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1)));
5 D: S" o3 g6 ^/ M8 y. i
X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
" `- _7 M! g: V1 f. _* b* X
subplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
: _3 B/ s3 b8 u; M6 ~
  n% W( n. m- b# K. z6 X# y& Y# yylabel('|X|'); title(['M = 10']);
5 {3 P8 x5 d  A
set(gca,'XTick',wtick,'YTick',magtick);

2 ?3 P0 S9 M  W; [6 T! _% M = 25
" ], R+ J# i, k9 }& l
M=25;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1)));
/ V  d- g& z! W
5 @. p/ v8 v& c$ V; gX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);

8 E8 B7 l$ \$ X' b" M! S* q. h. isubplot(2,2,2); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
% i" g) n$ h0 B$ b
title(['M = 25']);

set(gca,'XTick',wtick,'YTick',magtick);
" l* N) c( y5 n$ M8 s
% M = 50

M=50;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1)));
4 p$ A( ?# C" ?  U0 |
X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
# O# I! |3 u1 J; L* p& R
; f% a) ?5 @$ e5 H4 Ysubplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
* G+ r7 o7 Y; i+ n* S- k# D/ g6 y3 M" T
xlabel('\omega/\pi'); ylabel('|X|'); title(['M = 50']);

set(gca,'XTick',wtick,'YTick',magtick);
! R4 i0 T, ~* ?* z& m" e2 q# e
, o( ?. W' b* j% M = 101

" B- h2 D) Z: wM=101;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1)));
: f+ F8 R- V- n: G! w7 _7 K
  ?  G8 N# x: y) t3 cX=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]);
: g  r/ h3 l; {5 K0 c" U: H
xlabel('\omega/\pi'); title(['M = 101']);
( c2 f& C2 w9 f. @, t3 D) s& c
! X. N- o) Q# w& E: [set(gca,'XTick',wtick,'YTick',magtick);


( I0 N) ?# {; w) h' w哈明窗:
3 ]# o( @9 m( ~( W& H
6 q5 H. R% X" P+ e! L) k( Y% Hamming Window:
- }) g6 a; T$ W0 H
( K% |" x$ m6 i' V5 @clc; close all;
  \- y4 S, w# c! U
Hf_1=figure; set(Hf_1,'NumberTitle','off','Name','P0304d');
" g1 j! H8 I( y- @5 w2 v# N
w=linspace(-pi,pi,501); wtick=[-1:0.5:1]; magtick=[0:0.5:1.1];

% D- ]( a- s+ S% M = 10
0 J/ m! {/ Q* H0 x
M=10; n=0:M; x=(0.54-0.46*cos((2*pi*n)/(M-1)));

8 L2 j7 ]2 q3 ?% ~8 YX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
6 ~, b* e1 y8 d: V! h, A! P4 c6 Q
6 _! }" ?: z1 {0 K. F$ ysubplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
9 B4 b0 x) O$ e  Y, x# q
ylabel('|X|'); title(['M = 10']);

set(gca,'XTick',wtick,'YTick',magtick);
6 a- q' O" l6 u% p& @4 o
% M = 25
' \; O6 V: _- T
M=25; n=0:M; x=(0.54-0.46*cos((2*pi*n)/(M-1)));
- a! u. ]" L8 h
X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
: w' J/ I& Y: X; d8 q9 U
subplot(2,2,2); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
4 b- |# U2 N" H, M. O
3 q5 P/ Z* d- f4 S3 j! L% N  A$ |. Btitle(['M = 25']);

set(gca,'XTick',wtick,'YTick',magtick);
' f* g  f$ }2 E
% 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);

! p7 @, ?6 `+ _' M( xsubplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
' c4 G4 \5 x  J+ y7 C( K$ c' X9 f
2 s. D' B* N3 ]2 Y$ |, Hxlabel('\omega/\pi'); ylabel('|X|'); title(['M=50']);
( Y+ |; u6 v  G; f& `
* e& c: Z( \5 S! M2 _% gset(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);

! [- f+ @6 C1 }subplot(2,2,4); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

" N$ r1 m, H5 @* \2 I/ a  Rxlabel('\omega/\pi');title(['M=101']);

$ v- U' N; S8 S8 xset(gca,'XTick',wtick,'YTick',magtick);
) H! j' }* K1 I& p" M
# C8 D- R6 Q) a' I+ t) \

relchhiclty

谢谢分享