% Program plotting parallel rays focussed by refraction as they pass
% through a plano-parabolic lens. The x- and y-dimensions are not equal in
% order to enhance the refraction angles

clear; close all;

aconst = 1.6; % constant in x = ay^2

vid = VideoWriter('parabolaFocussingReflection.mp4','MPEG-4');
vid.FrameRate = 30;
vid.Quality = 100;
open(vid);
figure('ToolBar','none');
set(0,'defaultfigurecolor',[1 1 1]);
ycoord = (-1:0.0099999999:1)'; % range of y-values of incoming horizontal beam
xcoord = aconst.*ycoord.^2; % x-value for a given y-value on surface of parabola

xmax = xcoord(size(ycoord,1)) + 1;


for i = 10:size(ycoord,1)-10
    hold off;
    grad = (4*aconst*xcoord(i))^-0.5;
    theta = atan(grad);
    ldist = abs(ycoord(i))/tan(2*theta);
    xf = xcoord(i) - ldist;
    plot(xcoord,ycoord,'color',[0,0,0], 'LineWidth', 2.0); % plot parabola
    axis([-0.1 xmax -xmax/2-0.05 xmax/2+0.05]);
    axis square % square plot
    axis equal
    axis off
    hold on;
    xchord = [xf xcoord(i) xmax]; % x-coordinates of ray path between foci
    ychord = [0 ycoord(i) ycoord(i)]; % y-coordinates of ray path between foci
    plot(xchord,ychord, 'color',[1 0 0], 'Linewidth', 2.0); % reflection off upper half
    plot(xchord,-ychord, 'color',[1 0 0], 'Linewidth', 2.0); % reflection off lower half
    frame = getframe(gcf);
    writeVideo(vid,frame);
end
hold off;
% Output the movie as an mpg file
close(vid);
