MODULE Rect;

(* Procedures and predicates for rectangles. *)

IMPORT R2, Geometry, PS, Ellipse;

(* \section{Opposite Corner Representation} *)

(* The follow procedures represent a rectangle 
   oriented squarely with the page by two of its 
   opposite corners, "a" and "b". *)

PROC Draw(a, b) IS 
  IF 
    VAR ax, ay, bx, by IN 
      a = (ax, ay) AND b = (bx, by) -> 
        PS.MoveTo(a); 
        PS.LineTo((ax, by)); 
        PS.LineTo(b); 
        PS.LineTo((bx, ay)); 
        PS.Close() 
    END 
  FI 
END;

UI PointTool(Draw);

(* Append a closed rectangular path with opposite 
   corners "a" and "b" oriented squarely with the 
   page to the current path. *)

PRED On(p, a, b) IS 
  (E px, py, ax, ay, bx, by :: 
    p = (px, py) AND 
    a = (ax, ay) AND 
    b = (bx, by) AND 
    0 = 
      (px - ax) * (px - bx) * (py - ay) * 
      (py - by)) 
END;

UI PointTool(On);

(* The point "p" is on the rectangle with 
   opposite corners "a" and "b". 

   This constraint is actually true if "p" is on 
   any of the lines obtained by extending the 
   sides of the rectangle infinitely. However, if 
   the hint for "p" is inside the circle 
   circumscribing the rectangle and is not too 
   close to the center of the rectangle, the 
   solver will converge on a solution for "p" 
   that is actually on the rectangle. You can use 
   the "OnHint" function to compute a hint for 
   "p" that will converge on the actual rectangle 
   in most cases. *)

PRIVATE FUNC ecc = Ecc(wh) IS 
  (E w, h :: wh = (w, h) AND ecc = h / w) 
END;

PRIVATE FUNC p = OnCHint0(c, b, d) IS 
  (E wh = R2.Minus(b, c), maj :: 
    maj = R2.Plus(c, (CAR(wh), 0)) AND 
    p ~ (0.5, 0) REL (d, c) AND 
    Geometry.Colinear(d, p, c) AND 
    Ellipse.OnEcc(p, c, maj, Ecc(wh))) 
END;

/* See the comment for "OnCHint" below. This 
   procedure is here because it must preceed 
   "OnHint". */

FUNC p = OnHint(a, b, d) IS 
  (E c = Geometry.Mid(a, b) :: 
    p = OnCHint0(c, b, d)) 
END;

(* "p" is a good hint for a point on the 
   intersection of the rectangle with corners "a" 
   and "b" and the line through "d" and the 
   center of that rectangle. *)

FUNC res = Area(a, b) IS
  res = Geometry.DistX(a,b) * Geometry.DistY(a,b)
END;

(* "res" is the area of the rectangle with 
   opposite corners "a" and "b". *)

(* \section{Center-Corner Representation} *)

PRIVATE FUNC a = OtherCorner(c, b) IS
  a = (-1, 0) REL (c, b)
END;

/* "a" is the opposite corner of the rectangle
   from the point "b", where "c" is the center
   of the rectangle. */

(* The following procedures are like the ones 
   above, only they represent a rectangle by its 
   center "c" and one of its corners "b". *)

PROC DrawC(c, b) IS 
  Draw(OtherCorner(c, b), b)
END;

UI PointTool(DrawC);

(* Like "Draw" above, but the rectangle is 
   specified by its center "c" and one of its 
   corners "b". *)

PRED OnC(p, c, b) IS 
  On(p, OtherCorner(c, b), b)
END;

UI PointTool(OnC);

(* Like "On" above, but the rectangle is 
   specified by its center "c" and one of its 
   corners "b". The current point becomes "b". *)

FUNC p = OnCHint(c, b, d) IS 
  p = OnCHint0(c, b, d) 
END;

(* "p" is a good hint for a point on the 
   intersection of the line "cd" and the 
   rectangle with center "c" and corner "b". *)

