clear all;
%To use this code, just enter your input file prefix here - 
filename = '081506 mCherry Comparison';
%and the path to where your data should be saved.  If you don't specify an
%output location, it will  use the current directory.
% outputpath = '~/Desktop/Chassis-System/Data/PlateReader/DedicatedTranslation/';


[TimeData, LabelData] = platereaderread2MG(filename);

%%
%============
%set smoothing parameters
%============
method = 'fb';
if method == 'fb' | method == 'ff'
    cutofffreq = 0.15;
    neighbors = [];
elseif method == 'co'
    cutofffreq = 0.15;
    neighbors = 9;
elseif method == 'ma'
    cutofffreq = [];
    neighbors = 9;
end
plotunsmooth = 1;
plotsmooth = 0; 

Absbackgroundsamples = [30];
GFPbackgroundsamples = [27];
Mcherrybackgroundsamples = [26];

experimentsamples = [1:10];
names = ['RBS 2.0 - 1A3 in C86-Induced         ';...
         'RBS 2.0 - 1A3 in C86-Uninduced       ';...
         'RBS 0.3 - 1A3 in C86-Induced         ';...
         'RBS 0.3 - 1A3 in C86-Uninduced       ';...
         'RBS 8.0 - 1A3 in C86-Induced         ';...
         'RBS 8.0 - 1A3 in C86-Uninduced       ';...
         'RBS 4B - 1A3 in C86-Induced          ';...
         'RBS 4B - 1A3 in C86-Uninduced        ';...
         'RBS 2.8 - 1A3 in C86-Induced         ';...
         'RBS 2.8 - 1A3 in C86-Uninduced       ';...
         'Non RBS HK - 1A3 in C86-Induced      ';...
         'Non RBS HK - 1A3 in C86-Uninduced    ';...  
         'RBS 4B - 1A3 in MG1655-Induced       ';...
         'RBS 4B - 1A3 in MG1655-Uninduced     ';...
         'RBS 8.0 - 1A3 in MG1655-Induced      ';...
         'RBS 8.0 - 1A3 in MG1655-Uninduced    ';...
         'RBS 4B - 1A3 in X90-Induced          ';...
         'RBS 4B - 1A3 in X90-Uninduced        ';...
         'RBS 8.0 - 1A3 in X90-Induced         ';...
         'RBS 8.0 - 1A3 in X90-Uninduced       ';...
         'RBS 4B - 1A3 in X90(F-)-Induced      ';...
         'RBS 4B - 1A3 in X90(F-)-Uninduced    ';...
         'RBS 8.0 - 1A3 in X90(F-)-Induced     ';...
         'RBS 8.0 - 1A3 in X90(F-)-Uninduced   ';...
         'MG1655                               ';...
         'C86                                  ';...
         'X90                                  ';...    
         'EZ Media Plain                       ';...
         'EZ Media + Amp                       ';...
         'Ez Media + Amp + IPTG                '];
     wellnames = cellstr(names(experimentsamples,:));

%  All RBS's C86 Induced and Uninduced: 1:12,25:27
%  All RBS's C86 Induced: 1,3,5,7,9,11,25:27
%  All RBS's C86 Uninduced:  2,4,6,8,10,12,25:27
%  All RBS's MG1655 Induced and Uninduced: 13:16,25:27
%  All RBS's MG1655 Induced: 13,15,25:27
%  All RBS's MG1655 Uninduced:14,16,25:27
%  All RBS's X90 Induced and Uninduced: 17:20,25:27
%  All RBS's X90 Induced: 17,19,25:27
%  All RBS's X90 Uninduced:18,20,25:27
%  All RBS's X90(F-) Induced and Uninduced: 21:24,25:27
%  All RBS's X90(F-) Induced: 21,23,25:27
%  All RBS's X90(F-) Uninduced:22,24,25:27
%  RBS 4B Induced and Uninduced: 7,8,13,14,17,18,21,22,25:27
%  RBS 4B Induced: 7,13,17,21,25:27
%  RBS 4B Uninduced: 8,14,18,22,25:27
%  RBS 8 Induced and Uninduced: 5,6,15,16,19,20,23,24,25:27
%  RBS 8 Induced: 5,15,19,23,25:27
%  RBS 8 Uninduced: 6,16,20,24,25:27
%  Empty Strains: 25:27
%  Media Only: 28:30

