sigmaPre= (.1:.1:10);
sigma=sigmaPre.*sigmaPre;

% sigma=.01;
nu= -1./sigma;
w=nu.*sqrt(sigma/2);

w=15;
sigma=2*w^2;
Qofw = sqrt(pi)* exp(w.*w) .* erfc(w);
OnePlusCVsq  = .5 * Qofw.* ((1+2*w.*w).*Qofw -2*w) ./((1-w.*Qofw).*(1-w.*Qofw));
CV= sqrt(OnePlusCVsq-1);
mean=sqrt(2*sigma).*(1-w.*Qofw)./Qofw;
std= CV.*mean;
[mean std CV w]



sigmaPre= (.4:.1:10);
sigma=sigmaPre.*sigmaPre;
% sigma=.01;
nu= -1./sigma;
w=nu.*sqrt(sigma/2);
Qofw = sqrt(pi)* exp(w.*w) .* erfc(w);
OnePlusCVsq  = .5 * Qofw.* ((1+2*w.*w).*Qofw -2*w) ./((1-w.*Qofw).*(1-w.*Qofw));
CV= sqrt(OnePlusCVsq-1);
mean=sqrt(2*sigma).*(1-w.*Qofw)./Qofw;
std= CV.*mean;
[mean std CV w]


% Find nu Given sigma
% sigmaPre= (.4:.1:10);
% sigma=sigmaPre.*sigmaPre;
% sigma=.01;


clear meanVec;

sigmaVec=[.005 .02 .04 .08 .2 .4 .8 1 1.2 1.6 2 2.4 3 4 6 10 30 100 500 1000]


for isigma=1:(length(sigmaVec)-5)
    sigma=sigmaVec(isigma);
    Count=0;
    meanOld=0;
    for nu= -500.05:.01:100.05;
        Count=Count+1;
        w=nu.*sqrt(sigma/2);
        Qofw = sqrt(pi)* exp(w.*w) .* erfc(w);
        OnePlusCVsq  = .5 * Qofw.* ((1+2*w.*w).*Qofw -2*w) ./((1-w.*Qofw).*(1-w.*Qofw));
        CV= sqrt(OnePlusCVsq-1);
        mean=sqrt(2*sigma).*(1-w.*Qofw)./Qofw;
        std= CV.*mean;
        if (Count>1 & meanOld>1 & mean<=1), break; end
        meanOld=mean;
    end
    nuVec(isigma)=nu;
    meanVec(isigma)=mean;
    stdVec(isigma)=std;

    wVec(isigma)  = sqrt(2./sigma)*nu;
    cvVec(isigma) = std/mean;
    QofwVec(isigma) = sqrt(pi)* exp(w*w) .* erfc(w);
    logNVec(isigma) =  log(Qofw*sqrt(sigma/2));
 
    Entropy(isigma) = (Qofw+2*w-2*w*w*Qofw)/(2*Qofw) + logNVec(isigma);
         
end
 
 
 [nuVec'  sigmaVec' meanVec' stdVec' wVec' cvVec'  Entropy']
     nu        sigma     mean      std      ignore          cv

  -49.9900    0.0200    0.9998    0.1414   -7.0697    0.1414
  -24.9900    0.0400    0.9996    0.2000   -3.5341    0.2001
  -12.4900    0.0800    0.9994    0.2825   -1.7664    0.2826
   -4.9100    0.2000    0.9982    0.4287   -0.6944    0.4295
   -2.2500    0.4000    0.9994    0.5484   -0.3182    0.5487
   -0.7900    0.8000    0.9977    0.6596   -0.1117    0.6611
   -0.4800    1.0000    0.9995    0.6934   -0.0679    0.6937
   -0.2600    1.2000    0.9976    0.7184   -0.0368    0.7201
    0.0200    1.6000    0.9977    0.7567    0.0028    0.7585
    0.2000    2.0000    0.9947    0.7827    0.0283    0.7868
    0.3200    2.4000    0.9949    0.8038    0.0453    0.8079
    0.4400    3.0000    0.9976    0.8293    0.0622    0.8313
    0.5700    4.0000    0.9967    0.8567    0.0806    0.8595
    0.7100    6.0000    0.9925    0.8872    0.1004    0.8940
    0.8200   10.0000    0.9950    0.9225    0.1160    0.9271
    0.9400   30.0000    0.9956    0.9666    0.1329    0.9709
    0.9900  100.0000    0.9905    0.9811    0.1400    0.9906

    wVec
   -2.2356
   -1.1176
   -0.5586
   -0.2196
   -0.1006
   -0.0353
   -0.0215
   -0.0116
    0.0009
    0.0089
    0.0143
    0.0197
    0.0255
    0.0318
    0.0367
    0.0420
    0.0443
    0.0447
    0.0447

>> 

saveas(gcf,'Stim06Vec100.eps','eps')



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% ReadPlayBack.dat is Figure 1a
% First Figure


%  col 1             col 2         col4         col 5        ,     col 6
% time for observed ,  actual oberved   stim times  time for fit       fit

vf= figure;  set(vf,'name','Fig. 1a BM');
    plot(ReadPlayBack(2:736,1), ReadPlayBack(2:736,2) ,'LineWidth',2,'color',[0 0 .627]) ,hold
    plot(ReadPlayBack(2:2598,5),ReadPlayBack(2:2598,6),'LineWidth',2,'color',[1 0 1] )
    axis([0   75  -.4 2.3])
    h= gca;
    set(h,'FontName','Arial','FontSize',12,'FontWeight','bold','box','off')
 
    vw= legend(h,'measured','fit')
    legend('boxoff');
    set(vw,'Position',[.55 .764 .39 .114],'box','off')
   
for jj=1:20
    plot(ReadPlayBack(2:53,4),(-.3 +jj*.0025)  *ones(52,1),'.','MarkerSize',4);
end

xlabel('time(secs)','FontSize',16,'FontName','Arial','FontWeight','Bold')
ylabel('[DA] \muM','FontSize',16,'FontName','Arial','FontWeight','Bold')
scribe.textbox('String','Fig 1a','LineStyle','none','FontName','Arial', ...
       'Position',[.002,.008,.1,0.073],'FontWeight','Bold');

   hold
 delete('Fig1aBM.fig')
 saveas(gcf,'Fig1aBM.fig','fig')
 saveas(gcf,'Fig1aBM.eps','eps')
 saveas(gcf,'Fig1aBM.bmp','bmp')
 saveas(gcf,'Fig1aBM.tif','tif')
 saveas(gcf,'Fig1aBM.jpg','jpg')
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Saturday February 10th 2007  ;   Invert for the electrode

figure
fit = ReadPlayBack(2:2598,6);
fitFFT = fft(fit);
plot(abs(fitFFT),'.');
% Remove the electrode

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% OneHertzRegular.dat is Figure 1b

vf= figure;  set(vf,'name','Fig. 1b BM');
    plot(OneHertzRegular(1:2694,1), OneHertzRegular(1:2694,2) ,'LineWidth',2,'color',[0 0 .627]) ,hold
    plot(OneHertzRegular(1:2694,1),OneHertzRegular(1:2694,3),'LineWidth',2,'color',[1 0 1] )
%    axis([0   75  -.4 2.3])
    axis([0   43  -.4 1.3])
    h= gca;
    set(h,'FontName','Arial','LineWidth',2,'FontSize',12,'FontWeight','bold','box','off')

    vw= legend(h,'measured','fit')
    set(vw,'Position',[.15 .764 .39 .114],'box','off')

for jj=1:200
    plot(OneHertzRegular(1:4:124,4),(-.2 +jj*.00025)  *ones(124/4,1),'.','MarkerSize',1);
    plot(OneHertzRegular(1:4:124,4)+.005,(-.2 +jj*.00025)  *ones(124/4,1),'.','MarkerSize',1);
end

xlabel('time(secs)','FontSize',14,'FontName','Arial','FontWeight','Bold')
ylabel('[DA] \muM','FontSize',14,'FontName','Arial','FontWeight','Bold')
hold

scribe.textbox('String','Fig 1b','LineStyle','none','FontName','Arial', ...
       'Position',[.002,.008,.1,0.073],'FontWeight','Bold');


 delete('Fig1bBM.fig')
 saveas(gcf,'Fig1bBM.fig','fig')
 saveas(gcf,'Fig1bBM.eps','eps')
 saveas(gcf,'Fig1bBM.bmp','bmp')
 saveas(gcf,'Fig1bBM.tif','tif')
 saveas(gcf,'Fig1bBM.jpg','jpg')


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Total.dat is Figure 3
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

FigStrings =['Fig. 3a';'Fig. 3b';'Fig. 3c';'Fig. 3d'] ;
FigRootStr =['Fig3aBM';'Fig3bBM';'Fig3cBM';'Fig3dBM'] ;

mVec  =['1.00';'1.04';'1.08';'1.13'];
CVvec =['.000';'.264';'.440';'.879'];



for pp=1:4
    vf= figure;  set(vf,'name',FigStrings(pp,:));
    plot(Total(1:4290,1), Total(1:4290,pp+1) ,'LineWidth',2,'color',[0 0 .627]) ,hold
    axis([0   75  -.3 1.8])
    h= gca;
    set(h,'FontName','Arial','LineWidth',2,'FontSize',12,'FontWeight','bold','box','off')
 %   stims = 
    stims=find(Total(:,pp+5)==1);
    bb=length(stims);
    for jj=1:200
        plot(stims(:)/60,(-.2 +jj*.00025)       * ones(1,bb),'k.','MarkerSize',1);
        plot(stims(:)/60+.005,(-.2 +jj*.00025)  * ones(1,bb),'k.','MarkerSize',1);
    end

    xlabel('time(secs)','FontSize',14,'FontName','Arial','FontWeight','Bold')
    ylabel('[DA] \muM' ,'FontSize',14,'FontName','Arial','FontWeight','Bold')
    hold

    sa=scribe.textbox('String',FigStrings(pp,:),'LineStyle','none','FontName','Arial', ...
       'Position',[.002,.008,.1,0.073],'FontWeight','Bold');

    sb=scribe.textbox('String',strcat('m=',mVec(pp,:),'s, CV=',CVvec(pp,:)),'LineStyle','none','FontName','Arial', ...
       'Position',[.602,.708,.1,0.073],'FontWeight','Bold');

%    save(strcat(FigRootStr(pp,:),'.fig'))
%    saveas(gca,strcat(FigRootStr(pp,:),'.eps'),'eps')
end
scribe.textbox('String','Fig 3','LineStyle','none','FontName','Arial', ...
       'Position',[.008,-.018,.1,0.073],'FontWeight','Bold');


for pp=1:4
    IntegratedDopamine(pp) = sum(Total(1:4290,pp+1));
    MaxDopamine(pp) = max(Total(1:4290,pp+1));
end



IntegratedDopamine = (Total(4290,1)/4290)*IntegratedDopamine;


ID4=IntegratedDopamine(4);ID1=IntegratedDopamine(1);
MD4=MaxDopamine(4);MD1=MaxDopamine(1);
slope = (ID4-ID1)/(MD4 - MD1);
offs = (ID1 *MD4-MD1*ID4)/(MD4 - MD1);

ToPlot= [IntegratedDopamine' ,slope*MaxDopamine'+offs]

ToPlot1= [IntegratedDopamine' ,zeros(4,1)]
ToPlot2= [zeros(4,1) ,15*MaxDopamine']

y1l=13; y2l= (y1l-offs)/slope;
y1u=24; y2u= (y1u-offs)/slope;


vf= figure;   set(vf,'name','Fig. 4 BM');

    hl1= bar(ToPlot,'group'); hold
    axis([0  5   y1l y1u])
    ax1 = gca;
    set(ax1,'FontName','Arial','LineWidth',2,'FontSize',12,'FontWeight','bold','ylim',[11 25],'box','off')
    vw   =  legend(ax1,'Integrated Dopamine','Maximum Dopamine','box','off')
    set(vw,'box','off','position',[.2 .82 .39 .063])
    ylabel('[DA]\cdot time (\muM - secs)','FontSize',14,'FontName','Arial','FontWeight','Bold')

    ax2 = axes('Position',get(ax1,'Position'),  ... % 'XAxisLocation','top',
        'YAxisLocation','right', 'Color','none', 'XColor','k','YColor','k', ...
     'FontName','Arial','LineWidth',2,'FontSize',12,'FontWeight','bold','ylim',[y2l y2u],)
%    hl2 = bar(ToPlot2,'group');
%    set(ax2,'color',[0 0 .627],'Parent',ax2) ,
%    axis([0   43  -.4 1.3])
%    axis([0   43  -.4 1.3])

%xlabel('CV','FontSize',14,'FontName','Arial','FontWeight','Bold')
     ylabel('max Dopamine(\muM)','FontSize',14,'FontName','Arial','FontWeight','Bold')

     ss= .866 - (.866-.633)/2*(0:6);
   sc1 = scribe.line('x',[.88 .904],'y',[ss(1) ss(1)],'linewidth',3,'color',[1 1 1]) ;
   sc2 = scribe.line('x',[.88 .904],'y',[ss(2) ss(2)],'linewidth',3,'color',[1 1 1]) ;
   sc3 = scribe.line('x',[.88 .904],'y',[ss(3) ss(3)],'linewidth',3,'color',[1 1 1]) ;
   sc4 = scribe.line('x',[.88 .904],'y',[ss(5) ss(5)],'linewidth',3,'color',[1 1 1]) ;
   sc5 = scribe.line('x',[.88 .904],'y',[ss(6) ss(6)],'linewidth',3,'color',[1 1 1]) ;
   sc6 = scribe.line('x',[.88 .904],'y',[ss(7) ss(7)],'linewidth',3,'color',[1 1 1]) ;
 vx= legend(ax2,'max Dopamine')
 set(vx,'box','off','visible','on','Position',[.192 .752 .387 .063],'Color',[1 1 1]); %
 set(vx,'box','on');  set(vx,'box','off');%
 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% MaxDvsAmp.dat is Figure 4
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Col 1 CV	
Col 2 Poisson Parameters 		[DA] -secs						
Col 3 Regular Parameters
Col 4 
Col 5 
Col 6  Max DA
Col 7  2*(Poisson -J2*Max)
Col 8  Amp

Integrated Response and Maximum [DA] vs CV


vf= figure;   set(vf,'name','Fig. 4 BM');

    line(MaxDvsAmp(1:11,1), MaxDvsAmp(1:11,2) ,'LineStyle','--','LineWidth',2, ... 
         'Marker','square','MarkerSize',2,'color',[0 0 .627]) ,hold
    ax1 = gca;
    set(ax1,'FontName','Arial','LineWidth',2,'FontSize',12,'FontWeight','bold','ylim',[11 25])
    vw   =  legend(ax1,'Integrated Dopamine','box','off')
    set(vw,'box','off','position',[.2 .82 .39 .063])
    ylabel('[DA]\cdot time (\muM - secs)','FontSize',14,'FontName','Arial','FontWeight','Bold')

    ax2 = axes('Position',get(ax1,'Position'),  ... % 'XAxisLocation','top',
        'YAxisLocation','right', 'Color','none', 'XColor','k','YColor','k', ...
     'FontName','Arial','LineWidth',2,'FontSize',12,'FontWeight','bold','ylim',[0.5 2.5])
    hl2 = line(MaxDvsAmp(1:11,1), MaxDvsAmp(1:11,6) ,'LineWidth',2, ...
        'Marker','square','MarkerSize',2,'color',[0 0 .627],'Parent',ax2) ,
%    axis([0   43  -.4 1.3])
%    axis([0   43  -.4 1.3])

xlabel('CV','FontSize',14,'FontName','Arial','FontWeight','Bold')
ylabel('max Dopamine(\muM)','FontSize',14,'FontName','Arial','FontWeight','Bold')

 vx= legend(ax2,'max Dopamine')
 set(vx,'box','off','visible','on','Position',[.192 .752 .387 .063],'Color',[1 1 1]); %
 set(vx,'box','on');  set(vx,'box','off');%
 
 save('Fig4BM.fig')
 saveas(gca,'Fig4BM.eps','eps')
 saveas(gca,'Fig4BM.bmp','bmp')
 saveas(gca,'Fig4BM.tif','tif')
 saveas(gca,'Fig4BM.jpg','jpg')
 
 scribe.textbox('String','Fig 4','LineStyle','none','FontName','Arial', ...
       'Position',[.008,-.018,.1,0.073],'FontWeight','Bold');

% --------------------------------------------------------------   
%    Here is another Fig 4.   Unfinished; I eventually made a Fig 3e

Y = round(rand(5,3)*10);
subplot(2,2,1)
bar(Y,'group')
subplot(2,2,2) 
bar(Y,'stack')
title 'Stack'



vf= figure;   set(vf,'name','Fig. 4 BM');

    hl1= line(MaxDvsAmp(1:11,1), MaxDvsAmp(1:11,2) ,'LineStyle','--','LineWidth',2, ... 
         'Marker','square','MarkerSize',2,'color',[0 0 .627]) ,hold
    ax1 = gca;
    set(ax1,'FontName','Arial','LineWidth',2,'FontSize',12,'FontWeight','bold','ylim',[11 25])
    vw   =  legend(ax1,'Integrated Dopamine','box','off')
    set(vw,'box','off','position',[.2 .82 .39 .063])
    ylabel('[DA]\cdot time (\muM - secs)','FontSize',14,'FontName','Arial','FontWeight','Bold')

    ax2 = axes('Position',get(ax1,'Position'),  ... % 'XAxisLocation','top',
        'YAxisLocation','right', 'Color','none', 'XColor','k','YColor','k', ...
     'FontName','Arial','LineWidth',2,'FontSize',12,'FontWeight','bold','ylim',[0.5 2.5])
    hl2 = line(MaxDvsAmp(1:11,1), MaxDvsAmp(1:11,6) ,'LineWidth',2, ...
        'Marker','square','MarkerSize',2,'color',[0 0 .627],'Parent',ax2) ,
%    axis([0   43  -.4 1.3])
%    axis([0   43  -.4 1.3])

xlabel('CV','FontSize',14,'FontName','Arial','FontWeight','Bold')
ylabel('max Dopamine(\muM)','FontSize',14,'FontName','Arial','FontWeight','Bold')

 vx= legend(ax2,'max Dopamine')
 set(vx,'box','off','visible','on','Position',[.192 .752 .387 .063],'Color',[1 1 1]); %
 set(vx,'box','on');  set(vx,'box','off');%
 
 save('Fig4BM.fig')
 saveas(gca,'Fig4BM.eps','eps')
 saveas(gca,'Fig4BM.bmp','bmp')
 saveas(gca,'Fig4BM.tif','tif')
 saveas(gca,'Fig4BM.jpg','jpg')
 
 scribe.textbox('String','Fig 4','LineStyle','none','FontName','Arial', ...
       'Position',[.008,-.018,.1,0.073],'FontWeight','Bold');

 %  
%  set(ax1)
%     
% 
% hl1 = line(x1,y1,'Color','r');
% set(ax1,'XColor','r','YColor','r')
% 
% %%%%


% Next, create another axes at the same location as the first, 
%  placing the x-axis on top and the y-axis on the right. 
%    Set the axes Color to none to allow the first axes to be visible and color code the x- and y-axis to match the data. 





%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%           This is Figure 5

% Need to go to 

vf= figure;  set(vf,'name','Fig. 5a BM: Entropy');
plot(cvVec,Entropy,'LineWidth',2,'color',[0 0 .627]   )
     h= gca;
    set(h,'FontName','Arial','FontSize',12,'FontWeight','bold')
    axis([0   1  -1.5 1.05])
    
    
xlabel('CV','FontSize',14,'FontName','Arial','FontWeight','Bold')
ylabel('Entropy','FontSize',14,'FontName','Arial','FontWeight','Bold')
scribe.textbox('String','Fig 5','LineStyle','none','FontName','Arial', ...
       'Position',[.008,-.018,.1,0.073],'FontWeight','Bold');

 delete('Fig5BM.fig')    
 saveas(gcf,'Fig5BM.fig','fig')
 saveas(gcf,'Fig5BM.eps','eps')
 saveas(gcf,'Fig5BM.bmp','bmp')
 saveas(gcf,'Fig5BM.tif','tif')
 saveas(gcf,'Fig5BM.jpg','jpg')

    
    rcvVec=1:length(cvVec)
plot(cvVec,Entropy   - .5*log(sigmaVec(rcvVec)))
  plot(cvVec,Entropy   +.5*log(stdVec(rcvVec)))

  
   
 % ----- 

 vf= figure;  set(vf,'name','Fig. 5b BM: Entropy');

 cvForMaxDop(1,1:9) = MaxDvsAmp(2:10,1);
 IntegratedDop(1,1:9) =MaxDvsAmp(2:10,2);
 EntropyforMaxDop = interp1(cvVec,Entropy,cvForMaxDop)

EntropyforMaxDop(9) = Entropy(length(Entropy))

 plot(EntropyforMaxDop, IntegratedDop ,'LineWidth',2,'color',[0 0 .627]   )
    xlabel('Entropy','FontSize',14,'FontName','Arial','FontWeight','Bold')
    ylabel('Integrated Dopamine (\muM - secs)','FontSize',14,'FontName','Arial','FontWeight','Bold'); %[DA]\cdot time
     h= gca;
    set(h,'FontName','Arial','FontSize',12,'FontWeight','bold')
%    axis([0   1  -1.5 1.05])
    
      
scribe.textbox('String','Fig 5b','LineStyle','none','FontName','Arial', ...
       'Position',[.008,-.018,.1,0.073],'FontWeight','Bold');

 delete('Fig5BM.fig')    
 saveas(gcf,'Fig5BM.fig','fig')
 saveas(gcf,'Fig5BM.eps','eps')
 saveas(gcf,'Fig5BM.bmp','bmp')
 saveas(gcf,'Fig5BM.tif','tif')
 saveas(gcf,'Fig5BM.jpg','jpg')

    
    rcvVec=1:length(cvVec)
seplot(cvVec,Entropy   - .5*log(sigmaVec(rcvVec)))
  plot(cvVec,Entropy   +.5*log(stdVec(rcvVec)))