FUNC res = AreaC(c, b) IS
  res = Area(OtherCorner(c, b), b)
END;

(* "res" is the area of the rectangle with
   center "c" and corner "b". *)

(* \section{Pair representation} *)

(* The following routines encode a rectangle "r" 
   as a pair of points denoting two of its 
   opposite corners. 

   If "r = ((x1, y1), (x2, y2)" and we have: 

|    west = MIN(x1, x2)
|    east = MAX(x1, x2)
|    south = MIN(y1, y2)
|    north = MAX(y1, y2)

   then a point "p = (px, py)" is in the 
   rectangle "r" if: 

|    west <= px < east
|    south <= py < north

   We impose the restriction that if a rectangle 
   contains no points, then it must be equal to 
   "Rect.Empty". *)

CONST Empty = ((0, 0), (0, 0));

(* The empty rectangle. *)

PROC res := FromCorners(a, b) IS
  res := (a, b)
END;

(* "res" is the rectangle with opposite corners
   "a" and "b". *)

PROC res := FromCenterCorner(c,b) IS
  res := ((-1,0) REL (c,b), b)
END;

(* "res" is the rectangle with center "c" and
   corner "b". *)

PROC DrawR(r) IS 
  IF 
    VAR r1, r2 IN r = (r1, r2) -> Draw(r1, r2) END 
  FI 
END;

(* Like "Draw" above, but the rectangle "r" is 
   specified by a pair of points at opposite 
   corners of the rectangle. *)

PRED OnR(p, r) IS 
  (E r1, r2 :: r = (r1, r2) AND On(p, r1, r2)) 
END;

(* Like "On" above, but the rectangle "r" is 
   specified by a pair of points at opposite 
   corners of the rectangle. *)

FUNC p = OnRHint(r, d) IS 
  p = OnHint(CAR(r), CDR(r), d) 
END;

(* Like "OnHint" above, but the rectangle "r" is 
   specified by a pair of points at opposite 
   corners of the rectangle. *)

PROC p := North(r) IS 
  IF 
    VAR x1, y1, x2, y2 IN 
      r = ((x1, y1), (x2, y2)) -> 
        p := ((x1 + x2) / 2, MAX(y1, y2)) 
    END 
  FI 
END;

PROC p := South(r) IS 
  IF 
    VAR x1, y1, x2, y2 IN 
      r = ((x1, y1), (x2, y2)) -> 
        p := ((x1 + x2) / 2, MIN(y1, y2)) 
    END 
  FI 
END;

PROC p := East(r) IS 
  IF 
    VAR x1, y1, x2, y2 IN 
      r = ((x1, y1), (x2, y2)) -> 
        p := (MAX(x1, x2), (y1 + y2) / 2) 
    END 
  FI 
END;

PROC p := West(r) IS 
  IF 
    VAR x1, y1, x2, y2 IN 
      r = ((x1, y1), (x2, y2)) -> 
        p := (MIN(x1, x2), (y1 + y2) / 2) 
    END 
  FI 
END;

FUNC p = Center(r) IS 
  (E r1, r2 :: 
    r = (r1, r2) AND p = Geometry.Mid(r1, r2)) 
END;

PROC p := NW(r) IS 
  IF 
    VAR x1, y1, x2, y2 IN 
      r = ((x1, y1), (x2, y2)) -> 
        p := (MIN(x1, x2), MAX(y1, y2)) 
    END 
  FI 
END;

PROC p := NE(r) IS 
  IF 
    VAR x1, y1, x2, y2 IN 
      r = ((x1, y1), (x2, y2)) -> 
        p := (MAX(x1, x2), MAX(y1, y2)) 
    END 
  FI 
END;

PROC p := SE(r) IS 
  IF 
    VAR x1, y1, x2, y2 IN 
      r = ((x1, y1), (x2, y2)) -> 
        p := (MAX(x1, x2), MIN(y1, y2)) 
    END 
  FI 
END;

PROC p := SW(r) IS 
  IF 
    VAR x1, y1, x2, y2 IN 
      r = ((x1, y1), (x2, y2)) -> 
        p := (MIN(x1, x2), MIN(y1, y2)) 
    END 
  FI 
