×

经典功率谱估计

前端技术网 前端技术网 发表于2024-01-25 04:00:40 浏览2605 评论0

抢沙发发表评论

一、现代谱估计与经典谱估计的区别

1)经典谱估计以傅里叶变换为基础,分为直接法(即周期图法)和间接法,当样本数据很大时经典谱估计的效果是可以接受的,但是样本少时,此估计方法的效果往往不是很好,这是由于经典谱估计的天然缺陷造成的:经典谱估计认为除了样本数据以外的其他信号值全为零,这是不符合实际的,经典谱估计的分辨率低,不可避免的受到加窗的影响,而且它也不是真实的一致估计。

2)现代谱估计以模型为基础,利用采样数据建立模型,对数据进行外推,进而提高了谱估计的分辨率。但是现代谱估计外推的数据不能完全重构数据,只能用来估计功率谱。

经典功率谱估计

二者的区别,LZ应该可以自己归纳了~

二、信号功率谱怎么计算

用FFT求取信号频谱的实部和虚部,实部的平方价虚部的平方就是功率谱。周期性连续信号x(t)的频谱可表示为离散的非周期序列Xn,它的幅度频谱的功率谱平方│Xn│2所排成的序列,就被称之为该周期信号的“功率谱”。

对信号进行傅里叶变换,取sin部分为实部,cos部分为虚部,直接算实部和虚部的平方和,得到的就是频域功率谱的分布推荐使用matlab计算,因为一个函数FFT就可以算出来。

信号x(t)的功率谱密度计算方法:

1、先计算x(t)的傅立叶变换:X(jw),

2、取模:|X(jw)|,再平方:|X(jw)|^2,再除以样本长度:|X(jw)|^2/T

经典功率谱估计

3、就得到: x(t)的功率谱密度函数: Gxx(w)=|X(jw)|^2/T

扩展资料:

周期运动在功率谱中对应尖锋,混沌的特征是谱中出现"噪声背景"和宽锋。它是研究系统从分岔走向混沌的重要方法。在很多实际问题中(尤其是对非线性电路的研究)常常只给出观测到的离散的时间序列X1, X2, X3,Xn,那么如何从这些时间序列中提取前述的四种吸引子(零维不动点、一维极限环、二维环面、奇怪吸引子)的不同状态的信息。

可以运用数学上已经严格证明的结论,即拟合。我们将N个采样值加上周期条件Xn+i=Xi,则自关联函数(即离散卷积)为然后对Cj完成离散傅氏变换,计算傅氏系数。 Pk说明第k个频率分量对Xi的贡献,这就是功率谱的定义。

参考资料来源:百度百科-功率谱

三、从经典谱估计到现代谱估计

谐波分析最早可追溯到古代对时间的研究。18世纪伯努利(Bernoulli)、欧拉(Euler)和拉格朗日(Lagrange)等人对波动方程及其正弦解进行了研究,19世纪初叶,傅立叶(Fourier)证明了在有限时间段上定义的任何函数都可以用正弦和余弦分量的无限谐波的总和来表示。1898年舒斯特(Schuster)以傅立叶分析为基础来拟合待分析信号,研究太阳黑子数的周期变化,并提出了周期图的概念。1930年维纳(Wiener)发表了经典性论文《广义谐波分析》,对平稳随机过程的自相关函数和功率谱密度作了精确的定义,证明了二者之间存在着傅氏变换(以下简称傅氏变换)的关系,从而为功率谱分析奠定了坚实的统计学基础。由于1934年辛钦(Khintchine)也独立地证明了自相关函数和功率谱之间的傅氏变换关系,即维纳-辛钦(Wiener-Khintchine)定理。根据这个定理(详见第一章),平稳离散随机信号x(n)的自相关函数rxx(m)

rxx(m)=E[x(m+n)x*(n)](4-1)

与功率谱Pxx(ejω)之间构成一傅氏变换对,即

地球物理信息处理基础

若x(n)还是各态遍历性的,则其自相关函数可由它的一个采样时间序列用时间平均的方法求出,即

地球物理信息处理基础

在大多数应用中x(n)是实信号,于是上式可写成

地球物理信息处理基础

实际上,一般只能在时域观测到随机信号的有限个采样值(例如N个值),可表示为

xN(n)={x(0),x(1),…,x(N-1)}={x(n),n=0,1,…,N-1}

其自相关函数只能由这N个采样数据进行估计,常用有偏估计

地球物理信息处理基础

这是一种渐近一致估计,称之为采样自相关函数。

用采样自相关函数的傅氏变换作为功率谱的估计,这种方法是布莱克曼(R.Blackman)和杜基(J.Tukey)在1958年提出来的,称为功率谱估计的自相关法(简称BT法)。此方法需要先求出有限个观测数据估计自相关函数,然后再根据式(4-2)计算出功率谱。在快速傅氏变换(FFT)算法提出之前,这是一种最流行的功率谱估计方法。

