line source
+ − % WELCHDFT welch method with dft
+ − %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ − % WELCHDFT comput welch'averaged periodogram with dft
+ − %
+ − % CALL
+ − % Sf = welchdft(x,fs,f,Ns,olap,navs,order,win,psll)
+ − % Sf = welchdft(x,fs,f,Ns,olap,navs,order,win,[])
+ − % Sf = welchdft(x,fs,f,[],olap,navs,order,win,psll)
+ − % Sf = welchdft(x,fs,f,[],olap,navs,order,win,[])
+ − %
+ − %
+ − % INPUT
+ − %
+ − % - x, data series, Nx1 double
+ − % - fs, sampling frequency Hz, 1x1 double
+ − % - f, frequenci vector Hz, 1x1 double
+ − % - Ns, length of overlapping segments
+ − % - olap, overlap percentage
+ − % - navs, number of desired averages
+ − % - order, detrend order, -1,0,1,2,3,4,...
+ − % - win, window name. e.g 'BH92'
+ − % - psll, for Kaiser window only
+ − %
+ − % REFERENCES
+ − %
+ − % D. B. Percival and A. T. Walden, Spectral Analysis for Physical
+ − % Applications (Cambridge University Press, Cambridge, 1993) p 291.
+ − %
+ − % L Ferraioli 28-03-2011
+ − %
+ − % $Id: welchdft.m,v 1.1 2011/03/28 16:37:23 luigi Exp $
+ − %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ − function Sf = welchdft(x,fs,f,Ns,olap,navs,order,win,psll)
+ −
+ − N = numel(x);
+ − if isempty(olap)
+ − % compute shift factor
+ − n = floor((N-Ns)/(navs-1));
+ − else
+ − if olap>1
+ − olap = olap/100;
+ − end
+ − % compute Ns
+ − Ns = floor(N/(olap*(navs-1)+1));
+ − % compute shift factor
+ − n = olap*Ns;
+ − end
+ −
+ − Nf = numel(f);
+ − Sf = zeros(Nf,1);
+ −
+ − % run welch averages
+ − for ii=1:navs
+ − idx1 = n*(ii-1)+1;
+ − idx2 = n*(ii-1)+Ns;
+ − if idx2>N
+ − break
+ − end
+ − seg = x(idx1:idx2);
+ −
+ − % get estimate for a segment
+ − segxx = utils.math.computeDftPeriodogram(seg,fs,f,order,win,psll);
+ − % update sum
+ − Sf = Sf + segxx;
+ − end
+ − % complete average
+ − Sf = Sf./navs;
+ −
+ −
+ −
+ −
+ −
+ −
+ −
+ −
+ − end