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

uperrua

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

) {& y! i8 a  ]+ O

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

4 `% S* `. _2 p% e! N& k: Z0 E

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

6 _3 |9 P3 B3 ^
) _2 Z! y- C- G! C; \" d
( F) v# [: k0 }) @) I; R下面是几种窗函数归一化DTFT幅度的MATLAB程序:

附上DTFT函数（dtft.m）：

function [ X ] = dtft( x,n,w )

+ L- O. @7 ~/ C$ S9 y%   Computes Discrete-time Fourier Transform

) p, p* Q- G' M%   [X] = dtft(x,n,w)
4 G2 v9 z0 w  y( D9 H. D
%    X = DTFT values computed at w.frequencies
/ u! w4 |" N8 Y- `
6 B( K7 B$ g/ A9 _& H' V" ~+ c1 Q4 j%    x = finite duration sequence over n
+ I9 r- x% B9 ]  p9 S$ {% Q
%    n = sample position vector

& X3 m8 g& ~3 s! y9 |. N0 r%    w = frequency location vector

X = x*exp(-j*n'*w);
: j  w' W& I4 y, F- G' v, O
%
' t& v! k$ A2 p6 u; S5 F
% ^9 I8 z: P4 ~+ G3 p( J% N, Xend

# B  G' W/ U) k  F. U2 I. ^1 L( j

矩形窗:

%DTFT of a Rectangular Window, M=10,25,50,101
& v# v! ?8 d  Q: O4 d* j; `# n
7 o+ T  k9 B8 x2 u2 r  y# \" Y! Zclc; close all;

. O% U1 X! i! M4 c0 Z. z; ?& kHf_1=figure; set(Hf_1,'NumberTitle','off','Name','P0304a');

9 x# u2 x. @' Yw=linspace(-pi,pi,501); wtick=[-1:0.5:1]; magtick=[0:0.5:1.1];

2 a3 W3 }# D1 x% u% M=10

1 e) H% P2 H2 q/ W/ pM=10; n=0:M; x=ones(1,length(n));
+ T! P- ?9 i# z: p) O
2 Q. `) h: a3 u3 K  y4 t+ NX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
4 _) y9 Y% x: f; I3 b8 J
subplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

2 u: x  U* ?2 ]; r- Mylabel('|X|'); title(['M=10']);

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

7 H/ ^4 ^0 `7 V/ q7 M0 ~% M=25
4 o0 R: d# A" O
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]);

" Y. n6 g7 R4 R) U+ j% C% ytitle(['M=25']); set(gca,'XTick',wtick,'YTick',magtick);