1965年库利(Cooley)和杜基(Tukey)完善了著名的FFT算法,使计算傅氏变换的速度提高了两个数量级,运算量显著降低,这样DFT变换很快在各领域,特别是在工程实践中得到了广泛应用。由式(4-5)知,

为x(n)与x(-n)的卷积运算,因为

地球物理信息处理基础

若x(n)的傅氏变换为X(ejω),则x(-n)的傅氏变换等于X*(ejω)。对式(4-5)两端取傅氏变换,得到

地球物理信息处理基础

这表明:通过对随机数据直接进行离散傅氏变换,然后取其幅值的平方,再对多样本进行此种运算并取平均值作为功率谱的估计,即舒斯特的周期图,这种谱估计受到了人们的普遍重视,因为它不需要计算自相关函数,而直接计算功率谱。

周期图和自相关法以及它们的改进方法称为功率谱估计的经典方法,周期图和自相关法是经典功率谱估计的两个基本方法。由于FFT的出现,周期图和自相关法往往被结合起来使用,其步骤如下:

(1)对xN(n)补N个零,求

(2)由

作傅氏变换,得

,这时|m|≤M=N-1;

(3)对

加窗函数v(m),这时|m|≤M<<N-1,得

(4)利用

,求

的傅氏变换,即

地球物理信息处理基础

由周期图法得到的功率谱

,其估计方差并不随样本长度的增加而趋于零,

,出人意料的是,不管数据记录有多长,周期图和自相关法得到的估计都不是功率谱的良好估计。事实上,随着记录长度增加,这两种估计的随机起伏反而会更加严重!此外,它们存在着以下两个难以克服的固有缺点。

(1)频率分辨率(区分两个邻近频率分量的能力)不高。因为它们的频率分辨率(赫兹)反比于数据记录长度(秒)(即Δf=k/Tp=k/NT,k为常数,Tp=NT为数据的记录长度,T为采样周期),而实际应用中一般不可能获得很长的数据记录,即观察到的数据只能是有限个,而观察不到的数据被认为是0。这样,如果只有N个观测数据,而对于N以外的数据,信号仍有较强的相关性,那么估计出的功率谱就会出现很大的偏差。

(2)对于有限的观测数据,相当于将信号在时域内乘以矩形窗函数,因而在频域内则相当于使真正的功率谱与sinc函数进行卷积,由于sinc函数不同于δ函数,它有主瓣和旁瓣,这样使卷积后的功率谱不同于真正的功率谱。sinc函数的主瓣不是无限窄的,引起功率谱向附近频域扩展,造成谱的模糊,降低谱的分辨率;同时,由于sinc函数的旁瓣存在,导致能量向旁瓣中“泄漏”(称之为旁瓣泄漏),即引起频谱间的干扰,信号强的功率谱旁瓣影响信号弱的功率谱检测,严重时,会使主瓣产生很大失真,检测不出弱信号,或者把旁瓣误认为是信号,造成假信号。为了对经典功率谱估计进行改进,可以采用各种不同的窗函数,但其结果都是以增加主瓣宽度来换取旁瓣的压低,因此功率谱分辨率低是经典功率谱估计的致命缺点。

为了克服以上缺点,人们曾做过长期努力,提出了平均、加窗平滑等办法,在一定程度上改善了经典功率谱估计的性能。实践证明,对于长数据记录来说,以傅氏变换为基础的经典功率谱估计方法,的确是比较实用的。但是,经典方法始终无法根本解决频率分辨率和功率谱估计稳定性之间的矛盾,特别是在数据记录很短的情况下,这一矛盾显得尤为突出。这就促进了现代功率谱估计方法研究的展开。

现代功率谱估计方法主要是以随机过程(Stochastic Process)的参数模型(Parameter Model)为基础的,称之为参数模型方法。虽然说现代功率谱估计技术的研究和应用主要起始于60年代,但实际上,时间序列模型在非工程领域早已被采用,如Yule在1927年、Walker在1931年都曾使用过自回归模型预测描述经济的时间序列的发展趋势,而Prony则早在1795年就曾采用指数模型去拟合在气体化学实验中获得的数据。在统计学和数值分析领域中,人们也曾采用过模型方法。

