Home > NoiseTools > nt_idxx.m

# nt_idxx

## PURPOSE

nt_idxx(fname,iname,bsize,chunksize,nfft) - create an index file to summarize large data file

## SYNOPSIS

function [bstats,cstats,sstats]=nt_idxx(fname,iname,bsize,channels_to_keep,nfft,chunksize)

## DESCRIPTION

    nt_idxx(fname,iname,bsize,chunksize,nfft) - create an index file to summarize large data file

      fname: name of data file to index
      iname: name of index file to create [default fname with EXT = .idx]
      bsize: size of blocks over which to calculate stats [default: 100]
      channels_to_keep: ignore other channels
      nfft: fft size for psd [default: 1024]
      chunksize: size of chunks to read from disk [default: 500000]

     If bsize is a struct, the following fields are expected:
       bsize.basic:  blocksize to calculate basic statistics
       bsize.cov:  blocksize to calculate covariance [default: none]
       bsize.psd: blocksize to calculate psd [default: none]
     If bsize is numeric, it refers to basic statistics.

     Usage:
       nt_idx(fname,...): calculate index structs, store in index file

       [bstats,cstats,sstats]=nt_idx(fname,...): return index structs:
         bstats: basic stats (min, max, etc.)
         cstats: covariance
         sstats: psd

     NoiseTools

## CROSS-REFERENCE INFORMATION

This function calls:

- nt_double2int nt_double2int() - recode/decode double as integer to save space
- nt_greetings nt_greetings - display message the first time the toolbox is used
- nt_whoss size=nt_whoss - total Gbytes used by variables

This function is called by:

- nt_plotxx nt_plotxx(fname,bounds,chans) - plot using index file

