function P3()
% Tests Intersection(E1,E2)
% The examples...
E0 = struct('h',2,'k',3,'a',6,'b',4,'tau',-pi/3);
E(1) = struct('h',-1,'k',-3,'a',3,'b',1,'tau',-pi/10);
E(2) = struct('h',-1,'k',-3,'a',4,'b',3,'tau',-pi/5);
E(3) = struct('h',-1,'k',-3,'a',7,'b',5,'tau',-pi/5);
E(4) = struct('h',2,'k',10,'a',.9,'b',.8,'tau',0.1);
E(5) = struct('h',0,'k',9,'a',.9,'b',.8,'tau',0.1);
E(6) = struct('h',-2,'k',0,'a',2.9,'b',.2,'tau',-pi/4);
E(7) = struct('h',-2,'k',0,'a',4,'b',1.8,'tau',-pi/4);
E(8) = struct('h',-.5,'k',-1.5,'a',50,'b',4,'tau',-pi/3);

clc
for k=1:8
   close all
   figure
   hold on
   PlotEllipse(E0);
   PlotEllipse(E(k));
   axis([-10 8 -10 15])
   axis equal
   [x,y] = Intersection(E0,E(k));
   fprintf('\n\nExample %1d\n',k)
   if length(x)==0
       fprintf('No Intersection.')
       title('no intersection')
   elseif length(x)==1
       fprintf('x(1) = %15.12f  y(1) = %15.12f',x,y) 
       plot(x,y,'or')
   else
       fprintf('x(1) = %15.12f  y(1) = %15.12f\n',x(1),y(1))
       fprintf('x(2) = %15.12f  y(2) = %15.12f',x(2),y(2))
       plot(x,y,'or')
   end
   shg
   pause(2)
end

function PlotEllipse(E)
% Adds Ellipse E and its axes to the current figure window.
theta = linspace(0,2*pi);
h = E.h; k = E.k; a = E.a; b = E.b; tau = E.tau;
x0 = a*cos(theta); y0 = b*sin(theta);
cTau = cos(tau); sTau = sin(tau);
plot(h + cTau*x0 - sTau*y0 , k + sTau*x0 + cTau*y0,'k')
% Now display the foci and the major axis...
f = sqrt(a^2 - b^2);
plot([h+cTau*a h-cTau*a],[k+sTau*a k-sTau*a],'k')
plot([h+cTau*f h-cTau*f],[k+sTau*f k-sTau*f],'*k')