现代功率谱估计的提出主要是针对经典功率谱估计(周期图和自相关法)的分辨率和方差性能不好的问题而提出的。1967年Burg在地震学研究中受到线性预测滤波的启发,提出了最大熵谱估计方法,在提高分辨率方面作了最有意义的探索。1968年Parzen正式提出了自回归谱估计方法。1971年Van der Bos证明了一维最大熵谱估计与自回归谱估计等效。1972年出现的谱估计的Prony方法在数学上与自回归方法有某些类似。目前以自回归滑动平均模型为基础的谱估计已经比自回归模型谱估计具有更高的频率分辨率和更好的性能。1973年Pisarenko提出的谐波分解方法提供了可靠的频率估计方法。1981年Schmidt提出了谱估计的多信号分类(MUSIC)算法等。因此,现代功率谱分析主要有ARMA谱分析、最大似然、熵谱估计和特征分解四种方法。ARMA谱分析是一种建模方法,即通过平稳线性信号过程建立模型来估计功率谱密度;熵谱估计包括最大熵谱和最小交叉法;特征分解也叫特征构造法和子空间法,包括Pisarenko谐波分解法、Prony法、MUSIC法和ESPRIT法(用旋转不变技术估计参数方法)。

现代功率谱估计研究仍侧重于一维功率谱分析,而且大部分是建立在二阶矩(相关函数、方差、功率谱密度)基础上的。但由于功率谱密度是频率的实函数,缺少相位信息,因此,建立在高阶谱基础上的谱估计方**引起人们的注意,特别是双谱估计和三谱估计的研究受到了高度的重视。其它如多维谱估计、多通道谱估计等的研究也正在发展中。人们希望这些新方法能更多地在提取信息、估计相位和描述非线性等方面获得应用。

四、[matlab实现经典功率谱估计]matlab功率谱估计

1、直接法:

直接法又称周期图法,它是把随机序列x(n)的N个观测数据视为一能量有限的序列,直接计算x(n)的离散傅立叶变换,得X(k),然后再取其幅值的平方,并除以N,作为序列x(n)真实功率谱的估计。

Matlab代码示例:

clear;

Fs=1000;%采样频率

n=0:1/Fs:1;

%产生含有噪声的序列

xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n));

window=boxcar(length(xn));%矩形窗

nfft=1024;

[Pxx,f]=periodogram(xn,window,nfft,Fs);%直接法

plot(f,10*log10(Pxx));

2、间接法:

间接法先由序列x(n)估计出自相关函数R(n),然后对R(n)进行傅立叶变换,便得到x(n)的功率谱估计。

Matlab代码示例:

clear;

Fs=1000;%采样频率

n=0:1/Fs:1;

%产生含有噪声的序列

xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n));

nfft=1024;

cxn=xcorr(xn,"unbiased");%计算序列的自相关函数

CXk=fft(cxn,nfft);

Pxx=abs(CXk);

index=0:round(nfft/2-1);

k=index*Fs/nfft;

plot_Pxx=10*log10(Pxx(index+1));

plot(k,plot_Pxx);

3、改进的直接法:

对于直接法的功率谱估计,当数据长度N太大时,谱曲线起伏加剧,若N太小,谱的分辨率又不好,因此需要改进。

3.1、Bartlett法

Bartlett平均周期图的方法是将N点的有限长序列x(n)分段求周期图再平均。

Matlab代码示例:

clear;

Fs=1000;

n=0:1/Fs:1;

xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n));

nfft=1024;

window=boxcar(length(n));%矩形窗

noverlap=0;%数据无重叠

p=0.9;%置信概率

[Pxx,Pxxc]=psd(xn,nfft,Fs,window,noverlap,p);

index=0:round(nfft/2-1);

k=index*Fs/nfft;

plot_Pxx=10*log10(Pxx(index+1));

plot_Pxxc=10*log10(Pxxc(index+1));

figure(1)

plot(k,plot_Pxx);

pause;

figure(2)

plot(k,[plot_Pxx plot_Pxx-plot_Pxxc plot_Pxx+plot_Pxxc]);

3.2、Welch法

Welch法对Bartlett法进行了两方面的修正,一是选择适当的窗函数w(n),并再周期图计算前直接加进去,加窗的优点是无论什么样的窗函数均可使谱估计非负。二是在分段时,可使各段之间有重叠,这样会使方差减小。

Matlab代码示例:

clear;

Fs=1000;

n=0:1/Fs:1;

xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n));

nfft=1024;

window=boxcar(100);%矩形窗

window1=hamming(100);%海明窗

window2=blackman(100);%blackman窗

noverlap=20;%数据无重叠

range="half";%频率间隔为[0 Fs/2],只计算一半的频率

[Pxx,f]=pwelch(xn,window,noverlap,nfft,Fs,range);

[Pxx1,f]=pwelch(xn,window1,noverlap,nfft,Fs,range);

[Pxx2,f]=pwelch(xn,window2,noverlap,nfft,Fs,range);

plot_Pxx=10*log10(Pxx);

plot_Pxx1=10*log10(Pxx1);

plot_Pxx2=10*log10(Pxx2);

figure(1)

plot(f,plot_Pxx);

pause;

figure(2)

plot(f,plot_Pxx1);

pause;

figure(3)

plot(f,plot_Pxx2);

关于本次经典功率谱估计和现代谱估计与经典谱估计的区别的问题分享到这里就结束了,如果解决了您的问题,我们非常高兴。