本文同步发布在我的个人博客宅到没朋友,欢迎来玩。
1.前言
经典功率谱估计基于傅里叶变换的思想,典型代表为BT法和周期图法。
2.自相关函数
理论上求一个随机信号的自相关函数应该使用下面这个公式: R ( s , t ) = E [ X ( s ) x ( t ) ] R(s,t) = E[X(s)x(t)]\quad R(s,t)=E[X(s)x(t)]
但在实际应用中,我们只能得到一个随机信号有限长度的样本函数。
如果一个随机信号是均方遍历的,我们就可以用样本函数的时间自相关代替该随机信号的自相关函数,如下式:
(式1)
理论上只要样本函数无限长,两者就完全相等,但实际应用中我们只能得到有限长样本函数,所以用样本函数的时间自相关去近似随机信号的自相关函数是有误差的,但只要误差在容许范围内,这就是一个很好的近似。
如果m和N都比较大,上式运算量很大,如何减小运算量,实现自相关函数的快速计算?
这就是我之前文章讲过的,用 F F T FFT FFT,具体请移步还在按部就班的算信号自相关?FFT让你体验飞一般的感觉!
3.BT法进行功率谱估计
对上面所求得的 r ^ ( m ) \hat{r} \left( m \right) r^(m)做傅里叶变换,即: S B T ^ ( w ) = ∑ − M M r ^ ( m ) e j w m ∣ m ∣ ≤ N − 1 \hat{S_{BT}} \left( w \right) = \sum_{-M}^{M} \hat{r} \left( m \right )e^{jwm}\quad \left | m \right|\leq N-1 SBT^(w)=−M∑Mr^(m)ejwm∣m∣≤N−1
以上式为结果作为对真实功率谱 S ( w ) S \left( w \right) S(w)的估计,BT法也称为间接法。
注:一般情况下, M < N M<N M<N。
4.周期图法
周期图法又称为直接法,对式1两边直接求傅里叶变换,1式右边是两个信号卷积,时域相卷,频域相乘,所以我们可以得到周期图法表达式: S P R E ^ ( w ) = 1 N ∣ U N ( w ) ∣ 2 \hat{S_{PRE}} \left( w \right) = \frac{1}{N}\left | U_N \left( w \right ) \right |^{2} SPRE^(w)=N1∣UN(w)∣2
其中 U N ( w ) = ∑ n = 0 N − 1 u N ( n ) e − j w n U_N \left( w \right ) = \sum_{n=0}^{N-1}u_N \left( n \right )e^{-jwn} UN(w)=∑n=0N−1uN(n)e−jwn,这也是周期图法被叫做直接法的原因,它是直接观察数据的傅里叶变换求得的。
5.两者的联系与区别
BT法实质是在周期图法法的基础上加了一个矩形窗,即BT法是对周期图法法平滑,平滑使得BT法的方差小于周期图法,但分辨率下降。
6. M a t l a b Matlab Matlab仿真结果
图2 周期图法与BT法估计出的信号功率谱 (N=256,M=64)
可以看到周期图法发分辨率要明显好于BT法,但同时起伏也更大。
7. M a t l a b a Matlaba Matlaba代码
下载代码请移步BT法和周期图法。