% C86 Background Samples
% Absbackgroundsamples = [29];
% GFPbackgroundsamples = [26];
% Mcherrybackgroundsamples = [26];
% 
% X90 Background Samples
% Absbackgroundsamples = [29];
% GFPbackgroundsamples = [27];
% Mcherrybackgroundsamples = [27];
% Note: EZ (28) was contaminated, so I used EZ + Amp + IPTG (30). 

% Time Interval Used
% timepoints = [19:99];

replicates = 3;
samples = [experimentsamples,GFPbackgroundsamples,Mcherrybackgroundsamples,Absbackgroundsamples];
% Extract the data that we are interested in into three new variables.  By
% default, consider all timepoints and wells
% timepoints = [1:size(LabelData,1)-1];
wells = [1:size(LabelData,2)];
timepoints = [19:99];
% wells = [1:6];

rawGFP = LabelData(timepoints,wells,2);
rawMcherry = LabelData(timepoints,wells,3);
rawAbs = LabelData(timepoints,wells,1);
Time = TimeData(timepoints);

%If your data has a fixed number of replicates for each sample and control,
%and the replicates are in adjacent wells, use this code to reshape both 
%matrices so that the third dimension represents the different replicates
orderedAbs = reshape(rawAbs,[size(rawAbs,1),replicates,size(rawAbs,2)/replicates]);
orderedGFP = reshape(rawGFP,[size(rawGFP,1),replicates,size(rawGFP,2)/replicates]);
orderedMcherry = reshape(rawMcherry,[size(rawMcherry,1),replicates,size(rawMcherry,2)/replicates]);

%reorder the dimensions of the array so that the third dimension is the
%number of replicates.
orderedAbs = permute(orderedAbs,[1,3,2]);
orderedGFP = permute(orderedGFP,[1,3,2]);
orderedMcherry = permute(orderedMcherry,[1,3,2]);
orderedAbs = orderedAbs(:,samples,:);
orderedGFP = orderedGFP(:,samples,:);
orderedMcherry = orderedMcherry(:,samples,:);

%Background subtraction
%If you are subtracting one time course from another as a background
%subtraction, it would be better to wait until after the data has been
%filtered.  For now, we'll just subtract the mean of a media only well from
%the absorbance.
   
if ~isempty(Absbackgroundsamples)
    Absbackground = mean(orderedAbs(:,end,:),3);
    cellAbs = orderedAbs - repmat(Absbackground,[1,size(orderedAbs,2),size(orderedAbs,3)]);
   
else 
    cellAbs = orderedAbs;
    
end

cellOD = 3.113.*cellAbs-0.0158;

% cellOD = cellAbs;
%==============================
%Data Smoothing
%Put the data to be smoothed into a cell array so the function can handle
%arbitrary numbers of labels.
if ~isempty(method)
    timeseries{1} = cellOD;
    timeseries{2} = orderedGFP;
    timeseries{3} = orderedMcherry;
    [smoothtimeseries] = platereadersmooth(timeseries,method,cutofffreq,neighbors);
    %now get the smoothed data back out into individual arrays.
    fOD = smoothtimeseries{1};
    OD = cellOD;
    fGFP = smoothtimeseries{2};
    fMcherry = smoothtimeseries{3};
    GFP = orderedGFP;
    Mcherry = orderedMcherry;
else
    fOD = cellOD;
    fGFP = orderedGFP;
    fMcherry = orderedMcherry;
end


%===========================
%Background Subtraction
%You should probably be subtracting the fluorescence of cells not
%expressing GFP so the GFP signal is truly proportional to GFP.  This
%implementation assumes the OD of your samples is similar to that of the
%control at each time point.
% backgroundsample = size(fGFP,2);
if ~isempty(GFPbackgroundsamples)
    meanbackground = mean(fGFP(:,end-2,:),3);
