查看: 3798|回复: 10|关注: 0
关于功率谱的求法等问题
请问各位大神，网上有几种做功率谱的方法，例如直接法（periodogram）、间接法（自相关函数）、。。。等等，但是结果差别很大，
请问哪种方法是准确的？
还有功率谱和频谱的联系是什么？我看一种说法是 相差 根号2，不知所云。
关注者： 25
功率谱：信号先自相关再作FFT。
频& & 谱：信号直接作FFT。
1、 一个信号的频谱,只是这个信号从时域表示转变为频域表示,只是同一种信号的不同的表示方式而已,& &而功率谱是从能量的观点对信号进行的研究,其实频谱和功率谱的关系归根揭底还是信号和功率,能量等之间的关系。
2、 频谱是个很不严格的东西，常常指信号的Fourier变换，是一个时间平均（time average）概念；功率谱的概念是针对功率有限信号的(能量有限信号可用能量谱分析)，所表现的是单位频带内信号功率随频率的变换情况。保留频谱的幅度信息，但是丢掉了相位信息，所以频谱不同的信号其功率谱是可能相同的。
3、功率谱是随机过程的统计平均概念，平稳随机过程的功率谱是一个确定函数；而频谱是随机过程样本的Fourier变换，对于一个随机过程而言，频谱也是一个“随机过程”。（随机的频域序列）
4、功率概念和幅度概念的差别。此外，只能对宽平稳的各态历经的二阶矩过程谈功率谱，其存在性取决于二阶局是否存在并且二阶矩的Fourier变换收敛；而频谱的存在性仅仅取决于该随机过程的该样本的Fourier变换是否收敛。
1、功率谱可以从两方面来定义，一个是自相关函数的傅立叶变换，另一个是时域信号傅氏变换模平方然后除以时间长度。第一种定义就是常说的维纳辛钦定理，而第二种其实从能量谱密度来的。根据parseval定理，信号傅氏变换模平方被定义为能量谱，能量谱密度在时间上平均就得到了功率谱。
2、在频域分析信号分两种：
（1）.对确定性信号进行傅里叶变换，分析频谱信息。
（2）.随机信号的傅里叶信号不存在，转向研究它的功率谱。随机信号的功率谱和自相关函数是傅里叶变换对（即维纳辛钦定理）。
关注者： 25
（1）信号通常分为两类：能量信号和功率信号；
（2）一般来讲，能量信号其傅氏变换收敛（即存在），而功率信号傅氏变换通常不收敛，当然，若信号存在周期性，可引入特殊数学函数（Delta）表征傅氏变换的这种非收敛性；
（3）信号是信息的搭载工具，而信息与随机性紧密相关，所以实际信号多为随机信号，这类信号的特点是状态随机性随时间无限延伸，其样本能量无限。换句话说，随机信号（样本）大多属于功率信号而非能量信号，它并不存在傅氏变换，亦即不存在频谱；
（4）若撇开搭载信息的有用与否，随机信号又称随机过程，很多噪声属于特殊的随机过程，它们的某些统计特性具有平稳性，其均值和自相关函数具有平稳性。对于这样的随机过程，自相关函数蜕化为一维确定函数，前人证明该确定相关函数存在傅氏变换；
（5）能量信号频谱通常既含有幅度也含有相位信息；幅度谱的平方（二次量纲）又叫能量谱（密度），它描述了信号能量的频域分布；功率信号的功率谱（密度）描述了信号功率随频率的分布特点（密度：单位频率上的功率），业已证明，平稳信号功率谱密度恰好是其自相关函数的傅氏变换。对于非平稳信号，其自相关函数的时间平均（对时间积分，随时变性消失而再次退变成一维函数）与功率谱密度仍是傅氏变换对；
（6）实际中我们获得的往往仅仅是信号的一段支撑，此时即使信号为功率信号，截断之后其傅氏变换收敛，但此变换结果严格来讲不属于任何“谱”（进一步分析可知它是样本真实频谱的平滑：卷积谱）；
（7）对于（6）中所述变换若取其幅度平方，可作为平稳信号功率谱（密度）的近似，是为经典的“周期图法”；
（8）FFT是DFT的快速实现，DFT是DTFT的频域采样，DTFT是FT的频域延拓。人们不得已才利用DFT近似完成本属于FT的任务。若仅提FFT，是非常不专业的
关注者： 25
　　功率谱估计在现代信号处理中是一个很重要的课题，涉及的问题很多。在这里，结合matlab，做一个粗略介绍。功率谱估计可以分为经典谱估计方法与现代谱估计方法。经典谱估计中最简单的就是周期图法，又分为直接法与间接法。直接法先取N点数据的傅里叶变换（即频谱），然后取频谱与其共轭的乘积，就得到功率谱的估计；间接法先计算N点样本数据的自相关函数，然后取自相关函数的傅里叶变换，即得到功率谱的估计.都可以编程实现，很简单。在matlab中，周期图法可以用函数periodogram实现。
　　但是周期图法估计出的功率谱不够精细，分辨率比较低。因此需要对周期图法进行修正，可以将信号序列x(n)分为n个不相重叠的小段，分别用周期图法进行谱估计，然后将这n段数据估计的结果的平均值作为整段数据功率谱估计的结果。还可以将信号序列x(n)重叠分段，分别计算功率谱，再计算平均值作为整段数据的功率谱估计。
这2种称为分段平均周期图法，一般后者比前者效果好。加窗平均周期图法是对分段平均周期图法的改进，即在数据分段后，对每段数据加一个非矩形窗进行预处理，然后在按分段平均周期图法估计功率谱。相对于分段平均周期图法，加窗平均周期图法可以减小频率泄漏，增加频峰的宽度。welch法就是利用改进的平均周期图法估计估计随机信号的功率谱，它采用信号分段重叠，加窗，FFT等技术来计算功率谱。与周期图法比较，welch法可以改善估计谱曲线的光滑性，大大提高谱估计的分辨率。
　　现代谱估计主要针对经典谱估计分辨率低和方差性不好提出的，可以极大的提高估计的分辨率和平滑性。可以分为参数模型谱估计和非参数模型谱估计。参数模型谱估计有AR模型，MA模型，ARMA模型等；非参数模型谱估计有最小方差法和MUSIC法等。由于涉及的问题太多，这里不再详述，可以参考有关资料。
matlab中，现代谱估计的很多方法都可以实现。music方法用pmusic命令实现；pburg函数利用burg法实现功率谱估计；pyulear函数利用yule-walker算法实现功率谱估计等等。
　　另外，sptool工具箱也具有功率谱估计的功能。窗口化的操作界面很方便，而且有多种方法可以选择。
经典功率谱估计
直接法又称周期图法，它是把随机序列x(n)的N个观测数据视为一能量有限的序列，直接计算x(n)的离散傅立叶变换，得X(k)，然后再取其幅值的平方，并除以N，作为序列x(n)真实功率谱的估计。
间接法先由序列x(n)估计出自相关函数R(n)，然后对R(n)进行傅立叶变换，便得到x(n)的功率谱估计。
改进的直接法：
对于直接法的功率谱估计，当数据长度N太大时，谱曲线起伏加剧，若N太小，谱的分辨率又不好，因此需要改进。
Bartlett平均周期图的方法是将N点的有限长序列x(n)分段求周期图再平均。
功率谱：信号先自相关再作FFT。
频& & 谱：信号直接作FFT。
关于功率谱的求法，能否给我举一下例子（程序）？
我自己搜到的程序，原理跟你说的一样，但是结果差别很大，我就搞不清楚哪个方法是对的了
谁能帮我解释一下功率谱幅值的意义。 它和频谱的幅值之间有没有联系？
@ljelly，版主，我改写的程序，后三种方法的结果是相同了，不知道对不对？
我需要的东西是 “等带宽分析得到的功率谱”，以前办公室的做法是FFT化成单边谱后除以根2，得到功率谱有效值，我一直搞不清楚这种做法的目的。而且这跟我直接求功率谱的结果（幅值）也不同。
恳请版主指点迷津！！！！
Fs=1000;%采样频率
N=1000;%采样点数 Fs/N=t(采样总时间)
%频率分辨率=Fs/N需满足关注频率为频率分辨率的整数倍
n=0:1/Fs:(Fs/N-1/Fs);
xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+5*cos(2*pi*150*n)+rand(size(n));%数据
ff=n(1:N/2)*Fs;
subplot(421),plot(n,xn),title('时域图');
xlabel('时间');
ylabel('幅值');
%%fft
y=fft(xn,N);
y1=abs(y)*2/N;
plot_y=y1(1:(N/2));
%ploty2=plot_y/sqrt(2);
subplot(422),plot(ff,plot_y),title('频谱图（FFT）');
%stem(ff,plot_y,'r');
xlabel('频率/Hz');
ylabel('幅值');
%power spectrum（间接法）
%nfft=1024;%index=0:round(nfft/2-1);%k=index*Fs/
cxn=xcorr(xn,'coeff'); %计算序列的自相关函数
%xcorr(,'')中''内缺省or unbiased or biased or coeff
CXk=fft(cxn,N);
Pxx=abs(CXk);
plot_Pxx=10*log10(Pxx(1:N/2));
subplot(423),plot(ff,plot_Pxx),title('间接法（自相关函数法）功率谱');
xlabel('频率/Hz');
ylabel('功率谱密度');
%power spectrum（直接法)
window=boxcar(length(xn)); %矩形窗
%nfft=1024;
[Pxx,f]=periodogram(xn,window,N,Fs); %直接法
subplot(424),plot(ff,Pxx(1:N/2)),title('直接法（periodogram）功率谱');
xlabel('频率/Hz');
ylabel('功率谱密度');
%welch method
window=boxcar(length(xn)); %矩形窗
noverlap=20; %数据无重叠
range='half'; %频率间隔为[0 Fs/2]，只计算一半的频率
[Pxx,f]=pwelch(xn,window,noverlap,N,Fs,range);
%plot_Pxx=10*log10(Pxx);
subplot(425),plot(f,Pxx),title('welch法功率谱');
%Barlett method
window=boxcar(length(n)); %矩形窗
noverlap=0; %数据无重叠
p=0.9; %置信概率
[Pxx,Pxxc]=psd(xn,N,Fs,window,noverlap,p);
%plot_Pxx=10*log10(Pxx);&&
subplot(426),plot(Pxxc,Pxx/N*2),title('barlett法功率谱');
　　功率谱估计在现代信号处理中是一个很重要的课题，涉及的问题很多。在这里，结合matlab，做一个粗略介绍 ...
版主您好， 我需要的东西是“等带宽分析的功率谱”，我不是很清楚这个指什么，办公室以前的做法是将加速度时程做FFT然后转化到单边谱再除以根2。
这个与我直接求功率谱的结果（幅值）有很大不同，下边是我写的求功率谱的程序，后三种方法的结果是近似相等的，但与上边除根2的做法不一样。因为我要在频域计权，幅值对我很重要。
请问版主我这个问题如何解释？
Fs=1000;%采样频率
N=1000;%采样点数 Fs/N=t(采样总时间)
%频率分辨率=Fs/N需满足关注频率为频率分辨率的整数倍
n=0:1/Fs:(Fs/N-1/Fs);
xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+5*cos(2*pi*150*n)+rand(size(n));%数据
ff=n(1:N/2)*Fs;
subplot(421),plot(n,xn),title('时域图');
xlabel('时间');
ylabel('幅值');
%%fft
y=fft(xn,N);
y1=abs(y)*2/N;
plot_y=y1(1:(N/2));
%ploty2=plot_y/sqrt(2);
subplot(422),plot(ff,plot_y),title('频谱图（FFT）');
%stem(ff,plot_y,'r');
xlabel('频率/Hz');
ylabel('幅值');
%power spectrum（间接法）
%nfft=1024;%index=0:round(nfft/2-1);%k=index*Fs/
cxn=xcorr(xn,'coeff'); %计算序列的自相关函数
%xcorr(,'')中''内缺省or unbiased or biased or coeff
CXk=fft(cxn,N);
Pxx=abs(CXk);
plot_Pxx=10*log10(Pxx(1:N/2));
subplot(423),plot(ff,plot_Pxx),title('间接法（自相关函数法）功率谱');
xlabel('频率/Hz');
ylabel('功率谱密度');
%power spectrum（直接法)
window=boxcar(length(xn)); %矩形窗
%nfft=1024;
[Pxx,f]=periodogram(xn,window,N,Fs); %直接法
subplot(424),plot(ff,Pxx(1:N/2)),title('直接法（periodogram）功率谱');
xlabel('频率/Hz');
ylabel('功率谱密度');
%welch method
window=boxcar(length(xn)); %矩形窗
noverlap=20; %数据无重叠
range='half'; %频率间隔为[0 Fs/2]，只计算一半的频率
[Pxx,f]=pwelch(xn,window,noverlap,N,Fs,range);
%plot_Pxx=10*log10(Pxx);
subplot(425),plot(f,Pxx),title('welch法功率谱');
%Barlett method
window=boxcar(length(n)); %矩形窗
noverlap=0; %数据无重叠
p=0.9; %置信概率
[Pxx,Pxxc]=psd(xn,N,Fs,window,noverlap,p);
%plot_Pxx=10*log10(Pxx);&&
subplot(426),plot(Pxxc,Pxx/N*2),title('barlett法功率谱');
很遗憾没人回答。
我根据自己最近的摸索得出一些结论，列于下边，供其他人参考（注意：不一定对！！）
1、按照功率谱密度的定义（自相关函数）可以得出一个结论，用这个结果反推信号的有效值是对的，所以我认为这个方法的得到的功率谱幅度是标准的；
2、直接法及直接法的改进，网上有一些解释，目前没有搞清楚是怎么做的，希望楼下有人帮我解释一下。
关注者： 90
很遗憾没人回答。
我根据自己最近的摸索得出一些结论，列于下边，供其他人参考（注意：不一定对！！）
我记得我已经回答过你的问题了。
Powered by急急急，小白求助。跪求大大解释下出现的错误_matlab吧_百度贴吧
&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&签到排名：今日本吧第个签到，本吧因你更精彩，明天继续来努力！
本吧签到人数：0成为超级会员，使用一键签到本月漏签0次！成为超级会员，赠送8张补签卡连续签到：天&&累计签到：天超级会员单次开通12个月以上，赠送连续签到卡3张
关注：88,228贴子：
急急急，小白求助。跪求大大解释下出现的错误收藏
??? gm.m:&&&&&&&&& |Error: Expected a variable, function, or constant, found "end of line".
登录百度帐号推荐应用
为兴趣而生，贴吧更懂你。或离散时间傅里叶变换 DTFT_图文_百度文库
两大类热门资源免费畅读
续费一年阅读会员，立省24元！
离散时间傅里叶变换 DTFT
||文档简介
通信工程老师|
总评分4.6|
浏览量12373
&&维​纳​。​K​。​恩​格​尔​著​ ​数​字​信​号​处​理​-​-​-​使​用​M​A​T​L​A​B​ ​第​三​章
阅读已结束，如果下载本文需要使用
想免费下载本文？
下载文档到电脑，查找使用更方便
还剩63页未读，继续阅读
你可能喜欢查看: 11365|回复: 3|关注: 0
怎么用Matlab编写DTFT(离散时间傅立叶变换)?
如何用matlab画出一些经典序列的DTFT,我这边只有一个例子:&& n=-1:3;x=1:5;
k=0:500;w=(pi/500)*k;
X=x*(exp(-j*pi/500)).^(n'*k);
magX=abs(X);angX=angle(X);
subplot(2,1,1);plot(w/pi,magX);title('幅度响应');
ylabel('幅度');xlabel('以\pi为单位的频率');
subplot(2,1,2);plot(w/pi,angX);title('相位响应');
ylabel('相位/\pi');xlabel('以\pi为单位的频率');复制代码请各位帮忙,谢谢
fc1=2*pi*50;
t=1:1/fs:10;
L=length(t);
x=cos(fc1.*t)+2;
L=length(t);
NFFT=2^nextpow2(L);
freq_x=fft(x,NFFT)/L;
f = fs/2*linspace(0,1,NFFT/2);
plot(f,2*abs(freq_x(1:NFFT/2)));
谢谢楼上的,但能不能解释下上面的程序?
问题已经解决,正如我给的例子,可适用于任何序列的DTFT.
Powered by}