INTERFACE SMPoint;
(* 
   This module implements the sets of points on the sphere "S2". 
   A point is represented by its homogeneous coordinates in the
   projective 3-space. The coordinate "w", in the case of points in
   the sets "A" and "B" is left implicit since for every point
   in the sphere we have "w^2=x^2+y^2+z^2". Then,
   
      "S2" is the set of all cartesian points "(X,Y,Z)" on the
           sphere. In this case "X,Y,Z" are approximated by
           longreals. Note that "x^2+y^2+z^2=1".

      "A" is the set of points "(X,Y,Z)" on the sphere such that
          "X","Y", are "Z" are rational, i.e. the point "[w,x,y,z]"
          such that "w","x","y", and "z" are integer and
          "w^2=x^2+y^2+z^2".  In others words, "A" is the set of
          pytagorean quadruplets on the sphere.
          
      "B" is the set of points "[w,x,y,z]" on the sphere such that
          "x,y,z" are integers. In this case, "w=sqrt(x^2+y^2+z^2)"
          may not be integer. This set consists of every point on
          "S2" whose direction is rational.

      "C" consists of every point of the sphere such that there is a
          rational line "l" passing through it. These points are
          represented by the six Plucker coefficients of a line "l"
          such that "exit(l)=p". As we can prove, for every point "p"
          in "C \ A", there is only one rational line such that
          "exit(l)=p". For any point "p \in C", we assume a canonical
          representation "\stab(p)" such that, if "p \in C \ A" then
          "\stab(p)" is the unique rational line "l" such that
          "\exit(l) = p". On the other hand, if "p \in A" then
          "\stab(p) = O \join p".

      "Q2" is the set of points (in homogeneous coordinates) "[w,x,y]"
           of the projective plane "T^2" corresponding to the
           stereographic projection of the points of the set "A". Note
           that, for every point in this set the coordinates "w,x,y" 
           are integers.
          
      "T2" is the set of points "[w,x,y]" of the projective plane
           "T^2". The coordinates "w,x,y" are approximated by
           longreals.
          
   Created in Jul 19, 95 by Marcus Vinicius A. Andrade.
   Last edited in Jul., 28, 2000 Marcus Vinicius A. Andrade.            *)

IMPORT I3, I6, HI3, HLR3, LR3, Rd, Wr, Thread, Lex, FloatMode;

TYPE 
   AT =  RECORD a : I3.T END;
   BT =  RECORD b : I3.T END;
   CT =  RECORD c : I6.T END;
   S2 =  RECORD s : LR3.T END;
   Q2 =  RECORD q : I3.T END;
   T2 =  RECORD r : LR3.T END;

VAR 
  APNull : AT;        (* the tuple "[0,0,0]" used in invalid function return *)
  CPNull : CT;        (* the tuple "<<0,0,..0>>" used in invalid function return *)


PROCEDURE Stab(READONLY p: CT) : CT;
  (* 
    Return the canonical representation of the point "p" defined above *)

PROCEDURE IsApointNull(READONLY p: AT) : BOOLEAN;
  (*
    Returns TRUE if "p" is the tuple "[0,0,0]"   *)

PROCEDURE IsBpointNull(READONLY p: BT) : BOOLEAN;
  (*
    Returns TRUE if "p" is the tuple "[0,0,0]"   *)

PROCEDURE IsCpointNull(READONLY p: CT) : BOOLEAN;
  (*
    Returns TRUE if "p" is the tuple "<<0,0,0,0,0,0>>"   *)

PROCEDURE IsApoint(READONLY l: CT; VAR ap: AT) : BOOLEAN;
  (* 
    Returns TRUE if the CPoint represented by "l" corresponds to a 
    point of "A", that is, if "exit(l)" is an APoint. Besides, sets that 
    point to "ap". Otherwise, if "exit(l)" is not a point of "A", then 
    returns FALSE and sets "[0,0,0,0]" to "ap".                       *)

PROCEDURE IsCpoint(READONLY l: HI3.Line) : BOOLEAN;
  (*
    Checks if the line "l" defines a Cpoint, that is, if the line 
    with the coefficients "t" intersects S2.                         *)

PROCEDURE AToB(READONLY p: AT) : BT; 
  (*
    Returns the point in "B" corresponding to the point "p" in "A".   *)   

PROCEDURE AToC(READONLY p: AT) : CT; 
  (*
    Returns the point in "C" corresponding to the point "p" in "A".   *)   
 
PROCEDURE AToS(READONLY p: AT) : S2; 
  (*
    Returns the point in "S" corresponding to the point "p" in "A".   *)   

PROCEDURE BToC(READONLY p: BT) : CT; 
  (*
    Returns the point in "C" corresponding to the point "p" in "B".   *)   
 
PROCEDURE BToS(READONLY p: BT) : S2; 
  (*
    Returns the point in "S" correspondig to the point "p" in "B".    *)   
 
PROCEDURE CToS(READONLY p: CT) : S2; 
  (*
    Returns the point in "S" correspondig to the point "p" in "C".    *)   

PROCEDURE SToA(READONLY p: S2; eps: LONGREAL:=1.0d-6): AT;
  (*
    Returns a point "q" in "A" that is an approximation to the point 
    "p" in "S", that is, "|q - p| < eps". *)
    
PROCEDURE SToB(READONLY p: S2; eps: LONGREAL:=1.0d-6) : BT;
  (*
    Returns a point "q" in "B" that is an approximation to the point 
    "p" in ST, that is, "|q - p| < eps". *)

PROCEDURE SToC(READONLY p: S2; READONLY eps: LONGREAL:=1.0d-6) : CT;
  (*
    Returns a point in "C" that is an approximation to the point "p" 
    on the sphere. *)

PROCEDURE AStereo(READONLY p: AT): Q2;
  (*
    Returns the point in "Q" that is the stereographic projection of 
    the point "p" in "A". Since "w=sqrt(x^2+y^2+z^2)" is integer then 
    the result of this projection is a point in "Q". *) 

PROCEDURE BStereo(READONLY p: BT) : T2;
  (*
    Returns the point in "R" that is the stereographic projection of 
    the point "p" in "B". *) 
    
PROCEDURE CStereo(READONLY p: CT) : T2;
  (*
    Returns the point in "R" that is the stereographic projection of 
    the point "p" in "C". *) 

PROCEDURE SStereo(READONLY p: S2) : T2;
  (*
    Returns the point in "R" that is the stereographic projection of 
    the point "p" in "S". *) 

PROCEDURE InvAStereo(READONLY p: Q2) : AT;
  (*
    Returns a point in "A" that is an approximation to the inverse of 
    the stereographic projection of the point "p" in the 
    plane "R^2" with integers coordinates. *)

PROCEDURE InvBStereo(READONLY p: T2; READONLY eps: LONGREAL:=1.0d-6) : BT;
  (*
    Returns a point in "B" that is an approximation to the inverse of 
    the stereographic projection of the point "p" in the 
    plane "R^2". *)

PROCEDURE InvCStereo(READONLY q: T2; READONLY eps: LONGREAL:=1.0d-6) : CT;
  (*
    Returns a point in "C" that is an approximation to the inverse of 
    the stereographic projection of the point "p" in the 
    plane "R^2". *)

PROCEDURE InvSStereo(READONLY q: T2) : S2;
  (*
    Returns a point in "S" that is an approximation to the inverse of 
    the stereographic projection of the point "p" in the 
    plane "R^2". *)

PROCEDURE HomogenizeSPoint(READONLY p: S2) : HLR3.Point;
  (*
    Returns the homogeneous coordinates of the point corresponding to 
    p (a point on the sphere), i.e., returns "[1,p0,p1,p2]". *)

PROCEDURE ApointToHomogPoint(READONLY p: AT): HI3.Point;
  (* 
    Returns the four homogenous coordinates of the point "p" in "A". *)

PROCEDURE BpointToHomogPoint(READONLY p: BT; 
                             READONLY eps: LONGREAL:=1.0d-6) : HI3.Point;
  (* 
    Returns the four integer homogenous coordinates of the point "p" in "B"; 
    the w coordinate is approximated with precision "eps". *)

PROCEDURE CpointToHomogPoint(READONLY p: CT; 
                             READONLY eps: LONGREAL:=1.0d-6) : HI3.Point;
  (* 
    Returns the four integer homogenous coordinates of the point "p" in "C";
    each one of these coordinates is approximated with precision "eps". *)

PROCEDURE ApointFromHomogPoint(READONLY p: HI3.Point; 
                               READONLY eps: LONGREAL:=1.0d-6) : AT;
  (* 
    Returns the integer coordinates of the point in "A" corresponding to the 
    homogeneous point on the sphere "p".  These integer coordinates are 
    approximated with precision "eps". *)

PROCEDURE BpointFromHomogPoint(READONLY p: HI3.Point; 
                               READONLY eps: LONGREAL:=1.0d-6) : BT;
  (* 
    Returns the integer coordinates of the point in "B" corresponding to the
    homogeneous point "p". These integer coordinates are approximated with 
    precision "eps". *)

PROCEDURE CpointFromHomogPoint(READONLY p: HI3.Point; 
                               READONLY eps: LONGREAL:=1.0d-6) : CT;
  (* 
    Returns the integer coordinates of the point in "C" corresponding to the 
    homogeneous point "p".  These integer coordinates are approximated with 
    precision "eps". *)

PROCEDURE LineFromCpoint(READONLY p: CT) : HI3.Line;
  (* 
    Returns the line that defines the point "p" in "C". *)

PROCEDURE CpointFromLine(READONLY r: HI3.Line) : CT;
  (* 
    Given the line "r", returns the point of "C" defined by this line. *)

PROCEDURE ReadApoint(READONLY rd: Rd.T) : AT RAISES{Rd.Failure, Thread.Alerted, Lex.Error};
PROCEDURE ReadBpoint(READONLY rd: Rd.T) : BT RAISES{Rd.Failure, Thread.Alerted, Lex.Error};
PROCEDURE ReadBpoint2(READONLY rd: Rd.T) : BT RAISES{Rd.Failure, Thread.Alerted, Lex.Error};
PROCEDURE ReadCpoint(READONLY rd: Rd.T) : CT RAISES{Rd.Failure, Thread.Alerted, Lex.Error};
PROCEDURE ReadSpoint(READONLY rd: Rd.T) : S2 RAISES{Rd.Failure, Thread.Alerted, FloatMode.Trap, Lex.Error}; 
  (* 
    Reads an point of "A", "B", "C" and "S" on the given reader "rd". *)

PROCEDURE PrintApoint (wr: Wr.T; READONLY p: AT) RAISES{Thread.Alerted, Wr.Failure};
PROCEDURE PrintBpoint (wr: Wr.T; READONLY p: BT) RAISES{Thread.Alerted, Wr.Failure};
PROCEDURE PrintCpoint (wr: Wr.T; READONLY p: CT) RAISES{Thread.Alerted, Wr.Failure};
PROCEDURE PrintSpoint (wr: Wr.T; READONLY p: S2; ndig:CARDINAL:=15) RAISES{Thread.Alerted, Wr.Failure};
  (* 
    Prints a point of "A", "B", "C" or "S" on the writer "wr". *)

END SMPoint.