%     meanbackground = mean(fGFP(:,end-2,:),3);
    %Need to make this array the same size as the full array so we can subtract
    %arrays rather than using slow looping.
    meanbackgroundarray = repmat(meanbackground,[1,size(fGFP,2),size(fGFP,3)]);
    correctedGFP = fGFP - meanbackgroundarray;
    
    
    umeanbackground = mean(GFP(:,end-2,:),3);
%     umeanbackground = mean(GFP(:,end-2,:),3);
    %Need to make this array the same size as the full array so we can subtract
    %arrays rather than using slow looping.
    umeanbackgroundarray = repmat(umeanbackground,[1,size(GFP,2),size(GFP,3)]);
    ucorrectedGFP = GFP - umeanbackgroundarray;
    
else
    correctedGFP = fGFP;
    
    
    ucorrectedGFP = GFP;
    
end

if ~isempty(Mcherrybackgroundsamples)
    meanbackground2 = mean(fMcherry(:,end-1,:),3);
    %Need to make this array the same size as the full array so we can subtract
    %arrays rather than using slow looping.
    meanbackgroundarray2 = repmat(meanbackground2,[1,size(Mcherry,2),size(fMcherry,3)]);
    correctedMcherry = fMcherry - meanbackgroundarray2;
    
    
    umeanbackground2 = mean(Mcherry(:,end-1,:),3);
%     umeanbackground2 = mean(Mcherry(:,Mcherrybackgroundsamples,:),3);
    
    %Need to make this array the same size as the full array so we can subtract
    %arrays rather than using slow looping.
    umeanbackgroundarray2 = repmat(umeanbackground2,[1,size(Mcherry,2),size(Mcherry,3)]);
    ucorrectedMcherry = Mcherry - umeanbackgroundarray2;
    
else
    correctedMcherry = fMcherry;
    
    
    ucorrectedMcherry = Mcherry;
    
end

%===========================
%GFP and Mcherry Synthesis Rate Calculations

%find the amount of GFP per OD
GFPperOD = correctedGFP./fOD;
uGFPperOD = ucorrectedGFP./OD;
%Calculate the rate of accumulation of GFP per OD per unit time
dGFPperOD = diff(GFPperOD);
udGFPperOD = diff(uGFPperOD);
dt = mean(diff(Time));
dGFPperODdt = dGFPperOD./dt;
udGFPperODdt = udGFPperOD./dt;

%find the amount of Mcherry per OD
McherryperOD = correctedMcherry./fOD;
uMcherryperOD = ucorrectedMcherry./OD;
%Calculate the rate of accumulation of Mcherry per OD per unit time
dMcherryperOD = diff(McherryperOD);
udMcherryperOD = diff(uMcherryperOD);
dt = mean(diff(Time));
dMcherryperODdt = dMcherryperOD./dt;
udMcherryperODdt = udMcherryperOD./dt;

%Calculate the dilution rate, assuming exponential growth
dOD = diff(fOD);
udOD = diff(OD);

dODdt = dOD./dt;
udODdt = udOD./dt;

gammad = (1/log(2))*(dODdt./fOD(2:end,:,:));
ugammad = (1/log(2))*(udODdt./OD(2:end,:,:));

%calculate the degradation rate of GFP per OD
degradation = gammad.*GFPperOD(2:end,:,:);
udegradation = ugammad.*uGFPperOD(2:end,:,:);
%The synthesis rate is the sum of the accumulation rate and the degradation
%rate.
synthesis = dGFPperODdt + degradation;
usynthesis = udGFPperODdt + udegradation;

%calculate the degradation rate of Mcherry per OD
degradation2 = gammad.*McherryperOD(2:end,:,:);
udegradation2 = ugammad.*uMcherryperOD(2:end,:,:);
%The synthesis rate is the sum of the accumulation rate and the degradation
%rate.
synthesis2 = dMcherryperODdt + degradation2;
usynthesis2 = udMcherryperODdt + udegradation2;