$ k7 k, [, g  j! o+ _% M=50
+ a5 @$ F  ~; H
3 k& _9 ^6 b! fM=50; n=0:M; x=ones(1,length(n));
1 V" N8 k, s$ N0 }
& y2 D# p4 G$ H3 Y9 N( J% i, dX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
" C% i$ m& C# }' ^
subplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
0 }, s% W2 P, m9 D' p+ f
: j2 r8 S' K  G, O: I( dxlabel('\omega/\pi'); ylabel('|X|');
) a+ T" {2 x# [! Q1 [
6 D- \( K2 X$ e2 ttitle('M=50'); set(gca,'XTick',wtick,'YTick',magtick);
( l! R$ X4 G& E, A/ B( j& r; C: X$ O
# K! {; d5 H* u6 }9 r% M=101

M=101; n=0:M; x=ones(1,length(n));
' r$ u5 I& a: \6 B
X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
3 O9 ?5 ]# @* {) s( o# S
subplot(2,2,4); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

5 B2 J/ _! v* H4 wxlabel('\omega/\pi'); ylabel('|X|');
' d1 g# N3 Y9 X7 z8 p5 R: \1 M8 a. M
title(['M=101']);
& Q$ P+ p) o& o# u. X

三角窗:

% Triangular Window:
2 s8 |% {' \) c5 K
% DTFT of a Triangular Window,M = 10,25,50,101
: S& k( `# J' f
clc; close all;

Hf_1=figure; set(Hf_1,'NumberTitle','off','Name','P0304b');

$ f% W, k- \7 v. V" k$ Uw=linspace(-pi,pi,501); wtick=[-1:0.5:1]; magtick=[0:0.5:1.1];

- W' l  n0 d% K$ n! K; X5 `* |% M = 10

3 f- V  {1 o% o3 s& `0 IM=10; n=0:M; x=(1-(abs(M-1-(2*n))/(M+1)));
* `- H7 @* R/ N8 A
x=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
3 n! U  c7 D& ?4 a1 y- d
# I6 h8 s, {2 P# E5 f# q/ Dsubplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

ylabel('|X|'); title(['M = 10']);
5 ?5 Q4 c/ N( y$ Y9 D% r- n
" \: o" \& E: R+ Fset(gca,'XTick',wtick,'YTick',magtick);
* {# L, q1 w: B* ^& A0 [7 L2 Z
% M = 25
' Z* Y- f' |2 a9 h7 Y8 }; R* j
M=25; n=0:M; x=(1-(abs(M-1-(2*n))/(M+1)));
. X, d* \; G: G' e1 f* J
. F  v2 ~8 B( M2 N, k* pX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
. f4 S8 Z0 ?( b( H2 L" T  f
" e: l: I' D1 y" rsubplot(2,2,2); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

# @0 p# t/ @% u7 J- y5 u! stitle(['M = 25']); set(gca,'XTick',wtick,'YTick',magtick);

% M = 50

1 X. I2 h8 ^5 V# u' ]. D- e  _M=50; n=0:M; x=(1-(abs(M-1-(2*n))/(M+1)));
' I1 g+ l$ e% v" Y
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]);

* W; V4 p. V8 O1 X/ ~' s- L, Hxlabel('\omega/\pi'); ylabel('|X|'); title(['M = 50']);

, }! V- p/ j) e8 q: Zset(gca,'XTick',wtick,'YTick',magtick);

% M = 100

9 S; ~# f2 d9 Q3 X- L( EM=101;n=0:M; x=(1-(abs(M-1-(2*n))/(M+1)));

7 D3 K8 U! w; B# |& B. BX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
2 o( B1 r( T1 f; h  ~# B
: e- ~- [# B, e) i4 E% X# csubplot(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);
: k4 [0 z! y4 f' J


4 Y6 x2 Y8 R  @& o" ^: T海宁窗:
+ U# f* U& {$ P, P& z7 p% W
7 ~1 D" G' C2 Z: A6 [1 E& ^; Z) i9 v% Hann Window
, d2 o5 j# ]8 D8 W3 M: x2 {
) Z+ b6 n7 s) \/ `% DTFT of a Hann Window, M = 10,25,50,101

, ^1 ]* N; A( R- g" J+ }clc;close all;
2 T( e9 G! f% e0 X9 Q
+ k& V+ V& V1 RHf_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
. \5 g1 @2 p! ~" r
M=10;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1)));
: h' f& p4 t2 ~$ a: R
X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);

% \& \, Y3 }6 s0 b8 Msubplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
4 a- d& d$ e' o5 c
2 r5 s9 _8 [' P: Eylabel('|X|'); title(['M = 10']);

set(gca,'XTick',wtick,'YTick',magtick);
2 V& _$ [+ b. N7 h
% M = 25
1 F) E3 X" w2 A( U0 L* X; \, V9 G
, j. O3 \3 N' B/ G9 yM=25;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1)));
( {: q; _: i, {
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]);

, C) a$ j) W7 z% C( m( f: d title(['M = 25']);
! V. E6 {# {/ Y! d/ ^! r
% c* O# g7 U2 }; v, z  m6 dset(gca,'XTick',wtick,'YTick',magtick);

. Y2 P2 W2 i6 @) U7 i  Y% M = 50
$ E$ W2 i" Z* a) @; J* _- C# g
M=50;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1)));
8 k! a: f6 R; K+ J- V5 R- o7 e/ J* P
X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
: x, K5 a1 q& i  ?
) k; g% M0 l! S+ a) ]9 Osubplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

5 u( x( v. {  P1 h9 ~' lxlabel('\omega/\pi'); ylabel('|X|'); title(['M = 50']);

" L8 x$ m# K" yset(gca,'XTick',wtick,'YTick',magtick);
* R- a$ z9 c# ~0 x* Y! D0 S
4 @" Y6 T" k) O; y0 k: \' B% M = 101
0 `- p3 s6 O5 p& ]0 Q
5 K+ }! z" X/ W1 cM=101;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1)));
; V# E; g* S8 j. ?3 S( K
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]);
  K1 c' N0 d& u$ g! V
+ I/ z8 ^" s( J+ x  J5 S6 Pxlabel('\omega/\pi'); title(['M = 101']);

/ U0 @3 |1 q' a( Zset(gca,'XTick',wtick,'YTick',magtick);

  Z# @% j, C1 H! `0 {, ?2 V, r
" |0 I$ O. V& L! U! t3 b2 F* r哈明窗:
9 x6 P( _6 y( R; u
% Hamming Window:
& \7 x2 u! L* d4 `4 K
9 v8 T- i9 e3 ]" I$ C- hclc; 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];
+ ^6 x& S: |/ _+ Q7 }- z
) G9 h+ Q( w& t: t  i8 v% M = 10
' Q. w! S/ O) A! M% Q
M=10; n=0:M; x=(0.54-0.46*cos((2*pi*n)/(M-1)));

% m- j9 ?& R) |* Q7 o: p& [X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
& k" \/ m3 m2 |. S
subplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
) [% t! y& h7 _) y
8 M  }! f7 L( A+ Z  {9 R! o ylabel('|X|'); title(['M = 10']);

set(gca,'XTick',wtick,'YTick',magtick);
% {8 H0 h- d4 Z# B; k) }: x  R8 ?- {
% M = 25
3 v$ b! d* m& c6 C7 D4 R# R  [3 J
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);

) z  A# ^* ?8 P4 X5 [+ q. tsubplot(2,2,2); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

title(['M = 25']);
( [8 O5 _8 E+ E. ~
# I/ ^& k6 x4 Iset(gca,'XTick',wtick,'YTick',magtick);
) ]' g1 ]9 L6 ]1 C3 u
% M = 50
! Q% {* T' J) `8 z: |' }, e
. l; ^" w7 }, V9 H0 }: nM=50; n=0:M; x=(0.54-0.46*cos((2*pi*n)/(M-1)));

, x- i* |4 r" b, \* X" xX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
0 i+ O3 E+ m: P# d! u
$ l( M" P' A, t" lsubplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
$ O9 P* V# ~7 l$ c; l9 H6 A
4 t! M$ O0 g8 y/ N. |& yxlabel('\omega/\pi'); ylabel('|X|'); title(['M=50']);

set(gca,'XTick',wtick,'YTick',magtick);
' F  l# q9 ^0 |* C
% M = 101
' e( ]+ j: d5 D! R4 f
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);
  G/ W2 d* A2 F
9 I' V* j0 I3 N, Ssubplot(2,2,4); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

xlabel('\omega/\pi');title(['M=101']);

. m) v5 L% G" v6 O+ G; I- Sset(gca,'XTick',wtick,'YTick',magtick);

9 w" l4 {, @7 }# {7 h( N

relchhiclty

谢谢分享