END;

(* These procedures set "p" to the named 
   characteristic point on the boundry or at the 
   center of the rectangle "r". *)

PROC w := Width(r) IS 
  IF 
    VAR x1, y1, x2, y2 IN 
      r = ((x1, y1), (x2, y2)) -> 
        w := ABS(x2 - x1) 
    END 
  FI 
END;

(* "w" is the width of "r". *)

PROC h := Height(r) IS 
  IF 
    VAR x1, y1, x2, y2 IN 
      r = ((x1, y1), (x2, y2)) -> 
        h := ABS(y2 - y1) 
    END 
  FI 
END;

(* "h" is the height of "r". *)

PROC res := Size(r) IS 
  IF 
    VAR x1, y1, x2, y2 IN 
      r = ((x1, y1), (x2, y2)) -> 
        res := (ABS(x2 - x1), ABS(y2 - y1)) 
    END 
  FI 
END;

(* "res" is a pair "(w, h)", where "w" and "h" 
   are the width and height of "r", respectively. 
*)

PROC res := Plus(r, p) IS 
  IF 
    r = Empty -> res := Empty 
  | VAR r1, r2 IN 
      r = (r1, r2) -> 
        res := (R2.Plus(r1, p), R2.Plus(r2, p)) 
    END 
  FI 
END;

(* "res" is the rectangle obtained by translating 
   "r" by the vector "p". *)

PROC w, e, s, n := Edges(r) IS 
  IF 
    VAR r1x, r1y, r2x, r2y IN 
      r = ((r1x, r1y), (r2x, r2y)) -> 
        w := MIN(r1x, r2x); 
        e := MAX(r1x, r2x); 
        s := MIN(r1y, r2y); 
        n := MAX(r1y, r2y) 
    END 
  FI 
END;

(* Set "w" and "e" to the x-coordinates of the
   west and east edges of "r", and "s" and "n"
   to the y-coordinates of the south and north
   edges of "r". *)

PROC res := Inset(r, d) IS 
  VAR w, e, s, n IN 
    w, e, s, n := Edges(r); 
    w, s := w + d, s + d; 
    e, n := e - d, n - d; 
    IF 
      e < w OR n < s -> res := Empty 
    | res := ((w, s), (e, n)) 
    FI 
  END 
END;

(* "res" is the rectangle "r" inset by "d" 
   points. If "r" does not have width and height 
   at least "2*d", "res" is set to "Empty". *)

PROC res := Join(r, s) IS 
  IF 
    r = Empty -> res := s | s = Empty -> res := r 
  | VAR rW, rE, rS, rN, sW, sE, sS, sN IN 
      rW, rE, rS, rN := Edges(r); 
      sW, sE, sS, sN := Edges(s); 
      res := 
        ((MIN(rW, sW), MIN(rS, sS)), 
         (MAX(rE, sE), MAX(rN, sN))) 
    END 
  FI 
END;

(* "res" is the smallest rectangle containing 
   both "r" and "s". *)

PROC res := Meet(r, s) IS 
  IF 
    r = Empty OR s = Empty -> res := Empty 
  | VAR 
      rW, rE, rS, rN, sW, sE, sS, sN, resW, resE, 
      resS, resN 
    IN 
      rW, rE, rS, rN := Edges(r); 
      sW, sE, sS, sN := Edges(s); 
      resW := MAX(rW, sW); 
      resE := MIN(rE, sE); 
      resS := MAX(rS, sS); 
      resN := MIN(rN, sN); 
      IF 
        resW >= resE OR resS >= resN -> 
          res := Empty 
      | res := ((resW, resS), (resE, resN)) 
      FI 
    END 
  FI 
END;

(* "res" is the largest rectangle contained in 
   both "r" and "s". *)

FUNC res = AreaR(r) IS
  (E c1, c2 :: r = (c1, c2) AND
    res = Area(c1, c2))
END;

(* "res" is the area of the rectangle "r". *)