%At this point, we should average each of the replicates
meansynthesis = mean(synthesis,3);
umeansynthesis = mean(usynthesis,3);

meansynthesis2 = mean(synthesis2,3);
umeansynthesis2 = mean(usynthesis2,3);


unsmoothGFP = mean(orderedGFP,3);
unsmoothMcherry = mean(orderedMcherry,3);
unsmoothOD = mean(cellOD,3);
smoothGFP = mean(fGFP,3);
smoothMcherry = mean(fMcherry,3);
smoothOD = mean(fOD,3);


meanGFPperOD = mean(GFPperOD,3);
umeanGFPperOD = mean(uGFPperOD,3);

meanMcherryperOD = mean(McherryperOD,3);
umeanMcherryperOD = mean(uMcherryperOD,3);

meandGFPperODdt = mean(dGFPperODdt,3);
umeandGFPperODdt = mean(udGFPperODdt,3);

meandMcherryperODdt = mean(dMcherryperODdt,3);
umeandMcherryperODdt = mean(udMcherryperODdt,3);

meangammad = mean(gammad,3);
umeangammad = mean(ugammad,3);

meandegradation = mean(degradation,3);
umeandegradation = mean(udegradation,3);

meandegradation2 = mean(degradation2,3);
umeandegradation2 = mean(udegradation2,3);

%=======================
%Data Plotting

%Define some figure handles
raw = 1;
processed = 2;
surface = 4;
ODwindow = smoothOD > 0.1 & smoothOD < 0.3;
figure1columns = size(experimentsamples+1,2);
figure2columns = size(experimentsamples,2);
figure3columns = size(experimentsamples,2);
% timewindow = 10:97;
set(0,'DefaultAxesColorOrder',jet(figure1columns))
figure(raw);
whatiwant = [0 4*10^4]
PropertyValue1 = 7
PropertyValue2 = 5

subplot(1,3,1)
    hold off;
    if plotunsmooth
        semilogy(Time,unsmoothOD(:,1:figure1columns),'o','MarkerSize',5)
        hold on;
     
    end
    
    if plotsmooth
      semilogy(Time,smoothOD(:,1:figure1columns),'LineWidth',1.5); 
    end
    
    title('Raw and smoothed OD','FontSize',PropertyValue1)
    xlabel('Time (mins)','FontSize',PropertyValue1)
    ylabel('OD600','FontSize', PropertyValue1);
legendhandle = legend(wellnames,0);
set(legendhandle, 'FontSize', PropertyValue2)
set(gca, 'FontSize', PropertyValue1);

subplot(1,3,2)
    hold off;
    if plotunsmooth
        plot(smoothOD(:,1:figure1columns),unsmoothGFP(:,1:figure1columns),'o','MarkerSize',5)
        hold on;
    end
    
    if plotsmooth
        plot(smoothOD(:,1:figure1columns),smoothGFP(:,1:figure1columns),'LineWidth',1.5)
    end
    title('Raw and smoothed GFP','FontSize',PropertyValue1)
    xlabel('OD600','FontSize',PropertyValue1)
    ylabel('GFP counts', 'FontSize',PropertyValue1);
    xlim([0 1.5]);
    ylim(whatiwant);
set(gca, 'FontSize', PropertyValue1);
    
    subplot(1,3,3)
    hold off;
    if plotunsmooth
        plot(smoothOD(:,1:figure1columns),unsmoothMcherry(:,1:figure1columns),'o','MarkerSize',5)
        hold on;
    end
    
    if plotsmooth
        plot(smoothOD(:,1:figure1columns),smoothMcherry(:,1:figure1columns),'LineWidth',1.5)
    end
    title('Raw and smoothed mCherry','FontSize',PropertyValue1)
    xlabel('OD600','FontSize', PropertyValue1)
    ylabel('Mcherry counts', 'FontSize', PropertyValue1);
    xlim([0 1.5]);
    ylim(whatiwant);
set(gca, 'FontSize', PropertyValue1);
    figure(processed(1));

