  function ShowEccentric(xc,phi)
% Illustrates the Eccentric Model

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Using the subfunction Atan, here is the essential trigonometry 
% that locates the season start points...

c = cos(phi); s = sin(phi); 
tau = -c*xc + sqrt(1-(s*xc)^2);
Spring = [xc+tau*c  tau*s];
thetaSpring = MyAtan(Spring(2),Spring(1));

c = cos(phi+pi/2); s = sin(phi+pi/2); 
tau = -c*xc + sqrt(1-(s*xc)^2);
Summer = [xc+tau*c  tau*s];
thetaSummer = MyAtan(Summer(2),Summer(1));

c = cos(phi+pi); s = sin(phi+pi); 
tau = -c*xc + sqrt(1-(s*xc)^2);
Autumn = [xc+tau*c  tau*s];
thetaAutumn = MyAtan(Autumn(2),Autumn(1));

c = cos(phi+3*pi/2); s = sin(phi+3*pi/2); 
tau = -c*xc + sqrt(1-(s*xc)^2);
Winter = [xc+tau*c  tau*s];
thetaWinter = MyAtan(Winter(2),Winter(1));


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% The graphics...
close all
figure
set(gcf,'Position',[100 100 1200 800])
% Display the circle..
c = cos(linspace(0,2*pi,200));
s = sin(linspace(0,2*pi,200));
plot(c,s,'r',0,1.2,'.k',0,-1.2,'.k')
hold on
axis equal off
set(gcf,'color','k')
% Report the season angles...
text(-2.,1.,sprintf('thetaSpring    = %6.2f',180*thetaSpring/pi),'color','w','fontsize',18)
plot(-2.1,1,'.g','markersize',30)
plot(Spring(1),Spring(2),'.g','markersize',30)
text(-2.,.8,sprintf('thetaSummer = %6.2f',180*thetaSummer/pi),'color','w','fontsize',18)
plot(Summer(1),Summer(2),'.y','markersize',30)
plot(-2.1,.8,'.y','markersize',30)
text(-2,.6,sprintf('thetaAutumn  = %6.2f',180*thetaAutumn/pi),'color','w','fontsize',18)
plot(Autumn(1),Autumn(2),'.m','markersize',30)
plot(-2.1,.6,'.m','markersize',30)
text(-2,.4,sprintf('thetaWinter    = %6.2f',180*thetaWinter/pi),'color','w','fontsize',18)
plot(Winter(1),Winter(2),'.c','markersize',30)
plot(-2.1,.4,'.c','markersize',30)
% Report the values of xc and phi (in degrees)...
text(-2,-.2,sprintf('xc  = %6.3f',xc),'color','w','fontsize',18)
text(-2,-.4,sprintf('phi = %6.3f',phi*180/pi),'color','w','fontsize',18)
% Cross Hairs...
plot([Spring(1) Autumn(1)],[Spring(2) Autumn(2)],':w')
plot([Summer(1) Winter(1)],[Summer(2) Winter(2)],':w')
% The origin...
plot(0,0,'.w','markersize',30)
shg
% The simulation...
theta = 0;           % the angle of the ray
delta = pi/180;      % amount the ray is advanced each time step
nSteps = 800;        % the number of steps
for k=1:nSteps
    c = cos(k*delta);
    s = sin(k*delta);
    h1 = plot(c,s,'.w','markersize',30);
    h2 = plot([0 c],[0 s],'w','linewidth',2);
    pause(.03)
    shg
    delete(h1,h2) 
end

  function theta = MyAtan(y,x)
% Determines theta (0<=theta<2pi) so that (x,y)= (rcos(theta),rsin(theta))
% where r = sqrt(x^2+y^2).
theta = atan2(y,x);
if theta<0
    theta = 2*pi+theta;
end

    
