% This Matlab m file converts FDTD Solutions frequency 
% monitor field data into a CW movie.  
%
% Copyright 2013, Lumerical Solutions, Inc.


% load movie data from the .mat file
close all;
clear;
load full_E_field.mat;


% create phase vector
phi = linspace(0,4*2*pi,161); 	    % 4 optical cycles, 40 frames per cycle


% Rotate data and interpolate onto linearly spaced position
% vectors as requried by imagesc function
x2 = linspace(x(1),x(length(x)),length(x)*4);         % create linearly spaced position vectors
y2 = linspace(y(1),y(length(y)),length(y)*4);
[X,Y] = meshgrid(x,y);
[X2,Y2] = meshgrid(x2,y2);

Ex2 = flipud(Ex.');              % convert matrix order from FDTD form to matlab form
Ey2 = flipud(Ey.');              % These matrices will be used for the image plots of the fields
Ez2 = flipud(Ez.');
n2  = flipud(n.');
Ex2 = interp2(X,Y,Ex2,x2,y2');   % interpolate data onto linearly spaced position vectors, as required by imagesc
Ey2 = interp2(X,Y,Ey2,x2,y2');
Ez2 = interp2(X,Y,Ez2,x2,y2');
n2  = interp2(X,Y,n2,x2,y2');     




%%%%% create movie of Ey(t)
figure;
M = moviein(length(phi));	

clim = max(max(abs(Ey2)));
clim=[-clim,clim];                % set colorbar limits to max value of field

% calculate movie frames.  
for j=1:length(phi)
    hold off					% clear the frame
    imagesc(x2,y2,real(exp(-1i*phi(j))*Ey2),clim)	% plot field, multiplied by e^(i*phi)
    hold on 
    contour(x2,y2,real(n2),1,'k')			% plot index coutour
    h = findobj('Type','line');
    set(h,'LineWidth',2)
    axis off
    M(:,j) = getframe;				% record frame
end

movie(M);			% convert to avi movie format
v = VideoWriter('Ey_movie.avi','Uncompressed AVI');
open(v);
writeVideo(v,M);
close(v)



%%%%% create movie of |E(t)|^2
figure;
M = moviein(length(phi));	

clim = max(max( abs(Ex).^2+abs(Ey).^2+abs(Ez).^2 ));
clim = clim/2;                    % oversaturate by 2x
clim = [0,clim];                    % set colorbar limits to max value of field

% calculate movie frames.  
for j=1:length(phi)
    hold off					% clear the frame
    temp = (real(exp(-1i*phi(j))*Ex2)).^2 + (real(exp(-1i*phi(j))*Ey2)).^2 +(real(exp(-1i*phi(j))*Ez2)).^2;
    imagesc(x2,y2,temp,clim)	% plot field, multiplied by e^(i*phi)
    hold on 
    contour(x2,y2,real(n2),1,'k')			% plot index coutour
    h = findobj('Type','line');
    set(h,'LineWidth',2)
    axis off
    M(:,j) = getframe;				% record frame
end

movie(M);			% convert to avi movie format
v = VideoWriter('E2_movie.avi','Uncompressed AVI');
open(v);
writeVideo(v,M);
close(v)




%%%%% create movie of E(t) using vector plot
load vector_movie.mat;

[X,Y] = meshgrid(x,y);
X=X';
Y=Y';
Z=X*0;

figure;
%figure('units','normalized','outerposition',[0,0,1,1]);  % Make figure full screen
M = moviein(length(phi));	
v_scale=0.3e-6;

% calculate movie frames.  
for j=1:length(phi)
    hold off					% clear the frame 
    quiver3(X,Y,Z,real(exp(-1i*phi(j))*Ex)*v_scale,-real(exp(-1i*phi(j))*Ey)*v_scale,real(exp(-1i*phi(j))*Ez)*v_scale,0)
    axis([min(x),max(x),min(y),max(y),-1e-6,1e-6])
    zoom(4)
    M(:,j) = getframe;				% record frame
end

movie(M);			% convert to avi movie format
v = VideoWriter('Ey_movie_vector.avi','Uncompressed AVI');
open(v);
writeVideo(v,M);
close(v)