function [ractive rinactive bactive binactive output] = ANDgatetest()


n=20;
m=20;
rdurv=linspace(0,1000,n);
rintv=linspace(0,1,m);
[X,Y]=meshgrid(rdurv,rintv);

ractive=zeros(n,m);
rinactive=zeros(n,m);
bactive=zeros(n,m);
binactive=zeros(n,m);
output=zeros(n,m);

for j=1:n
    dur=rdurv(j);
    for k=1:m
        rint=rintv(k);
        %display(['----------------------------------------------------']);
        %display(['--------- Emission Duration ' num2str(dur) ' seconds ---------']);
        %display(['----------------------------------------------------']);


        % 1 min 60
        % 1 h 3600
        % 1 d 86400 
        global p
        p=struct();
        %red light
        p.rstart=200;
        p.rdur=dur;
        p.rint=rint*10;
        %blue light
        p.bstart=200;
        p.bdur=dur/5;
        p.bint=rint/10;
        [t,y]=ANDGate();


        %% tRNA
        tRNAthresh=40;
        rstart=find(y(:,14)>tRNAthresh,1);
        rtend=t(rstart + find(y(rstart:end,14)<tRNAthresh,1)-1);
        if size(rstart-p.rstart,1)>0
            ractive(j,k)=rstart-p.rstart;
        else
            ractive(j,k)=-20;
        end
        if size(rtend-p.rstart-p.rdur,1)>0
            rinactive(j,k)=rtend-p.rstart-p.rdur;
        else
            rinactive(j,k)=-20;
        end
        %display(['REDLIGHT Time for activation: ' num2str(rtstart-p.rstart) ' s - Time for deactivation: ' num2str(rtend-p.rstart-p.rdur) ' s '])


        %% T7RNAP

        T7RNAPmrnathresh=1;
        bstart=find(y(:,16)>T7RNAPmrnathresh,1);
        btstart=t(bstart);
        btend=t(bstart + find(y(bstart:end,16)<T7RNAPmrnathresh,1)-1);
        if size(bstart-p.bstart,1)>0
            bactive(j,k)=btstart-p.bstart;
        else
            bactive(j,k)=-20;
        end
        if size(btend-p.bstart-p.bdur,1)>0
            binactive(j,k)=btend-p.bstart-p.bdur;
        else
            binactive(j,k)=-20;
        end
        %display(['BLUELIGHT Time for activation: ' num2str(btstart-p.rstart) ' s - Time for deactivation: ' num2str(btend-p.rstart-p.rdur) ' s '])

        %% Output
        output(j,k)=y(end,20);
        %display(['Final output of the color is ' num2str(output) ' nM'])
        
end
end

figure(1)
clf
surf(X',Y',ractive)
title('Activation Time Red Light')
ylabel('light intensity')
xlabel('exposure time (s)')
zlabel('activation time')

figure(2)
clf
surf(X',Y',rinactive)
title('Deactivation Time Red Light')
ylabel('light intensity')
xlabel('exposure time (s)')
zlabel('deactivation time')

figure(3)
clf
surf(X',Y',bactive)
title('Activation Time Blue Light')
ylabel('light intensity')
xlabel('exposure time (s)')
zlabel('activation time')

figure(4)
clf
surf(X',Y',binactive)
title('Deactivation Time Blue Light')
ylabel('light intensity')
xlabel('exposure time (s)')
zlabel('deactivation time')

figure(5)
clf
surf(X',Y',output)
title('Output of dye ')
ylabel('light intensity')
xlabel('exposure time (s)')
zlabel('dye output (nM)')

end