INTERFACE Triang;

(* 
  Mostly-triangulated meshes for automatic topology visualization
  
  Created 1993 by Rober M. Rosi
*)

IMPORT Oct, Random;
IMPORT LR3;

(* === ELEMENTS === *)
  
TYPE 
  Arc = Oct.Arc;  
  (* A "Triang.Arc" is an "Oct.Arc" whose "edge" field is a "Triang.Edge". *)

  Edge <: PublicEdge;
  PublicEdge = Oct.Edge OBJECT
      exists: BOOLEAN := TRUE;       (* FALSE for immaterial edges. *)
      spring: BOOLEAN := TRUE;       (* TRUE if edge must be smoothed *)
      style: EdgeStyle := 0;          (* Edge rendering style *)
    METHODS
      init(): Edge;
        (* 
          Initializes "self" as an isolated nonloop edge. *)
    END;
  EdgeStyle = CARDINAL;
  
TYPE
  Node <: PublicNode;
  PublicNode = OBJECT 
      num: CARDINAL;    (* Vertex or face number *)
    END; 
  
  Vertex <: PublicVertex;
  PublicVertex = Node OBJECT 
      exists: BOOLEAN := TRUE;        (* FALSE for ghost vertices. *)
      fixed: BOOLEAN := FALSE;        (* TRUE if position is fixed *)
      style: VertexStyle := 0;        (* Vertex renderning style *)
    END;
  VertexStyle = CARDINAL;
  
  Face <: PublicFace;  
  PublicFace = Node OBJECT
      exists: BOOLEAN := TRUE;        (* FALSE for ghost faces *)
      patch: PatchNum := 0;           (* Patch number for area balancing. *)
      style: FaceStyle := 0;          (* Face rendering style *)
    END;
  PatchNum = CARDINAL;
  FaceStyle = CARDINAL;

(* === ELEMENT CREATION === *)

PROCEDURE MakeEdge(): Arc;
  (* Creates a new unattached edge with distinct endpoints. *)
  
PROCEDURE MakeVertex(): Vertex;
  (* Creates an unattached vertex record *)

PROCEDURE MakeFace(): Face;
  (* Creates an unattached face record *)

(* === ARC PROPERTIES === *)

PROCEDURE Org(a: Arc): Node;
  (* The origin of arc "a". *)

PROCEDURE SetOrg(a: Arc; n: Node);
  (* Sets the origin of arc "a" to be "n". *)

PROCEDURE SetAllOrgs(a: Arc; n: Node);
  (* Does "SetOrg(t,n)" for all arcs "t" with same origin as "a". *)

PROCEDURE Left(a: Arc): Node;
  (* The left node of "a". *)

PROCEDURE SetLeft(a: Arc; n: Node);
  (* Sets the left face node of "a" to "n".  Same as "SetOrg(Tor(a),n)". *)

PROCEDURE SetAllLefts(a: Arc; n: Node);
  (* Does SetLeft(t,n) for all arcs "t" with same left face as "a". *)

PROCEDURE OrgV(a: Arc): Vertex;
  (* The origin of arc "a", narrowed to "Vertex". *)

PROCEDURE LeftF(a: Arc): Face;
  (* The left node of "a", narrowed to "Face". *)

(* === CONSTRUCTION TOOLS === *)

PROCEDURE Connect(a, b: Arc): Arc;
  (*
    Adds a new edge from Org(a) to Dest(b). The new edge
    will be spliced in as Onext(a) and Dnext(b).
    Returns the new edge. *)

PROCEDURE AddSpokes(a: Arc): Arc;
  (*
    Adds a new vertex in the middle of Left(a), connected
    by new edges ("spokes") to all vertices of Left(a).  Returns the 
    spoke Onext(a). *)

PROCEDURE MakeGrid(nx, ny: CARDINAL): ARRAY [0..7] OF Arc;
  (*
    Builds a grid of "nx" by "ny" square cells,
    each divided into four triangles. Returns the corner
    edges, numbered as shown below:
    
                  [6]               [5]
                +-->--+-----+-----+--<--+
                |\   /|\   /|\   /|\   /|
                | \ / | \ / | \ / | \ / |
             [7]v  x  |  x  |  x  |  x  v[4]
                | / \ | / \ | / \ | / \ |
                |/   \|/   \|/   \|/   \|
                +-----+-----+-----+-----+
                |\   /|\   /|\   /|\   /|
                | \ / | \ / | \ / | \ / |
                |  x  |  x  |  x  |  x  |
                | / \ | / \ | / \ | / \ |
                |/   \|/   \|/   \|/   \|
                +-----+-----+-----+-----+
                |\   /|\   /|\   /|\   /|
                | \ / | \ / | \ / | \ / |
             [0]^  x  |  x  |  x  |  x  ^[3]
                | / \ | / \ | / \ | / \ |
                |/   \|/   \|/   \|/   \|
                +-->--+-----+-----+--<--+
                  [1]               [2]
      
    Even-numbered arcs lie on the front side, odd-numbered arcs
    lie on the back side.
  *)
    
