Different output from pwelch and my own implementation of welch

3 次查看（过去 30 天）
Jan Slavotinek
Jan Slavotinek 2017-2-3
回答： Jeremy 2017-2-10
Hi!
I am currently trying to understand and implement how to estimate power spectrum as i need it for further analysis. So far, for analysing of my data i have used pwelch function, which make a nice output. Now i need to implement somehow my own estimate using only fft. So i have studied how welch method works (divide data into more segments (possibly overlaping, using hanning window, and do fft for each and then average the result).
So i have written my own function where i try to do estimate. Problem is, thht output is different from what produce pwelch function.
My welch:
function [freq, ampl] = func_FFT2_welch(s, Fs, parts, OverlapFactor)
FrameLength = 2^(nextpow2(length(s)/parts) -1);
Overlap = int32(OverlapFactor * FrameLength);
Shift = FrameLength - Overlap;
part = s(1 : FrameLength);
s(1: Shift) = [];
tmpFFT = abs(fft(part .* hanning(FrameLength)))/FrameLength;
ampl = 2*tmpFFT(1 : int32(FrameLength/2+1));
cnt = 1;
while(length(s) >= FrameLength)
part = s(1 : FrameLength);
s(1: Shift) = [];
tmpFFT = abs(fft(part .* hanning(FrameLength)))/FrameLength;
ampl = ampl + 2*tmpFFT(1 : int32(FrameLength/2+1));
cnt = cnt +1;
end
ampl = ampl / cnt;
ampl = ampl';
freq = Fs/2*linspace(0,1, int32(FrameLength/2+1));
end
See plot i have enclosed. Blue plot is pwelch output (desired). Red is irrelevant to my question. Yellow is simple fft over whole data array (it has a lot of noise, shape is fuzzy) Purple is output of my implementation of welch (whole data array divided into n frames = 2^(nextpow2(length(s)/parts) -1) long, each windowed and overlap 30%). Noise cancelation is good, but problem is, that the overall shape of spectrum is changed.
for the future use i will need to find a linear approximation of some parts of the spectrum. I dont like the way how my function changed the slope of the spectrum.
Note that result is in dB ... Also note that in this case i have called my function with parts = 50 and OverlapFactor = 0.3;
i call pwelch with same overlap and with various "frame length". This only changes variance, but a slope of a the spectrum is unchanged by this.

请先登录，再回答此问题。

回答（1 个）

Jeremy
Jeremy 2017-2-10
At first glance the blue line looks like the square of the purple line. Looking at your code confirms this, pwelch is calculating a PSD which is based on each FFT result x the conjugate of the FFT result so it is a quantity squared result.

类别

在 Help Center 和 File Exchange 中查找有关 Spectral Measurements 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by