Converting spectrum.pwelch//psd to equivalent pwelch call
5 次查看(过去 30 天)
显示 更早的评论
I have successfully converted a number of my code blocks from the older method to pwelch. Each of those explicitely called out a window type and segment length.
I can not convert the simplest case because I do not understand what the defaults for spectrum.pwelch & psd are. I do not think I am specifying the correct segment length or window to match the defaults in the original code. I have tried a number of window and segment lengths to try and get close to original psdEst spectrum ... but not even close thus far.
Original Code
h=spectrum.welch;
psdEst = psd(h, timeSeries, 'FS',44100); timeSeries is a 20000 sample sinewave type signal
Replacement Code
psdEst = pwelch(timeSeries,[],[],[],44100);
1 个评论
Mathieu NOE
2024-9-5
Each of those explicitely called out a window type and segment length.
I don't see that neither in the original code or replacement code that you posted . The old and new pwelch functions may have different default settings so we have to use more parameters to make sure you are still having matching results
回答(1 个)
William Rose
2024-9-5
As you can see, the use of spectrum.pwelch is not recommended, and the use of pwelch is recommended.
1 个评论
William Rose
2024-9-5
I like your choice of psdEst for the output, because it is a good reminder that when we compute a spectrum, we are computing an estimate of the true but unknown spectrum.
You can observe default values of spectrum.pwelch by omitting the semicolon when defining h:
h=spectrum.welch
Create a signal: two sinusoids plus noise.
fs=44100; % sampling rate (Hz)
t=0:1/fs:0.1; % time vector (0.1 sec)
f1=7789; f2=10103; % frequencies (Hz)
x=cos(2*pi*f1*t)+sin(2*pi*f2*t)+randn(size(t));
Find spectrum by two methods:
psdEst1=psd(h,x,'FS',fs);
[psdEst2,f2]=pwelch(x,64,32,[],fs);
Plot the two spectra:
figure;
subplot(211), plot(psdEst1.Frequencies,psdEst1.Data,'-b.')
grid on; ylabel('PSD'); title('spectrum.welch')
subplot(212), plot(f2,psdEst2,'-r.')
grid on; xlabel('Frequency (Hz)'); ylabel('PSD'); title('pwelch')
The spectra look the same. I chose the values for pwelch() arguments window and overlap, based on the values displayed when I defined h.
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Spectral Estimation 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!