% Cartoon of Frank-van-der-Merwe thin-film growth being monitored by XRR
close all; clear

vid = VideoWriter('filmGrowthFvdMcartoon.mp4','MPEG-4');
vid.Quality = 100;
vid.FrameRate = 60;
open(vid);
figure(1)
set(gcf,'units','pixels','position',[0 0 1920 1080],'ToolBar','none');
set(0,'defaultfigurecolor',[1 1 1]);
set(gca,'linewidth',7);

lp = [-0.3 -0.5 0.8]; % Lighting direction 
ec = [0.25 0.28 0.5]; % Phil blue
fc = [1.0 0.83 0]; % Gold
dr = [0.5 0.0 0]; % Dark red

pos1 = [0 0 1 1]; 
LL = 1;
noCells = 40; % Number of cells

% Cylinder for detector
theta = 0:pi/180:2*pi;
circlex = sin(theta);
circlez = cos(theta);
[a,b,c] = cylinder(0.1,500);

[V,F] = platonic_solid(2,sqrt(3.00)); % Cube, side length = 1
V(:,3) = 0.14*V(:,3)-0.14; % Substrate, top surface (+/- 1)^2 @ z = 0

[V2,F2] = platonic_solid(2,sqrt(3)/noCells); % Unit cell of deposited material, edge length = 2/noCells
V2(:,3) = V2(:,3)+1/noCells; % Shift upwards by half a u.c. so bottom surface @ z = 0

