rmeddis@23: function UTIL_showMAP (options) rmeddis@24: % UTIL_showMAP produces summaries of the output from MAP's mmost recent run rmeddis@24: % All MAP outputs are stored in global variables and UTIL_showMAP rmeddis@24: % simply assumes that they are in place. rmeddis@24: % rmeddis@23: % options rmeddis@24: % options.printModelParameters=1; % print model parameters rmeddis@24: % options.showModelOutput=1; % plot all stages output rmeddis@24: % options.printFiringRates=1; % mean activity at all stages rmeddis@24: % options.showACF=1; % SACF (probabilities only) rmeddis@24: % options.showEfferent=1; % plot of efferent activity rmeddis@24: % options.surfProbability=0; % HSR (probability) surf plot rmeddis@24: % options.fileName=[]; % parameter filename for plot title rmeddis@23: rmeddis@23: dbstop if warning rmeddis@23: rmeddis@23: global dt ANdt saveAN_spikesOrProbability savedBFlist saveMAPparamsName... rmeddis@23: savedInputSignal TMoutput OMEoutput ARattenuation ... rmeddis@23: DRNLoutput IHC_cilia_output IHCrestingCiliaCond IHCrestingV... rmeddis@23: IHCoutput ANprobRateOutput ANoutput savePavailable tauCas ... rmeddis@23: CNoutput ICoutput ICmembraneOutput ICfiberTypeRates MOCattenuation rmeddis@23: global OMEParams DRNLParams IHC_cilia_RPParams IHCpreSynapseParams rmeddis@23: global AN_IHCsynapseParams MacGregorParams MacGregorMultiParams rmeddis@23: rmeddis@23: rmeddis@23: restorePath=path; rmeddis@23: addpath ( ['..' filesep 'utilities'], ['..' filesep 'parameterStore']) rmeddis@23: rmeddis@23: if nargin<1 rmeddis@23: options=[]; rmeddis@23: end rmeddis@23: % defaults (plot staged outputs and print rates only) rmeddis@24: if ~isfield(options,'printModelParameters') rmeddis@24: options.printModelParameters=0; end rmeddis@23: if ~isfield(options,'showModelOutput'),options.showModelOutput=1;end rmeddis@23: if ~isfield(options,'printFiringRates'),options.printFiringRates=1;end rmeddis@23: if ~isfield(options,'showACF'),options.showACF=0;end rmeddis@23: if ~isfield(options,'showEfferent'),options.showEfferent=0;end rmeddis@23: if ~isfield(options,'surfProbability'),options.surfProbability=0;end rmeddis@23: if ~isfield(options,'fileName'),options.fileName=[];end rmeddis@23: rmeddis@23: rmeddis@23: if options.printModelParameters rmeddis@23: % Read parameters from MAPparams<***> file in 'parameterStore' folder rmeddis@23: % and print out all parameters rmeddis@23: cmd=['MAPparams' saveMAPparamsName ... rmeddis@23: '(-1, 1/dt, 1);']; rmeddis@23: eval(cmd); rmeddis@23: end rmeddis@23: rmeddis@23: if options.printFiringRates rmeddis@23: %% print summary firing rates rmeddis@23: fprintf('\n\n') rmeddis@23: disp('summary') rmeddis@23: disp(['AR: ' num2str(min(ARattenuation))]) rmeddis@23: disp(['MOC: ' num2str(min(min(MOCattenuation)))]) rmeddis@23: nANfiberTypes=length(tauCas); rmeddis@23: if strcmp(saveAN_spikesOrProbability, 'spikes') rmeddis@23: nANfibers=size(ANoutput,1); rmeddis@23: nHSRfibers=nANfibers/nANfiberTypes; rmeddis@23: duration=size(TMoutput,2)*dt; rmeddis@23: disp(['AN: ' num2str(sum(sum(ANoutput(end-nHSRfibers+1:end,:)))/... rmeddis@23: (nHSRfibers*duration))]) rmeddis@23: rmeddis@23: nCNneurons=size(CNoutput,1); rmeddis@23: nHSRCNneuronss=nCNneurons/nANfiberTypes; rmeddis@23: disp(['CN: ' num2str(sum(sum(CNoutput(end-nHSRCNneuronss+1:end,:)))... rmeddis@23: /(nHSRCNneuronss*duration))]) rmeddis@23: disp(['IC: ' num2str(sum(sum(ICoutput))/duration)]) rmeddis@23: % disp(['IC by type: ' num2str(mean(ICfiberTypeRates,2)')]) rmeddis@23: else rmeddis@23: disp(['AN: ' num2str(mean(mean(ANprobRateOutput)))]) rmeddis@23: [PSTH pointsPerBin]= UTIL_makePSTH(ANprobRateOutput, dt, 0.001); rmeddis@23: disp(['max max AN: ' num2str(max(max(... rmeddis@23: PSTH/pointsPerBin )))]) rmeddis@23: end rmeddis@23: end rmeddis@23: rmeddis@23: rmeddis@23: %% figure (99) summarises main model output rmeddis@23: if options.showModelOutput rmeddis@23: plotInstructions.figureNo=99; rmeddis@23: signalRMS=mean(savedInputSignal.^2)^0.5; rmeddis@23: signalRMSdb=20*log10(signalRMS/20e-6); rmeddis@23: rmeddis@23: % plot signal (1) rmeddis@23: plotInstructions.displaydt=dt; rmeddis@23: plotInstructions.numPlots=6; rmeddis@23: plotInstructions.subPlotNo=1; rmeddis@23: plotInstructions.title=... rmeddis@23: ['stimulus: ' num2str(signalRMSdb, '%4.0f') ' dB SPL']; rmeddis@23: r=size(savedInputSignal,1); rmeddis@23: if r==1, savedInputSignal=savedInputSignal'; end rmeddis@23: UTIL_plotMatrix(savedInputSignal', plotInstructions); rmeddis@23: rmeddis@23: % stapes (2) rmeddis@23: plotInstructions.subPlotNo=2; rmeddis@23: plotInstructions.title= ['stapes displacement']; rmeddis@23: UTIL_plotMatrix(OMEoutput, plotInstructions); rmeddis@23: rmeddis@23: % DRNL (3) rmeddis@23: plotInstructions.subPlotNo=3; rmeddis@23: plotInstructions.title= ['BM displacement']; rmeddis@23: plotInstructions.yValues= savedBFlist; rmeddis@23: UTIL_plotMatrix(DRNLoutput, plotInstructions); rmeddis@23: rmeddis@23: switch saveAN_spikesOrProbability rmeddis@23: case 'spikes' rmeddis@23: % AN (4) rmeddis@23: plotInstructions.displaydt=ANdt; rmeddis@23: plotInstructions.title='AN'; rmeddis@23: plotInstructions.subPlotNo=4; rmeddis@23: plotInstructions.yLabel='BF'; rmeddis@23: plotInstructions.yValues= savedBFlist; rmeddis@23: plotInstructions.rasterDotSize=1; rmeddis@23: if length(tauCas)==2 rmeddis@23: plotInstructions.plotDivider=1; rmeddis@23: else rmeddis@23: plotInstructions.plotDivider=0; rmeddis@23: end rmeddis@23: if sum(sum(ANoutput))<100 rmeddis@23: plotInstructions.rasterDotSize=3; rmeddis@23: end rmeddis@23: UTIL_plotMatrix(ANoutput, plotInstructions); rmeddis@23: rmeddis@23: % CN (5) rmeddis@23: plotInstructions.displaydt=ANdt; rmeddis@23: plotInstructions.subPlotNo=5; rmeddis@23: plotInstructions.title='CN spikes'; rmeddis@23: if sum(sum(CNoutput))<100 rmeddis@23: plotInstructions.rasterDotSize=3; rmeddis@23: end rmeddis@23: UTIL_plotMatrix(CNoutput, plotInstructions); rmeddis@23: rmeddis@23: % IC (6) rmeddis@23: plotInstructions.displaydt=ANdt; rmeddis@23: plotInstructions.subPlotNo=6; rmeddis@23: plotInstructions.title='IC'; rmeddis@23: if size(ICoutput,1)>3 rmeddis@23: if sum(sum(ICoutput))<100 rmeddis@23: plotInstructions.rasterDotSize=3; rmeddis@23: end rmeddis@23: UTIL_plotMatrix(ICoutput, plotInstructions); rmeddis@23: else rmeddis@23: plotInstructions.title='IC (HSR) membrane potential'; rmeddis@23: plotInstructions.displaydt=dt; rmeddis@23: plotInstructions.yLabel='V'; rmeddis@23: plotInstructions.zValuesRange= [-.1 0]; rmeddis@23: UTIL_plotMatrix(ICmembraneOutput, plotInstructions); rmeddis@23: end rmeddis@23: rmeddis@23: otherwise % probability (4-6) rmeddis@23: plotInstructions.displaydt=dt; rmeddis@23: plotInstructions.numPlots=2; rmeddis@23: plotInstructions.subPlotNo=2; rmeddis@23: plotInstructions.yLabel='BF'; rmeddis@23: if nANfiberTypes>1, rmeddis@23: plotInstructions.yLabel='LSR HSR'; rmeddis@23: plotInstructions.plotDivider=1; rmeddis@23: end rmeddis@23: plotInstructions.title='AN - spike probability'; rmeddis@23: UTIL_plotMatrix(ANprobRateOutput, plotInstructions); rmeddis@23: end rmeddis@23: end rmeddis@23: rmeddis@23: if options.surfProbability rmeddis@23: %% surface plot of probability rmeddis@23: figure(97), clf rmeddis@23: % select only HSR fibers at the bottom of the matrix rmeddis@23: ANprobRateOutput= ANprobRateOutput(end-length(savedBFlist)+1:end,:); rmeddis@23: PSTHbinWidth=0.001; rmeddis@23: PSTH=UTIL_PSTHmakerb(ANprobRateOutput, ANdt, PSTHbinWidth); rmeddis@23: [nY nX]=size(PSTH); rmeddis@23: time=PSTHbinWidth*(1:nX); rmeddis@23: surf(time, savedBFlist, PSTH) rmeddis@23: shading interp rmeddis@23: set(gca, 'yScale','log') rmeddis@23: xlim([0 max(time)]) rmeddis@23: ylim([0 max(savedBFlist)]) rmeddis@23: zlim([0 1000]) rmeddis@23: xlabel('time (s)') rmeddis@23: ylabel('best frequency (Hz)') rmeddis@23: zlabel('spike rate') rmeddis@23: view([-20 60]) rmeddis@23: % view([0 90]) rmeddis@23: title ([options.fileName ': ' num2str(signalRMSdb,'% 3.0f') ' dB']) rmeddis@23: end rmeddis@23: rmeddis@23: rmeddis@23: %% plot efferent control values as dB rmeddis@23: if options.showEfferent rmeddis@23: plotInstructions=[]; rmeddis@23: plotInstructions.figureNo=98; rmeddis@23: figure(98), clf rmeddis@23: plotInstructions.displaydt=dt; rmeddis@23: plotInstructions.numPlots=2; rmeddis@23: plotInstructions.subPlotNo=1; rmeddis@23: plotInstructions.zValuesRange=[ -25 0]; rmeddis@23: plotInstructions.title= ['AR strength. Signal level= ' ... rmeddis@23: num2str(signalRMSdb,'%4.0f') ' dB SPL']; rmeddis@23: UTIL_plotMatrix(20*log10(ARattenuation), plotInstructions); rmeddis@23: rmeddis@23: plotInstructions.subPlotNo=2; rmeddis@23: plotInstructions.yValues= savedBFlist; rmeddis@23: plotInstructions.yLabel= 'BF'; rmeddis@23: plotInstructions.title= ['MOC strength']; rmeddis@23: plotInstructions.zValuesRange=[ -25 0]; rmeddis@23: subplot(2,1,2) rmeddis@23: % imagesc(MOCattenuation) rmeddis@23: UTIL_plotMatrix(20*log10(MOCattenuation), plotInstructions); rmeddis@23: colorbar rmeddis@23: end rmeddis@23: rmeddis@23: %% ACF plot if required rmeddis@23: if options.showACF rmeddis@23: tic rmeddis@23: method.dt=dt; rmeddis@23: method.segmentNo=1; rmeddis@23: method.nonlinCF=savedBFlist; rmeddis@23: rmeddis@23: minPitch= 80; maxPitch= 4000; numPitches=100; % specify lags rmeddis@23: pitches=10.^ linspace(log10(minPitch), log10(maxPitch),numPitches); rmeddis@23: pitches=fliplr(pitches); rmeddis@23: filteredSACFParams.lags=1./pitches; % autocorrelation lags vector rmeddis@23: filteredSACFParams.acfTau= .003; % time constant of running ACF rmeddis@23: filteredSACFParams.lambda= 0.12; % slower filter to smooth ACF rmeddis@23: filteredSACFParams.lambda= 0.01; % slower filter to smooth ACF rmeddis@23: rmeddis@23: filteredSACFParams.plotACFs=0; % special plot (see code) rmeddis@23: filteredSACFParams.plotFilteredSACF=0; % 0 plots unfiltered ACFs rmeddis@23: filteredSACFParams.plotMoviePauses=.3; % special plot (see code) rmeddis@23: rmeddis@23: filteredSACFParams.usePressnitzer=0; % attenuates ACF at long lags rmeddis@23: filteredSACFParams.lagsProcedure= 'useAllLags'; rmeddis@23: % filteredSACFParams.lagsProcedure= 'useBernsteinLagWeights'; rmeddis@23: % filteredSACFParams.lagsProcedure= 'omitShortLags'; rmeddis@23: filteredSACFParams.criterionForOmittingLags=3; rmeddis@23: filteredSACFParams.plotACFsInterval=200; rmeddis@23: rmeddis@23: if filteredSACFParams.plotACFs rmeddis@23: % plot original waveform on ACF plot rmeddis@23: figure(13), clf rmeddis@23: subplot(4,1,1) rmeddis@23: t=dt*(1:length(savedInputSignal)); rmeddis@23: plot(t,savedInputSignal) rmeddis@23: xlim([0 t(end)]) rmeddis@23: title(['stimulus: ' num2str(signalRMSdb, '%4.0f') ' dB SPL']); rmeddis@23: end rmeddis@23: rmeddis@23: % plot original waveform on summary/smoothed ACF plot rmeddis@23: figure(96), clf rmeddis@23: subplot(2,1,1) rmeddis@23: t=dt*(1:length(savedInputSignal)); rmeddis@23: plot(t,savedInputSignal) rmeddis@23: xlim([0 t(end)]) rmeddis@23: title(['stimulus: ' num2str(signalRMSdb, '%4.0f') ' dB SPL']); rmeddis@23: rmeddis@23: rmeddis@23: % compute ACF rmeddis@23: switch saveAN_spikesOrProbability rmeddis@23: case 'probability' rmeddis@23: inputToACF=ANprobRateOutput.^0.5; rmeddis@23: otherwise rmeddis@23: inputToACF=ANoutput; rmeddis@23: end rmeddis@23: rmeddis@23: disp ('computing ACF...') rmeddis@23: [P, BFlist, sacf, boundaryValue] = ... rmeddis@23: filteredSACF(inputToACF, method, filteredSACFParams); rmeddis@23: disp(' ACF done.') rmeddis@23: rmeddis@23: % SACF rmeddis@23: subplot(2,1,2) rmeddis@23: imagesc(P) rmeddis@23: ylabel('periodicities (Hz)') rmeddis@23: xlabel('time (s)') rmeddis@23: title(['running smoothed (root) SACF. ' saveAN_spikesOrProbability ' input']) rmeddis@23: pt=[1 get(gca,'ytick')]; % force top xtick to show rmeddis@23: set(gca,'ytick',pt) rmeddis@23: set(gca,'ytickLabel', round(pitches(pt))) rmeddis@23: tt=get(gca,'xtick'); rmeddis@23: set(gca,'xtickLabel', round(100*t(tt))/100) rmeddis@23: end rmeddis@23: rmeddis@23: path(restorePath)