subplot(3,2,1)
    hold off;
    if plotunsmooth
        plot(smoothOD(:,1:figure2columns),umeanGFPperOD(:,1:figure2columns),'o','MarkerSize',5)
        hold on;
    end
    
    if plotsmooth
        plot(smoothOD(:,1:figure2columns),meanGFPperOD(:,1:figure2columns),'LineWidth',1.5)
    end
    title('GFP per OD');
%ylim([-2000, 150000]);
    
subplot(3,2,2)
    hold off; 
    if plotunsmooth
        plot(smoothOD(2:end,1:figure2columns),umeandGFPperODdt(:,1:figure2columns),'o','MarkerSize',5)
        hold on;
      
    end
    if plotsmooth
            plot(smoothOD(2:end,1:figure2columns),meandGFPperODdt(:,1:figure2columns),'LineWidth',1.5)
    end
    title('Accumulation rate of GFP per OD');
     ylim([-200, 500]);

    subplot(3,2,3)
    hold off;
    if plotunsmooth
        plot(smoothOD(2:end,1:figure2columns),umeangammad(:,1:figure2columns).*60,'o','MarkerSize',5)
        hold on;
    end
    if plotsmooth
        plot(smoothOD(2:end,1:figure2columns),meangammad(:,1:figure2columns).*60,'LineWidth',1.5)
    end
    title('Dilution rate (doublings/hour)');
    ylim([-1,2]);

    subplot(3,2,4)
    hold off;
    if plotunsmooth
        plot(smoothOD(2:end,1:figure2columns),umeandegradation(:,1:figure2columns),'o','MarkerSize',5)
        hold on;
    end
    if plotsmooth
        plot(smoothOD(2:end,1:figure2columns),meandegradation(:,1:figure2columns),'LineWidth',1.5)
    end
   title('Degradation rate');
%ylim([-200,500])
    subplot(3,2,5)
    hold off;
    if plotunsmooth
        plot(smoothOD(2:end,1:figure2columns),umeansynthesis(:,1:figure2columns),'o','MarkerSize',5)
        hold on;
    end
    if plotsmooth
        plot(smoothOD(2:end,1:figure2columns),meansynthesis(:,1:figure2columns),'LineWidth',1.5)
    end
    title('Synthesis rate');
    %ylim([-200, 700]);
legend(wellnames,0);

%Figure 3 - Mcherry Plots

figure(3)

subplot(3,2,1)
    hold off;
    if plotunsmooth
        plot(smoothOD(:,1:figure3columns),umeanMcherryperOD(:,1:figure3columns),'o','MarkerSize',5)
        hold on;
    end
    
    if plotsmooth
        plot(smoothOD(:,1:figure3columns),meanMcherryperOD(:,1:figure3columns),'LineWidth',1.5)
    end
    title('Mcherry per OD');
%ylim([-2000, 150000]);
    
    subplot(3,2,2)
    hold off; 
    if plotunsmooth
        plot(smoothOD(2:end,1:figure3columns),umeandMcherryperODdt(:,1:figure3columns),'o','MarkerSize',5)
        hold on;
       ylim([-100, 200]);
    end
    if plotsmooth
            plot(smoothOD(2:end,1:figure3columns),meandMcherryperODdt(:,1:figure3columns),'LineWidth',1.5)
    end
    
    title('Accumulation rate of Mcherry per OD');
       
subplot(3,2,3)
    hold off;
    if plotunsmooth
        plot(smoothOD(2:end,1:figure3columns),umeangammad(:,1:figure3columns).*60,'o','MarkerSize',5)
        hold on;
    end
    if plotsmooth
        plot(smoothOD(2:end,1:figure3columns),meangammad(:,1:figure3columns).*60,'LineWidth',1.5)
    end
    title('Dilution rate (doublings/hour)');
    ylim([-1,2]);

    subplot(3,2,4)
    hold off;
    if plotunsmooth
        plot(smoothOD(2:end,1:figure3columns),umeandegradation2(:,1:figure3columns),'o','MarkerSize',5)
        hold on;
    end
    if plotsmooth
        plot(smoothOD(2:end,1:figure3columns),meandegradation2(:,1:figure3columns),'LineWidth',1.5)
    end
   title('Mcherry Degradation rate');
