基于matlab的经典功率谱估计(有源程序)
- 格式:doc
- 大小:23.00 KB
- 文档页数:3
随机信号利用经典谱估计法估计信号的功率谱作业综述:给出一段信号“asd.wav”,利用经典谱估计法的原理,通过不同的谱估计方法,求出信号的功率谱密度函数。
采用MATLAB语言,利用MATLAB语言强大的数据处理和数据可视化能力,通过GUI的对话框模板,使操作更为简便!在一个GUI界面中,同时呈现出不同方法产生出的功率谱。
这里给出了几种不同的方法:BT法,周期图法,平均法以及Welch法。
把几种不同方法所得到的功率谱都呈现在一个界面中,便于对几种不同方法得到的功率谱作对比。
一.题目要求给出一段信号及采样率,利用经典谱估计法估计出信号的功率谱。
二.基本原理及方法经典谱估计的方法,实质上依赖于传统的傅里叶变换法。
它是将数据工作区外的未知数据假设为零,相当于数据加窗,主要方法有BT法,周期图法,平均法以及Welch法。
1. BT法(Blackman-Tukey)● 理论基础:(1)随机序列的维纳-辛钦定理由于随机序列{X(n)}的自相关函数Rx(m)=E[X(n)X(n+m)]定义在离散点m 上,设取样间隔为,则可将随机序列的自相关函数用连续时间函数表示为等式两边取傅里叶变换,则随机序列的功率谱密度(2)谱估计BT法是先估计自相关函数Rx(m)(|m|=0,1,2…,N-1),然后再经过离散傅里叶变换求的功率谱密度的估值。
即其中可有式得到。
2. 周期图法● 理论基础:周期图法是根据各态历经随机过程功率谱的定义来进行谱估计的。
在前面我们已知,各态历经的连续随机过程的功率谱密度满足式中是连续随机过程第i个样本的截取函数的频谱。
对应在随机序列中则有由于随机序列中观测数据仅在的点上存在,则的N点离散傅里叶变换为:因此有随机信号的观测数据的功率谱估计值(称“周期图”)如下:由于上式中的离散傅里叶变换可以用快速傅里叶变换计算,因此就可以估计出功率谱。
3.平均法:● 理论基础:平均法可视为周期图法的改进。
周期图经过平均后会使它的方差减少,达到一致估计的目的,有一个定理:如果是不相关的随机变量,且都有个均值及其方差,则可以证明它们的算术平均的均值为。
功率谱估计 matlab
在MATLAB中进行功率谱密度估计可以使用多种方法,其中最常
用的是基于信号处理工具箱中的函数。
功率谱密度估计是一种用于
分析信号频谱特性的方法,它可以帮助我们了解信号中不同频率成
分的能量分布情况。
在MATLAB中,可以使用periodogram函数来对信号进行功率谱
密度估计。
该函数可以接受原始信号作为输入,并返回频率和对应
的功率谱密度估计值。
另一个常用的函数是pwelch,它可以对信号
进行Welch方法的功率谱估计,该方法是一种常用的频谱估计方法,可以减小估计值的方差。
除了这些内置函数,MATLAB还提供了其他一些工具和函数用于
功率谱密度估计,比如spectrogram函数用于计算信号的短时功率
谱密度估计,cpsd函数用于计算信号的交叉功率谱密度估计等。
在进行功率谱密度估计时,需要注意选择合适的窗函数、重叠
比例等参数,以保证估计结果的准确性和可靠性。
此外,还需要考
虑信号长度、采样频率等因素对功率谱密度估计的影响。
总之,在MATLAB中进行功率谱密度估计有多种方法和工具可供选择,需要根据具体的应用场景和要求来选择合适的方法和函数进行使用。
希望这些信息能对你有所帮助。
在matlab中,功率谱估计是信号处理和频谱分析中常用的一种方法。
通过对信号的频谱特性进行估计,可以有效地分析信号的功率分布情况,从而为信号处理和系统设计提供重要的参考信息。
在matlab中,提供了多种功率谱估计的函数,以下将对其中几种常用的函数进行介绍和分析。
1. periodogram函数periodogram函数是matlab中用于估计信号功率谱密度的函数之一。
它基于傅里叶变换将离散时间信号转换成频域信号,然后计算频域信号的功率谱密度。
其调用格式为:[Pxx, F] = periodogram(x,window,nfft,fs)其中,x为输入的离散时间信号,window为窗函数,nfft为离散傅里叶变换的点数,fs为信号的采样频率。
periodogram函数返回的Pxx 为功率谱密度估计值,F为对应的频率。
2. pwelch函数pwelch函数也是用于估计功率谱密度的函数,它采用了Welch方法,通过对信号进行分段处理,然后对各段信号进行傅里叶变换,并对各段功率谱密度进行平均。
其调用格式为:[Pxx, F] = pwelch(x,window,noverlap,nfft,fs)其中,x为输入的离散时间信号,window为窗函数,noverlap为相邻分段的重叠点数,nfft为离散傅里叶变换的点数,fs为信号的采样频率。
pwelch函数返回的Pxx为功率谱密度估计值,F为对应的频率。
3. cpsd函数cpsd函数用于估计信号的交叉功率谱密度,即两个信号之间的频谱特性。
其调用格式为:[Pxy, F] = cpsd(x,y,window,noverlap,nfft,fs)其中,x和y为输入的两个离散时间信号,window为窗函数,noverlap为相邻分段的重叠点数,nfft为离散傅里叶变换的点数,fs为信号的采样频率。
cpsd函数返回的Pxy为交叉功率谱密度估计值,F为对应的频率。
4. mscohere函数mscohere函数用于估计信号的相干函数,即两个信号之间的相关性。
功率谱估计 matlab
在MATLAB中,可以使用多种方法来进行功率谱密度(PSD)的估计。
以下是一些常用的方法:
1. 通过信号处理工具箱中的函数进行估计:
MATLAB的信号处理工具箱提供了一些内置函数来进行功率谱密度估计,比如pwelch()和periodogram()函数。
这些函数可以直接对信号进行处理并估计其功率谱密度。
2. 基于频谱估计的方法:
在MATLAB中,你可以使用基于频谱估计的方法来进行功率谱密度估计,比如传统的傅里叶变换、Welch方法、Bartlett方法、Blackman-Tukey方法等。
这些方法可以通过MATLAB中的相关函数来实现,比如fft()函数用于傅里叶变换,pwelch()函数用于Welch 方法估计等。
3. 使用自相关函数:
自相关函数可以用于估计信号的功率谱密度。
在MATLAB中,你
可以使用xcorr()函数来计算信号的自相关函数,然后对自相关函
数进行傅里叶变换来得到功率谱密度估计。
4. 基于模型的方法:
MATLAB中还提供了一些基于模型的方法来进行功率谱密度估计,比如Yule-Walker方法、Maximum Entropy方法等。
你可以使用相
应的函数来实现这些方法,比如pyulear()函数用于Yule-Walker
方法估计。
总的来说,MATLAB提供了丰富的工具和函数来进行功率谱密度
的估计,你可以根据具体的需求和信号特性选择合适的方法来进行
估计。
希望这些信息能够帮助到你。
功率谱估计案例 matlab在MATLAB中进行功率谱估计有许多不同的方法和工具。
其中,常用的方法包括周期图法(periodogram method)、Welch方法、Bartlett方法、Blackman-Tukey方法、自回归模型(autoregressive model)和傅里叶变换法等。
这些方法可以用于估计信号的功率谱密度,进而分析信号的频谱特性。
以周期图法为例,MATLAB提供了periodogram函数来实现功率谱估计。
用户可以直接输入信号数据并指定采样频率,函数将返回频率和对应的功率谱估计结果。
使用periodogram函数可以轻松地对信号进行功率谱分析,并可视化频谱特性。
另外,MATLAB还提供了pwelch函数来实现Welch方法,该方法可以对信号进行分段处理并计算每个段的功率谱估计,最后将结果进行平均以得到最终的功率谱密度估计。
这种方法可以降低估计的方差,更适用于非平稳信号的功率谱分析。
除了内置函数外,MATLAB还提供了丰富的工具箱,如信号处理工具箱(Signal Processing Toolbox)和控制系统工具箱(Control System Toolbox),这些工具箱中包含了更多高级的功率谱估计方法和工具,用户可以根据具体需求选择合适的方法进行功率谱分析。
在实际应用中,用户还可以结合MATLAB中的数据处理和可视化功能,对功率谱估计结果进行进一步分析和展示。
通过MATLAB强大的编程功能,用户可以灵活地定制功率谱估计的流程,并将分析结果以图表或报告的形式输出,从而更好地理解信号的频谱特性。
综上所述,MATLAB提供了丰富的功率谱估计方法和工具,用户可以根据具体需求选择合适的方法进行功率谱分析,并结合MATLAB 的数据处理和可视化功能进行全面的信号频谱特性分析。
matlab 功率谱计算在MATLAB中,可以使用多种方法来计算信号的功率谱。
下面我将从多个角度介绍几种常用的方法。
方法一,使用fft函数计算功率谱。
1. 首先,将信号进行零均值化,即减去信号的均值。
2. 然后,使用fft函数对零均值化后的信号进行傅里叶变换,得到频域表示。
3. 对频域表示进行平方运算,得到每个频率分量的幅度平方。
4. 最后,对幅度平方进行归一化处理,即除以信号长度和采样频率的乘积,得到功率谱密度。
示例代码如下:matlab.% 假设信号为x,采样频率为Fs.x = % 输入信号。
Fs = % 采样频率。
% 零均值化。
x = x mean(x);% 计算功率谱。
N = length(x); % 信号长度。
X = fft(x); % 傅里叶变换。
Pxx = (abs(X).^2)/(NFs); % 幅度平方归一化。
% 绘制功率谱图。
f = (0:N-1)(Fs/N); % 频率轴。
plot(f, 10log10(Pxx));xlabel('频率 (Hz)');ylabel('功率谱密度 (dB/Hz)');方法二,使用pwelch函数计算功率谱。
MATLAB还提供了pwelch函数,可以更方便地计算信号的功率谱密度估计。
pwelch函数使用了Welch方法,可以自动进行分段加窗、重叠和平均处理,得到更准确的功率谱估计结果。
示例代码如下:matlab.% 假设信号为x,采样频率为Fs.x = % 输入信号。
Fs = % 采样频率。
% 计算功率谱。
[Pxx, f] = pwelch(x, [], [], [], Fs);% 绘制功率谱图。
plot(f, 10log10(Pxx));xlabel('频率 (Hz)');ylabel('功率谱密度 (dB/Hz)');以上是两种常用的计算信号功率谱的方法,你可以根据实际需求选择适合的方法进行计算。
基于Matlab功率谱密度估计方法摘要在实际情况下, 许多平稳信号无法导出数学表达式, 要准确获取这些信号的功率谱密度存在一定的困难。
根据维纳-辛钦 (Wiener Khintchine)定理,提出一种基于Matlab编程实现这类信号的功率谱密度的估计方法。
通过仿真实验表明, 该方法简单易行,准确性较高。
关键词:平稳信号;功率谱密度;估计方法;Estimation Method for Power Spectral Density of StationarySignalsAbstractIn practice m any mathematical expressions can no t be derived from stationary signals as there are: some difficulties in getting the power spectral density of the signals According to Wiener Khintchine theorem, the paper suggested an estimation method for power spectral density of the signals based on Matlab programming Simulation results show that the method is simple and comparatively accurate.iim. Keywords:Stationary signal; Powerspectral density; estuation method; 前言信号的功率谱密度在通信系统的设计、信号传输、信号分析及信号处理等方面有很重要的作用。
频谱分析是信号在频率域上的重要手段,他反映了信号的频率成分以及分布情况。
信号频谱估计是信号分析的重要手段,也是信号综合的前提在实际情况下, 许多信号很难导出闭合的数学表达式, 例如密勒码 ( Miler code)、双极性不归零码 ( NRZI) 等信号至今在时域中还没有数学表达式, 虽然Hecht& Guida在1969年导出了密勒码的功率谱密度的表达式并做出了图形,但推导过程及得出的表达式非常复杂。
实验功率谱估计实验目的:1、掌握最大熵谱估计的基本原理。
2、了解最终预测误差(FPE)准则。
3、掌握周期图谱估计的基本原理。
4、掌握传统谱估计中直接法与间接法之间的关系。
5、复习快速傅里叶变换与离散傅里叶变换之间关系。
实验内容:1、设两正弦信号的归一化频率分别为0.175和0.20,用最大熵法编程计算信噪比S/N=30dB、N=32点时该信号的最大熵谱估计结果。
2、用周期图法编程计算上述信号的谱估计结果。
程序示例:1、最大熵谱估计clc;N=32;SNR=30;fs=1;t=1:N;t=t/fs;y=sin(2*pi*0.175*t)+sin(2*pi*0.20*t);x = awgn(y,SNR);M=1;P(M)=0;Rx(M)=0;for n=1:NP(M)=P(M)+(abs(x(n)))^2;ef(1,n)=x(n);eb(1,n)=x(n);endP(M)=P(M)/N;Rx(M)=P(M);M=2;A=0;D=0;for n=M:NA=A+ef(M-1,n)*eb(M-1,n-1);D=D+(abs(ef(M-1,n)))^2+(abs(eb(M-1,n-1)))^2; endxishu=-2*A/D;a(M-1,M-1)=-2*A/D;P(M)=P(M-1)*(1-(abs(xishu))^2);FPE(M-1)=P(M)*(N+M)/(N-M);TH=FPE(M-1);for n=M:Nef(M,n)=ef(M-1,n)+xishu*eb(M-1,n-1);eb(M,n)=eb(M-1,n-1)+xishu*ef(M-1,n);endM=M+1;A=0;D=0;for n=M:NA=A+ef(M-1,n)*eb(M-1,n-1);D=D+(abs(ef(M-1,n)))^2+(abs(eb(M-1,n-1)))^2;endxishu=-2*A/D;a(M-1,M-1)=-2*A/D;P(M)=P(M-1)*(1-(abs(xishu))^2);FPE(M-1)=P(M)*(N+M)/(N-M);for m=1:M-2a(M-1,m)=a(M-2,m)+xishu*a(M-2,M-1-m);endwhile FPE(M-1)<THTH=FPE(M-1);for n=M:Nef(M,n)=ef(M-1,n)+xishu*eb(M-1,n-1);eb(M,n)=eb(M-1,n-1)+xishu*ef(M-1,n);endM=M+1;A=0;D=0;for n=M:NA=A+ef(M-1,n)*eb(M-1,n-1);D=D+(abs(ef(M-1,n)))^2+(abs(eb(M-1,n-1)))^2;endxishu=-2*A/D;a(M-1,M-1)=-2*A/D;P(M)=P(M-1)*(1-(abs(xishu))^2);FPE(M-1)=P(M)*(N+M)/(N-M);for m=1:M-2a(M-1,m)=a(M-2,m)+xishu*a(M-2,M-1-m);endendT=1/fs;sum1=0;f=0.01:0.01:0.5;for m=1:M-1;sum1=sum1+a(M-1,m)*exp(-j*2*pi*m*f*T);ends1=(abs(1+sum1)).^2;s=P(M)*T./s1;plot(f,10*log10(s),'k');xlabel('f/fs');ylabel('功率谱/dB');2、周期图谱估计clc;clear;N=32;SNR=30;fs=1;t=1:N;t=t/fs;y=sin(2*pi*0.175*t)+sin(2*pi*0.20*t);x = awgn(y,SNR);sum1=0;f=0.05:0.01:0.5;for m=1:Nsum1=sum1+x(m)*exp(-j*2*pi*m*f);ends=(abs(sum1)).^2/N;plot(f,10*log10(s),'k');xlabel('f/fs');ylabel('功率谱/dB');实验结果:1、最大熵法估计结果:2、周期图法估计结果:。
第25卷第2期 湖北民族学院学报(自然科学版) Vol .25 No .22007年6月 Journal of Hubei I nstitute f or Nati onalities (Natural Science Editi on ) Jun .2007收稿日期:2006-12-11.基金项目:国家民委自然科学基金项目资助(MZZ04004).作者简介:刘嵩(1979-),男,硕士,主要从事模式识别与Matlab 的研究.基于Matlab 的W elch 法功率谱估计刘 嵩1,罗 敏2(1.湖北民族学院信息工程学院,湖北恩施445000;2.湖北民族学院理学院,湖北恩施445000)摘要:在经典谱估计领域,由W elch 提出的修正周期图法获得了有效的应用,在多数情况下,其频率分辨率、估计方差能满足实际应用.运用M atlab 软件进行W elch 算法仿真与分析,从而得出一些有价值的结论.关键词:功率谱估计;周期图;窗函数;分辨率中图分类号:TP311.13文献标识码:A 文章编号:1008-8423(2007)02-0215-02功率谱估计(PSE )是用有限长的数据来估计信号的功率谱,它对于认识一个随机信号来讲是重要的,是数字信号处理的重要研究内容之一.而MAT LAB 是一个强大的数值计算软件,程序设计自由度大,程序的可移植性好并具有强大的图形处理功能,因此可以利用MAT LAB 来实现功率谱估计.本文介绍了功率谱估计W elch 算法的MAT LAB 实现.1 W elch 算法目前经典谱估计算法有Barlett 法、W elch 法、Nattall 法等,W elch 算法是由W elch 提出的修正周期图法,是经典谱估计中获得有效应用一种算法[1].W elch 算法谱估计采取数据分段加窗处理再求平均的办法,先分别求出每段的谱估计,然后进行总平均[2].根据概率统计理论,若将原长度为N 的数据分成L 段,每段长度取M =N /L,且各段数据互为独立,则估计的方差将只有原来不分段的1/L,达到一致估计的目的.其算法如式(1):G LM (ω)=6L l =1G M ,l (ω)/L (1) 但若L 增加M 减小,则分辨率下降.相反,若L 减小M 增加,虽然估计的偏差减小,但估计的方差增大.所以在实际中必须兼顾分辨率与方差的要求适当选取L 和M 的值.在分段时为了减小因分段数增加给分辨率带来的影响,采取各段数据有一定重叠的方法[3].2 W elch 算法的M atlab 实现2.1 仿真实现W elch 法谱估计流程如图1.图1 W elch 法谱估计流程图Fig .1 Fl owchart of s pectru m esti m ati on of W elch method在Matlab 中利用p sd 函数与pwelch 函数都可以实现W elch 法的功率谱估计,其方法是一样的,只是部分参数设置有所不同.利用Matlab 产生一随机信号,然后用pwelch 函数进行谱估计.算法的主要源程序如下: nfft=1024;window=har mm ing(100);noverlap=20;[Pxx,f]=pwelch(xn,window,noverlap,nfft,Fs);p l ot_Pxx=103l og10(Pxx);p l ot(f,p l ot_pxx)title(’W elch法功率谱估计(ha mm ing)’)图2 不同窗函数的功率谱估计Fig.2 The PSE of different window functi ons2.2 结果分析上面参考程序选用了Har mm ing窗,还可以类似的将其修改为Chebyshev窗、Rectangular窗,得到仿真波形见图2.由图2分析可知,在其他条件一致的情况下,Ha mm ing窗的噪声水平较低,而Chebyshev窗的分辨率较好,Reclangular窗由于对信号的突然截断及无法避免的吉布斯现象,频谱泄露严重,噪声水平较高.所以窗函数的选择一定要针对不同的信号和不同的处理目的加以选用才能收到良好效果[4].3 结论利用Matlab实现W elch算法,从而完成功率谱估计,不仅实现方便,形象直观,还有利于分析窗函数对功率谱估计的影响.如何选取最佳窗函数,提高频谱分辨率,是课题组正致力研究的问题,上述方法和结论对功率谱估计有一定的参考价值.参考文献:[1] 飞思科技.MAT LAB7辅助信号处理技术应用[M].北京:电子工业出版社,2005.[2] 常建平,李海林.随机信号分析[M].北京:科学出版社,2006.[3] 金连文,韦岗.现代数字信号处理简明教程[M].北京:清华大学出版社,2004.[4] 姚武川,姚天任.经典谱估计方法的MAT LAB分析[J].华中理工大学学报,2000,28(4):45-48.Ana lysis of W elch M ethod i n PSE Ba sed on M a tl abL I U Song1,LUO M in2(1.School of I nfor mati on Engineering,Hubei I nstitute for Nati onalities,Enshi445000,China;2.School of Science,Hubei I nstitute f or Nati onalities,Enshi445000,China)Abstract:I n the field of classical chart,the a mended peri odic table put for ward byW elch is p r oved t o be effective in app licati on.Under many conditi ons,its frequency distincti on and esti m ate variance can satisfy the p ractical ap2 p licati on.This paper gives an examp le t o si m ulate the m with Matlab,analyses the perf or mances of these methods and dra ws s ome useful conclusi ons.Key words:power s pectrum esti m ati on;peri odic table;window functi on;res oluti on612 湖北民族学院学报(自然科学版) 第25卷。
[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),并再周期图计算前直接加进去,加窗的优点是无论什么样的窗函数均可使谱估计非负。
经典功率谱估计
直接法:
直接法又称周期图法,它是把随机序列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));
间接法:
间接法先由序列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);
改进的直接法:
对于直接法的功率谱估计,当数据长度N太大时,谱曲线起伏加剧,若N太小,谱的分辨率又不好,因此需要改进。
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]);
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);。