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

uperrua

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


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

如果被测信号是随机或者未知的，或者是一般使用者对窗函数不大了解，要求也不是特别高时，可以选择汉宁窗，因为它的泄漏、波动都较小，并且选择性也较高。但在用于校准时选用平顶窗较好，因为它的通带波动非常小，幅度误差也较小。
( G9 u' R- \3 T

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

$ k0 W$ q2 D% ~8 }- h. N3 S

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

6 L1 w5 i0 a% A9 h8 r) L; @
1 x$ W' D# _1 W4 i6 ^5 k
& V6 A9 u9 U: d& O下面是几种窗函数归一化DTFT幅度的MATLAB程序:
7 l0 C' _/ O" L4 t
附上DTFT函数（dtft.m）：
: i: D3 _9 H* V( d( o' g
$ f! s" Q' k2 [function [ X ] = dtft( x,n,w )

%   Computes Discrete-time Fourier Transform
5 }8 ?3 H3 g! q& M
%   [X] = dtft(x,n,w)
5 i1 e4 E. m: D7 O, R6 P
%    X = DTFT values computed at w.frequencies

# F- o) n8 |2 T" Z. z4 Y% L8 A4 y%    x = finite duration sequence over n
% H- I& q, b5 ?- _4 n
' u+ D: S, @' g8 E%    n = sample position vector
$ R* p& Y% W& X# {1 d) _
%    w = frequency location vector

X = x*exp(-j*n'*w);

8 j$ K& K9 |, S  C! W  |6 H+ I%
- z0 ]/ s% H  h& q! ^
end


' Y  R  g1 P9 O
1 F6 H+ o8 @( J  J  Z. H矩形窗:

%DTFT of a Rectangular Window, M=10,25,50,101

( B1 x: U) O" O6 }) ]' `; Vclc; close all;

) s1 C; Q1 |" {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];
- s4 A- F' [1 r! T$ c3 B" h
- H; \5 N  ~: e3 U% M=10

: t( P- I. T7 V3 \8 p, D8 m: JM=10; n=0:M; x=ones(1,length(n));

X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
+ ^2 j7 k3 f- G
subplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
* f% N& ]) i1 v* }
ylabel('|X|'); title(['M=10']);
8 D2 `# A/ M8 D& a' r! A
' X* m, d( J7 Hset(gca,'XTick',wtick,'YTick',magtick);

% M=25

, [4 v- k+ O5 x# `* l4 K. F* ]M=25; n=0:M; x=ones(1,length(n));
+ Z* x% C* U' q2 n) x5 @1 t+ T3 k
, d  p1 x" ?8 L5 I9 c2 l# {8 uX=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]);
  A& ?" u8 W4 U  `: a; A
. Y/ v+ b5 V  u1 n# q, f& r) ytitle(['M=25']); set(gca,'XTick',wtick,'YTick',magtick);

% M=50

