load('Figure4S1A_Data'); % same as Figure1_Data.mat
       
%% SYNCH to birth or G1/S
chArea=2; chGem=1; ctrl=1; sns=4; 
for nn=1:height(fulltracks)
    L=size(fulltracks.FullTracks{nn,1},1); 
    clear GemRiseVec
    SynchToStart{nn}=ones(L,300)*NaN;
    SynchTo_G1S{nn}=ones(L,600)*NaN;
    SynchTo_G1S_GR{nn}=ones(L,600)*NaN;
    SynchTo_G1S_GR_SL{nn}=ones(L,600)*NaN;
    SynchToStart_G1{nn}=false(L,300);
    for ii=1:L
        kk=ii;
        xdata=fulltracks.FullTracks{nn,1}(kk).Feat.Ti-fulltracks.FullTracks{nn,1}(kk).Feat.Ti(1)+1;
        ydata=fulltracks.FullTracks{nn,1}(kk).Feat.Area;
        zdata=fulltracks.FullTracks{nn,1}(kk).GrowthRate.GR;
        zzdata=fulltracks.FullTracks{nn,1}(kk).GrowthRateSL;
        xdata=xdata(4:end-4); ydata=ydata(4:end-4); 
        gemrise=round(fulltracks.FullTracks{nn,1}(kk).Gemrise);
        if length(zdata)>length(ydata) && gemrise>5 && size(ydata,1)>gemrise
            zdata=zdata(1:length(ydata));
            GemRiseVec(ii)=gemrise;
            L2=length(ydata); L3=length(zzdata);
            SynchTo_G1S_GR{nn}(ii,(302-gemrise):300)=zdata(1:(gemrise-1)); % GR from birth to 300 (G1/S)
            SynchTo_G1S_GR{nn}(ii,301:(301+L2-gemrise))=zdata(gemrise:end); % GR from 301 (G1/S) to end
            SynchTo_G1S_GR_SL{nn}(ii,(302-gemrise):300)=zzdata(1:(gemrise-1)); % GR from birth to 300 (G1/S)
            SynchTo_G1S_GR_SL{nn}(ii,301:(301+L3-gemrise))=zzdata(gemrise:end); % GR from 301 (G1/S) to end
        end
        if gemrise>5 && size(ydata,1)>gemrise
            GemRiseVec(ii)=gemrise;
            L2=length(ydata);
            SynchToStart_G1{nn}(ii,1:gemrise)=true(1,gemrise); % logical: 1 for G1, 0 for after
            SynchToStart{nn}(ii,1:L2)=ydata'; % nuclear area from birth forward
            SynchTo_G1S{nn}(ii,(302-gemrise):300)=ydata(1:(gemrise-1)); % area from birth to 300 (G1/S)
            SynchTo_G1S{nn}(ii,301:(301+L2-gemrise))=ydata(gemrise:end); % area from 301 (G1/S) to end
            SynchToStart1=SynchToStart{nn};
            SynchToStart2=SynchToStart{nn};
            SynchToStart1(~SynchToStart_G1{nn})=NaN; % make cells not in G1 into nan
            SynchToStart2(SynchToStart_G1{nn})=NaN; % make cells in G1 into nan
        end
    end  
end
%% plot
tt=[];
for kk=1:size(SynchTo_G1S{ctrl},1)
    temp=SynchTo_G1S{ctrl,1}(kk,:);
    sz=temp(290:300);
    tt=[tt;sz];
end
tt=sort(tt);
prctSMctrl=prctile(tt(:,1),10);
prctBGctrl=prctile(tt(:,1),90);

sma=[]; big=[];
for kk=1:size(SynchTo_G1S{ctrl},1)
    temp=SynchTo_G1S{ctrl,1}(kk,:);
    sz=temp(290);
    if sz<=prctSMctrl
        sma=[sma;temp];
    elseif sz>=prctBGctrl
        big=[big;temp];
    else
        continue
    end
end
% r=randperm(size(sma,1),37);
% rr=randperm(size(big,1),37); 
figure; hold on
for kk=1:36
    ss = kk; bb = kk; %ss=r(kk); bb=rr(kk);
    plot(1:600,smoothdata(sma(ss,:),'movmean'),'b-','linewidth',1.5);
    plot(1:600,smoothdata(big(bb,:),'movmean'),'r-','linewidth',1.5);
end
axis([175 400 350 1200]);line([300 300],get(gca,'ylim'),'color','k'); 
set(gca,'ytick',400:100:1200,'xtick',[180 220 260 300 340 380],'xticklabel',-25:10:20,'fontsize',14);
xlabel('Time (hrs from G1/S)'); ylabel('Nuclear Area (pix.)')