% % Subplot of 3D cartoon of film growth being monitored by XRR
subplot('position',pos1);
% Draw substrate
ps = patch('Faces',F,'Vertices',V,'FaceColor',ec,...
    'FaceAlpha',1,'EdgeColor','none','FaceLighting','flat',...
    'DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0); 
hold on
% Draw detector 
detWall = surf(a,b,c/4+1.7,'FaceColor',[0.7 0.7 0.7],'LineStyle','none',...
    'FaceLighting','flat','DiffuseStrength',0.5);
rotate(detWall,[0 1 0],90-atand(0.25/1.5),[0,0,2/noCells]);
detEnd = fill3(0.1*circlex,0.1*circlez,0*circlez+1.95,[0.5 0.5 0.5],'LineStyle','none',...
    'FaceLighting','flat','DiffuseStrength',0.5,'SpecularStrength',0.5);
rotate(detEnd,[0 1 0],90-atand(0.25/1.5),[0,0,2/noCells]);


% Draw incident and reflected x-ray beams
kinArrow = mArrow3([-1.5 0 0.25],[0 0 0],'color',fc, ...
    'stemWidth',0.01,'tipWidth',0.025,'facealpha',1);
koutArrow = mArrow3([0 0 0],[1.5 0 0.25],'color',fc, ...
    'stemWidth',0.01,'tipWidth',0.025,'facealpha',0.5);

set(gca, 'Projection','perspective');
light('Position',lp,'Style','infinite');
axis equal;
axis off;
view(35,35);
xlim([-LL*1.5 LL*2]);
ylim([-LL LL]);
zlim([-0.25 0.6]);

ii = 0; % Controls when a frame is recorded
frame = getframe(gcf);
writeVideo(vid,frame);
S(1) = 0;
for jj = 1:2 % Two monolayers of growth
    kk = 0; % Runs from 0 to noCells^2 for each ML and describes the ML coverage
    fillArray = zeros(noCells);
    while (any(fillArray(:) == 0)) % Incomplete monolayer
        x1 = randi([1 noCells],1); % Random integer x-coordinate between 1 and noCells
        y1 = randi([1 noCells],1); % Random integer y-coordinate between 1 and noCells
        % Limit 3x3 box neighborhood to within boundaries of substrate
        if (x1 == 1)
            xlo = x1;
            xhi = x1 + 2;
        else
            xlo = x1 - 1;
            xhi = x1 + 1;
        end
        if (x1 == noCells)
            xhi = x1;
            xlo = x1 - 2;
        end
        if (y1 == 1)
            ylo = y1;
            yhi = y1 + 2;
        else
            ylo = y1 - 1;
            yhi = y1 + 1;
        end
        if (y1 == noCells)
            yhi = y1;
            ylo = y1 - 2;
        end
        % Neighborhood of position (x1,y1)
        neighborBox = fillArray(xlo:xhi,ylo:yhi);
        if (sum(neighborBox(:)) < (xhi - xlo + 1)*(yhi - ylo + 1)) % At least one unoccupied site in 3 x 3 neighborhood
            % Position list of unoccupied sites within 3 x 3 neighborhood 
            [row,column] = find(neighborBox==0);
            rowSize = size(row,1); % Number of unoccupied sites within 3 x 3 neighborhood 
            selectedUCidx = randi([1 rowSize]); % Pick one of unoccupied sites at random 
            % Shift 3 x 3 region by one u.c. if selected site is on the
            % edge of the substrate 
            if (x1 == 1)
                x1 = x1 + row(selectedUCidx)-1;
            elseif (x1 == noCells)
                x1 = x1 + row(selectedUCidx)-3;
            else
                x1 = x1 + row(selectedUCidx)-2;
            end
            if (y1 == 1)
                y1 = y1 + column(selectedUCidx)-1;
            elseif (y1 == noCells)
                y1 = y1 + column(selectedUCidx)-3;
            else
                y1 = y1 + column(selectedUCidx)-2;
            end
            kk = kk + 1; % Increment number of occupied sites by one
            V3(:,1) = -1 + 1/noCells + V2(:,1) + 2*(x1-1)/noCells;
            x1 = round(V3(1,1)*noCells/2+(1+noCells/2)); % x-coordinate of u.c. in units of u.c.
            V3(:,2) = -1 + 1/noCells + V2(:,2) + 2*(y1-1)/noCells;
            y1 = round(V3(1,2)*noCells/2+(1+noCells/2)); % y-coordinate of u.c. in units of u.c.
            V3(:,3) = V2(:,3)+(jj-1)*2/noCells; % z-coordinate of u.c. in units of u.c.
            % Draw u.c. in position (x1,y1) 
            ps2 = patch('Faces',F2,'Vertices',V3,'FaceColor',dr*(1+(jj-1)/8),...
                'FaceAlpha',1,'EdgeColor','none','FaceLighting','flat',...
                'DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);
            fillArray(x1,y1) = 1; % Flag this site as occupied            
        end
        
        delete(kinArrow);
        delete(koutArrow);
        
        arrowHt = kk*(2/noCells)/noCells^2 + (jj-1)*2/noCells;
        kinArrow = mArrow3([-1.5 0 0.25+arrowHt],[0 0 0+arrowHt],'color',fc, ...
            'stemWidth',0.01,'tipWidth',0.025,'facealpha',1);
        roughNow = 1.4*(1 - abs((noCells^2/2)-kk)/(noCells^2/2));
        koutAlpha = exp(-roughNow.^2);        
        koutArrow = mArrow3([0 0 0+arrowHt],[1.5 0 0.25+arrowHt],'color',fc, ...
            'stemWidth',0.01,'tipWidth',0.025,'facealpha',koutAlpha);
        ii = ii + 1; % Number of attempts within an incomplete monolayer to complete that monolayer
        % Add frame to video every second attempt to fill a u.c. in order
        % to speed things up, even at 60 Hz frame rate 
        if (mod(ii,2) == 0) 
            frame = getframe(gcf);
            writeVideo(vid,frame);
            S(ii/2) = kk;
        end
    end
end
frame = getframe(gcf);
writeVideo(vid,frame);

noFrames = floor(ii/2);

% _________________________________________________________________________
% Plot out of reflectivity and monolayer coverage 

close all
clear vid
figure(2)
set(gcf,'ToolBar','none');

vid = VideoWriter('filmGrowthFvdMxrrSignal.mp4','MPEG-4');
vid.Quality = 100;
vid.FrameRate = 60;
open(vid);
set(0,'defaultfigurecolor',[1 1 1]);
set(gca,'linewidth',7); 
% Equation describing roughness, which is greatest when occupancy of ML is
% 50% 
rough = 1.4*(1 - abs((noCells^2/2)-S)/(noCells^2/2));
refl = exp(-rough.^2);

for jj = 1:noFrames
    newplot
    hold off
    yyaxis left 
    % Plot reflectivity
    plot(refl(1:jj)/max(refl(:)),'r','LineWidth', 2.0)
    set(gca,'fontsize',18);
    set(gca,'TickLength',[0.02, 1]);
    set(gca,'linewidth',2);
    xlabel('Time [arb. units]');
    set(gca,'YTickLabel',[]);
    set(gca,'XTickLabel',[]);
    set(gca,'xtick',[])
    set(gca,'ytick',[])
    xlim([0 noFrames]);
    ylim([0 1.05]);
    set(gca, 'YColor','k')
    ylabel('Reflectivity [arb. units]','Color','r');

    yyaxis right 
    % Plot monolayer coverage 
    plot(S(1:jj)/max(S(:)),'b','LineWidth', 2.0)
    set(gca,'fontsize',18);
    set(gca,'TickLength',[0.02, 1]);
    set(gca,'linewidth',2);
    set(gca,'YTickLabel',[])
    set(gca,'XTickLabel',[])
    set(gca,'xtick',[])
    set(gca,'ytick',[])
    xlim([0 noFrames]);
    ylim([0 1.05]);
    set(gca, 'YColor','k')
    ylabel('Monolayer coverage','Color','b');
    frame = getframe(gcf);
    writeVideo(vid,frame);
end
close(vid)

