当前位置: 首页 > >

信号处理课件第11章 经典谱估计

发布时间:

第十一章 经典谱估计
11.1 概述 11.2 自相关函数的估计 11.3 经典谱估计的基本方法 11.4 经典谱估计的质量 11.5 经典谱估计的改进 11.6 经典谱估计算法比较 11.7 短时傅里叶变换

11.1 概述
请抓住并搞清楚如下四个问题:
? ? ? ?

功率谱为什么要估计? 如何估计? 如何评价估计质量? 如不理想,如何改进?

*稳随机信号功率谱的两个定义:
集总*均

两者等效

求极限运算

求均值运算

随机信号的 单个样本

*稳信号 X ( n )
PX ( e
j?

单一样本

x(n, i) ? x(n)
j?

截短
)

xM (n )
) PM ( e
j?

)

Px ( e

可将 x M ( n ) 看作能量信号,因此,可对它作 傅立叶变换,并得到功率谱:

x M ( n ) 的功率谱 PM ( e j? ) 和单个样本的 问题 :

功率谱

Px ( e
PX ( e

j?

)
)

有何关系?和整个随机信号的

功率谱

j?

有何关系?

1. 求极限:

2. 求均值:

单一样本的功率谱不能收敛到所有样本的 功率谱,因此必须有求均值运算,此即如 下定义的来历:

各态遍历信号也是如此。

证明:

双求和变成 单求和:
证明了两个公 式等效。所以 自相关函数是 集总自相关。

功率谱的两个定义都要求:样本无穷多,时 间无限长,即需要集总*均。

实际工作中,我们往往能得到的是:
1. 单一的样本; 2. 单一样本的有限长数据; 问题:如何用这单一样本的有限长数据去估

计原随机信号真实的自相关函数和功
率谱

11.2 自相关函数估计
目的:自身估计的需要;
功率谱估计的需要 定义: 集总自相关

时间自相关

估计方法: 实际求出的 自相关函数 从估计方法上看,实际上是把随机信号“视 为”单样本有限长的确定性信号。问题是:

*似质量如何
Estimation
Estimate Estimator(估计子)

自相关函数估计的质量:

单个样本

偏差
? rx ( m ) ?

1 N

N ?1? m

?

x(n) x(n ? m )

n?0

1. 偏差

估计方法

来自定义

所有样本

所以:

含义

渐*无偏估计

对固定的N,此结论 给出了m的选取原则

那儿来的 三角窗?
在数据上加矩形窗,长度为 N ,该矩形窗函数 的自相关函数正是三角窗!注意矩形窗加在数 据上,三角窗加在相关函数上,体现在估计的 自相关函数的均值上。

2.方差

来自定义 包含两项 前面结果

方差

四阶统计量!

由:

零均值高 斯分布

最后 导出

有:
N ? ?, N ? ?, ? b ia[ r ( m )] ? 0 ? var[ r ( m )] ? 0

渐*一 致估计

3.自相关函数的计算
已知单个样本的 N 点数据 估计
? rx ( m )

两个方法:
(1) 直接按定义:

利用

? ? rx ( m ) ? r x ( ? m )

? rx ( m ), m ? ? ( N ? 1) ~ ( N ? 1)

最大长度

(2) 利用FFT:
Step1: 将 x N ( n ) 补 N 个零得 x 2 N ( n ) ; Step2: 对 x 2 N ( n ) 做FFT,得 X 2 N ( k ) ;

Step3: 对 X 2 N ( k ) 求幅*方,得 X 2 N ( k ) ;
Step4: 由 X 2 N ( k ) 得
2

2

1 N

X 2 N (k )

2



对其作IFFT,得 思考: 和

。 有何关系

自相关函数的另一个估计方法(估计子):

很容易证明: 是 的无偏估计, 但方差性能不好。在一些谱估计的方法中, 有时用到该公式。

要求:很好掌握自相关函数 的估计方法及估计性质。

11.3 经典谱估计
问题的提出:对随机信号 X ( n ) ,我们往往 只能得到它的: 1. 单一的样本 x ( n , i ) ? x ( n ) ; 并且仅是 2. 单一样本的有限长数据;

如何用这 N 数据去估计原随机信号真实的功

率谱

Px ( e

j?

)

经典谱估计中有两个基本的方法:
1.周期图(Periodogram)法:

X

N

(k ) ?

?

N ?1

x ( n )W N

nk

,

n?0

? (k ) ? 1 X PP E R N

N

(k )

2

思路:对 x N ( n ) 做DTFT(DFT),得到频谱;对
该频谱求幅*方,再除以N,即得到“周期图” 功率谱,以此作为对真谱的估计。

2.自相关(Blackman-Tukey BT法)法:
Step1 Step2

因为先要估计自相关函数,所以 又称间接法。与此相对应,周期 图法又称直接法。

3.直接法和间接法的关系:

需要考虑两种情况:
(一)
M ? N ?1 M ? N ?1

(二)

数据的范围 自相关函数 的范围

(一) M ? N ? 1 比较用两种方法的估计出的离散谱:

2N 点的谱,把所能估计出的自相关函数都使 用上了,而估计自相关函数时,把 N 点数据 也全都使用上了。

? (k ) ? 1 X PP E R N

N

(k )

2

? 2 N (k ) ? 1 X (k ) PP E R 2N N

2

对 x N ( n ) 补N 个 零,做DFT,得到

IFFT

结论: 在 M ? N ? 1 时,直接法和间接法 估计的结果是一样的。 使用间接法时,往往取 M ? N ? 1 , 这时二者是不一样的。因此,直接法可看 作是间接法的特例。

思考:x N ( n ) 不补零, 即:
? (k ) ? 1 X PP E R N
N

(k )

2

N点离散谱

如何和

相等?

N点离散谱

(二) M ? N ? 1 相当于只用了部分自相关函数

所以:
加在自相关函数上。目的是将 其截短。第二次加窗。

直接法和间接法之间的关系

11.4经典谱估计的质量
也分两种情况讨论
均值 方差 (一)
M ? N ?1 M ? N ?1 M ? N ?1

主要考察的是

无偏估计

一致估计 周期图和 自相关法 是等效的, 统一考虑

1. 偏差
估计值的均值

自相关 函数估 计的性 质

于是有:

的真实功率谱; 的频谱谱;

的频谱谱;
三角窗; 注意: 三角窗频谱恒为正

由于 最后有:

如何理解这一结果

因为:

所以:

M ? N ?1

周期图和自相关法 都是渐*无偏估计

2. 方差

定义:
又遇到四阶矩问题,直接求解困难。

思路:
(1)假定
(2)求 在

是高斯零均值的随机过程;
处的协方差 :

(3)令

,则
? ? ? co v[ P (? 1 ), P (? 1 )] ? var[ P (? )]

有关方差公式的推导不作要求。主要 是掌握结论,并用来说明问题。

求解的关键

推导的结果:方差

(1)N

? ?



经典功率谱估计不是一致估计

解释:
D0 (? )

??

?B 2

B 2

?

?

D 0 (? ? ? )

D 0 (? ? ? )

B 2

?? ?? ?

B 2

??

? ? ??

? ??

?

?

推导的结果:协方差

(2)

若 假定 有

的主瓣宽度为

; 内:

在主瓣外为零;

那么,在频率范围

D0 (? )

??

?B 2

B 2

?

?

D 0 (? 2 ? ? )

D 0 (? 1 ? ? )

?1 ? ? 2 ? B

??

? ? ?? 2

? ? ?1

?

?

D 0 (? 1 ? ? )

D0 (? ? ? 2 )

?1 ? ? 2 ? B
??

? ? ?1

? ? ?2

?

?



处, 在 处不相关;

说明:随机变量

后果:使估计出的谱曲线起伏加剧;

原因:功率谱的定义中即要求极限,又要求均
值;而实际的估计方法,仅靠单次实现 的有限长,无极限、又无均值运算,因 此产生上述问题。

设想:增大数据长度,效果如何

增大, 的主瓣( B )将变窄,因 此,引起不相关的区域进一步增多,从而引 起谱曲线的更加起伏,实际上是方差变大。
通常,增加 ,会提高谱的分辨率,对经 典谱估计来说,增加 固然会有利于提高 分辨率,但谱曲线的起伏令使用者难以接受, 这是经典谱估计的一个致命缺点。 分辨率和方差(体现在曲线起伏上),是 经典谱估计中的一对矛盾。

N=16 对白 噪声 在不 同长 度情 况下 估计 出的 谱曲 线:
10 0 -10 -20 -30 -40 0 0.25 0.5 -20 0 0 10

N =32

-10

0.25

0.5

10 0 -10 -20 -30 -40 0 0.25 0.5

10 0 -10 -20 -30 -40 0 0.25 0.5

N=64

N =128

经典谱估计质量的讨论:
(二)

? :加在估计的自相关函数上 rx ( m )



周期图谱估计和自相关法的谱估计 不再一样!

1. 偏差

谁的主瓣比较宽
N ? ? ? W (? ) ? ? (? )

假定1: 是慢 变谱,在 的主瓣内*似为 一个常数 假定2
1 2?

? ? V (? ) d ?
?

?

? v (0 ) ? 1

窗函数的一般要求

也是渐*无偏估计!

2.方差
考虑特殊情况, 为白噪序列,其 功率谱应为常数,即 P (? ) ? ?
2

时 对白噪声功率 谱估计的方差 时 对白噪声功率 谱估计的方差

两种情况下估计的方差之比:

: 方差改 进之比

取哈明窗:

1. 在 加上 后,估计的谱 的偏差劣于 M=N-1 时估计的谱,而方 差优于 M=N-1 时估计的谱; 2. 上加窗 以后,估计谱 方差的改进体现在两个方面:
(1)

(2)在 因为B变大,

的范围上, 估计的谱

曲线变得 *滑些

不相关的点变少。

3. 方差的减小是以牺牲分辨率为代价的!
原主瓣宽,取决于 现主瓣宽,取决于

若分辨率能满足要求,则这样做 是有意义的,即既保证了分辨率,又 使估计出的谱较为*滑。

11.5 直接法估计的改进
任务:改进
对 估计的性能;

目标:主要是改进方差的性能

方法:*滑与*均;
1. *滑(Smoothing) 用 对 的加窗来实现

*滑

2. *均(Average) 理论依据: L个独立同分布随即变量和的 分布,方差减小 倍,即:

将一个较长的信号分成若干段,对每一段求功 率谱,每一段的功率谱都是随机变量,然后* 均之。类似相干*均,用以弥补经典谱估计中 缺少的求均值运算。注意:信号应是*稳的, 且每一段的统计特性基本一样。

(1) Bartlett*均 将 分成 段,每段 点,即

?

?

思考:如何确定 M 或者 L?

每一段谱

*均后谱

*均后估计出的功率谱的性能如何?
在数据上加了数据窗 宽度是

结果,在自相关函数上引入了窗函数


类似

的自相关;
引入的

统计性能分析: (1)偏差增大,分辨率进一步下降;
(2)方差减小,但到不了 倍

思考:如何确定M 或者L?

(2) Welch*均

特点:交叠分段

若重叠一半,段数

变大

:不一定是矩形窗,如Hamming窗

归一化因子, 保证无偏估计
? E PP E R (? ) ? P (? )

?

?

Welch *均是常用的经典谱估计方 法,MATLAB中有相应的命令

假定1: 是慢 变谱,在 的主瓣内*似为 一个常数)

