rmeddis@34: function testANprob(targetFrequency,BFlist, levels, ...
rmeddis@29:     paramsName, paramChanges)
rmeddis@38: % testANprob generates rate/level functions for AN and brainstem units.
rmeddis@38: %  also other information like PSTHs, MOC efferent activity levels.
rmeddis@38: % A 'probability' model is used.
rmeddis@38: % e.g.
rmeddis@34: % testANprob(1000,1000, -10:10:80, 'Normal')
rmeddis@34: 
rmeddis@29: global IHC_VResp_VivoParams  IHC_cilia_RPParams IHCpreSynapseParams
rmeddis@29: global AN_IHCsynapseParams
rmeddis@32: global ANprobRateOutput dt ANtauCas
rmeddis@29: global ARattenuation MOCattenuation
rmeddis@29: 
rmeddis@29: AN_spikesOrProbability='probability';
rmeddis@29: 
rmeddis@29: dbstop if error
rmeddis@29: addpath (['..' filesep 'MAP'], ['..' filesep 'utilities'], ...
rmeddis@29:     ['..' filesep 'parameterStore'],  ['..' filesep 'wavFileStore'],...
rmeddis@29:     ['..' filesep 'testPrograms'])
rmeddis@29: 
rmeddis@35: if nargin<5, paramChanges=[]; end
rmeddis@35: if nargin<4, paramsName='Normal'; end
rmeddis@35: if nargin<3, levels=-10:10:80; end
rmeddis@35: if nargin==0, targetFrequency=1000; BFlist=1000; end
rmeddis@29: 
rmeddis@29: nLevels=length(levels);
rmeddis@29: 
rmeddis@29: toneDuration=.2;
rmeddis@29: rampDuration=0.002;
rmeddis@29: silenceDuration=.02;
rmeddis@29: localPSTHbinwidth=0.001;
rmeddis@29: 
rmeddis@29: % Use only the first frequency in the GUI targetFrequency box to defineBF
rmeddis@29: % targetFrequency=stimulusParameters.targetFrequency(1);
rmeddis@29: % BFlist=targetFrequency;
rmeddis@29: 
rmeddis@29: AN_HSRonset=zeros(nLevels,1);
rmeddis@29: AN_HSRsaturated=zeros(nLevels,1);
rmeddis@29: AN_LSRonset=zeros(nLevels,1);
rmeddis@29: AN_LSRsaturated=zeros(nLevels,1);
rmeddis@29: 
rmeddis@29: AR=zeros(nLevels,1);
rmeddis@29: MOC=zeros(nLevels,1);
rmeddis@29: 
rmeddis@29: figure(15), clf
rmeddis@34: set(gcf,'position',[980   356   401   321])
rmeddis@29: drawnow
rmeddis@29: 
rmeddis@29: %% guarantee that the sample rate is at least 10 times the frequency
rmeddis@29: sampleRate=50000;
rmeddis@29: while sampleRate< 10* targetFrequency
rmeddis@29:     sampleRate=sampleRate+10000;
rmeddis@29: end
rmeddis@29: 
rmeddis@29: %% adjust sample rate so that the pure tone stimulus has an integer
rmeddis@29: %% numver of epochs in a period
rmeddis@29: dt=1/sampleRate;
rmeddis@29: period=1/targetFrequency;
rmeddis@29: 
rmeddis@29: %% main computational loop (vary level)
rmeddis@29: levelNo=0;
rmeddis@29: for leveldB=levels
rmeddis@29:     levelNo=levelNo+1;
rmeddis@29: 
rmeddis@29:     fprintf('%4.0f\t', leveldB)
rmeddis@29:     amp=28e-6*10^(leveldB/20);
rmeddis@29: 
rmeddis@29:     time=dt:dt:toneDuration;
rmeddis@29:     rampTime=dt:dt:rampDuration;
rmeddis@29:     ramp=[0.5*(1+cos(2*pi*rampTime/(2*rampDuration)+pi)) ...
rmeddis@29:         ones(1,length(time)-length(rampTime))];
rmeddis@29:     ramp=ramp.*fliplr(ramp);
rmeddis@29: 
rmeddis@29:     silence=zeros(1,round(silenceDuration/dt));
rmeddis@29: 
rmeddis@29:     % create signal (leveldB/ targetFrequency)
rmeddis@29:     inputSignal=amp*sin(2*pi*targetFrequency'*time);
rmeddis@29:     inputSignal= ramp.*inputSignal;
rmeddis@29:     inputSignal=[silence inputSignal];
rmeddis@29: 
rmeddis@29:     %% run the model
rmeddis@29:     showPlotsAndDetails=0;
rmeddis@29: 
rmeddis@29: 
rmeddis@29:     MAP1_14(inputSignal, 1/dt, BFlist, ...
rmeddis@29:         paramsName, AN_spikesOrProbability, paramChanges);
rmeddis@29: 
rmeddis@32:     nTaus=length(ANtauCas);
rmeddis@29: 
rmeddis@29:     %LSR (same as HSR if no LSR fibers present)
rmeddis@29:     [nANFibers nTimePoints]=size(ANprobRateOutput);
rmeddis@29: 
rmeddis@29:     numLSRfibers=1;
rmeddis@29:     numHSRfibers=numLSRfibers;
rmeddis@29: 
rmeddis@29:     LSRspikes=ANprobRateOutput(1:numLSRfibers,:);
rmeddis@29:     PSTH=UTIL_PSTHmaker(LSRspikes, dt, localPSTHbinwidth);
rmeddis@29:     PSTHLSR=PSTH/(localPSTHbinwidth/dt);  % across fibers rates
rmeddis@29:     PSTHtime=localPSTHbinwidth:localPSTHbinwidth:...
rmeddis@29:         localPSTHbinwidth*length(PSTH);
rmeddis@29:     AN_LSRonset(levelNo)= max(PSTHLSR); % peak in 5 ms window
rmeddis@29:     AN_LSRsaturated(levelNo)= mean(PSTHLSR(round(length(PSTH)/2):end));
rmeddis@29: 
rmeddis@29:     % HSR
rmeddis@29:     HSRspikes= ANprobRateOutput(end- numHSRfibers+1:end, :);
rmeddis@29:     PSTH=UTIL_PSTHmaker(HSRspikes, dt, localPSTHbinwidth);
rmeddis@29:     PSTH=PSTH/(localPSTHbinwidth/dt); % sum across fibers (HSR only)
rmeddis@29:     AN_HSRonset(levelNo)= max(PSTH);
rmeddis@29:     AN_HSRsaturated(levelNo)= mean(PSTH(round(length(PSTH)/2): end));
rmeddis@29: 
rmeddis@29:     figure(15), subplot(2,2,4)
rmeddis@35:     hold off, bar(PSTHtime,PSTH, 'k')
rmeddis@29:     hold on,  bar(PSTHtime,PSTHLSR,'r')
rmeddis@29:     ylim([0 1000])
rmeddis@29:     xlim([0 length(PSTH)*localPSTHbinwidth])
rmeddis@29:     set(gcf,'name',[num2str(BFlist), ' Hz: ' num2str(leveldB) ' dB']);
rmeddis@29: 
rmeddis@29:     AR(levelNo)=min(ARattenuation);
rmeddis@29:     MOC(levelNo)=min(MOCattenuation(length(MOCattenuation)/2:end));
rmeddis@29: 
rmeddis@29: 
rmeddis@29:     figure(15), subplot(2,2,3)
rmeddis@34:     plot(20*log10(MOC), 'k'), hold on
rmeddis@34:     plot(20*log10(AR), 'r'),  hold off
rmeddis@34:     title(' MOC/AR'), ylabel('dB attenuation')
rmeddis@29:     ylim([-30 0])
rmeddis@29: 
rmeddis@34: end % level
rmeddis@29: 
rmeddis@29: figure(15), subplot(2,2,3)
rmeddis@34: plot(levels,20*log10(MOC), 'k'), hold on
rmeddis@34: plot(levels,20*log10(AR), 'r'),  hold off
rmeddis@34: title(' MOC/AR'), ylabel('dB attenuation')
rmeddis@29: ylim([-30 0])
rmeddis@29: xlim([0 max(levels)])
rmeddis@29: 
rmeddis@29: fprintf('\n')
rmeddis@29: toneDuration=2;
rmeddis@29: rampDuration=0.004;
rmeddis@29: silenceDuration=.02;
rmeddis@29: nRepeats=200;   % no. of AN fibers
rmeddis@29: fprintf('toneDuration  %6.3f\n', toneDuration)
rmeddis@29: fprintf(' %6.0f  AN fibers (repeats)\n', nRepeats)
rmeddis@29: fprintf('levels')
rmeddis@29: fprintf('%6.2f\t', levels)
rmeddis@29: fprintf('\n')
rmeddis@29: 
rmeddis@29: 
rmeddis@29: % ---------------------------------------------------- display parameters
rmeddis@29: 
rmeddis@29: 
rmeddis@29: nRows=2; nCols=2;
rmeddis@29: 
rmeddis@29: % AN rate - level ONSET functions
rmeddis@29: subplot(nRows,nCols,1)
rmeddis@29: plot(levels,AN_LSRonset,'ro'), hold on
rmeddis@38: plot(levels,AN_HSRonset,'ko', 'MarkerEdgeColor','k', 'markerFaceColor','k'), hold off
rmeddis@29: ylim([0 1000]),  xlim([min(levels) max(levels)])
rmeddis@29: ttl=['tauCa= ' num2str(IHCpreSynapseParams.tauCa)];
rmeddis@29: title( ttl)
rmeddis@29: xlabel('level dB SPL'), ylabel('peak rate (sp/s)'), grid on
rmeddis@34: text(0, 800, 'AN onset', 'fontsize', 14)
rmeddis@29: 
rmeddis@29: % AN rate - level ADAPTED function
rmeddis@29: subplot(nRows,nCols,2)
rmeddis@29: plot(levels,AN_LSRsaturated, 'ro'), hold on
rmeddis@38: plot(levels,AN_HSRsaturated, 'ko', 'MarkerEdgeColor','k', 'markerFaceColor','k'), hold off
rmeddis@35: maxYlim=340;
rmeddis@35: ylim([0 maxYlim])
rmeddis@29: set(gca,'ytick',0:50:300)
rmeddis@29: xlim([min(levels) max(levels)])
rmeddis@29: set(gca,'xtick',[levels(1):20:levels(end)])
rmeddis@29: %     grid on
rmeddis@29: ttl=[   'spont=' num2str(mean(AN_HSRsaturated(1,:)),'%4.0f')...
rmeddis@29:     '  sat=' num2str(mean(AN_HSRsaturated(end,1)),'%4.0f')];
rmeddis@29: title( ttl)
rmeddis@29: xlabel('level dB SPL'), ylabel ('adapted rate (sp/s)')
rmeddis@35: text(0, maxYlim-50, 'AN adapted', 'fontsize', 14), grid on
rmeddis@29: 
rmeddis@29: allData=[ levels'  AN_HSRonset AN_HSRsaturated...
rmeddis@29:     AN_LSRonset AN_LSRsaturated ];
rmeddis@29: fprintf('\n levels \tANHSR Onset \tANHSR adapted\tANLSR Onset \tANLSR adapted\tCNHSR\tCNLSR\tICHSR  \tICLSR \n');
rmeddis@29: UTIL_printTabTable(round(allData))
rmeddis@29: 
rmeddis@29: 
rmeddis@29: UTIL_showStruct(IHC_cilia_RPParams, 'IHC_cilia_RPParams')
rmeddis@29: UTIL_showStruct(IHCpreSynapseParams, 'IHCpreSynapseParams')
rmeddis@29: UTIL_showStruct(AN_IHCsynapseParams, 'AN_IHCsynapseParams')
rmeddis@29: 
rmeddis@29: 
rmeddis@29: allData=[ levels'  AN_HSRonset AN_HSRsaturated...
rmeddis@29:     AN_LSRonset AN_LSRsaturated ];
rmeddis@29: fprintf('\n levels \tANHSR Onset \tANHSR adapted\tANLSR Onset \tANLSR adapted\tCNHSR\tCNLSR\tICHSR  \tICLSR \n');
rmeddis@29: UTIL_printTabTable(round(allData))
rmeddis@29: