Upload files to ''

code/Fig_1/Fig_1.m

@@ -0,0 +1,40 @@
+load('data_fig_1.db','-mat');
+figure;
+
+%% theta frequency
+frwin=3;
+for a=1:2
+    rx=size(data.unwrapped_phase_right{a},1);
+    freq_right{a}=([data.unwrapped_phase_right{a}(:,1+frwin:end) nan(rx,frwin)]-[nan(rx,frwin) data.unwrapped_phase_right{a}(:,1:end-frwin)])/2/pi/0.02/frwin/2;
+    rx=size(data.unwrapped_phase_left{a},1);
+    freq_left{a}=([data.unwrapped_phase_left{a}(:,1+frwin:end) nan(rx,frwin)]-[nan(rx,frwin) data.unwrapped_phase_left{a}(:,1:end-frwin)])/2/pi/0.02/frwin/2;    
+end
+
+%% acceleration
+win=10;    %gfr=fr;
+
+
+    rx=size(data.x_right,1);
+    acc_right=-([nan(rx,win) data.x_right(:,1:end-win)]+[data.x_right(:,1+win:end) nan(rx,win)]-2*data.x_right)/0.02/0.02/win/win;
+    rx=size(data.x_left,1);    
+    acc_left=([nan(rx,win) data.x_left(:,1:end-win)]+[data.x_left(:,1+win:end) nan(rx,win)]-2*data.x_left)/0.02/0.02/win/win;    
+
+
+ fct=0; 
+th=8; axlim=[-120  220]; 
+subplot(2,1,1)
+
+gacc=acc_right;
+gfr=freq_right;
+do_plots;
+
+set(AX(2),'XTick',[]);
+
+subplot(2,1,2)
+    gfr=freq_left; 
+    gacc=acc_left;    
+
+do_plots;
+set(AX(2),'XTick',[]);
+
+

code/Fig_1/Fig_1.m~

@@ -0,0 +1,48 @@
+load('data_fig_1.db','-mat');
+
+
+subplot(2,1,1)
+
+%% theta frequency
+frwin=3;
+for a=1:3
+    rx=size(data.unwrapped_phase_right{a},1);
+    freq_right{a}=([data.unwrapped_phase_right{a}(:,1+frwin:end) nan(rx,frwin)]-[nan(rx,frwin) data.unwrapped_phase_right{a}(:,1:end-frwin)])/2/pi/0.02/frwin/2;
+    rx=size(data.unwrapped_phase_left{a},1);
+    freq_left{a}=([data.unwrapped_phase_left{a}(:,1+frwin:end) nan(rx,frwin)]-[nan(rx,frwin) data.unwrapped_phase_left{a}(:,1:end-frwin)])/2/pi/0.02/frwin/2;    
+end
+
+%% acceleration
+win=10;    %gfr=fr;
+
+
+    rx=size(x,1);
+    acc_right=-([nan(rx,win) x(:,1:end-win)]+[x(:,1+win:end) nan(rx,win)]-2*x)/0.02/0.02/win/win;
+
+
+axlim=0; fct=0; th=0; 
+th=ths(r); axlim=[-120  220]; %fct=0.0127;
+new_figfig
+% set(AX(2),'XTickLabel',[]);
+% set(AX(1),'XTickLabel',[]);
+
+set(AX(2),'XTick',[]);
+% L=legend('acceleration','CA1','MEC','Location', 'NorthEast');
+% set(L,'Color','none')
+% set(L,'box','off')
+% set(L,'FontSize',10,'FontName','Arial')
+subplot(2,1,2)
+% win=5;
+if r==3
+    gfr=data.unwrapped_phase_left; gacc=lacc;
+    rx=size(lx,1);
+    gacc=([nan(rx,win) lx(:,1:end-win)]+[lx(:,1+win:end) nan(rx,win)]-2*lx)/0.02/0.02/win/win;    
+else
+    gfr=fr; gacc=acc; 
+    rx=size(x,1);
+    gacc=-([nan(rx,win) x(:,1:end-win)]+[x(:,1+win:end) nan(rx,win)]-2*x)/0.02/0.02/win/win;    
+end
+% axlim=[-160 260];
+new_figfig
+set(AX(2),'XTick',[]);
+

code/Fig_1/data_fig_1.db

@@ -0,0 +1 @@
+/annex/objects/MD5-s23588531--4e064a672abd8adcfb19ca2e5e57e4f9

code/Fig_1/do_plots.m

@@ -0,0 +1,62 @@
+
+% selx=423:624;
+selx=451:601;
+bx=(-500:500)/50;
+col1=.1*[1 1 1];%[.0 .6 .6];
+col2=[255 175 3]/255;%[77 255 255]/255;%[255 195 77]/255;%[225 54 255]/255;%[255 110 238]/255;[0 .63 .54];
+col3=[255 0 0]/255;%[106 77 255]/255;%[77 255 106]/255;%[.60 .65 .0];%[77 255 106]/255;%[0 1 1]*.9;
+xxx=gacc(:,selx);
+mxxx=nanmean(xxx);
+sxxx=nanstd(xxx)/2;
+hold on;
+
+
+[AX,H1,H2]=plotyy(bx(selx),nanmean(gacc(:,selx)),bx(selx),[nanmean(gfr{1}(:,selx))' nanmean(gfr{2}(:,selx))']);
+uistack(AX(1), 'top')
+
+set(AX(1),'Color','none')
+set(AX(1),'FontName','Arial');
+set(AX(2),'FontName','Arial');
+if th==0
+    th=nanmean(nanmean(gfr{2}(:,[423:450 600:624])));
+    th=.5*round(th/.5);
+end
+if fct==0
+    fct= (nanmean([max(nanmean(gfr{2}(:,selx))) max(nanmean(gfr{1}(:,selx)))])-th)/max(mxxx);
+end
+set(H1,'LineWidth',1);set(H2,'LineWidth',2);
+set(H1,'LineStyle','--','Color',col1);
+ set(H2(1),'Color',col2)
+set(H2(2),'Color',col3)
+if axlim==0
+    axlim=1.01*[1.15*min(mxxx-sxxx) max(mxxx+sxxx)];
+end
+axs=-1000:100:1000;
+set(AX(1),'YLim',axlim);
+set(AX(1),'YLim',axlim); %axlim=get(AX(2),'YLim');
+set(AX(1),'YTick',axs);
+set(AX(1),'XTick',-3:1:3);
+set(AX(2),'YLim',axlim*fct+th-0.0);%[6.332 10]);
+set(AX(2),'YTick',6:11);
+set(AX(1),'XLim',[-1 2]);
+set(AX(2),'XLim',[-1 2]);
+set(AX(1),'box','off');
+set(AX(1),'YColor','k');
+set(AX(2),'YColor',col3);
+AX(1)=reg_plots(AX(1));
+AX(2)=reg_plots(AX(2));
+
+
+
+AX_annotations = copyobj(AX(1),get(AX(1),'Parent'));
+A=area(AX_annotations,bx(selx),[mxxx-sxxx; 2*sxxx]'); 
+set(A(1),'Visible','off');
+set(A(2),'FaceColor',col1+(1-col1)*.8);
+set(A(2),'EdgeColor','none');
+uistack(AX_annotations, 'bottom');
+set(AX_annotations,'Color',[1 1 1]);
+set(AX_annotations,'XTick',[],'YTick',[]);
+dsp=diff(axlim)/200*120;
+ set(AX(1),'dataaspectratio',[1 dsp 1]);
+ set(AX(2),'dataaspectratio',[1 dsp*fct 1]);
+ set(AX_annotations,'dataaspectratio',[1 dsp 1]);

code/Fig_1/reg_plots.m

@@ -0,0 +1,3 @@
+function ax=reg_plots(ax)
+set(ax,'FontSize',13);
+set(ax,'FontName','Arial');

code/Fig_1b/Fig_1c.m

@@ -0,0 +1,60 @@
+load('data_fig_1c.db','-mat');
+
+ecp=bsxfun(@rdivide,bsxfun(@minus,eeghp.ecL.p,nanmean(eeghp.ecL.p)),nanmean(amp.ecL.p));
+ecm=bsxfun(@rdivide,bsxfun(@minus,eeghp.ecL.m,nanmean(eeghp.ecL.m)),nanmean(amp.ecL.m));
+
+egp=bsxfun(@rdivide,bsxfun(@minus,eeg.ecL.p,nanmean(eeg.ecL.p)),nanmean(amp.ecL.p));
+egm=bsxfun(@rdivide,bsxfun(@minus,eeg.ecL.m,nanmean(eeg.ecL.m)),nanmean(amp.ecL.m));
+%%
+ks=20;
+sel=4801:14400;
+aux=ecp(sel,:);
+[idx,Cp]=kmeans(aux',ks);
+aux=ecm(sel,:);
+[idxm,Cm]=kmeans(aux',ks);
+%%
+close all;
+f=figure;
+pp=zeros(19201,ks);
+pm=pp;
+for j=1:ks
+   pp(:,j)=mean(ecp(:,idx==j),2);
+   pm(:,j)=mean(ecm(:,idxm==j),2);
+end
+
+C=corr(pp(sel(2401:4800),:),pm(sel(2401:4800),:));
+C(1:ks+1:end)=-1;
+nsp=8;
+fp=zeros(1,nsp);
+fm=fp;
+phm=fp;
+for j=1:nsp
+    [mC,rw]=max(C);
+    [~,cl]=max(mC);
+    rw=rw(cl);
+  
+
+fp(j)=rw;
+fm(j)=cl;
+   
+   C(:,cl)=-1;
+   C(rw,:)=-1;
+   phm(j)=angle(mean(exp(1i*phase.ecL.m(9601,bsxfun(@eq,idxm,fm(j)))),2));
+end
+[~,sindex]=sort(phm,'descend');
+for k=1:nsp
+    j=sindex(k);
+    subplot(nsp,1,k)
+    rw=fp(j); cl=fm(j);
+    disp([C(rw,cl),corr(pp(sel(2401:4800),rw),pm(sel(2401:4800),cl)),sum(idx==rw),sum(idxm==cl)])
+
+    auxp=mean(ecp(:,idx==fp(j)),2); 
+    plot((-9600:9600)/4800,auxp,'c','LineWidth',2);
+   hold on;
+   auxm=mean(ecm(:,idxm==fm(j)),2);
+   plot((-9600:9600)/4800,auxm,'b','LineWidth',2);
+   axis([-.5 .5 1.4*[-1 1]]);   
+   set(gca,'XTick',-0.5:.5:.5,'XTickLabel',[]);
+   golden(.3)
+end
+set(gca,'XTick',-0.5:.5:.5,'XTickLabel',-0.5:.5:.5,'FontName','Arial');

code/Fig_1b/data_fig_1c.db

@@ -0,0 +1 @@
+/annex/objects/MD5-s1715354492--3acf785d1a3a26cefc5fae02b184747d

code/Fig_2/Fig_2.m

@@ -0,0 +1,225 @@
+load('data_fig_2.db','-mat');
+
+res={};
+bins=(1:300)*.002;
+res.allcells.positive=0*bins;
+res.allcells.negative=0*bins;
+res.allcells.shuff=zeros(1000,numel(bins));
+typ='allcells';
+delta_r=5.0666667;
+delta_l=30.4;
+if 1
+
+    for i=1:numel(alldata.ts)
+        
+        
+        data=alldata.pos{i};
+        ts=alldata.ts{i};
+        
+        %%
+        nlaps=min(numel(data.log.end_r),numel(data.log.end_l));
+        
+        right=(data.log.start_r(1:nlaps)+delta_r)';
+        left=(data.log.start_l(1:nlaps)+delta_l)';  
+
+        dwin=.5;
+
+        %% speed
+        selfast=is_inside(ts',data.log.start_r'+1,data.log.start_r'+4) | is_inside(ts',data.log.start_l'+31.5,data.log.start_l'+34.5);
+        selslow=is_inside(ts',data.log.start_l'+5,data.log.start_l'+25) | is_inside(ts',data.log.start_r'+8,data.log.start_r'+28);
+        res.(typ).fast(i,:)=hist(autoc2(ts(selfast)),bins)/6;    
+        res.(typ).slow(i,:)=hist(autoc2(ts(selslow)),bins)/40;    
+        
+        %% acceleration
+
+        res.(typ).negative=res.(typ).negative + hist(autoc2(ts(is_inside(ts',right,right+dwin))),bins);    
+        res.(typ).positive=res.(typ).positive + hist(autoc2(ts(is_inside(ts',left,left+dwin))),bins);                       
+        
+        for j=1:1000
+            shuff_fast=data.log.start_r'+1+3*rand(1,nlaps);
+            shuff_slow=data.log.start_l'+5+20*rand(1,nlaps);          
+            shuff=[shuff_fast shuff_slow];
+            shuff=shuff(randperm(numel(shuff),numel(shuff)/2));
+            res.(typ).shuff(j,:)=res.(typ).shuff(j,:)+hist(autoc2(ts(is_inside(ts',shuff,shuff+dwin))),bins);
+        end       
+    
+    end
+end
+
+%%
+%%
+close all
+
+bins=(1:150)*.004;
+
+
+%%
+sel=bins>0.01 & bins<0.3;
+types={'allcells'};
+nshuf=1000;
+sname='shuff_speed.mat';
+if exist(sname,'file')
+    load(sname,'-mat')
+else
+    sh_f1={};sh_f2={};
+    for ty=1:numel(types)
+        typ=types{ty};
+    sh_f1.(typ)=zeros(1,nshuf);
+    sh_f2.(typ)=zeros(1,nshuf);      
+for j=1:nshuf
+    if mod(j,100)==0
+        disp(j);
+    end
+    sel=rand(1,size(res.(typ).fast,1))>.5;
+    
+    aux1=mean([reshape(res.(typ).fast(sel,:),[],150); reshape(res.(typ).slow(~sel,:),[],150)]);
+    aux2=mean([reshape(res.(typ).fast(~sel,:),[],150); reshape(res.(typ).slow(sel,:),[],150)]);
+    aux1(end)=0;
+    aux2(end)=0;
+    aux1=aux1/max(aux1);    
+    aux2=aux2/max(aux2);    
+    c1=do_fit(bins,aux1);
+    c2=do_fit(bins,aux2);
+    
+    sh_f1.(typ)(j)=c1.f;
+    sh_f2.(typ)(j)=c2.f;
+ 
+end
+
+    end
+save(sname,'sh_f1','sh_f2','-mat');
+end
+%% speed
+
+for ty=numel(types)
+    typ=types{ty};
+    f=figure;
+    hold on;
+    aux=mean(reshape(res.(typ).fast,[],150))/8;
+    aux(end)=0;
+    mxaux=max(aux);
+    aux=aux/max(aux);
+ffast=do_plot(bins,aux,[0 1 0],[.5 1 .5]);    
+    auxs=mean(reshape(res.(typ).slow,[],150))/40;
+    auxs(end)=0;
+    mxauxs=max(auxs);
+
+    auxs=auxs/max(auxs);
+    fslow=do_plot(bins,auxs,[1 0 0],[1 .5*[1 1]]); 
+
+subplot(2,2,1)
+axis([0 .3 0.3 1.2]);
+gl=1.5;
+
+set(gca,'YTick',0:.5:2,'FontName','Arial');
+subplot(2,2,2)
+axis([0 .3 .15*[-1 1] ]);
+
+set(gca,'YTick',-.6:.1:.6,'FontName','Arial');
+
+
+subplot(2,2,3)
+fbins=7.5:.01:8.5;
+h=hist([sh_f1.(typ) sh_f2.(typ)],fbins);
+b=bar(fbins,h,'LineStyle','none','FaceColor',.5*[1 2 1],'BarWidth',1.0);
+mx=ceil(2*max(h));
+hold on;
+legend off;
+line((fslow)*[1 1],[0 mx],'LineWidth',2,'Color',[1 1 1]*.5);
+line(ffast*[1 1],[0 mx],'LineWidth',2,'Color',[1 .5 .5]);
+axis([7.6 8.0 0 380]);
+axis([7.5 9.5 0 380]);
+
+set(gca,'YTick',0:100:300)
+
+subplot(224)
+fbins=7.5:.01:8.5;
+h=hist([sh_f1.(typ) sh_f2.(typ)],fbins);
+b=bar(fbins,h,'LineStyle','none','FaceColor',.5*[1 2 1],'BarWidth',1.0);
+mx=ceil(2*max(h));
+hold on;
+legend off;
+line((fslow)*[1 1],[0 mx],'LineWidth',2,'Color',[1 1 1]*.5);
+line(ffast*[1 1],[0 mx],'LineWidth',2,'Color',[1 .5 .5]);
+axis([7.6 8.0 0 300]);
+
+set(gca,'YTick',0:100:300)
+
+
+end
+
+%% Acceleration
+gs=0;
+figure;
+for ty=numel(types)
+    typ=types{ty};
+
+    hold on;
+    aux=gauss_boxcar(mean(reshape(res.(typ).negative,2,150)),gs);
+    aux=aux/max(aux);
+fm=do_plot(bins,aux,[0 0 0],.5*[1 1 1]);    
+    aux=gauss_boxcar(mean(reshape(res.(typ).shuff(1,:),2,150)),gs);
+    aux=aux/max(aux);
+    fc=do_plot(bins,aux,[0 1 0],[.5 1 .5]);    
+    aux=gauss_boxcar(mean(reshape(res.(typ).positive,2,150)),gs);
+    aux=aux/max(aux);
+    fp=do_plot(bins,aux,[1 0 0],[1 .5*[1 1]]); 
+end
+subplot(2,2,1)
+axis([0 .3 0.3 1.2]);
+
+set(gca,'YTick',0:.5:2,'FontName','Arial');
+subplot(2,2,2)
+axis([0 .3 .15*[-1 1] ]);
+
+set(gca,'YTick',-.6:.1:.6,'FontName','Arial');
+
+
+sname='shuff_acc.mat';
+if exist(sname,'file')
+    load(sname,'-mat')
+else
+    
+for j=1:nshuf
+    if mod(j,100)==0
+        disp(j);
+    end
+
+    acc_f(j)=all_fit(res.(typ).shuff(j,:));
+
+end
+
+save(sname,'acc_f','-mat');
+end
+%%
+subplot(2,2,3)
+fbins=7:.05:9;
+h=hist(acc_f,fbins);
+h=h;%/sum(h);
+b=bar(fbins,h,'LineStyle','none','FaceColor',.5*[1 2 1],'BarWidth',1.0);
+mx=ceil(2.1*max(h));
+gg=fit(fbins',h','gauss1');
+hold on;
+p1=plot(gg);
+set(p1,'Color',[.5 1 .5]-.5);
+legend off;
+% line(fc*[1 1],[0 mx],'LineWidth',2,'Color',[0 1 0]);
+line(fm*[1 1],[0 mx],'LineWidth',2,'Color',[1 1 1]*.5);
+line(fp*[1 1],[0 mx],'LineWidth',2,'Color',[1 .5 .5]);
+axis([7.5 9.5 0 280]);
+set(gca,'YTick',0:100:300);
+p_n=(1-erf((fm-gg.b1)/gg.c1)); 
+p_p=(1-erf((fp-gg.b1)/gg.c1));
+
+
+title(['neg: ' num2str(1-erf(abs(fm-gg.b1)/gg.c1)) '   pos: ' num2str(p_p)]);
+
+
+
+
+
+%%
+
+
+
+

code/Fig_2/Fig_2.m~

@@ -0,0 +1,196 @@
+load('data_fig_2.db','-mat');
+
+res={};
+bins=(1:300)*.002;
+res.allcells.positive=0*bins;
+res.allcells.negative=0*bins;
+res.allcells.shuff=zeros(1000,numel(bins));
+typ='allcells';
+delta_r=5.0666667;
+delta_l=30.4;
+if 1
+
+    for i=1:numel(alldata.ts)
+        
+        
+        data=alldata.pos{i};
+        ts=alldata.ts{i};
+        
+        %%
+        nlaps=min(numel(data.log.end_r),numel(data.log.end_l));
+        
+        right=(data.log.start_r(1:nlaps)+delta_r)';
+        left=(data.log.start_l(1:nlaps)+delta_l)';  
+
+        dwin=.5;
+
+        %% speed
+        selfast=is_inside(ts',data.log.start_r'+1,data.log.start_r'+4) | is_inside(ts',data.log.start_l'+31.5,data.log.start_l'+34.5);
+        selslow=is_inside(ts',data.log.start_l'+5,data.log.start_l'+25) | is_inside(ts',data.log.start_r'+8,data.log.start_r'+28);
+        res.(typ).fast(i,:)=hist(autoc2(ts(selfast)),bins)/6;    
+        res.(typ).slow(i,:)=hist(autoc2(ts(selslow)),bins)/40;    
+        
+        %% acceleration
+
+        res.(typ).negative=res.(typ).negative + hist(autoc2(ts(is_inside(ts',right,right+dwin))),bins);    
+        res.(typ).positive=res.(typ).positive + hist(autoc2(ts(is_inside(ts',left,left+dwin))),bins);                       
+        
+        for j=1:1000
+            shuff_fast=data.log.start_r'+1+3*rand(1,nlaps);
+            shuff_slow=data.log.start_l'+5+20*rand(1,nlaps);          
+            shuff=[shuff_fast shuff_slow];
+            shuff=shuff(randperm(numel(shuff),numel(shuff)/2));
+            res.(typ).shuff(j,:)=res.(typ).shuff(j,:)+hist(autoc2(ts(is_inside(ts',shuff,shuff+dwin))),bins);
+        end       
+    
+    end
+    save('partial_res.db','res','-mat');
+end
+
+%%
+%%
+close all
+
+bins=(1:150)*.004;
+
+
+%%
+sel=bins>0.01 & bins<0.3;
+types={'allcells'};
+nshuf=1000;
+sname='shuff_speed.mat';
+if 0%exist(sname,'file')
+    load(sname,'-mat')
+else
+    sh_f1={};sh_f2={};
+    for ty=1:numel(types)
+        typ=types{ty};
+    sh_f1.(typ)=zeros(1,nshuf);
+    sh_f2.(typ)=zeros(1,nshuf);      
+%     aux_all=[.5*(res.(typ).fast(:,1:2:end)+res.(typ).fast(:,2:2:end)); .5*(res.(typ).slow(:,1:2:end)+res.(typ).slow(:,2:2:end)); ]
+for j=1:nshuf
+    if mod(j,100)==0
+        disp(j);
+    end
+    sel=rand(1,size(res.(typ).fast,1))>.5;
+    
+    aux1=mean([reshape(res.(typ).fast(sel,:),[],150); reshape(res.(typ).slow(~sel,:),[],150)]);
+    aux2=mean([reshape(res.(typ).fast(~sel,:),[],150); reshape(res.(typ).slow(sel,:),[],150)]);
+    aux1(end)=0;
+    aux2(end)=0;
+    aux1=aux1/max(aux1);    
+    aux2=aux2/max(aux2);    
+    c1=do_fit(bins,aux1);
+    c2=do_fit(bins,aux2);
+    
+    sh_f1.(typ)(j)=c1.f;
+    sh_f2.(typ)(j)=c2.f;
+ 
+end
+
+    end
+save(sname,'sh_f1','sh_f2','-mat');
+end
+%% speed
+
+
+
+
+for ty=numel(types)
+    typ=types{ty};
+    f=figure;
+    hold on;
+    aux=mean(reshape(res.(typ).fast,[],150))/8;
+    aux(end)=0;
+    mxaux=max(aux);
+    aux=aux/max(aux);
+ffast=do_plot(bins,aux,[0 1 0],[.5 1 .5]);    
+    auxs=mean(reshape(res.(typ).slow,[],150))/40;
+    auxs(end)=0;
+    mxauxs=max(auxs);
+
+    auxs=auxs/max(auxs);
+    fslow=do_plot2(bins,auxs,[1 0 0],[1 .5*[1 1]]); 
+
+subplot(2,2,1)
+axis([0 .3 0.3 1.2]);
+gl=1.5;
+%golden(gl);
+set(gca,'YTick',0:.5:2,'FontName','Arial');
+subplot(2,2,2)
+axis([0 .3 .15*[-1 1] ]);
+%golden(gl);
+set(gca,'YTick',-.6:.1:.6,'FontName','Arial');
+
+
+subplot(2,2,3)
+fbins=7.5:.01:8.5;
+h=hist([sh_f1.(typ) sh_f2.(typ)],fbins);
+b=bar(fbins,h,'LineStyle','none','FaceColor',.5*[1 2 1],'BarWidth',1.0);
+mx=ceil(2*max(h));
+hold on;
+legend off;
+line((fslow)*[1 1],[0 mx],'LineWidth',2,'Color',[1 1 1]*.5);
+line(ffast*[1 1],[0 mx],'LineWidth',2,'Color',[1 .5 .5]);
+axis([7.6 8.0 0 380]);
+axis([7.5 9.5 0 380]);
+
+set(gca,'YTick',0:100:300)
+
+subplot(224)
+fbins=7.5:.01:8.5;
+h=hist([sh_f1.(typ) sh_f2.(typ)],fbins);
+b=bar(fbins,h,'LineStyle','none','FaceColor',.5*[1 2 1],'BarWidth',1.0);
+mx=ceil(2*max(h));
+hold on;
+legend off;
+line((fslow)*[1 1],[0 mx],'LineWidth',2,'Color',[1 1 1]*.5);
+line(ffast*[1 1],[0 mx],'LineWidth',2,'Color',[1 .5 .5]);
+axis([7.6 8.0 0 300]);
+
+set(gca,'YTick',0:100:300)
+
+
+end
+
+%% Acceleration
+gs=0;
+for ty=numel(types)
+    typ=types{ty};
+%     f=figure;
+%     plot(bins,gauss_boxcar(sum(res.(typ).right.pre,1),gs),'Color',[1 .5 .5]);
+    hold on;
+    aux=gauss_boxcar(mean(reshape(res.(typ).right,2,150)),gs);
+    aux=aux/max(aux);
+fm=myplot33(bins,aux,[0 0 0],.5*[1 1 1]);    
+%     plot(bins,gauss_boxcar(sum(res.(typ).left.pre,1),gs),'Color',.5*[1 1 1]);
+
+
+
+%%
+
+
+
+   
+%     plot(bins,aux/max(aux),'Color',[1 0 0]);
+    aux=gauss_boxcar(mean(reshape(res.(typ).control(1,:),2,150)),gs);
+    aux=aux/max(aux);
+    fc=myplot3(bins,aux,[0 1 0],[.5 1 .5]);    
+%     plot(bins,aux/max(aux),'Color',[0 1 0]);
+
+    aux=gauss_boxcar(mean(reshape(res.(typ).left,2,150)),gs);
+    aux=aux/max(aux);
+    fp=myplot3(bins,aux,[1 0 0],[1 .5*[1 1]]); 
+end
+subplot(2,2,1)
+axis([0 .3 0.3 1.2]);
+gl=1.5;
+golden(gl);
+set(gca,'YTick',0:.5:2,'FontName','Arial');
+subplot(2,2,2)
+axis([0 .3 .15*[-1 1] ]);
+golden(gl);
+set(gca,'YTick',-.6:.1:.6,'FontName','Arial');
+
+
+

code/Fig_2/all_fit.m

@@ -0,0 +1,8 @@
+function f=all_fit(aux)
+autoc_step=.002;
+autoc_bins=1:300; 
+bins=(autoc_bins)*autoc_step*2; bins=bins(1:150);
+    aux=mean(reshape(aux,[],150));
+    aux=aux/max(aux);    
+    c=do_fit(bins,aux);
+    f=c.f;

code/Fig_2/autoc2.m

@@ -0,0 +1,19 @@
+function ht=autoc2(ts,mn,mx)%,autoc_bins,autoc_step)
+    if ~exist('mn','var')
+        mn=0;
+        mx=1;
+    end
+
+    overl=5000;
+    ksel=1:min(numel(ts),overl); %memory overload
+    ht=[];
+    if numel(ts)>1
+
+        while ~isempty(ksel)
+            haux=bsxfun(@minus,ts(ksel),ts(ksel)');
+            ht=[ht; haux(haux>mn & haux<mx)];
+            ksel=ksel(end)+1:min(numel(ts),ksel(end)+overl);
+        end
+    end
+
+end

code/Fig_2/data_fig_2.db

@@ -0,0 +1 @@
+/annex/objects/MD5-s1008483337--77108f3ab8132a9e1122cb1ff7829911

code/Fig_2/do_fit.m

@@ -0,0 +1,14 @@
+function [c,g]=do_fit(bins,aux)
+
+sel=bins>0.015 & bins<0.3;
+
+
+fo = fitoptions('Method','NonlinearLeastSquares',...
+               'Lower',[0.2,0,0,0,6,-1,-1],...
+               'Upper',[2,20,1,Inf,12,1,1],...
+               'StartPoint',[1,3,.2,10,8,0,0]);
+% ft = fittype('a-b*x+c*cos(2*pi*f*x)*exp(-d*x)','options',fo);
+ft = fittype('a*exp(-b*x)+c*cos(2*pi*f*x)*exp(-d*x)+g+h*x','options',fo);
+
+
+[c, g] = fit(bins(sel)',aux(sel)',ft);

code/Fig_2/do_plot.m

@@ -0,0 +1,45 @@
+function ff=do_plot(bins,aux,col,col2)    
+% aux=aux/max(aux);
+% 
+% sel=bins>0.015 & bins<0.3;
+% 
+% 
+% fo = fitoptions('Method','NonlinearLeastSquares',...
+%                'Lower',[0.2,0,0,0,6,-1,-1],...
+%                'Upper',[2,20,1,Inf,12,1,1],...
+%                'StartPoint',[1,3,.2,10,8,0,0]);
+% % ft = fittype('a-b*x+c*cos(2*pi*f*x)*exp(-d*x)','options',fo);
+% ft = fittype('a*exp(-b*x)+c*cos(2*pi*f*x)*exp(-d*x)+g+h*x','options',fo);
+
+% aux=aux/max(aux);
+[curve2, ~] = do_fit(bins,aux);
+
+aux=aux/feval(curve2,0);
+[curve2, gof2] = do_fit(bins,aux);
+
+
+subplot(2,2,1)
+p1=plot(curve2);
+set(p1,'Color',col2,'LineWidth',2);
+hold on;
+plot(bins,aux,'Color',col);
+plot(bins,curve2.a*exp(-curve2.b*bins)+curve2.g+curve2.h*bins,'c');
+
+legend off;
+axis([0 .3 0 1.2])
+
+subplot(2,2,2)
+
+plot(bins,curve2.c.*cos(2*pi*curve2.f.*bins).*exp(-curve2.d.*bins),'Color',col2,'LineWidth',2);
+hold on;
+plot(bins,aux-feval(curve2,bins)'+curve2.c.*cos(2*pi*curve2.f.*bins).*exp(-curve2.d.*bins),'Color',col);
+
+% disp(curve2.f);
+% curve2
+% gof2
+disp(['freq: ' trunc(curve2.f,3) '  R^{2}: ' trunc(gof2.rsquare,3)]);
+ff=curve2.f;
+axis([0 .3 -.2 .2])
+
+
+% plot(bins,.2+.9*exp(-6*bins).*(1+.1*cos(2*pi*bins*9)),'b');

code/Fig_3/Fig_3.m

@@ -0,0 +1,149 @@
+close all;
+% load('N:\kropff\analyisis\figs\fig2a\all open data_q95.db','-mat');
+% load('N:\kropff\analyisis\figs\fig2a\all car data_q95.db','-mat');
+% % % load('all open data.db','-mat');
+% % % load('all car data.db','-mat');
+% db_all.car=db_car;
+% db_all.open=db_open;
+% clear db_car; clear db_open;
+%     addpath('/gpfs/work/kropff/common_functions/');
+
+dats={'ec', 'hc'};
+% Haccum={};
+w1=0; w2=3;
+jet2=interp1(1:64,jet,(1+w1:(64-w1-w2-1)/(64-1):64-w2));
+
+% load(['all ' dat ' data.db'],'-mat')
+dang=pi/50;
+dacc=20;%15;
+accbins=0:dacc:200; %350
+angbins=-pi:dang:pi; angbins(end)=[]; angbins=angbins+dang/2;  
+cols=[ [1 0 0]; [255 175 3]/255];
+
+
+for d=1:2%1:2
+% dat='open';
+% dat='car';
+dat=dats{d};
+
+% close all;
+
+
+
+ 
+% close all;
+histlim=5*dacc;
+
+
+
+angbins=-pi:dang:pi; angbins(end)=[]; angbins=angbins+dang/2;  
+Ha=zeros(numel(angbins),numel(accbins)); Oa=Ha;
+Hva=Ha;
+
+myfl=['data_fig_3.db'];
+
+nr=numel(unique(who('-file',myfl)));
+L={};
+L.Fr=cell(1,-1+2*numel(accbins));
+L.BFr=cell(-1+2*numel(accbins),7);
+% L.Occ=zeros(1,21);
+% L.Fr=L.Occ; L.BFr=L.Occ;
+for i=[1 3:nr] %2 hippocampal drive not in hippocampus
+    rat=['rat' num2str(i)];
+    db_all.(dat)=load(myfl,rat,'-mat');
+    db_all.(dat)=db_all.(dat).(rat);
+
+    sel=true(size(db_all.(dat).accel));
+    if isempty(sel)
+        warning([num2str(i) 'empty']);
+        continue;
+    end
+    
+
+    sel=sel&db_all.(dat).t_dist>5;
+
+    sel1=conv(double(sel), ones(25,1)/25,'same');
+    sel=sel1>.99;
+    selx=abs(cos(db_all.(dat).anga-db_all.(dat).angv))>.9;
+    xxx=db_all.(dat).accel.*sign(cos(db_all.(dat).anga-db_all.(dat).angv));
+
+    frw=50*5;
+%     [Hs,Occs,Hv]=freq_hist(db_all.(dat).accel(sel),mod(db_all.(dat).anga(sel)-db_all.(dat).angv(sel)+pi,2*pi)-pi,db_all.(dat).freq(sel)-fr_baseline(sel),accbins,angbins);
+    [Hs,Occs,Hv]=freq_hist(db_all.(dat).accel(sel),mod(db_all.(dat).anga(sel)-db_all.(dat).angv(sel)+pi,2*pi)-pi,db_all.(dat).freq(sel),accbins,angbins);
+
+    % Linear results
+    
+    rxxx=round(xxx/dacc);
+%     L.Occ=arrayfun(@(x) sum(rxxx(sel & selx)==x),-10:10);
+    ws=50*[1 2 5 10 20 50 200];
+    aend=accbins(end)/dacc;
+    for w=1:numel(ws)
+        fr_baseline=nanconv(db_all.(dat).freq,ones(ws(w),1)/ws(w)); %10 s window
+        
+        for x=-aend:aend
+            xj=aend+1+x;
+            if w==1
+                L.Fr{xj}=[L.Fr{xj}; db_all.(dat).freq(rxxx==x & sel & selx)];
+            end
+            L.BFr{xj,w}=[L.BFr{xj,w}; fr_baseline(rxxx==x & sel & selx)];        
+        end    
+    end
+    
+    Ha=Ha+Hs.*Occs; Oa=Oa+Occs; Hva=Hva+Hv.*Occs;
+    Haccum{d}=Ha;
+    Oaccum{d}=Oa;
+    Laccum{d}=L;
+
+end
+
+end
+%%
+close all
+for d=1:2
+angbins(101)=angbins(1)+2*pi;
+
+f=figure;
+hold off;
+H=Haccum{d}; O=Oaccum{d};
+H(O>0)=H(O>0)./O(O>0);
+histlim=max(O(:))/10000;
+
+H(O<histlim)=nan;
+Hs=gauss_smoothing_2D(repmat(H,3,1),[1.5 3]);
+h1=size(H,1);
+Hs=Hs(h1+(1:h1),:);
+
+polar3d(flipud(Hs'),angbins(1),angbins(end),accbins(1),accbins(size(Hs,2)),1,'surf'); shading flat;view(2); 
+axis equal; view([-90,90]); axis tight; 
+hold on;
+
+colormap(jet2);
+
+cc=colorbar;
+set(cc,'YTick',7.5:.5:10);
+name=['fig2a_' dat];
+plot3([0 100],200*[1 1],[15 15],'k')
+caxis([8.1 9.7])
+
+
+
+name=['colormap_' dats{d}];
+end
+%%
+f=figure;
+accbins2=[fliplr(-accbins(2:end)) accbins];
+for d=2
+    L=Laccum{d};
+    sel=cellfun(@numel,L.Fr);
+    sel=sel>max(sel)/10000;
+    mn=cellfun(@mean,L.Fr);
+    eom=cellfun(@std,L.Fr)./sqrt(cellfun(@numel,L.Fr));
+    errorbar(accbins2(sel),mn(sel),eom(sel),'Color',cols(d,:));
+    hold on;
+    plot(accbins2(sel),mn(sel),'Color',cols(d,:),'LineWidth',3);
+    dat_acc=arrayfun(@(x) accbins2(x)*ones(size(L.Fr{x})),1:numel(L.Fr),'UniformOutput',0);
+ 
+end
+axis([-220 240 7.95 10]);
+set(gca,'YTick',6:.5:12)
+

code/Fig_3/Fig_3.m~

@@ -0,0 +1,180 @@
+close all;
+% load('N:\kropff\analyisis\figs\fig2a\all open data_q95.db','-mat');
+% load('N:\kropff\analyisis\figs\fig2a\all car data_q95.db','-mat');
+% % % load('all open data.db','-mat');
+% % % load('all car data.db','-mat');
+% db_all.car=db_car;
+% db_all.open=db_open;
+% clear db_car; clear db_open;
+%     addpath('/gpfs/work/kropff/common_functions/');
+
+dats={'ec', 'hc'};
+% Haccum={};
+w1=0; w2=3;
+jet2=interp1(1:64,jet,(1+w1:(64-w1-w2-1)/(64-1):64-w2));
+
+% load(['all ' dat ' data.db'],'-mat')
+dang=pi/50;
+dacc=20;%15;
+accbins=0:dacc:200; %350
+angbins=-pi:dang:pi; angbins(end)=[]; angbins=angbins+dang/2;  
+cols=[ [1 0 0]; [255 175 3]/255];
+
+
+for d=1:2%1:2
+% dat='open';
+% dat='car';
+dat=dats{d};
+
+% close all;
+
+
+
+ 
+% close all;
+histlim=5*dacc;
+
+
+
+angbins=-pi:dang:pi; angbins(end)=[]; angbins=angbins+dang/2;  
+Ha=zeros(numel(angbins),numel(accbins)); Oa=Ha;
+Hva=Ha;
+
+myfl=['data_fig_3.db'];
+
+nr=numel(unique(who('-file',myfl)));
+L={};
+L.Fr=cell(1,-1+2*numel(accbins));
+L.BFr=cell(-1+2*numel(accbins),7);
+% L.Occ=zeros(1,21);
+% L.Fr=L.Occ; L.BFr=L.Occ;
+for i=[1 3:nr] %2 hippocampal drive not in hippocampus
+    rat=['rat' num2str(i)];
+    db_all.(dat)=load(myfl,rat,'-mat');
+    db_all.(dat)=db_all.(dat).(rat);
+
+    sel=true(size(db_all.(dat).accel));
+    if isempty(sel)
+        warning([num2str(i) 'empty']);
+        continue;
+    end
+    
+
+    sel=sel&db_all.(dat).t_dist>5;
+
+    sel1=conv(double(sel), ones(25,1)/25,'same');
+    sel=sel1>.99;
+    selx=abs(cos(db_all.(dat).anga-db_all.(dat).angv))>.9;
+    xxx=db_all.(dat).accel.*sign(cos(db_all.(dat).anga-db_all.(dat).angv));
+
+    frw=50*5;
+%     [Hs,Occs,Hv]=freq_hist(db_all.(dat).accel(sel),mod(db_all.(dat).anga(sel)-db_all.(dat).angv(sel)+pi,2*pi)-pi,db_all.(dat).freq(sel)-fr_baseline(sel),accbins,angbins);
+    [Hs,Occs,Hv]=freq_hist(db_all.(dat).accel(sel),mod(db_all.(dat).anga(sel)-db_all.(dat).angv(sel)+pi,2*pi)-pi,db_all.(dat).freq(sel),accbins,angbins);
+
+    % Linear results
+    
+    rxxx=round(xxx/dacc);
+%     L.Occ=arrayfun(@(x) sum(rxxx(sel & selx)==x),-10:10);
+    ws=50*[1 2 5 10 20 50 200];
+    aend=accbins(end)/dacc;
+    for w=1:numel(ws)
+        fr_baseline=nanconv(db_all.(dat).freq,ones(ws(w),1)/ws(w)); %10 s window
+        
+        for x=-aend:aend
+            xj=aend+1+x;
+            if w==1
+                L.Fr{xj}=[L.Fr{xj}; db_all.(dat).freq(rxxx==x & sel & selx)];
+            end
+            L.BFr{xj,w}=[L.BFr{xj,w}; fr_baseline(rxxx==x & sel & selx)];        
+        end    
+    end
+    
+    Ha=Ha+Hs.*Occs; Oa=Oa+Occs; Hva=Hva+Hv.*Occs;
+    Haccum{d}=Ha;
+    Oaccum{d}=Oa;
+    Laccum{d}=L;
+
+end
+
+end
+%%
+close all
+for d=1:2
+angbins(101)=angbins(1)+2*pi;
+
+f=figure;
+hold off;
+H=Haccum{d}; O=Oaccum{d};
+H(O>0)=H(O>0)./O(O>0);
+histlim=max(O(:))/10000;
+
+H(O<histlim)=nan;
+Hs=gauss_smoothing_2D(repmat(H,3,1),[1.5 3]);
+h1=size(H,1);
+Hs=Hs(h1+(1:h1),:);
+
+polar3d(flipud(Hs'),angbins(1),angbins(end),accbins(1),accbins(size(Hs,2)),1,'surf'); shading flat;view(2); 
+axis equal; view([-90,90]); axis tight; 
+hold on;
+
+colormap(jet2);
+
+cc=colorbar;
+set(cc,'YTick',7.5:.5:10);
+name=['fig2a_' dat];
+plot3([0 100],200*[1 1],[15 15],'k')
+caxis([8.1 9.7])
+
+
+
+name=['colormap_' dats{d}];
+end
+%%
+f=figure;
+accbins2=[fliplr(-accbins(2:end)) accbins];
+for d=2
+    L=Laccum{d};
+    sel=cellfun(@numel,L.Fr);
+    sel=sel>max(sel)/10000;
+    mn=cellfun(@mean,L.Fr);
+    eom=cellfun(@std,L.Fr)./sqrt(cellfun(@numel,L.Fr));
+    errorbar(accbins2(sel),mn(sel),eom(sel),'Color',cols(d,:));
+    hold on;
+    plot(accbins2(sel),mn(sel),'Color',cols(d,:),'LineWidth',3);
+    dat_acc=arrayfun(@(x) accbins2(x)*ones(size(L.Fr{x})),1:numel(L.Fr),'UniformOutput',0);
+ 
+end
+axis([-220 240 7.95 10]);
+set(gca,'YTick',6:.5:12)
+
+%%
+f=figure;
+accbins2=[fliplr(-accbins(2:end)) accbins];
+for d=2
+    subplot(1,2,d)
+    for j=[3 ]
+        
+    L=Laccum{d};
+    L.Fr2=arrayfun(@(x) L.Fr{x}-L.BFr{x,j},1:numel(L.Fr),'UniformOutput',0);
+    sel=cellfun(@numel,L.Fr2);
+    sel=sel>max(sel)/10000 & abs(accbins2)<210;
+    mn=cellfun(@mean,L.Fr2);
+    eom=cellfun(@std,L.Fr2)./sqrt(cellfun(@numel,L.Fr2));
+    errorbar(accbins2(sel),mn(sel),eom(sel),'Color',cols(d,:)*(j-1)/6);
+    hold on;
+    plot(accbins2(sel),mn(sel),'Color',cols(d,:)*(j-1)/6,'LineWidth',3);
+    end
+
+    
+    axis([-220 240 -.25 2]);
+    golden(2.5)
+set(gca,'YTick',-2:.5:12)
+
+    
+    dat_acc=arrayfun(@(x) accbins2(x)*ones(size(L.Fr2{x})),1:numel(L.Fr2),'UniformOutput',0);    
+%     [cr_neg,p_neg]=corr(cell2mat(dat_acc(sel & accbins2<0)'),cell2mat(L.Fr2(sel & accbins2<0)'))
+%     [cr_neg,p_neg]=corr(accbins2(sel & accbins2<0)',mn(sel & accbins2<0)');
+    lneg_bas = fitlm(accbins2(sel & accbins2<0)',mn(sel & accbins2<0)','linear')
+    lpos_bas = fitlm(accbins2(sel & accbins2>0)',mn(sel & accbins2>0)','linear')
+    
+end

code/Fig_3/data_fig_3.db

@@ -0,0 +1 @@
+/annex/objects/MD5-s1831968181--c6324a1e14a92575b78e9b2d687066c3

code/Fig_3/gauss_smoothing_2D.m

@@ -0,0 +1,12 @@
+function map=gauss_smoothing_2D(map, sigma)
+% filt=fspecial('gaussian',[ceil(3*sigma(1)) ceil(3*sigma(2))]);%,sigma(1));
+l1=ceil(3*sigma(1));
+l2=ceil(3*sigma(2));
+[X,Y]=meshgrid(-l1:l1,-l2:l2);
+filt=exp(-.5*((X/sigma(1)).^2+(Y/sigma(2)).^2));
+filt=filt/sum(filt(:));
+[m1,m2]=size(map);
+
+% map=repmat(map,3,3); % angular data can be tricky
+map=nanconv(map,filt,'same','edge','2d' );
+% map=map(m1+(1:m1),m2+(1:m2));

code/Fig_3d/Fig_3d.m

@@ -0,0 +1,116 @@
+close all;
+res3={};
+fcol=figure;
+    acclim=50; 
+    savename=['intrinsic_results_ns_res' num2str(acclim) '.mat'];
+    
+    
+
+    load(savename,'-mat');
+    types=fields(res2);
+    for ty=1:numel(types)
+    typ=types{ty};
+
+ f1=figure
+
+    gs=1;
+    res=res2;
+
+
+        
+    autoc_step=.002;
+    autoc_bins=1:300; 
+    bins=(autoc_bins)*autoc_step;
+    phbins_norm=(pi/60:pi/60:5*pi)/2/pi;
+
+%%
+
+    hold on;
+    aux=nanmean(res.(typ).right);
+    aux=aux/max(aux);
+    fm=do_plot(bins,aux,[0 0 0],.5*[1 1 1]);    
+
+   
+    aux=nanmean(res.(typ).left);
+    aux=aux/max(aux);
+    fc=do_plot(bins,aux,[1 0 0],[1 .5 .5]);    
+    
+    aux=res.(typ).control2(1,:);
+    aux=aux/max(aux);
+    fneut=do_plot(bins,aux,[0 122 255]/255,[51 226 217]/255);     
+
+    subplot(2,2,1)
+    set(gca,'XTick',0:.1:.6);
+    axis([0 .3 0 1]);
+    set(gca,'YTick',0:.5:2,'FontName','Arial');
+    ss=subplot(2,2,2)
+    set(gca,'XTick',0:.1:.6);
+    axis([0.0 .2 0.25*[-1 1] ]);
+    set(gca,'YTick',-.6:.2:.6,'FontName','Arial');
+    title(typ);
+
+    %%
+        fsh=res.(typ).fcontrol2(:);
+    
+    %%
+        subplot(2,2,4)
+    fbins=6:.02:9;
+        h=hist(fsh,fbins);
+
+        h=h/sum(h);
+        bar(fbins,h,'EdgeColor','none','FaceColor',[51 226 217]/255,'BarWidth',.8);
+
+        gg=fit(fbins',h(:),'gauss1');
+        hold on;
+        p1=plot(gg);
+        set(p1,'Color',[0 122 255]/255);
+        legend off;
+        mx=ceil(1.1*max(max(h),gg.a1));
+        line(fm*[1 1],[0 mx],'LineWidth',2,'Color',[1 1 1]*0);
+        line(fc*[1 1],[0 mx],'LineWidth',2,'Color',[1 0 0]);
+        axis([7.3 9 0 0.3]);
+        set(gca,'XTick',7:.5:10,'YTick',0:0.1:1);
+
+
+%%
+if ty<numel(types)
+figure(fcol)
+ subplot(2,3,ty)
+        h=hist(fsh,fbins);
+ 
+
+        h=h/sum(h);
+        bar(fbins,h,'EdgeColor','none','FaceColor',[51 226 217]/255,'BarWidth',1.0);
+
+        gg=fit(fbins',h(:),'gauss1');
+        hold on;
+        p1=plot(gg);
+        set(p1,'Color',[0 122 255]/255);
+        legend off;
+        mx=ceil(1.1*max(max(h),gg.a1));
+        line(fm*[1 1],[0 mx],'LineWidth',1,'Color',[1 1 1]*0);
+        line(fc*[1 1],[0 mx],'LineWidth',1,'Color',[1 0 0]);
+        axis([7.3 9 0 0.3]);
+
+        set(gca,'XTick',7:.5:9,'XTickLabel',{7,'',8,'',9});
+else
+
+    ylim([0 .5]);
+    xlim([7.9 9])
+        subplot(2,2,3)
+
+
+                
+                subplot(2,2,2)
+                axis([.09 .15 -.05 .08])
+             
+                
+
+end
+  end
+ 
+  
+    
+
+
+%%

code/Fig_3d/Fig_3d.m~

@@ -0,0 +1,120 @@
+close all;
+res3={};
+fcol=figure;
+    acclim=50; 
+    savename=['intrinsic_results_ns_res' num2str(acclim) '.mat'];
+    
+    
+
+    load(savename,'-mat');
+    types=fields(res2);
+    for ty=1:numel(types)
+    typ=types{ty};
+%         vs={'pn','pm','cfm','cfc','cshm','cshst','cshp99','cshp95','cshp5','cshp1'};
+%         for j=1:numel(vs)
+%             eval([vs{j} '=nan*acclist;']);
+%         end
+ f1=figure
+
+    gs=1;
+    res=res2;
+
+
+        
+    autoc_step=.002;
+    autoc_bins=1:300; 
+    bins=(autoc_bins)*autoc_step;
+    phbins_norm=(pi/60:pi/60:5*pi)/2/pi;
+
+%%
+
+    hold on;
+    aux=nanmean(res.(typ).right);
+    aux=aux/max(aux);
+    fm=do_plot(bins,aux,[0 0 0],.5*[1 1 1]);    
+
+   
+    aux=nanmean(res.(typ).left);
+    aux=aux/max(aux);
+    fc=do_plot(bins,aux,[1 0 0],[1 .5 .5]);    
+    
+    aux=res.(typ).control2(1,:);
+    aux=aux/max(aux);
+    fneut=do_plot(bins,aux,[0 122 255]/255,[51 226 217]/255);     
+
+    subplot(2,2,1)
+    set(gca,'XTick',0:.1:.6);
+    axis([0 .3 0 1]);
+    set(gca,'YTick',0:.5:2,'FontName','Arial');
+    ss=subplot(2,2,2)
+    set(gca,'XTick',0:.1:.6);
+    axis([0.0 .2 0.25*[-1 1] ]);
+    set(gca,'YTick',-.6:.2:.6,'FontName','Arial');
+    title(typ);
+
+    %%
+        fsh=res.(typ).fcontrol2(:);
+    
+    %%
+        subplot(2,2,4)
+    fbins=6:.02:9;
+        h=hist(fsh,fbins);
+
+        h=h/sum(h);
+        bar(fbins,h,'EdgeColor','none','FaceColor',[51 226 217]/255,'BarWidth',.8);
+
+        gg=fit(fbins',h(:),'gauss1');
+        hold on;
+        p1=plot(gg);
+        set(p1,'Color',[0 122 255]/255);%[.5 1 .5]-.5);
+        legend off;
+        mx=ceil(1.1*max(max(h),gg.a1));
+        line(fm*[1 1],[0 mx],'LineWidth',2,'Color',[1 1 1]*0);
+        line(fc*[1 1],[0 mx],'LineWidth',2,'Color',[1 0 0]);
+        axis([7.3 9 0 0.3]);
+        set(gca,'XTick',7:.5:10,'YTick',0:0.1:1);
+
+
+%%
+if ty<numel(types)
+figure(fcol)
+ subplot(2,3,ty)
+        h=hist(fsh,fbins);
+ 
+
+        h=h/sum(h);
+        bar(fbins,h,'EdgeColor','none','FaceColor',[51 226 217]/255,'BarWidth',1.0);
+
+        gg=fit(fbins',h(:),'gauss1');
+        hold on;
+        p1=plot(gg);
+        set(p1,'Color',[0 122 255]/255);
+        legend off;
+        mx=ceil(1.1*max(max(h),gg.a1));
+        line(fm*[1 1],[0 mx],'LineWidth',1,'Color',[1 1 1]*0);
+        line(fc*[1 1],[0 mx],'LineWidth',1,'Color',[1 0 0]);
+        axis([7.3 9 0 0.3]);
+
+        set(gca,'XTick',7:.5:9,'XTickLabel',{7,'',8,'',9});%,'YTick',0:0.5:1);
+
+else
+
+    ylim([0 .5]);
+    xlim([7.9 9])
+        subplot(2,2,3)
+
+
+                
+                subplot(2,2,2)
+                axis([.09 .15 -.05 .08])
+             
+                
+
+end
+  end
+ 
+  
+    
+
+
+%%

code/Fig_3d/do_fit.m

@@ -0,0 +1,18 @@
+function [c,g]=do_fit(bins,aux)
+c={};
+c.f=nan;
+g=nan;
+
+if ~any(isnan(aux))
+sel=bins>0.015 & bins<0.3;
+
+
+fo = fitoptions('Method','NonlinearLeastSquares',...
+               'Lower',[0.2,0,0,0,6,-1,-1],...
+               'Upper',[2,20,1,Inf,12,1,1],...
+               'StartPoint',[1,3,.2,1,8,0,0]);
+ft = fittype('a*exp(-b*x)+c*cos(2*pi*f*x)*exp(-d*x)+g+h*x','options',fo);
+
+
+[c, g] = fit(bins(sel)',aux(sel)',ft);
+end

code/Fig_3d/do_plot.m

@@ -0,0 +1,29 @@
+function ff=do_plot(bins,aux,col,col2)    
+
+[curve2, ~] = do_fit(bins,aux);
+
+aux=aux/feval(curve2,0);
+[curve2, gof2] = do_fit(bins,aux);
+
+
+subplot(2,2,1)
+p1=plot(curve2);
+set(p1,'Color',col2,'LineWidth',2);
+hold on;
+plot(bins,aux,'Color',col);
+
+legend off;
+axis([0 .3 0 1.2])
+
+subplot(2,2,2)
+
+plot(bins,curve2.c.*cos(2*pi*curve2.f.*bins).*exp(-curve2.d.*bins),'Color',col2,'LineWidth',2);
+hold on;
+plot(bins,aux-feval(curve2,bins)'+curve2.c.*cos(2*pi*curve2.f.*bins).*exp(-curve2.d.*bins),'Color',col);
+
+
+disp(['freq: ' trunc(curve2.f,3) '  R^{2}: ' trunc(gof2.rsquare,3)]);
+ff=curve2.f;
+axis([0 .3 -.2 .2])
+
+

code/Fig_5/Fig_5.m

@@ -0,0 +1,44 @@
+function Fig_5()
+x=-50:1:50; % 50 samples per second
+close all;
+f2=figure;
+speed=.92*exp(-(x/15).^2);
+acc=get_acceleration(x,0,speed)/5;
+acc0=get_acceleration0(x,0,acc);
+th=smooth_theta2(acc,1-exp(-.1));
+ff=90; fg=.8;
+
+subplot(2,3,1)
+plot(x,[speed*fg; acc*ff*fg*.8; acc*ff]');
+
+set(gca,'YTick',-1:1,'XTick',-100:50:100);
+axis([-50 50 -1 1])
+axis square
+% hold on;
+
+subplot(2,3,2)
+plot(x,[speed*fg; acc*ff*fg*.8; acc0*ff]');
+set(gca,'YTick',-1:1,'XTick',-100:50:100);
+axis([-50 50 -1 1])
+axis square
+
+subplot(2,3,3)
+plot(x,[speed*fg; acc*ff*fg*.8; th*ff]');
+
+set(gca,'YTick',-1:1,'XTick',-100:50:100);
+axis([-50 50 -1 1])
+axis square
+[corr(speed',acc') corr(acc',acc');
+corr(speed',acc0') corr(acc',acc0');
+corr(speed',th') corr(acc',th')]
+
+
+
+
+function acc=get_acceleration(x,x0,speed)
+x=x-x0;
+acc=[0 speed(3:end)-speed(1:end-2) 0]/2;
+
+
+function acc=get_acceleration0(x,x0,acc)
+acc(acc<0)=0;

code/Fig_5/Fig_5.m~

@@ -0,0 +1,44 @@
+function Fig_5()
+x=-50:1:50;
+close all;
+f2=figure;
+speed=.92*exp(-(x/15).^2);
+acc=get_acceleration(x,0,speed)/5;
+acc0=get_acceleration0(x,0,acc);
+th=smooth_theta2(acc,1-exp(-.1));%r(i)*(1+decay)>th(i-1)
+ff=90; fg=.8;
+
+subplot(2,3,1)
+plot(x,[speed*fg; acc*ff*fg*.8; acc*ff]');
+
+set(gca,'YTick',-1:1,'XTick',-100:50:100);
+axis([-50 50 -1 1])
+axis square
+% hold on;
+
+subplot(2,3,2)
+plot(x,[speed*fg; acc*ff*fg*.8; acc0*ff]');
+set(gca,'YTick',-1:1,'XTick',-100:50:100);
+axis([-50 50 -1 1])
+axis square
+
+subplot(2,3,3)
+plot(x,[speed*fg; acc*ff*fg*.8; th*ff]');
+
+set(gca,'YTick',-1:1,'XTick',-100:50:100);
+axis([-50 50 -1 1])
+axis square
+[corr(speed',acc') corr(acc',acc');
+corr(speed',acc0') corr(acc',acc0');
+corr(speed',th') corr(acc',th')]
+
+
+
+
+function acc=get_acceleration(x,x0,speed)
+x=x-x0;
+acc=[0 speed(3:end)-speed(1:end-2) 0]/2;
+
+
+function acc=get_acceleration0(x,x0,acc)
+acc(acc<0)=0;

code/Fig_5/smooth_theta2.m

@@ -0,0 +1,12 @@
+function th=smooth_theta2(acc0,decay)
+
+racc=acc0;
+racc(acc0<0)=0;
+th=racc(1)*ones(size(racc));
+for i=2:numel(racc)
+    if racc(i)*(1+decay)>th(i-1) 
+        th(i)=racc(i); 
+    else        
+        th(i)=th(i-1)*(1-decay);    
+    end
+end

code/Fig_6/Fig_6.m

@@ -0,0 +1,67 @@
+load('data_fig_6.db','-mat');
+close all;
+
+
+%%
+f=figure;
+subplot(2,1,1)
+sel=1:numel(data.v);
+speeds=data.v(sel)';
+accs=[0 speeds(3:end)-speeds(1:end-2) 0]/.04;
+ths=smooth_theta2(accs,1-exp(-.1)); 
+
+raccs=accs;
+raccs(raccs<0)=0;
+
+plot(sel/50,(speeds)/std(speeds),'r');
+hold on;
+
+
+plot(sel/50,accs/std(accs),'b');
+plot(sel/50,ths/std(accs),'g');
+xlim([0 10])
+
+
+ths=8+ths/100;
+da=20; abins=-5:5;
+dv=5; vbins=0:12;
+rv=round(speeds/dv);
+ra=round(accs/da);
+vfr=arrayfun(@(x) mean(ths(rv==x)),vbins);
+afr=arrayfun(@(x) mean(ths(ra==x)),abins);
+m8=vfr(1);
+[th_f5,th_s5]=mybp(ths-m8,0.5);
+[th_f1,th_s1]=mybp(ths-m8,0.1);
+[th_f3,th_s3]=mybp(ths-m8,0.3);
+
+vfr_f=arrayfun(@(x) mean(th_f5(rv==x))+m8,vbins);
+vfr_s=arrayfun(@(x) mean(th_s5(rv==x))+m8,vbins);
+
+afr_f=arrayfun(@(x) mean(th_f1(ra==x))+m8,abins);
+afr_s=arrayfun(@(x) mean(th_s1(ra==x))+m8,abins);
+
+
+subplot(2,2,3)
+plot(vbins*dv,[vfr; vfr_s; vfr_f]');
+axis([0 50 7.8 9.1])
+golden(2.5)
+
+subplot(2,2,4)
+plot(abins*da,[afr; afr_s; afr_f]');
+axis([-120 140 7.8 9.1])
+set(gca,'YTick',6:.5:12);
+
+
+
+[corr(speeds',ths') corr(accs',ths');
+corr(speeds',th_f5') corr(accs',th_f1');
+corr(speeds',th_s5') corr(accs',th_s1')]
+drawnow;
+
+
+
+
+function [th_f,th_s]=mybp(ths,f1)
+    th_f=real(bandpass_and_hilbert(ths',f1,100,50)');
+    th_s=real(bandpass_and_hilbert(ths',0,f1,50)');
+end

code/Fig_6/Fig_6.m~

@@ -0,0 +1,79 @@
+% f1=figure;
+mypath='/mnt/Vol2/kropff/analyisis/figs/theta_new/fig2abc-openfield theta/clean_test/';
+files=dir([mypath 'peak_matr_final_f2_*.db']);
+
+close all;
+for ifile=14%1:numel(files) % 14 10
+p=load([mypath files(ifile).name],'-mat');%3 6 8
+data=p.out.data{1};
+
+%%
+% close all;
+f=figure;
+subplot(2,1,1)
+sel=1:numel(data.v);
+speeds=data.v(sel)';
+accs=[0 speeds(3:end)-speeds(1:end-2) 0]/.04;
+% ths=smooth_theta(accs,0.1); %0.06
+ths=smooth_theta2(accs,1-exp(-.1)); %0.06
+
+raccs=accs;
+raccs(raccs<0)=0;
+
+plot(sel/50,(speeds)/std(speeds),'r');
+hold on;
+
+
+plot(sel/50,accs/std(accs),'b');
+plot(sel/50,ths/std(accs),'g');
+xlim([0 10])
+
+
+%
+% sel2=true(size(speeds));
+ths=8+ths/100;
+da=20; abins=-5:5;
+dv=5; vbins=0:12;
+rv=round(speeds/dv);
+ra=round(accs/da);
+vfr=arrayfun(@(x) mean(ths(rv==x)),vbins);
+afr=arrayfun(@(x) mean(ths(ra==x)),abins);
+m8=vfr(1);%mean(ths);
+[th_f5,th_s5]=mybp(ths-m8,0.5);
+[th_f1,th_s1]=mybp(ths-m8,0.1);
+[th_f3,th_s3]=mybp(ths-m8,0.3);
+
+vfr_f=arrayfun(@(x) mean(th_f5(rv==x))+m8,vbins);
+vfr_s=arrayfun(@(x) mean(th_s5(rv==x))+m8,vbins);
+
+afr_f=arrayfun(@(x) mean(th_f1(ra==x))+m8,abins);
+afr_s=arrayfun(@(x) mean(th_s1(ra==x))+m8,abins);
+
+
+subplot(2,2,3)
+plot(vbins*dv,[vfr; vfr_s; vfr_f]');
+axis([0 50 7.8 9.1])
+golden(2.5)
+
+subplot(2,2,4)
+plot(abins*da,[afr; afr_s; afr_f]');
+axis([-120 140 7.8 9.1])
+set(gca,'YTick',6:.5:12);
+
+
+
+[corr(speeds',ths') corr(accs',ths');
+corr(speeds',th_f5') corr(accs',th_f1');
+corr(speeds',th_s5') corr(accs',th_s1')]
+drawnow;
+end
+
+
+
+
+function [th_f,th_s]=mybp(ths,f1)
+    th_f=real(bandpass_and_hilbert(ths',f1,100,50)');
+    th_f=abs(th_f).*cos(angle(th_f));
+    th_s=bandpass_and_hilbert(ths',0,f1,50)';
+    th_s=abs(th_s).*cos(angle(th_s));
+end

code/Fig_6/bandpass_and_hilbert.m

@@ -0,0 +1,65 @@
+%%% copyright 2012, Emilio Kropff
+%%% output:  first dim = number of frequencies
+%%%         second dim = number of trials
+%%%          third dim = number of timestamps
+
+
+
+function h=bandpass_and_hilbert(x,flow,fhigh,sampling)
+    l=size(x,1); 
+
+    NFFT= 2^nextpow2(l); % power of 2 is a faster fft 
+
+    x_fft=fft(x,NFFT);
+
+    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+    % Hilbert Transform, taken from the MATLAB function hilbert(). h is
+    % the 'filter' applied on the fft of the signal. The filter
+    % consists in multiplying the fft by i*sign(freq). But the output
+    % is hilbert(x)=x-i*transform(x), which results in a total filter that multyplies
+    % by 2 positive frequencies and by 0 negative frequencies (and by 1 the
+    % two central frequencies).
+    % http://en.wikipedia.org/wiki/Hilbert_transform#Relationship_with_
+    % the_Fourier_transform
+    hilbert_filter  = zeros(NFFT,1); 
+    % NFFT is even (it is a power of 2) so there are two central frequencies
+    hilbert_filter([1 NFFT/2+1]) = 1; 
+    hilbert_filter(2:NFFT/2) = 2;   %positive freqs (upper part of the circle)        
+
+    
+    fft_freqs=sampling/2*(linspace(0,1,NFFT/2+1))'; %This includes only the upper part of the circle
+    
+    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+    % Bandpass Filter
+    bandpass_filter=bsxfun(@gt,fft_freqs,flow) & bsxfun(@lt,fft_freqs,fhigh);
+    bandpass_filter=[bandpass_filter; bandpass_filter(end-1:-1:2,:)]; 
+
+    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+    % Finally, perform the bandpass and hilbert transform in one step
+    % (blocked)
+    
+%     nfreqmax=floor(.5*10^8/numel(x_fft));
+    nfreqmax=floor(.1*10^8/numel(x_fft));
+    fstart=1; fend=min(numel(flow),fstart-1+nfreqmax);
+    h=NaN(l,numel(flow));
+    while fstart<=fend
+        haux=ifft(bsxfun(@times,x_fft,bsxfun(@times,bandpass_filter(:,fstart:fend),hilbert_filter)));%combined_filter(:,ones(1,car.nlaps)));
+        h(:,fstart:fend)=haux(1:l,:);
+        fstart=fend+1;
+        fend=min(numel(flow),fstart-1+nfreqmax);
+    end
+
+%     h=squeeze(h);
+%     if nargout>1
+%         bpx=zeros(nf,nx,l);
+%         for f=1:nf
+%             bandpass_filter=(fft_freqs>flow).*(fft_freqs<fhigh);
+%             bandpass_filter=[bandpass_filter; bandpass_filter(end-1:-1:2)];             
+%             hilbert_transf=ifft(bsxfun(@times,x_fft,bandpass_filter));
+%             bpx(f,:,:)=hilbert_transf(1:l,:);
+%         end
+%         bpx=squeeze(bpx);
+%     end    
+
+
+ 

code/Fig_6/data_fig_6.db

@@ -0,0 +1 @@
+/annex/objects/MD5-s12799163--20e6d1d287bbcaba6b3df7730dec5ee5

code/Fig_6/smooth_theta2.m

@@ -0,0 +1,12 @@
+function th=smooth_theta2(acc0,decay)
+
+racc=acc0;
+racc(acc0<0)=0;
+th=racc(1)*ones(size(racc));
+for i=2:numel(racc)
+    if racc(i)*(1+decay)>th(i-1) 
+        th(i)=racc(i); 
+    else        
+        th(i)=th(i-1)*(1-decay);    
+    end
+end