% Cartoon movie of electric-field vectors induced by a moving charge 
% This version moves the electron in the plane of the field grid 

clear; close all;
vid = VideoWriter('electricField2.mp4','MPEG-4');
vid.Quality = 100;
frmRate = 30;
vid.FrameRate = frmRate;
open(vid);
figure('units','pixels','position',[0 0 3840 2160],'ToolBar','none');
set(0,'defaultfigurecolor',[1 1 1]); % white background = [1 1 1]
set(gca,'linewidth',7);
lp = [2 1 1];
myBlue = [0.5 0.55 0.88];
myBlue2 = [0.61 0.67 0.92];
pos1 = [-0.1 -0.2 1.2 1.5];

[xx,yy,zz] = sphere(70); 
r = 0.5; % Radius of electron
ampl = 1.25;
lambda = 16; % Wavelength in terms of grid dimensions
% Want c to assume a value so that waves move one wavelength/second
c = lambda/frmRate; % c = f lambda; f = 1s = frmRate frames

xMx = 28;
yMx = 28;
x = linspace(-xMx,xMx,2*xMx+1);
y = linspace(-yMx,yMx,2*yMx+1);
[X,Y] = meshgrid(x,y);
Z = 0.0*X;

subplot('position',pos1);
for frm = 0:29 % Time in units of 1/frmRate
    newplot
    hold on
    plot3(X,Y,Z,'.y')
    axis equal
    axis off
    set(gca, 'Projection','perspective');
    set(gca,'View',[44,40]);
    set(gcf,'Color','k');
    campos([43 41 25]);
    xlim([-xMx-0.1 xMx+0.1]);
    ylim([-yMx-0.1 yMx+0.1]);
    zlim([-4*ampl 4*ampl]);

    accy = sin(2*pi*c*frm/lambda); % Acceleration of electron in y-direction
    posy = -sin(2*pi*c*frm/lambda); % Position of electron in y-direction
    R = (X.^2 + (Y-posy).^2).^0.5; % Matrix of distances of each element in the meshgrid from present position of electron 
    for ii = 1:2*xMx+1
        for jj = 1:2*yMx+1
            if (X(jj,ii) ~= 0) || (Y(jj,ii) ~= 0)
                % Value of position of electron as detected by field point 
                % after field propagates out from centre
                posyHere = -sin(2*pi*c*(frm - abs(R(jj,ii))/c)/lambda); 
                alpha = atan2((Y(jj,ii)-posyHere),X(jj,ii)); % Angle of grid elements relative to electron
                accperp = -ampl*posyHere*cos(alpha); % Perpendicular component of acceleration 
                Efield = accperp/R(jj,ii);
                arrowStartx = X(jj,ii);
                arrowStarty = Y(jj,ii);
                arrowStartz = Z(jj,ii);
                arrowEndx = arrowStartx - accperp*sin(alpha);
                arrowEndy = arrowStarty + accperp*cos(alpha);
                arrowEndz = 0; 
                % Determine whether arrow has nonzero length because
                % mArrow3 totally craps out if it is zero. Took me forever
                % to debug!! 
                delxa = arrowEndx - arrowStartx; 
                delya = arrowEndy - arrowStarty; 
                delza = arrowEndz - arrowStartz; 
                arrowLen = (delxa^2 + delya^2 + delza^2)^0.5; 
                if (arrowLen > 0)
                    h = mArrow3([arrowStartx arrowStarty arrowStartz],...
                        [arrowEndx arrowEndy arrowEndz],...
                        'color',myBlue2,'stemWidth',0.07,'tipWidth',0.1,...
                        'facealpha',(abs(Efield)/ampl)^0.35);
                end 
            end
        end
    end 
    % Oscillating electron 
    surf(xx*r,yy*r+posy,zz*r,'FaceColor',myBlue,'EdgeColor','none');
    light('Position',lp,'Style','infinite');

    frame = getframe(gcf);
    writeVideo(vid,frame);
    hold off
end
% Output the movie as an mpg file
close(vid);