## SOURCE CODE

    0001 function [bstats,cstats,sstats]=nt_idxx(fname,iname,bsize,channels_to_keep,nfft,chunksize)
    0002 %nt_idxx(fname,iname,bsize,chunksize,nfft) - create an index file to summarize large data file
    0003 %
    0004 %  fname: name of data file to index
    0005 %  iname: name of index file to create [default fname with EXT = .idx]
    0006 %  bsize: size of blocks over which to calculate stats [default: 100]
    0007 %  channels_to_keep: ignore other channels
    0008 %  nfft: fft size for psd [default: 1024]
    0009 %  chunksize: size of chunks to read from disk [default: 500000]
    0010 %
    0011 % If bsize is a struct, the following fields are expected:
    0012 %   bsize.basic:  blocksize to calculate basic statistics
    0013 %   bsize.cov:  blocksize to calculate covariance [default: none]
    0014 %   bsize.psd: blocksize to calculate psd [default: none]
    0015 % If bsize is numeric, it refers to basic statistics.
    0016 %
    0017 % Usage:
    0018 %   nt_idx(fname,...): calculate index structs, store in index file
    0019 %
    0020 %   [bstats,cstats,sstats]=nt_idx(fname,...): return index structs:
    0021 %     bstats: basic stats (min, max, etc.)
    0022 %     cstats: covariance
    0023 %     sstats: psd
    0024 %
    0025 % NoiseTools
    0026 nt_greetings;
    0027
    0028
    0029
    0030 assert(nargin>0, '!');
    0031 if nargin<2 ; iname=[]; end
    0032 if nargin<3 || isempty(bsize)
    0033     bsize.basic=100;
    0034 end
    0035 if nargin<4; channels_to_keep=[]; end
    0036 if nargin<5 || isempty(nfft); nfft=1024; end
    0037 if nargin<6 || isempty(chunksize); chunksize=500000; end
    0038
    0039 if isnumeric(bsize); tmp=bsize; bsize=[]; bsize.basic=tmp; end
    0040
    0041 % check for FieldTrip
    0042 try, ft_version; catch, disp('You must download FieldTrip from http://www.fieldtriptoolbox.org'); return; end
    0043
    0044 % use separate structs to make it easy to read just one kind of stats from file
    0045 bstats=[]; % index structure for basic stats
    0046 cstats=[]; % index structure for covariance
    0047 sstats=[]; % index structure for spectrogram
    0048 bstats.fname=fname;
    0049
    0050 % read header
    0051 h=ft_read_header(fname);
    0052 bstats.header=h;
    0053 bstats.sr=h.Fs;
    0054 bstats.nsamples=h.nSamples;
    0055 bstats.label=h.label;
    0056 bstats.nchans=h.nChans;
    0057
    0058 if isempty(channels_to_keep); channels_to_keep=1:bstats.nchans; end
    0059 if any(channels_to_keep>bstats.nchans); error('!'); end
    0060 bstats.channels_to_keep=channels_to_keep;
    0061 bstats.nchans=numel(channels_to_keep);
    0062
    0063 % allocate basic stats arrays:
    0064 nbasic=ceil(bstats.nsamples/bsize.basic); % total number of blocs for basic stats
    0065 bstats.min=zeros(nbasic,bstats.nchans);
    0066 bstats.max=zeros(nbasic,bstats.nchans);
    0067 bstats.mean=zeros(nbasic,bstats.nchans);
    0068 bstats.rms=zeros(nbasic,bstats.nchans);
    0069 bstats.card=zeros(nbasic,1,'uint32');
    0070
    0071 chunksize=floor(chunksize/bsize.basic)*bsize.basic;
    0072
    0073 % allocate covariance array
    0074 if isfield(bsize,'cov')
    0075     tmp=log2(bsize.cov/bsize.basic);
    0076     assert(tmp==round(tmp), ...
    0077         'bsize.cov should be power of 2 times bsize.basic');
    0078     ncov=ceil(bstats.nsamples/bsize.cov);
    0079     cstats.cov=zeros(ncov,bstats.nchans,bstats.nchans);
    0080     cstats.card=zeros(ncov,1,'uint32');
    0081     chunksize=floor(chunksize/bsize.cov)*bsize.cov;
    0082 end
    0083
    0084 % allocate psd array
    0085 if isfield(bsize,'psd');
    0086     tmp=log2(bsize.psd/bsize.basic);
    0087     assert(tmp==round(tmp), ...
    0088         'bsize.psd should be power of 2 times bsize.basic');
    0089     npsd=ceil(bstats.nsamples/bsize.psd);
    0090     sstats.psd=zeros(npsd,bstats.nchans,nfft/2+1);
    0091     sstats.card=zeros(npsd,1,'uint32');
    0092     sstats.nfft=nfft;
    0093     chunksize=floor(chunksize/bsize.psd)*bsize.psd;
    0094 end
    0095
    0096
    0097 foffset=0;
    0098 boffset=0;
    0099 coffset=0;
    0100 soffset=0;
    0101
    0102 while true
    0103
    0104     %if file_offset>=i.nsamples; break; end
    0105
    0106     % read chunk from disk
    0107     begsample=foffset+1;
    0108     endsample=min(foffset+chunksize,bstats.nsamples);
    0109     x=ft_read_data(fname, 'begsample',begsample,'endsample',endsample);
    0110     x=x'; % --> time X channels
    0111     x=x(:,channels_to_keep);
    0112
    0113     % fold chunk into blocks
    0114     n=floor(size(x,1)/bsize.basic); % number of blocks in this chunk
    0115     xb=x(1:n*bsize.basic,:);
    0116     xb=reshape(xb,[bsize.basic,n,bstats.nchans]);
    0117     bstats.min(boffset+(1:n),:)=min(xb);
    0118     bstats.max(boffset+(1:n),:)=max(xb);
    0119     bstats.mean(boffset+(1:n),:)=mean(xb);
    0120     bstats.rms(boffset+(1:n),:)=sqrt(mean(xb.^2));
    0121     bstats.card(boffset+(1:n),:)=bsize.basic;
    0122     boffset=boffset+n;
    0123
    0124     % extra bit at end of file?
    0125     if size(x,1)>n*bsize.basic
    0126         tmp=x(n*bsize.basic+1,:);
    0127         bstats.min(boffset+1,:)=min(tmp);
    0128         bstats.max(boffset+1,:)=max(tmp);
    0129         bstats.mean(boffset+1,:)=mean(tmp);
    0130         bstats.rms(boffset+1,:)=sqrt(mean(tmp.^2));
    0131         bstats.card(boffset+1,:)=size(tmp,1);
    0132     end
    0133
    0134     foffset=foffset+n*bsize.basic;
    0135
    0136     if isfield(bstats, 'cov')
    0137         n=floor(size(x,1)/bsize.cov); % number of blocks
    0138         xb=x(1:n*bsize.cov,:);
    0139         xb=reshape(xb,[bsize.cov, n, bstats.nchans]);
    0140         for iBlock=1:n
    0141             tmp=squeeze(xb(:,iBlock,:));
    0142             cstats.cov(coffset+iBlock,:,:) = tmp'*tmp;
    0143             cstats.cardcov(coffset+iBlock,:)=bsize.cov;
    0144         end
    0145         coffset=coffset+size(xb,2);
    0146         if size(x,1)>n*bsize.cov
    0147             tmp=x(n*bsize.cov+1:end,:);
    0148             cstats.cov(coffset+1,:,:)=tmp'*tmp;
    0149             cstats.cardcov(coffset+1,:)=size(tmp,1);
    0150         end
    0151     end
    0152
    0153     if isfield(bstats, 'psd')
    0154         n=floor(size(x,1)/bsize.psd); % number of blocks
    0155         xb=x(1:n*bsize.psd,:);
    0156         xb=reshape(xb,[bsize.psd, n, bstats.nchans]);
    0157         for iBlock=1:n
    0158             tmp=squeeze(xb(:,iBlock,:));
    0159             sstats.psd(soffset+iBlock,:,:) = pwelch(tmp, nfft, 'power')';
    0160             sstats.cardpsd(soffset+iBlock,:,:)=bsize.psd;
    0161         end
    0162         soffset=soffset+size(xb,2);
    0163         if size(x,1)>n*bsize.psd
    0164             tmp=x(n*bsize.psd+1:end,:);
    0165             sstats.psd(soffset+1,:,:) = pwelch(tmp, nfft, 'power')';
    0166             sstats.cardpsd(soffset+1,:)=size(tmp,1);
    0167         end
    0168     end
    0169
    0170     nt_whoss
    0171     disp([foffset boffset, coffset]);
    0172
    0173     if endsample>=bstats.nsamples; break; end;
    0174 end
    0175
    0176 if ~nargout
    0177     if isempty(iname)
    0178         [FILEPATH,NAME,EXT]=fileparts(fname);
    0179     else
    0180         [FILEPATH,NAME,EXT]=fileparts(iname);
    0181     end
    0182     iname=[FILEPATH,filesep,NAME,'.idxx'];
    0183     bstats.min=nt_double2int(bstats.min);
    0184     bstats.max=nt_double2int(bstats.max);
    0185     bstats.mean=nt_double2int(bstats.mean);
    0186     bstats.rms=nt_double2int(bstats.rms);
    0187     save(iname, 'bstats', 'cstats', 'sstats');
    0188     clear bstats cstats sstats;
    0189 end
    0190
    0191
    0192
    0193