PROCEDURE GridCell(a: Arc; READONLY x, y: CARDINAL): Arc;
  (*
    Returns the grid cell that lies "x" columns east and "y" rows west
    of the cell whose bottom arc is "a".
    
    More precisely, given an east-pointing arc "a" in the grid, 
    not on the north border, returns another east-pointing
    arc whose origin lies "x" columns east and "y" rows north of Org(a).
  *)
  
TYPE GridVisitProc = PROCEDURE (e: Arc; x, y: CARDINAL);
  
PROCEDURE EnumGridVertices(
    a: Arc; 
    nx, ny: CARDINAL;
    visit: GridVisitProc;
  );
  (*
    Assumes "a" is the east-pointing arc out of the 
    southwest corner of a triangulated grid.
    Calls "visit(e, x, y)" for each vertex of the grid,
    where "e" is one arc out of the vertex, and "x,y"
    are its indices.  
    
    The vertices are numbered from [0,0] to [2*nx, 2*ny];
    both indices are even for cell corners, and both are
    odd for cell centers. 
  *)
    
PROCEDURE SetQuarterGridStyles(
    c: Arc; 
    order: CARDINAL;
    vertexStyle: VertexStyle;
    edgeStyle: EdgeStyle;
    joinStyle: VertexStyle;
    faceStyle: FaceStyle;
    facePatch: PatchNum;
  );
  (*
    Given an unglued grid of "order" by "order" cells, such as
    produced by makeGrid, sets the style and the logical patch numbers 
    in all vertices, edges, and triangles of one quarter of that grid.
    The quarter-grid id the one adjacent to the corner arc "c" 
    
    The "vertexStyle" applies only to the grid corner Org(c);
    the styles of other vertices are left unchanged.

    The "edgeStyle" applies only to the grid edges in half of the
    patch's main diagonal, starting with Oprev(c).
    The "joinStyle" applies to the grid vertices joining those edges.
    
    The "faceStyle" properties applies only to the traingles 
    in one quarter of the grid, between the the grid's side that 
    starts with "c" and the two diagonals of the patch.

    The properties of all other vertices, edges,and faces of the
    patch are left unchanged.
    *)

PROCEDURE Glue(a, b: Arc; n: CARDINAL): Arc;
  (* 
    Identifies "n" edges from Left(a), starting with "a" and
    turning counterclockwise, with "n" edges from Right(b), starting with
    "b" and turning clockwise.  
    
    The edges in the "b" chain are removed; the edges in the "a"
    chain are left in the structure.  Returns the last edge 
    of the "a" chain.  *)

(* === GLOBAL PROCEDURES === *)

PROCEDURE NumberVertices(a: Arc): CARDINAL;
  (* Enumerates all (primal) vertices reachable from "a" by chains
     of "Sym" and "Onext", and assigns them distinct serial numbers 
     from 0 up.  Returns the number of vertices found. *)
    