%ylim([-200,500])
    
    subplot(3,2,5)
    hold off;
    if plotunsmooth
        plot(smoothOD(2:end,1:figure3columns),umeansynthesis2(:,1:figure3columns),'o','MarkerSize',5)
        hold on;
    end
    if plotsmooth
        plot(smoothOD(2:end,1:figure3columns),meansynthesis2(:,1:figure3columns),'LineWidth',1.5)
    end
    title('Mcherry Synthesis rate');
    %ylim([-200, 700]);
 legend(wellnames,0);

%  figure(4)
%     subplot(1,1,1)
%     if plotunsmooth
%         plot3(Time(1:98),umeandMcherryperODdt(:,1:figure3columns),umeandGFPperODdt(:,1:figure2columns),'.','MarkerSize',20)
%     end
%     if plotsmooth
%         plot3(Time(1:98),meandMcherryperODdt(:,1:figure3columns),meandGFPperODdt(:,1:figure2columns),'LineWidth',1.5)
%     end
%     title('mCherry vs. GFP Levels')
%     xlabel('Time','FontSize', PropertyValue1)
%     ylabel('mCherry Accumulation Rate', 'FontSize', PropertyValue1);
%     zlabel('GFP Acccumulation Rate', 'FontSize', PropertyValue1);
%     axis ij
%     grid on
%     ylim([-20, 120])
%     zlim([-20, 20])
%     legend(wellnames,0);
%     legendhandle = legend(wellnames,0);
%     set(legendhandle, 'FontSize', PropertyValue2)
    
    figure(4)
    subplot(1,1,1)
    if plotunsmooth
        plot3(Time(:),umeanMcherryperOD(:,1:figure3columns),umeanGFPperOD(:,1:figure2columns),'LineWidth',1.5)
%         plot3(Time(:),umeanMcherryperOD(:,1:figure3columns),umeanGFPperOD(:,1:figure2columns),'.','MarkerSize',20)    
    end
    if plotsmooth
        plot3(Time(:),meanMcherryperOD(:,1:figure3columns),meanGFPperOD(:,1:figure2columns),'LineWidth',1.5)
    end
    title('mCherry vs. GFP Levels')
    xlabel('Time','FontSize', PropertyValue1)
    ylabel('mCherry per OD', 'FontSize', PropertyValue1);
    zlabel('GFP per OD', 'FontSize', PropertyValue1);
    axis ij
    grid on
    xlim([0, 1200])
    ylim([0, 40000])
    zlim([-3000, 20000])
    legend(wellnames,0);
    legendhandle = legend(wellnames,0);
    set(legendhandle, 'FontSize', PropertyValue2)
    
%     figure(4)
%     subplot(1,1,1)
%     if plotunsmooth
%         mesh(umeanMcherryperOD(:,1),umeanGFPperOD(:,1),'.','MarkerSize',20)
%     end
%     if plotsmooth
%         mesh(meandMcherryperODdt(:,1),meandGFPperODdt(:,1),'LineWidth',1.5)
%     end
%     title('mCherry vs. GFP Levels')
%     xlabel('Time','FontSize', PropertyValue1)
%     ylabel('mCherry per OD', 'FontSize', PropertyValue1);
%     zlabel('GFP per OD', 'FontSize', PropertyValue1);
%     axis ij
%     grid on
%     xlim([400, 1200])
%     ylim([0, 20000])
%     zlim([-3000, 20000])
%     legend(wellnames,0);
%     legendhandle = legend(wellnames,0);
%     set(legendhandle, 'FontSize', PropertyValue2)

% figure(4)
% subplot(1,1,1)
%     if plotunsmooth
%     
%     tri = delaunay(Time,umeanMcherryperOD(:,1))
% 
%     trisurf(tri,Time,umeanMcherryperOD(:,1),zeros(size(Time)))
% end
%  
   

