diff testPrograms/testRP2.m @ 38:c2204b18f4a2 tip

End nov big change

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/testPrograms/testRP2.m	Mon Nov 28 13:34:28 2011 +0000
@@ -0,0 +1,295 @@
+function testRP2(MAPparamsName,paramChanges)
+% testIHC evaluates IHC I/O function
+% multiple BF can be used but only one is easier to interpret.
+% e.g. testRP(1000,'Normal',{});
+
+global experiment method inputStimulusParams
+global stimulusParameters IHC_VResp_VivoParams IHC_cilia_RPParams
+savePath=path;
+addpath (['..' filesep 'utilities'],['..' filesep 'MAP'])
+dbstop if error
+
+figure(4), clf,
+set (gcf, 'name', ['IHC'])
+set(gcf,'position',[613   354   360   322])
+drawColors='rgbkmcy';
+drawnow
+
+if nargin<3
+    paramChanges=[];
+end
+if nargin<2
+    MAPparamsName='Normal';
+end
+if nargin<3
+    BF=800;
+end
+
+levels=-20:10:100;
+nLevels=length(levels);
+toneDuration=.05;
+silenceDuration=.01;
+sampleRate=50000;
+dt=1/sampleRate;
+
+allIHC_RP_peak=[];
+allIHC_RP_dc=[];
+
+%% Ruggero
+%%Ruggero data
+RuggeroData=[
+    0	2.00E-10;
+    10	5.00E-10;
+    20	1.50E-09;
+    30	2.50E-09;
+    40	5.30E-09;
+    50	1.00E-08;
+    60	1.70E-08;
+    70	2.50E-08;
+    80	4.00E-08;
+    90	6.00E-08;
+    100	1.50E-07;
+    110	3.00E-07;
+    ];
+
+
+BF=10000;
+targetFrequency=BF;
+
+IHC_RP_peak=zeros(nLevels,1);
+IHC_RP_min=zeros(nLevels,1);
+IHC_RP_dc=zeros(nLevels,1);
+
+time=dt:dt:toneDuration;
+
+rampDuration=0.004;
+rampTime=dt:dt:rampDuration;
+ramp=[0.5*(1+cos(2*pi*rampTime/(2*rampDuration)+pi)) ...
+    ones(1,length(time)-length(rampTime))];
+ramp=ramp.*fliplr(ramp);
+
+silence=zeros(1,round(silenceDuration/dt));
+
+toneStartptr=length(silence)+1;
+toneMidptr=toneStartptr+round(toneDuration/(2*dt)) -1;
+toneEndptr=toneStartptr+round(toneDuration/dt) -1;
+
+levelNo=0;
+for leveldB=levels
+    levelNo=levelNo+1;
+    % replicate at all levels
+    amp=28e-6*10^(leveldB/20);
+    
+    %% create signal (leveldB/ frequency)
+    inputSignal=amp*sin(2*pi*targetFrequency'*time);
+    inputSignal= ramp.*inputSignal;
+    inputSignal=[silence inputSignal silence];
+    inputStimulusParams.sampleRate=1/dt;
+    %         global IHC_ciliaParams
+    
+    %% disable efferent for fast processing
+    
+    %% run the model
+    global  DRNLoutput IHC_cilia_output IHCrestingCiliaCond IHCrestingV...
+        IHCoutput
+    AN_spikesOrProbability='probability';
+    
+    MAP1_14(inputSignal, sampleRate, BF, ...
+        MAPparamsName, AN_spikesOrProbability, paramChanges);
+    
+    % DRNL
+    DRNLoutput=DRNLoutput;
+    DRNL_peak(levelNo,1)=max(DRNLoutput(toneMidptr:toneEndptr));
+    DRNL_min(levelNo,1)=min(DRNLoutput(toneMidptr:toneEndptr));
+    DRNL_dc(levelNo,1)=mean(DRNLoutput(toneMidptr:toneEndptr));
+    
+end
+%% plot DRNL
+subplot(2,2,1)
+referenceDisp=10e-9;
+semilogy(levels,DRNL_peak, 'linewidth',2), hold on
+semilogy(RuggeroData(:,1),RuggeroData(:,2),'o')
+title(['BM: Ruggero ' num2str(BF) ' Hz'])
+ylabel ('displacement(m)'), xlabel('dB SPL')
+xlim([min(levels) max(levels)]), ylim([1e-10 1e-7])
+grid on
+
+
+%% Dallos
+BF=800;
+targetFrequency=BF;
+
+IHC_RP_peak=zeros(nLevels,1);
+IHC_RP_min=zeros(nLevels,1);
+IHC_RP_dc=zeros(nLevels,1);
+
+time=dt:dt:toneDuration;
+
+rampDuration=0.004;
+rampTime=dt:dt:rampDuration;
+ramp=[0.5*(1+cos(2*pi*rampTime/(2*rampDuration)+pi)) ...
+    ones(1,length(time)-length(rampTime))];
+ramp=ramp.*fliplr(ramp);
+
+silence=zeros(1,round(silenceDuration/dt));
+
+toneStartptr=length(silence)+1;
+toneMidptr=toneStartptr+round(toneDuration/(2*dt)) -1;
+toneEndptr=toneStartptr+round(toneDuration/dt) -1;
+
+levelNo=0;
+for leveldB=levels
+    levelNo=levelNo+1;
+    % replicate at all levels
+    amp=28e-6*10^(leveldB/20);
+    
+    %% create signal (leveldB/ frequency)
+    inputSignal=amp*sin(2*pi*targetFrequency'*time);
+    inputSignal= ramp.*inputSignal;
+    inputSignal=[silence inputSignal silence];
+    inputStimulusParams.sampleRate=1/dt;
+    %         global IHC_ciliaParams
+    
+    %% disable efferent for fast processing
+    method.DRNLSave=1;
+    method.IHC_cilia_RPSave=1;
+    method.IHCpreSynapseSave=1;
+    method.IHC_cilia_RPSave=1;
+    method.segmentDuration=-1;
+    moduleSequence=1:4;
+    
+    %% run the model
+    global  DRNLoutput IHC_cilia_output IHCrestingCiliaCond IHCrestingV...
+        IHCoutput
+    AN_spikesOrProbability='probability';
+    
+    MAP1_14(inputSignal, sampleRate, BF, ...
+        MAPparamsName, AN_spikesOrProbability, paramChanges);
+    
+    % DRNL
+    DRNLoutput=DRNLoutput;
+    DRNL_peak(levelNo,1)=max(DRNLoutput(toneMidptr:toneEndptr));
+    DRNL_min(levelNo,1)=min(DRNLoutput(toneMidptr:toneEndptr));
+    DRNL_dc(levelNo,1)=mean(DRNLoutput(toneMidptr:toneEndptr));
+    
+    % cilia
+    IHC_ciliaData=IHC_cilia_output;
+    IHC_ciliaData=IHC_ciliaData;
+    IHC_cilia_peak(levelNo,1)=...
+        max(IHC_ciliaData(toneMidptr:toneEndptr));
+    IHC_cilia_min(levelNo,1)=...
+        min(IHC_ciliaData(toneMidptr:toneEndptr));
+    IHC_cilia_dc(levelNo,1)=...
+        mean(IHC_ciliaData(toneMidptr:toneEndptr));
+    
+    % RP
+    IHC_RPData=IHCoutput;
+    IHC_RPData=IHC_RPData;
+    IHC_RP_peak(levelNo,1)=...
+        max(IHC_RPData(toneMidptr:toneEndptr));
+    IHC_RP_min(levelNo,1)=...
+        min(IHC_RPData(toneMidptr:toneEndptr));
+    IHC_RP_dc(levelNo,1)=...
+        mean(IHC_RPData(toneMidptr:toneEndptr));
+    
+end
+
+%% Dallos and Harris data
+%% plot receptor potentials
+figure(4)
+subplot(2,2,3)
+% RP I/O function min and max
+restingRP=IHC_RP_peak(1);
+toPlot= [fliplr(IHC_RP_min(:,1)') IHC_RP_peak(:,1)'];
+microPa=   28e-6*10.^(levels/20);
+microPa=[-fliplr(microPa) microPa];
+plot(microPa,toPlot, 'linewidth',2)
+
+dallosx=[-0.9	-0.1	-0.001	0.001	0.01	0.9];
+dallosy=[-8	-7.8	-6.5	11	16.5	22]/1000 + restingRP;
+subplot(2,2,3)
+hold on, plot(dallosx,dallosy, 'o')
+plot([-1 1], [restingRP restingRP], 'r')
+title(' Dallos(86) data at 800 Hz')
+ylabel ('receptor potential(V)'), xlabel('Pa')
+ylim([-0.08 -0.02]), xlim([-1 1])
+grid on
+
+
+%% Patuzzi
+BF=7000;
+targetFrequency=BF;
+
+IHC_RP_peak=zeros(nLevels,1);
+IHC_RP_min=zeros(nLevels,1);
+IHC_RP_dc=zeros(nLevels,1);
+
+time=dt:dt:toneDuration;
+
+rampDuration=0.004;
+rampTime=dt:dt:rampDuration;
+ramp=[0.5*(1+cos(2*pi*rampTime/(2*rampDuration)+pi)) ...
+    ones(1,length(time)-length(rampTime))];
+ramp=ramp.*fliplr(ramp);
+
+silence=zeros(1,round(silenceDuration/dt));
+
+toneStartptr=length(silence)+1;
+toneMidptr=toneStartptr+round(toneDuration/(2*dt)) -1;
+toneEndptr=toneStartptr+round(toneDuration/dt) -1;
+
+levelNo=0;
+for leveldB=levels
+    levelNo=levelNo+1;
+    % replicate at all levels
+    amp=28e-6*10^(leveldB/20);
+    
+    %% create signal (leveldB/ frequency)
+    inputSignal=amp*sin(2*pi*targetFrequency'*time);
+    inputSignal= ramp.*inputSignal;
+    inputSignal=[silence inputSignal silence];
+    inputStimulusParams.sampleRate=1/dt;
+    %         global IHC_ciliaParams
+    
+    %% run the model
+    global  DRNLoutput IHC_cilia_output IHCrestingCiliaCond IHCrestingV...
+        IHCoutput
+    AN_spikesOrProbability='probability';
+    
+    MAP1_14(inputSignal, sampleRate, BF, ...
+        MAPparamsName, AN_spikesOrProbability, paramChanges);
+    
+    % RP
+    IHC_RPData=IHCoutput;
+    IHC_RP_peak(levelNo,1)=...
+        max(IHC_RPData(toneMidptr:toneEndptr));
+    IHC_RP_min(levelNo,1)=...
+        min(IHC_RPData(toneMidptr:toneEndptr));
+    IHC_RP_dc(levelNo,1)=...
+        mean(IHC_RPData(toneMidptr:toneEndptr));
+end
+
+
+%% RP I/O function min and max
+figure(4)
+subplot(2,2,4)
+restingRP=IHC_RP_peak(1);
+peakRP=max(IHC_RP_peak);
+plot(levels, IHC_RP_peak, 'linewidth',2)
+hold on
+plot(levels, IHC_RP_dc, ':', 'linewidth',2)
+plot([min(levels) max(levels)], [restingRP restingRP], 'r')
+xlim([min(levels) max(levels)])
+% animal data
+sndLevel=[5	15	25	35	45	55	65	75];
+RPanimal=restingRP+[0.5	2	4.6	5.8	6.4	7.2	8	10.2]/1000;
+% could be misleading when restingRP changes
+RPanimal=-0.060+[0.5	2	4.6	5.8	6.4	7.2	8	10.2]/1000;
+hold on, plot(sndLevel,RPanimal,'o')
+
+grid on
+title(' 7 kHz Patuzzi')
+ylabel ('RP(V) peak and DC'), xlabel('dB SPL')
+ylim([-0.07 -0.04])
+path(savePath);
+disp(paramChanges)