nt_deboing

y=nt_deboing(x,events) - fit, remove ringing associated with events

  y: clean data

  x: data to clean
  events: list of event samples

 NoiseTools.

0001 function y=nt_deboing(x,events)
0002 %y=nt_deboing(x,events) - fit, remove ringing associated with events
0003 %
0004 %  y: clean data
0005 %
0006 %  x: data to clean
0007 %  events: list of event samples
0008 %
0009 % NoiseTools.
0010 nt_greetings;
0011 
0012 if nargin<2; error('!'); end
0013 
0014 if size(x,2)>1; 
0015     error('x should be column vector'); 
0016 end
0017 
0018 ORDER=10; % order of polynomial trend
0019 NSAMPLES=100; % number of samples over which to estimate impulse response
0020 EXTRA=50; % samples before stimulus to anchor trend
0021 NNUM=8;NDEN=8; % number of filter coeffs
0022 THRESH=3; % threshold for robust detrending
0023 
0024 % remove events too close to beginning or end
0025 events(find(events<EXTRA))=[];
0026 events(find(events>size(x,1)-NSAMPLES))=[];
0027 
0028 y=x;
0029 for iEvent=1:numel(events)
0030     event=events(iEvent);
0031     
0032     % select portion to fit filter response, remove polynomial trend
0033     event_response=x(event-EXTRA:event+NSAMPLES);
0034     event_response=nt_detrend(event_response,ORDER,[],[],THRESH); 
0035     event_response=event_response(EXTRA+1:end);
0036        
0037     event_response=event_response(2:end); % not sure why...
0038     
0039     % estimate filter parameters
0040     event_response=[event_response;zeros(size(event_response))]; % helps ensure stable filter
0041     [B,A]=stmcb(event_response,NNUM,NDEN); 
0042     
0043     % estimate filter response to event
0044     model=filter(B,A,[1;zeros(NSAMPLES,1)]);
0045     y(event+(1:size(model,1)))=x(event+(1:size(model,1))) - model;
0046     %figure(1); clf; plot(nt_demean([event_response(1:size(model,1)),model,x(event+(1:size(model,1)))])); pause
0047 end
0048 
0049 if nargout==0
0050     % plot, don't return
0051     figure(1); clf;
0052     plot([x,y]);
0053     
0054     clear y
0055 end
0056 
0057