pwelch vs fft, plotting with 10*log10

57 views (last 30 days)
Ganeshkumar M
Ganeshkumar M on 6 Sep 2018
Commented: Honglei Chen on 11 Sep 2018
Hi Guys! been at this for some time and really hope someone can explain.
1. I have 2 signals x and y. Though i get a similar shape of the plot, the power amplitude differs when I use fft or pwelch. Would there be any reason for this?
2. When plotting, why do we need to plot with 10*log10(pxx) instead of only using pxx? Even when I use either ways of plotting, I cant seem to emulate the yrange in the example plot (attached). The closest plot seems to be plot(fxx, 10*log10(pxx), fyy, 10*log10(pyy)) but the yrange is different.
3. When we perform subsequent calculations, such as Granger causality, should we use the 10*log10(pxx) value or just the pxx value?
Thank you so much!
code as below:
fs = 120;
N = 500;
T = 0:1/fs:N/fs;
win_len = N;
noverlap = 0; % no overlap since discrete 97 data
nfft = 500;
trials = 1000;
x = randn(trials,N); y = randn(trials,N); r1 = 0.9; r2 = 0.7;
for trial = 1:trials
for t = 4:N
x(trial,t) = 2*r1*cos(40*(1/fs)*2*pi)*x(trial,t-1) - r1^2*x(trial,t-2) + randn(1);
y(trial,t) = -0.356*x(trial,t-1) + 2*r2*cos(10*(1/fs)*2*pi)*y(trial,t-1) + 0.7136*x(trial,t-2) - r2^2*y(trial,t-2) -0.356*x(trial,t-3) + randn(1);
end
end
xA = reshape(x',1,size(x,1)*size(x,2));
yA = reshape(y',1,size(y,1)*size(y,2));
[pxx,fxx] = pwelch(xA,hanning(win_len),noverlap, nfft, fs);
[pyy,fyy] = pwelch(yA,hanning(win_len),noverlap, nfft, fs);
[pxy,fxy] = cpsd(xA,yA,hanning(win_len),noverlap,nfft, fs);
px = fft_manual(x,win_len,fs);
py = fft_manual(y,win_len,fs);
plot(fxx, 10*log10(pxx), fyy, 10*log10(pyy))
figure; plot(fxx,pxx,fyy,pyy)
figure; plot(fxx,px,fyy,py)
  1 Comment
Ganeshkumar M
Ganeshkumar M on 6 Sep 2018
Edited: Ganeshkumar M on 6 Sep 2018
Here is the paper im trying to replicate where x1 = x and x2 = y:
Stokes, P. A., & Purdon, P. L. (2018). Correction for Stokes and Purdon, A study of problems encountered in Granger causality analysis from a neuroscience perspective. Proceedings of the National Academy of Sciences, 115(29), E6964–E6964. http://doi.org/10.1073/pnas.1809324115

Sign in to comment.

Answers (1)

Honglei Chen
Honglei Chen on 6 Sep 2018
1. pwelch uses Welch method, which involves windowing and averaging on top of fft, that's why they are different. The basic shape is probably similar but FFT probably gives you the best resolution while pwelch provides smoother spectrum
2. 10*log10 is the conversion to dB. It's a monotonic function so again the relative relation in the spectrum is preserved. The purpose of it is to compress the dynamic range so both big and small features can be easily seen.
3. For the downstream processing, in general you should use pxx, not the dB value.
HTH
  2 Comments
Ganeshkumar M
Ganeshkumar M on 11 Sep 2018
Thanks for your answer :) May ask for a bit more clarification?
1. I performed fft and then averaged it over the number of trials. I guess this is the same as the pwelch? The resolution seems to be the same but the peak amplitude is still higher when using the averaged fft.
2. I understand how 10*log10 transforms the graph. It follows the same shape of the example. but the magnitude is totally different. Why is this the case?
Thanks for the clarification so far!
Honglei Chen
Honglei Chen on 11 Sep 2018
I don't know how you averaged it among several trials, but one thing pwelch does is to segment the data and average among different segments. In addition, fft gives you the spectrum while pwelch by default outputs spectral density.
Are you talking about the y-axis range or the value itself? It may be related to the fact that they are not really the same physical quantity.

Sign in to comment.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!