nt_bias_fft

[c0,c1]=nt_bias_fft(x,freq,nfft) - covariance with and w/o filter bias

 x: data 
 freq: row vector of normalized frequencies to keep (wrt sr)
 nfft: fft size

 The filter has zeros at all frequencies except those immediately inferior
 or superior to values in vector freq.
 
 If freq has two rows, keep frequencies between corresponding values on
 first and second row.

 NoiseTools

0001 function [c0,c1]=nt_bias_fft(x,freq,nfft)
0002 %[c0,c1]=nt_bias_fft(x,freq,nfft) - covariance with and w/o filter bias
0003 %
0004 % x: data
0005 % freq: row vector of normalized frequencies to keep (wrt sr)
0006 % nfft: fft size
0007 %
0008 % The filter has zeros at all frequencies except those immediately inferior
0009 % or superior to values in vector freq.
0010 %
0011 % If freq has two rows, keep frequencies between corresponding values on
0012 % first and second row.
0013 %
0014 % NoiseTools
0015 
0016 if max(freq(:))>0.5; error('frequencies should be <= 0.5'); end
0017 if nfft>size(x,1); error('nfft too large'); end
0018 
0019 filt=zeros(floor(nfft/2)+1,1);
0020 
0021 if size(freq,1)==1
0022     for k=1:size(freq,2)
0023         idx=round(freq(1,k)*nfft+0.5);
0024         filt(idx)=1;
0025     end
0026 elseif size(freq,1)==2
0027     for k=1:size(freq,2)
0028         idx=round(freq(1,k)*nfft+0.5) : round(freq(2,k)*nfft+0.5);
0029         filt(idx)=1;
0030     end
0031 else
0032     error('freq should have one or two rows');
0033 end
0034 
0035 filt=[filt;flipud(filt(2:end-1))];
0036 
0037 %plot(filt); pause
0038 
0039 w=hanning(nfft);
0040 
0041 %plot(filt); return
0042 
0043 [m,n,o]=size(x);
0044 c0=nt_cov(x);
0045 c1=zeros(size(c0));
0046 for j=1:o
0047     nframes=ceil((m-nfft/2)/(nfft/2));
0048     for k=1:nframes
0049         idx=(k-1)*nfft/2;
0050         idx=min(idx,m-nfft);
0051         z=x(idx+1:idx+nfft,:,j);
0052         Z=fft(nt_vecmult(z,w));
0053         Z=nt_vecmult(Z,filt);
0054         c1=c1+real(Z'*Z);
0055     end
0056 end
0057 
0058 %[todss,fromdss,ratio,pwr]=dss0(c0,c1);