假定2 1 ? ? ?W 2?
?

2

(? ) d ? ? 1

Welch*均法的方差比Barttlett方法有明显 的减小,而偏差几乎没有减小 3. Nottall 法:*滑与*均相结合

三种改进方法:

11.6 总结与比较

请掌握如下的方法: 白噪 声1 白噪 声2

由自 己指定
H (z)

H(z)

H(z)

两个输 出都是 随机信 号

令:
则: 得到
y (n)

构成一 复信号

的功率谱;

在 y ( n ) 的基础上再加上四个复正弦,归 一化频率分别是:

调整

,可以得到不同的信噪比,本例取

这样, 的真实功率谱可得到,并可画 出。我们可以此作为比较各种算法的依据。
Px (? ) ? Py (? ) ?

? 2? A
k ?1

4

2 k

? (? ? ? k )

实际工作中,对信号 总取有限长,如 ,由这128点去“求”功 率谱,得到的当然是估计值。

0

0

-20

-20

-40

-40

-60 -0.5 0

-0.25

0 (a)

0.25

0.5

-60 -0.5 0

-0.25

0 (b)

0.25

0.5

-20

-20

-40

-40

-60 -0.5

-0.25

0 (c)

0.25

0.5

-60 -0.5

-0.25

0 (d)

0.25