/ g& w( s3 o' K% v& }: C% ZM=50; n=0:M; x=ones(1,length(n));
4 k! b9 Q7 b' a2 C* w2 ?- ^
$ R/ U- s) C+ p. YX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
+ [4 f: _# |# w; y  T% d, X7 r
0 c3 ?  c0 h3 M# L) ^, Asubplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

+ D, v( k. T8 q4 n! N3 P' Axlabel('\omega/\pi'); ylabel('|X|');
# G& [# i2 j, y, e: _6 _7 \
title('M=50'); set(gca,'XTick',wtick,'YTick',magtick);

% M=101
: h) d3 \/ N0 \7 j% U2 I
& _# }4 {. A' D* L1 iM=101; n=0:M; x=ones(1,length(n));

X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
$ q' f' D& x/ s  s0 e1 v. a
subplot(2,2,4); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
1 C7 ?0 ?- d: P9 d, S
xlabel('\omega/\pi'); ylabel('|X|');
# ?, J5 H  p& r3 r: r" m5 P8 P1 p
title(['M=101']);
9 J, A% w0 X1 K5 B( E7 ~8 Q
  z- c% \4 Z. E* |( b9 y
; c. |7 {# O& r; v& a! h+ _: L1 V, L三角窗:
* T0 N( G! D# L4 B6 S7 h1 V8 ^
1 d5 }. z  p1 q7 b; |# D" k% Triangular Window:

% DTFT of a Triangular Window,M = 10,25,50,101

2 q; q- Y( ?+ ]; S& M! v5 N# `; gclc; 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];
) |1 c9 \+ M& L+ |5 r" }
% M = 10

+ o2 e+ ^/ x0 f: j# M% sM=10; n=0:M; x=(1-(abs(M-1-(2*n))/(M+1)));
* P6 J/ m! Q5 e# j/ `8 v1 ]1 W
7 Z8 X' r3 J5 z% fx=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);

/ ~$ F8 o5 N8 P8 a- \subplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

ylabel('|X|'); title(['M = 10']);

$ }2 \: I# N4 k# f" l1 T+ Zset(gca,'XTick',wtick,'YTick',magtick);
9 K& W1 a3 G( E' ^2 p
; a1 A! |4 d, u$ ?' M+ K8 s% M = 25

5 |7 j( D" M( l2 T9 gM=25; n=0:M; x=(1-(abs(M-1-(2*n))/(M+1)));

6 b* R$ |9 G& ?X=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);

" h  W3 W0 T9 _# K& o* s; Jsubplot(2,2,2); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
; W$ H$ _3 ?) c! [, Y8 I
title(['M = 25']); set(gca,'XTick',wtick,'YTick',magtick);
2 j' ~" {; d- b+ h7 L, G6 w* R6 Y
+ S  ~' C, q1 l& B9 t- D% M = 50
. L/ l  K+ i% N1 Z3 J
M=50; n=0:M; x=(1-(abs(M-1-(2*n))/(M+1)));

9 ?1 J# [: Z# m7 y* ]+ f# dX=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]);

- H7 z1 {* @  C) _xlabel('\omega/\pi'); ylabel('|X|'); title(['M = 50']);
2 d1 J( D- r2 J3 R2 D5 \7 ]  }7 o
set(gca,'XTick',wtick,'YTick',magtick);
8 R7 a) y" f. V$ C3 ~0 x" g
/ K& l+ F7 ?' H% M = 100

M=101;n=0:M; x=(1-(abs(M-1-(2*n))/(M+1)));
* K$ W( j! w( k% s6 J7 p6 y
0 N% V/ [% s! v+ g; mX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);

5 P1 b0 U# Q5 L$ K) fsubplot(2,2,4); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

. \! S) t7 [+ o  T3 Z1 zxlabel('\omega/\pi'); title(['M = 101']);
& l, q% l3 r  O) p
9 v6 x/ j* t* [5 S6 j$ ~2 f, O4 uset(gca,'XTick',wtick,'YTick',magtick);

# O  Z0 {7 o- p# _% O# u

' U& D2 P3 |9 k海宁窗:

; ]3 ]# K% i5 i1 h9 Z- ^* o  K% `% Hann Window

% DTFT of a Hann Window, M = 10,25,50,101

0 W: [( x5 I4 A& @, Yclc;close all;
5 g5 g1 A8 C, |, L( T
Hf_1 = figure; set(Hf_1,'NumberTitle','off','Name','P0304c');
( ]7 J" ^0 Z$ S8 T! X+ w' U4 E6 A
w=linspace(-pi,pi,501); wtick=[-1:0.5:1]; magtick=[0:0.5:1.1];

% M = 10
$ E4 \9 ?4 v) `# T# ?: `/ }
M=10;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1)));
" i8 m0 `, x4 n8 K" l. @5 B/ l  \
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]);

& n( r) O1 E) P+ [% F1 {) @ylabel('|X|'); title(['M = 10']);

set(gca,'XTick',wtick,'YTick',magtick);
: d6 @; u. }1 O7 M' u$ q- p
% M = 25
* u0 g, {4 T/ Y% b, c
M=25;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1)));
$ m* n* z+ l5 w
& M/ \" T. ^" p* K( zX=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]);
- f: s: s6 V/ z- }9 b
title(['M = 25']);

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

