求信号的频率和幅值修正

yxiankun

我分别用psd和pwelch两种方法对幅值为10，频率为15的正弦信号进行分析
为什么用两种方法得到的幅值会相差这么多，请高人指点。
正确的幅值是不是应该为100啊，也就是10的平方？
还有为什么幅值会随窗长的改变而改变啊
怎样才能求出准确的功率谱幅值啊？
谢谢！
clc,clear,close all;
fs=100;
t=0:1/fs:15;
A=10*sin(15*2*pi*t);
N=512;
window=hanning(N);
noverlap=N/2;
figure
plot(t,A)
[P1,f1]=psd(A,length(A),fs,window,noverlap,'none');
figure
plot(f1,P1)
[P2,f2]=pwelch(A,window,noverlap,length(A),100);
figure
plot(f2,P2)

yxiankun

大家帮我啊

hahaer

楼主的方法我认为非常的不妥，至少在功率谱估计的时候你选择FFT点数为采样数据的长度，那么将会有对应长度的0被补上，这样你得不到准确的功率谱估计结果。
其次，PSD方法基本上是已经淘汰的一种方法了，注意看下matlab中关于PSD方法的一个帮助，有个 NOTE 1: To express the result of PSD, Pxx, in units of Power per Hertz multiply Pxx by 1/Fs [1].
而且，PSD对于实数求得的是单边谱，所以对应结果要乘以2，如果非要用你的方法，可以按如下方式进行：
clc,clear,close all;
fs=200;
t=0:1/fs:15;
A=10*sin(15*2*pi*t);
N=512;
NFFT = 1024;
w = window(@hann,N);
noverlap=20;
figure
plot(t,A)
[P1,P2,f1]=psd(A,NFFT,fs,w,noverlap,'none');
figure
P1 = (P1./fs*2);
plot(f1,P1)

[P2,f2]=pwelch(A,w,noverlap,NFFT,fs);
figure
plot(f2,P2)

hahaer

      Fs = 1000;   t = 0:1/Fs:.296;
       x = cos(2*pi*t*200)+randn(size(t));  
       h = spectrum.welch;                  % Create a Welch spectral estimator.
       psd(h,x,'Fs',Fs);                    % Calculate and plot the one-sided PSD.
       hpsd = psd(h,x,'ConfLevel',p);    % PSD with confidence level
       figure,plot(hpsd)

这个是Matlab自带的帮助，可以通过设定不同的方法：
----------------------------
    Power spectral density (PSD)
    ----------------------------
    burg - Burg                             periodogram - Periodogram
    cov  - Covariance                       welch       - Welch
    mcov - Modified Covariance              yulear      - Yule-Walker AR
    mtm  - Thomson multitaper method (MTM)

进行功率谱估计

hahaer

还有一个重要的概念就是：功率谱峰值不一定对应的是该频点的频谱幅值的平方，这点与能量谱需要区分开来，功率谱的实际含义应该是沿频点积分后得到的是信号的功率，也就是功率谱下包含的面积吧，我是这样理解。

yxiankun

谢谢hahaer的帮助。
采用你和程序后两者画出的图是一样的。
另外P1 = (P1./fs*2)中乘2是不是因为采用了Hanning窗后的修正系数啊？
如果是的话，那么在Pwelch中为什么不用乘2了啊

hahaer

另外P1 = (P1./fs*2)中乘2是不是因为采用了Hanning窗后的修正系数啊？
如果是的话，那么在Pwelch中为什么不用乘2了啊
第一个psd直接求的单边谱没有进行幅度折叠叠加
Pwelch好像是加过了，你看下帮助文件吧，其实你最开始的程序最主要的问题还是NFFT的设置问题，其他没有什么问题的

hahaer

乘2不是窗函数的修正系数，有点困了，眼睛疼，睡觉去了，窗函数的修正好像有过一篇文章，现在记不得在哪里了，哪天找到发给你吧

yxiankun

放入hahaer了