0.5

(a)真实谱;(b)周期图;(c)Welch*均,四 段,无迭合,Hamming窗;(d)同c, 但迭合16点

(e)BT法,M=32;(f)BT法,M=16
0 0

-20

-20

-40

-40

-60 -0.5

-0.25

0 (e)

0.25

0.5

-60 -0.5

-0.25

0 (f)

0.25

0.5

经典功率谱估计的特点:
1. 物理概念明确,可用FFT快速算法。所以

是大众化的谱估计方法; 2. 对周期图,分辨率受到
对自相关法,分辨率受到

的限制;
的限制;

3. 方差性能不好,不是一致估计, N 增 大时谱曲线反而起伏加剧;

4. 改进方法是“*滑”与“*均”,改进的 目
的是减小方差,但牺牲了分辨率; 5. 注意窗函数的作用与影响:

加在数据上的窗函数:
各个窗 函数的 产生加在自相关函数上的延迟窗: 作用及 影响是 什么?

11.7 短时傅里叶变换
*稳信号:均值、方差及均方都不随时间变 化,自相关函数仅和两个观察时间的差有关, 和观察的具体位置无关 ; 非*稳信号:均值、方差都随时间变化,自相 关函数也和观察的时间位置有关,信号的频率 也随时间而变化,如语音、脑电及其他含有较 多突变分量的信号 。其一阶、二阶统计量和 功率谱的估计显然不能简单地使用*稳信号的 估计方法,必须考虑其时变因素。