% M = 50
* |8 ~; F4 B  R, u
( o4 V6 q) ]$ ^, MM=50;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1)));

: e  Y/ c4 k- SX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);

- k% s5 x4 {. I. c3 x- s$ [subplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
( r0 A" R6 u5 S/ O; u
- e  d9 o7 o2 g* p; v( [$ S6 Mxlabel('\omega/\pi'); ylabel('|X|'); title(['M = 50']);
$ e2 Y/ ~  T" \) a3 S! z# W- b
9 h2 P1 U, Q9 S% d. ?' @set(gca,'XTick',wtick,'YTick',magtick);
8 B4 f* ]8 h. t+ R2 _* t
% M = 101
% W' ?+ \& q/ X; h0 Q
M=101;n=0:M; x=0.5*(1-cos((2*pi*n)/(M-1)));

3 N7 F7 t0 J$ e( ~& @: X' |% tX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);

, e0 B+ C! d7 U1 h: Nsubplot(2,2,4); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
* ^; [& B5 f9 e3 ]. E
* w( J# F$ q  C  H* hxlabel('\omega/\pi'); title(['M = 101']);

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


2 r& c3 A8 P5 p: m哈明窗:

% Hamming Window:
; G( i/ c  L% x! S9 f! S
& I. z! A# ]0 k( a8 M, [clc; close all;
2 G/ W& R) D) `$ t2 s
0 ?2 `# n7 `! i# n. m: fHf_1=figure; set(Hf_1,'NumberTitle','off','Name','P0304d');
; w" d* _: H: M# I& l
, h: J/ l+ X+ M) Sw=linspace(-pi,pi,501); wtick=[-1:0.5:1]; magtick=[0:0.5:1.1];
0 S: j6 Q8 ^2 F7 C5 j9 O
% M = 10
% E4 e/ e3 L+ [  w
M=10; n=0:M; x=(0.54-0.46*cos((2*pi*n)/(M-1)));
0 _! e6 K3 [7 N5 x
! d+ u5 ^) g! `5 Y0 e3 T3 RX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
' E5 L0 S5 Q5 L% e
2 m* R8 Q/ Q6 w6 O. t0 B' r+ Ksubplot(2,2,1); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

ylabel('|X|'); title(['M = 10']);
0 e0 r$ F1 H8 C9 n
set(gca,'XTick',wtick,'YTick',magtick);

3 ^. w, |# f/ F$ C; ?4 |% M = 25
8 ?9 t+ i. Y$ N1 K/ C9 k  n
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);
# b2 s! N7 I1 C' t
3 W$ |% x, X$ ]+ g# b( |2 Qsubplot(2,2,2); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

3 P! }7 _) D. ?/ G9 ?title(['M = 25']);

set(gca,'XTick',wtick,'YTick',magtick);
) m% Q( J* T8 ~. r( c; A+ e
% M = 50

* e3 W0 ?& n1 I* X5 CM=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);
- R" w3 M, s+ S  t7 ?* S0 I/ s  j
& d7 v' o2 c9 u9 bsubplot(2,2,3); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);

xlabel('\omega/\pi'); ylabel('|X|'); title(['M=50']);

) o0 ~, m" y2 }7 X9 U) Pset(gca,'XTick',wtick,'YTick',magtick);
1 y3 i! l8 @9 x1 d% V
% M = 101
6 r1 W3 U1 \$ {% m6 K6 Q+ a
M=101; n=0:M; x=(0.54-0.46*cos((2*pi*n)/(M-1)));
9 v# a. o* C! ^( O
7 K& Q/ A9 u& v  R7 D, VX=dtft(x,n,w); magX=abs(X); magX=magX/max(magX);
4 f& R& @# Q" [# r) L3 K' C
& q) P$ ~& M; g1 ~# _2 Gsubplot(2,2,4); plot(w/pi,magX,'LineWidth',1.5); axis([-1 1 0 1.1]);
' t' O6 B0 p! f8 i1 ^& t
/ V# u3 q$ A* pxlabel('\omega/\pi');title(['M=101']);

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



7 k: H0 C) ^! _) O

relchhiclty

谢谢分享