MODULE Ellipse;

(* Procedures and predicates for 
   ellipses with arbitrary 
   orientations. *)

IMPORT Geometry, PS, Rel;

(* ECCENTRICITY REPRESENTATION *)

(* The following routines define an 
   ellipse "[a, b, e]" by its center 
   "a", the tip of its major axis 
   "b", and its eccentricity "e" 
   (the real-valued ratio of the 
   length of its minor axis to the 
   length its major axis). This 
   representation is somewhat more 
   efficient than the one used later 
   in this module, but cannot be 
   clicked in. *)

PRIVATE PROC DrawQuarterEcc(a, r) IS 
  VAR b IN 
    b := PS.CurrentPoint(); 
    IF 
      VAR 
        c = (0, r) REL (a, b), 
        d 
      = (1, r * 0.5523) REL (a, b), 
        e = (0.5523, r) REL (a, b) 
      IN 
        PS.CurveTo(d, e, c) 
      END 
    FI 
  END 
END;

/* Add a quarter ellipse with 
   eccentricity "r" starting at the 
   current point about the center 
   "a" to the current path. Positive 
   values of "r" create the curve in 
   a counter-clockwise direction 
   from the current point around 
   "a", while negative values of "r" 
   create the curve in a clockwise 
   direction around "a". */

PRIVATE PROC DrawSemiEcc(a, e) IS 
  DrawQuarterEcc(a, e); 
  DrawQuarterEcc(a, 1 / e) 
END;

/* Add half an ellipse of 
   eccentricity "r" about "a" 
   starting from the current point. 
   Positive values of "r" create the 
   curve in a counter-clockwise 
   direction from the current point 
   around "a", while negative values 
   of "r" create the curve in a 
   clockwise direction around "a". */

PROC DrawEcc(a, b, e) IS 
  PS.MoveTo(b); 
  DrawSemiEcc(a, -e); 
  DrawSemiEcc(a, -e); 
  PS.Close() 
END;

PROC DrawEccCC(a, b, e) IS 
  PS.MoveTo(b); 
  DrawSemiEcc(a, e); 
  DrawSemiEcc(a, e); 
  PS.Close() 
END;

(* Add a clockwise and counter- 
   clockwise ellipse "[a, b, e]", 
   respectively, to the current 
   path. The current point becomes 
   "b". *)

PRED OnEcc(p, a, b, e) IS 
  (E q = Rel.Inv(p, a, b), 
     x = CAR(q), y = CDR(q) :: 
    x * x + (y * y) / (e * e) = 1) 
END;

(* The point "p" is on the ellipse 
   "[a, b, e]". *)

(* THREE-POINT REPRESENTATION *)

(* The following routines define an 
   ellipse "[a, b, c]" by its center 
   "a", a point "b" at the tip of 
   one of its axes, and a point "c" 
   whose distance from "a" is the 
   length of the axis perpendicular 
   to "ab". *)

PROC Draw(a, b, c) IS 
  VAR ab, ac IN 
    ab := Geometry.Dist(a, b); 
    ac := Geometry.Dist(a, c); 
    DrawEcc(a, b, ac / ab) 
  END 
END;

PROC DrawCC(a, b, c) IS 
  VAR ab, ac IN 
    ab := Geometry.Dist(a, b); 
    ac := Geometry.Dist(a, c); 
    DrawEccCC(a, b, ac / ab) 
  END 
END;

UI PointTool(Draw);

UI PointTool(DrawCC);

(* Add a clockwise and counter- 
   clockwise ellipse "[a, b, c]", 
   respectively, to the current 
   path. The current point becomes 
   "b". *)

PROC DrawSemi(a, c) IS 
  VAR b, ab, ac, e IN 
    b := PS.CurrentPoint(); 
    ab := Geometry.Dist(a, b); 
    ac := Geometry.Dist(a, c); 
    e := ac / ab; 
    DrawSemiEcc(a, -e) 
  END 
END;

PROC DrawSemiCC(a, c) IS 
  VAR b, ab, ac, e IN 
    b := PS.CurrentPoint(); 
    ab := Geometry.Dist(a, b); 
    ac := Geometry.Dist(a, c); 
    e := ac / ab; 
    DrawSemiEcc(a, e) 
  END 
END;

UI PointTool(DrawSemi);

UI PointTool(DrawSemiCC);

(* Add a clockwise or 
   counter-clockwise semi-ellipse to 
   the current path starting at the 
   current point. The semi-ellipse 
   has center "a", the length of its 
   major axis is the distance from 
   "a" to the current point, and the 
   length of its minor axis is the 
   distance from "a" to "c". *)

PRED On(p, a, b, c) IS 
  (E ab, ac :: 
    ab = Geometry.Dist(a, b) AND 
    ac = Geometry.Dist(a, c) AND 
    OnEcc(p, a, b, ac / ab)) 
END;

UI PointTool(On);

(* The point "p" is on the ellipse 
   "[a, b, c]". *)