view utilities/UTIL_CAPgenerator.m @ 38:c2204b18f4a2 tip

End nov big change

function wholeNerveCAP  = UTIL_CAPgenerator...
    (ANresponse, dt, BFlist, numANfibers, plotCAP)
% 
% Generates a compound action potential by convolving an impulse repsonse,
% as defined by Mark Chertoff (2004, JASA), with the response of the 
% auditory nerve.
%
%
% -e(-k*time)*SIN(2*PI()*f*time)
%
% mu(t) = e^-kt * sin(omega*t)
% omega = 2 * pi * f
%
%
% Robert T. Ferry
% 01st May 2008
%
%

nChannels=length(BFlist);
[r nSpikeEpochs]=size(ANresponse);

wholeNerveCAP       = [];
channelCAP          = [];
e                   = exp(1);
k                   = 1000;
impulseDuration     = 0.01;
impulseFrequency    = 750;
impulseTime         = dt:dt:impulseDuration;
impulseResponse     = -e.^(-k*impulseTime).*sin(2*pi*impulseFrequency*impulseTime);
impulseResponse=impulseResponse-mean(impulseResponse);

% WholeNerveCAP
ANoutput = sum(ANresponse, 1);
convolvedWholeNerveCAP = conv(ANoutput, impulseResponse(1,:));
% truncate
convolvedWholeNerveCAP=convolvedWholeNerveCAP(1:nSpikeEpochs);

% apply measurement time constant
sampleRate=1/dt;
upperFreq=sampleRate/4;
lowPassCutOff=40;
wholeNerveCAP = UTIL_Butterworth(convolvedWholeNerveCAP, dt, lowPassCutOff, upperFreq, 1);
% or do not do this
% wholeNerveCAP = convolvedWholeNerveCAP;

% Plot output?

CAPtime = dt:dt:dt*length(wholeNerveCAP);

if (plotCAP == 1)
    figure(9)

    subplot(3,1,1), plot(impulseTime, impulseResponse)
    title('Impulse response')
    xlabel('Time (s)'), ylabel('Amplitude (Pa)')
    xlim([0 max(impulseTime)]), ylim([-inf inf])
	
    subplot(3,1,2), plot(CAPtime, wholeNerveCAP)
    title(['AN CAP (whole-nerve)  '   num2str(length(BFlist)) ' channels' num2str(numANfibers) ' fibers/ch'])
    xlabel('Time (s)'), ylabel('Amplitude (Pa)')
    xlim([0 max(CAPtime)]), ylim([-inf inf])   
	
end