Home > NoiseTools > nt_bias_fft.m

# nt_bias_fft

## PURPOSE

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

## SYNOPSIS

function [c0,c1]=nt_bias_fft(x,freq,nfft)

## DESCRIPTION

    [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

## CROSS-REFERENCE INFORMATION

This function calls:

- nt_cov [c,tw]=nt_cov(x,shifts,w) - time shift covariance
- nt_vecmult y=nt_vecmult(x,v) - multiply all rows or columns of matrix by vector

This function is called by:

- nt_zapline [y,yy]=nt_zapline(x,fline,nremove,p,plotflag) - remove power line artifact

## SOURCE CODE

    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(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);