TYPE
  Topology = RECORD
      NV: INTEGER; (* Number of vertices *)
      NF: CARDINAL; (* Number of faces *)       
      NE: CARDINAL; (* Number of edges *)

      vertex: REF ARRAY OF Vertex;  (* The vertices, indexed by "num" *)
      face : REF ARRAY OF Face;     (* The faces, indexed by "num" *)
      edge: REF ARRAY OF Arc;       (* One arc for each edge *)
      out: REF ARRAY OF Arc;        (* One arc for each vertex where:
                                                Org(out[v]) = vertex[v] *)
      side: REF ARRAY OF Arc;       (* One arc for each face where: 
                                                 Org(side[f]) = face[f] *)
    END;
    (*
      A triangulation's topology in array form. *)

TYPE
  AdjacencyMatrix = ARRAY OF ARRAY OF BOOLEAN;

PROCEDURE MakeTopology(a: Arc): Topology;
  (*
    Builds the vertex/edge tables for the given Triang structure.
    Assumes "NumberVertices" has been called. *)
    
PROCEDURE MakeAdjacencyMatrix(READONLY top: Topology): REF AdjacencyMatrix;
  (*
    Builds the adjacency matrix for the topology "top". *)
    
PROCEDURE GetVariableVertices(READONLY top: Topology; VAR vr: ARRAY OF BOOLEAN);
  (*
    Sets "vr[v] := TRUE" for every vertex "v" that has "exists=TRUE" 
    and "fixed=FALSE". *)
    
PROCEDURE MaxVertexStyle(READONLY v: ARRAY OF Vertex): VertexStyle;
  (*
    Returns the maximum "v[i].style" for all "i". *)

PROCEDURE MaxEdgeStyle(READONLY e: ARRAY OF Edge): EdgeStyle;
  (*
    Returns the maximum "e[i].style" for all "i". *)

PROCEDURE MaxFaceStyle(READONLY f: ARRAY OF Face): FaceStyle;
  (*
    Returns the maximum "f[i].style" for all "i". *)

PROCEDURE MaxPatchNum(READONLY f: ARRAY OF Face): PatchNum;
  (*
    Returns the maximum "f[i].patch" for all "i". *)

(* === GEOMETRIC TOOLS === *)

TYPE
  Coords = ARRAY OF LR3.T; 
  (*
    An assignment of the vertices of some triangulation 
    to points of "R^3". *)

PROCEDURE InitCoords(
    coins: Random.T; 
    VAR c: Coords; 
    radius: REAL := 1.0
  );
  (*
    Fills c with random coordinates in the range
    [-radius __ +radius]. *)

PROCEDURE Barycenter(READONLY top: Topology; READONLY c: Coords): LR3.T;
  (*
    Returns the barycenter of all existing vertices. *)
    
PROCEDURE Displace(READONLY top: Topology; d: LR3.T; VAR c: Coords);
  (*
    Displaces all existing vertices by "d". *)

PROCEDURE Scale(READONLY top: Topology; s: LONGREAL; VAR c: Coords);
  (*
    Scales the coordinates of all existing vertices by "s". *)

PROCEDURE MeanVertexNorm(READONLY top: Topology; READONLY c: Coords): LONGREAL;
  (*
    The average distance of existing vertices from the origin, in the 
    root-mean-square sense; that is,
    "sqrt(sum(norm(c[v])^2, v IN VExist))/|VExist|)". *)
    
PROCEDURE MeanEdgeLength(READONLY top: Topology; READONLY c: Coords): LONGREAL;
  (*
    The average length OF existing edges, in the root-mean-square sense; 
    that is "sqrt(sum(dist(c[org(e)], c[dst(e)])^2, e IN EExist))/|EExist|)". *)

PROCEDURE NormalizeVertexNorms(
    READONLY top: Topology; 
    VAR c: Coords; 
    norm: LONGREAL := 1.0d0;
  );
  (*
    Shifts and scales all existing vertices so that they have
    barycenter (0,0,0) and mean square distance from the origin
    equal to "norm". *)

PROCEDURE NormalizeEdgeLengths(
    READONLY top: Topology; 
    VAR c: Coords;
    length: LONGREAL := 1.0d0;
  );
  (*
    Shifts and scales all existing vertices so that they have
    barycenter (0,0,0), and the mean square length of existing edges
    is "length". *)

PROCEDURE FaceCross(a: Arc; READONLY c: Coords): LR3.T;
  (*
    A vector perpendicular to face "Left(a)", whose length is
    twice the area of the face. Returns the zero vector if
    the face has "exists=FALSE". *)

PROCEDURE FaceNormal(a: Arc; READONLY c: Coords): LR3.T;
  (*
    Normal of face "Left(a)"; same as "LR3.Dir(FaceCross(a, c))".
    Returns an arbitrary unit vector if the face has zero area. *)
    
PROCEDURE FaceBarycenter(a: Arc; READONLY c: Coords): LR3.T;
  (*
    The barycenter of face "LeftF(a)".  Assumes all
    the vertices of the face exist. *)

PROCEDURE VertexCross(a: Arc; READONLY c: Coords): LR3.T;
  (*
    A vector approximately orthogonal to the 
    surface at Org(a), whose length is proportional the
    area of the polygon defined by the existing neighbors of "Org(a)". *)
  
PROCEDURE VertexNormal(a: Arc; READONLY c: Coords): LR3.T;
  (*
    Estimated normal at Org(a), considering only neighbors
    that exist; same as "LR3.Dir(VertexCross(a, c))".
    Returns an arbitrary unit vector if "VertexCross(a, c)" is zero. *)
    
PROCEDURE NeighborBarycenter(a: Arc; READONLY c: Coords): LR3.T;
  (*
    Barycenter of the existing neighbors of "Org(a)". *)
  
PROCEDURE ComputeAllVertexNormals(
    READONLY top: Topology; 
    READONLY c: Coords;
  ): REF ARRAY OF LR3.T;
  (*
    Returns a vector with the result of VertexNormal
    applied to each vertex of "top". *)
  
END Triang.