Reference level P0 in psd() logarithm?
显示 更早的评论
Hi,
I have written code to implement Yule-Walker algorithm myself and compared it against psd() result In my code, I directly take the log10 of of the squared DTFT, and found a constant difference around 100 that my value is above psd() value.
According to http://en.wikipedia.org/wiki/Decibel the definition of power decibel requires a reference level P0, and I therefore suspect it is this P0 that lowered psd() value. Could anyone tell me the value of this P0?
Bob
采纳的回答
Wayne King
2012-1-14
1 个投票
Hi Bob, In PSD estimates, you can take this to be 1.
The use of the logarithm in PSD calculations is different than the use of the logarithm in measurements like sound pressure level where the reference level is critical.
The reason you typically take the log in PSD estimation is for variance stabilization. The variance in a PSD estimate like the periodogram is frequency-dependent. Taking the log() removes that frequency dependence and makes the PSD estimate much easier to interpret.
Have you looked at aryule() in the Signal Processing Toolbox?
Basically the Yule-Walker spectral estimators:
- use aryule() to obtain the AR coefficients
- return the error variance along with the AR coefficients from aryule
- use freqz() to estimate the frequency response
- use e*abs(h).^2 as the basis for PSD estimate where e is the prediction error from aryule() and h is the frequency response from freqz()
- scale appropriately with the sampling frequency (if in Hz) or (2*pi) if in radians/sample.
更多回答（1 个）
Wayne King
2012-1-14
You're welcome Bob, just to reiterate what I stated above, if you try the following:
x = randn(1000,1);
[a,e] = aryule(x,4);
[h,w] = freqz(1,a,512);
Pxx = e*abs(h).^2;
Pxx(2:end-1) = 2*Pxx(2:end-1);
Pxx = Pxx./(2*pi);
plot(w./pi,10*log10(Pxx),'k'); grid on;
set(gca,'xlim',[0 1]);
figure;
pyulear(x,4,512);
You'll see what I mean about the scaling. aryule() obtains the biased estimate of the autocorrelation and then solves the Yule-Walker equations using levinson().
类别
在 帮助中心 和 File Exchange 中查找有关 Linear Predictive Coding 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