方法:分段,每一小段可看作是*稳的。

概念:

x(t ) ? L ( R)
2

其STFT定义为:
*

STFTx (t , ?) ? ? x(? ) gt ,? (? ) d? ? ? x(? ) g (? ? t )e
* ? j??

d?

?? x(? ), g (? ? t )e

j??

?

式中

g t , ? (? ) ? g (? ? t ) e
|| g (? ) ||? 1

j? ?

并且窗函数应取对称函数。

x(τ )
x(? ) g (? ? t1 ) x(? ) g (? ? t2 )
x(? ) g (? ? t3 )

0

t1

t2

t3

τ

Ω
FT FT FT

0

t1

t2

t3

t

概念:

“谱图(spectrogram)”
2

STFTx (t , ?) ?

? x(? ) g (? ? t )e

? j??

d?

2

? S x (t , ?)

谱图是恒正的,且是实的。

由于
所以

|| g (? ) ||? 1

?? S

x

(t , ?)dtd? ? E x

谱图是信号能 量的分布。

考虑 x(t ) 是随机信号的一个样本,谱图可实现 信号功率谱的估计。注意,它们是 (t , ?) 的函数

将信号 x(t ) 变换为一个二维函数 为信号的联合时频分析:
*

(t , ?)

的方法称

STFTx (t , ?) ? ? x(? ) gt ,? (? )d?
S x (t , ?) ?| STFTx (t , ?) |
W x (t , ? ) ?
2

STFT 谱图
?
2 )e
? j? ?

? x (t ?
*

?
2

) x (t ?
*

d?

Wigner分布
C x (t , ? ) ? 1 2?

???

x ( u ? ?2 ) x ( u ? ?2 ) g (? , ? ) e

? j ( ? t ? ? ? ? u? )

dud? d?

Cohen类分布

例1

信号x(n)由三个不同频率的正弦首尾相接 所组成,即
? sin (? 1 n ), ? x ( n ) ? ? sin (? 2 n ), ? sin (? n ), 3 ?
1

0 ? n ? N1 ? 1 N1 ? n ? N 2 ? 1 N2 ? n ? N ?1
Signal in time

Real part

0 -1

N ? N 2 ? N1
Energy spectral density

Linear scale

|STFT| , Lh=48, Nf=192, lin. scale, contour, Thld=5%

2

Frequency [Hz]
159517975 0

? 3 ? ? 2 ? ?1

0.4

0.3

0.2

0.1

0

50

100

150 200 Time [s]

250

300

350

例2

线性频率调制信号(chirp)
x ( n ) ? ex p ( j? n ) ? ex p ( jn? n )
2

其频率 与时间 成正比

非*稳

例2

两个chirp信号,一个频率随时间增长, 一个频率随时间减小,求它们和的谱图。

与本章内容有关的MATLAB文件
1. pwelch.m
本文件用Welch*均法估计一个信号的功率 谱,其基本调用格式是: [Px, F] = pwelch(x, Nfft, Fs, window, Noverlap) 式中 x 是随机信号,Fs是抽样频率,Nfft是对x作 FFT 时 的 长 度 , window 是 选 用 的 窗 函 数 , Noverlap是估计x的功率谱时每一段叠合的长度。 缺省时,Nfft = 256 , noverlap = 0, window = Hanning(Nfft), Fs = 2。输出的Px是估计出的功率 谱,按上述调用格式给出的是幅*方值,F是频 率轴坐标。

2. spectrum.m
功能和pwelch.m类似,可用Welch*均法来估 计一个信号的自功率谱,还可用于估计两个信 号的互功率谱。 3.specgram.m 估计信号的谱图,但实际上估计的是其短时傅 里叶变换。该文件主要针对非*稳信号,当然 也可用于*稳信号,甚至确定性信号。




友情链接: