12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152 |
- close all; clear all; clc;
- T=150*1e-3; %[s] duration of post-stimulus time window
- Fs=1e3; %[Hz] sampling frequency
- Ncs=1; % Number of channels
- Nts=T*Fs; % Number of timepoints
- yoff=0.25; % Y-offset
- ts=linspace(0,T,Nts); %[s] time points after stimulus
- ft = @(p,t)(p(5)*exp(p(1)*p(3)+0.5*p(2)^2*p(3)^2-p(3)*t)...
- .*normcdf(t,p(1)+p(2)^2*p(3),p(2))+p(4)*normcdf(t,p(1),p(2))+p(6));
- % Generate some synthetic data with random latency ranging 20÷70 [ms]
- mu=zeros(1,Ncs); lt=zeros(1,Ncs); xt=zeros(Ncs,Nts);
- for ch=1:Ncs
- mu(ch)=(20+randi(50,1));
- xt(ch,:)=feval(ft,[mu(ch) 10 3 0.1 6 yoff],ts*1e3);
- lt(ch)=ts(find(xt(ch,:)>=(yoff+1/3*max(xt(ch,:)-yoff)),1,'first'));
- end
- yt=xt+.01*randn(Ncs,Nts);
- % Get Latencies
- [latTimes0, lsqCurveTimes0, lsqFittedCurve0, lsqFittedParam0, lsqCurveFitErr0] = computeLatencies_offset(xt,ts*1e3,10);
- [latTimes, lsqCurveTimes, lsqFittedCurve, lsqFittedParam, lsqCurveFitErr] = computeLatencies_offset(yt,ts*1e3,10);
- %% Plot results
- % multi-line plot offset
- ploffs0=0;%repmat((1:Ncs)'*max(yt(:)),1,Nts);
- ploffs=0;%repmat((1:Ncs)'*max(yt(:)),1,length(lsqCurveTimes));
- figure(); hold on
- plot(-inf,-inf,'b-'); plot(-inf,-inf,'k-'); plot(-inf,-inf,'r-');
- plot(-inf,-inf,'b*'); plot(-inf,-inf,'ko'); plot(-inf,-inf,'ro');
- plot(-inf,-inf,'bx'); plot(-inf,-inf,'kx'); plot(-inf,-inf,'rx');
- plot(ts*1e3,xt+ploffs0,'b-'); hold on;
- plot(ts*1e3,yt+ploffs0,'k-'); hold on;
- plot(lsqCurveTimes,lsqFittedCurve+ploffs,'r-'); hold on;
- plot(lt*1e3,ploffs0(:,1),'b*');
- plot(latTimes0,ploffs0(:,1),'ko');
- plot(latTimes,ploffs0(:,1),'ro');
- plot(mu,ploffs0(:,1),'bx');
- plot(lsqFittedParam0(:,1),ploffs0(:,1),'kx');
- plot(lsqFittedParam(:,1),ploffs0(:,1),'rx');
- xlabel('time (ms)');
- ylabel('channels'); %set(gca,'YTick',[]);
- hl=legend('x(t) sample signal','y(t)=x(t)+noise','f(t) fitted function','true latency',...
- 'latency fit(x)','latency fit(y)','true mean','mean fit(x)','mean fit(y)');
- hl.EdgeColor='none'; hl.Location='southeast';
|