% figure(surface);
% timeaxis = Time(2:end);
% samplestoplot = [1:4:29];
% [xx,yy] = meshgrid(timeaxis,1:length(samplestoplot));
% mesh(xx,yy,meansynthesis(:,samplestoplot)');
% xlabel('Time (mins)')
% ylabel('Plotted Sample #');
%

indices = fOD>0.1 & fOD<0.3;

numberofsamplesinrange = sum(indices);
correspondingsynthesisrates = synthesis.*indices(2:end,:,:);
meansynthesisrates = sum(correspondingsynthesisrates)./numberofsamplesinrange;

synthesisratesminusmean = synthesis-repmat(meansynthesisrates,[size(synthesis,1),1,1]);
synthesisratesminusmean = synthesisratesminusmean.*indices(2:end,:,:);
synthesisratesminusmeansquared = synthesisratesminusmean.^2;
varsynthesisrates = sum(synthesisratesminusmeansquared)./numberofsamplesinrange;
sdsynthesisrates = varsynthesisrates.^0.5;
sesynthesisrates = sdsynthesisrates./sqrt(numberofsamplesinrange);
samplemeansynthesisrates = squeeze(mean(meansynthesisrates,3));
samplestsynthesisrates = squeeze(mean(sdsynthesisrates,3));
samplesesynthesisrates = squeeze(mean(sesynthesisrates,3));

% figure


% [b,a] = butter(2,0.125);
% label = 1;
% d = diff(LabelData(timepoints,wells,label));
% fd = filter(b,a,d);
% f = filter(b,a,LabelData(timepoints,wells,label) - mean(LabelData(1,wells,label)));
% f = f + mean(LabelData(1,wells,label));
% s = smooth(LabelData(timepoints,wells(1),label),13,'moving');
% noise = LabelData(timepoints,wells(1),label) - s;
% 
% 
% df = diff(f);
% % wells = [48];
% fftplot = 3;
% figure(fftplot)
% hold off;
% subplot(1,3,1)
% hold off;
% PlotFFT(noise(:,1),1/mean(diff(TimeData)));
% % PlotFFT(diff(LabelData(timepoints,wells(1),label)),1/mean(diff(TimeData)));
% hold on;
% % PlotFFTred(LabelData(timepoints,wells(1),label)-mean(LabelData(timepoints
% % ,wells(2),label)),1/mean(diff(TimeData)));
% 
% subplot(1,3,2)
% hold off;
% plot(TimeData(timepoints(1:end-1)),d,'.-b');
% hold on;
% plot(TimeData(timepoints(1:end-1)),df,'.-r');
% plot(TimeData(timepoints(1:end-1)),fd,'.-g');
% subplot(1,3,3)
% plot(TimeData(timepoints),LabelData(timepoints,wells,label),'.-b');
% hold on;
% plot(TimeData(timepoints),f,'.-r');
%%
% plot(fdifferentials2,'-g');
% fdifferentials(50:end) = 0;
% differentials2 = ifft(fdifferentials,512);
% differentials2 = abs(differentials2(1:143));


% sabsorbance = smooth(LabelData(:,wells,1),17,'rloess');
% sfluorescence = smooth(LabelData(:,wells,2),17,'rloess');
% snormalized = sfluorescence./sabsorbance;
% sdifferentials = diff(snormalized)./mean(diff(TimeData));


% s2absfft = fft(LabelData(:,wells,1),512);
% s2absfft(255:end) = 0;
% s2absorbance = ifft(s2absfft,512);
% s2absorbance = s2absorbance(1:length(TimeData));
% s2flufft = fft(LabelData(:,wells,2),512);
% s2flufft(255:end) = 0;
% s2fluorescence = ifft(s2flufft,512);
% s2fluorescence = s2fluorescence(1:length(TimeData));
% 
% s2normalized = s2fluorescence./s2absorbance;
% s2differentials = diff(s2normalized)./mean(diff(TimeData));

% platereaderplotting2(TimeData,LabelData(:,wells,1),LabelData(:,wells,2),normalized,differentials)
% % figure
% % platereaderplotting2(TimeData,sabsorbance,sfluorescence,snormalized,sdifferentials)
% figure
% platereaderplotting2(TimeData,LabelData(:,wells,1),LabelData(:,wells,2),normalized,fdifferentials)

% platereaderplotting2(TimeData,s2absorbance(1:length(TimeData)),s2fluorescence,s2normalized,s2differentials)

%legend(LabelNames(wells,:),'Location','EastOutside');


%[dadt] = calculategrowthrate(time, absorbance)
%%
% figure
% transformed = fft(LabelData(:,wells,1),512);
% P = transformed.*conj(transformed)/512;
% f = 0.2*(0:256)/512;
% plot(f,P(1:257))
% %%
% diffs = fft(differentials,512);
% P3 = diffs.*conj(diffs)/512;
% %figure
% plot(f,P3(1:257))
% [b,a] = butter(2,0.15)
% filtered = filter(b,a,differentials);
% 
% range = 132:138;
% diffs(range) = 0;
% diffs(512-range) = 0;
% rev = ifft(diffs,512);
% figure
% hold on;
% plot(TimeData,rev(1:length(TimeData)),'-g')
% %plot(TimeData(1:end-1),differentials,'-r')
% plot(TimeData(1:end-1),filtered,'-b')

% windowSize = 5;
% manAbs = filter(ones(1,windowSize)/windowSize,1,nAbs(:,wells)-repmat(nAbs(1,wells),length(timepoints),1));
% manAbs = manAbs+repmat(nAbs(1,wells),length(timepoints),1);
% 
% manGFP = filter(ones(1,windowSize)/windowSize,1,nGFP(:,wells)-repmat(nGFP(1,wells),length(timepoints),1));
% manGFP = manGFP+repmat(nGFP(1,wells),length(timepoints),1);
% 
% sAbs(:,1) = smooth(nAbs(:,wells(1)),5,'moving');
% sAbs(:,2) = smooth(nAbs(:,wells(2)),5,'moving');
% 
% sGFP(:,1) = smooth(nGFP(:,wells(1)),5,'moving');
% sGFP(:,2) = smooth(nGFP(:,wells(2)),5,'moving');
% 
% figure(5);
% 
% subplot(1,2,1)
% hold off;
% plot(Time,nGFP(:,wells),'-b');
% hold on;
% plot(Time,fGFP,'-r');
% plot(Time,manGFP,'-g');
% plot(Time,sGFP,'-k');
% 
% subplot(1,2,2)
% hold off;
% plot(Time,nAbs(:,wells),'-b');
% hold on;
% plot(Time,fAbs,'-r');
% plot(Time,manAbs,'-g');
% plot(Time,sAbs,'-k');

% nGFP = LabelData(timepoints,:,2) - repmat(LabelData(1,:,2),length(timepoints),1);


% [b,a] = butter(2,0.15);
% fGFP = filter(b,a,nGFP(:,wells)-repmat(nGFP(1,wells),length(timepoints),1));
% fGFP = fGFP+repmat(nGFP(1,wells),length(timepoints),1);
% 
% [b2,a2] = butter(2,0.3);
% fAbs = filter(b2,a2,nAbs(:,wells)-repmat(nAbs(1,wells),length(timepoints),1));
% fAbs = fAbs+repmat(nAbs(1,wells),length(timepoints),1);

% noise = diff(nAbs(:,wells));
% figure(10)
% plot(Time(1:end-1),noise);
% figure(11)
% PlotFFT(noise(:,1),1/mean(diff(TimeData)));
% hold on;
% PlotFFTred(noise(:,2),1/mean(diff(TimeData)));
% PlotFFT(noise(:,3),1/mean(diff(TimeData)));
% PlotFFTred(noise(:,4),1/mean(diff(TimeData)));
% PlotFFT(noise(:,5),1/mean(diff(TimeData)));
% PlotFFTred(noise(:,6),1/mean(diff(TimeData)));