function theta = mars_angle( d )
%
%   Find the angle (in degrees) at which Mars appears in the sky 
%   on day d of the year 2000.

  Pe = 365.256;     %  period of Earth orbit (days)
  Pm = 686.980;     %  period of Mars orbit (days)
  d0 = 187;         %  day of 2000 on which Mars crosses angle 0
  Re = 1.496e8;     %  radius of Earth orbit (km)
  Rm = 2.279e8;     %  radius of Mars orbit (km)

  phi_e = 2*pi*d/Pe;
  phi_m = 2*pi*(d-d0)/Pm;

  xe = Re*cos( phi_e );
  ye = Re*sin( phi_e );
  xm = Rm*cos( phi_m );
  ym = Rm*sin( phi_m );

  
  
  % disp('xm - xe = ')
  % disp(xm - xe)
  % disp('ym - ye = ')
  % disp(ym - ye)
  
  theta = atan( ( ym - ye )/( xm - xe ) );
  % theta = atan2( ( ym - ye ), ( xm - xe ) );

 theta = 180*theta/pi;   %  convert from radians to degrees
