(* From ReadTopology *)

      (* Read edge data: *)
      EVAL Mis.ReadCommentsJS(rd, '|');
      FOR j := 0 TO top.NE-1 DO 
        Lex.Skip(rd);
        WITH 
          ne = Lex.Int(rd),      (* index to edge *) 
          e = top.edge[ne]
        DO
          <* ASSERT ne = j *>
          e.num := ne;
          e.pa := Map3D.ReadHandle(rd,map);
        END;
      END;
      Lex.Skip(rd);

      (* Read face data: *)
      EVAL Mis.ReadCommentsJS(rd, '|');
      FOR j := 0 TO top.NF-1 DO 
        Lex.Skip(rd);
        WITH 
          nf = Lex.Int(rd),      (* index to face *) 
          f = top.face[nf]
        DO
          <* ASSERT nf = j *>
          f.num := nf;
          f.pa := Map3D.ReadHandle(rd,map);
        END;
      END;
      Lex.Skip(rd);

======================================================================

PROCEDURE DeleteFaceEdge (a : Handle) =
  BEGIN
    SpliceFacets(a, prevF(a));
    SpliceFacets(Flip(a), prevF(Flip(a))); 
  END DeleteFaceEdge;

From Topology:
      
      vertex: REF ARRAY OF Vertex;
      edge: REF ARRAY OF Edge;
      face: REF ARRAY OF Face; 
      cell: REF ARRAY OF Cell;   

      der: CARDINAL; (* Edge Ring Degree *)
      bdr: CARDINAL; (* Where (0) indicates without boundary *)
                     (*       (1) indicates with    boundary *)
                     (*   and (2) indicates that cells are octahedra *)

From MakeTopology:

    top.cell := NEW(REF ARRAY OF Cell, top.NP);
    top.region     := NEW(REF ARRAY OF Handle, top.NP);
    top.der        := Map3D.DegreeEdgeRing(top.face[0].pa);
    top.bdr        := bdr;


    FOR i:= 0 TO top.NFE-1 DO
      VAR
        c: Handle := top.fe[i];
        v,p: Node;
        vi,pi: CARDINAL;
      BEGIN
        FOR k:= 0 TO 1 DO
          v := Org(c); 
          vi := v.num;
          top.vertex[vi] := v;
          top.out[vi] := c;
          p := Pneg(c);
          IF p # NIL THEN
            pi := p.num;
            top.cell[pi] := p;
            top.region[pi] := RevE(Dual(c);
            <* ASSERT Pneg(c) = Org(top.region[pi]) *>
          END;
          c := Flip(c); 
        END
      END
    END;

  PROCEDURE Der(READONLY top: TOPOLOGY): CARDINAL =
    (* Maximum number of sides in any face *)
    VAR an: Handle; der: CARDINAL := 0; d: CARDINAL;
    BEGIN
      FOR i := 0 TO top.NF-1 DO
        WITH a = top.fref[i] DO
          an := a; d := 0;
	  REPEAT INC(d); an := nextE(an) UNITL an = a
	END;
	IF d > der THEN der := d END;
      END;
      RETURN der
    END Der;

  PROCEDURE Bdr(READONLY top: Topology): CARDINAL =
    (* 1 if the map has a free border (i.e. any cell is NIL or non-existent *)
    BEGIN
      FOR i := 0 TO top.NF-1 DO
        WITH a = top.fref[i] DO
	  IF (NOT Exists(Pneg(a))) OR (NOT Exists(Ppos(a))) THEN RETURN 1 END
	END
      END;
      RETURN 0
    END Bdr;

      WriteCommentsJS(tp, "\nEdge data:\n", '|');
      FOR i := 0 TO top.NE-1 DO
        WITH er = top.eref[i] DO 
          Wr.PutText(tp, Fmt.Pad(Fmt.Int(i), eWidth)); 
          Wr.PutText(tp, "  ");
          Map3D.PrintHandle(tp, er, eWidth+1, TRUE);  
        END
      END;

      WriteCommentsJS(tp, "\nFace data:\n", '|');
      FOR i := 0 TO top.NF-1 DO
        WITH fr = top.fref[i] DO 
          Wr.PutText(tp, Fmt.Pad(Fmt.Int(i), fWidth)); 
          Wr.PutText(tp, "  ");
          Map3D.PrintHandle(tp, fr, fWidth+1, TRUE);  
        END
      END;

        Wr.PutText(tp, "der        ");
        Wr.PutText(tp, Fmt.Pad(Fmt.Int(Der(top)),m) & "\n");
        Wr.PutText(tp, "bdr        ");
        Wr.PutText(tp, Fmt.Pad(Fmt.Int(Bdr(top)),m) & "\n");
 
PROCEDURE NumberNeighborVertex(neigh: REF ARRAY OF Vertex): CARDINAL;
  (*
     Return the number neighbor vertex of OrgV(a). *)



PROCEDURE Neighbors(a: Handle;  READONLY top: Topology): REF ARRAY OF Vertex =
  CONST
    IniStackSize = 1024;
  VAR
    vstack := NEW(REF ARRAY OF Vertex, IniStackSize);
    nstack: CARDINAL := 0;
    v : CARDINAL;

  BEGIN
    v := OrgV(a).num;
    FOR i := 0 TO top.NE-1 DO
      WITH
        ei = top.edge[i],
	ei0 = ei.vertex[0].num,
	ei1 = ei.vertex[1].num
      DO
	IF (v = ei0) OR (v = ei1)
	  THEN
            IF v = ei0
	      THEN
	        vstack[nstack] := top.vertex[ei1];
                INC(nstack);
	      ELSE
	        vstack[nstack] := top.vertex[ei0];
                INC(nstack);
            END
        END
      END
    END;          
      
    WITH
      r = NEW(REF ARRAY OF Vertex, nstack)
    DO
      r^ := SUBARRAY(vstack^, 0, nstack);
      RETURN r;
    END;
  END Neighbors;
  


PROCEDURE WriteTable(
    name: TEXT; 
    READONLY top: Topology;
    comments: TEXT := " ";
    debug: BOOLEAN := FALSE;
    ); 
  (* This procedure create tables of vertex/edges/faces/cells for one
     triangulation. 
  *)


PROCEDURE NumberVertices(a: Handle): CARDINAL;
  (* Enumerates all (primal) vertices reachable from "a" by chains of
     "Onext", "Onext_1" and "Flip", and assigns them distinct serial
     numbers from 0 up.  Returns the number of vertices found.Use the
     procedure travesal "EnumVertices" of library "libm3triang".
  *)

(*  
PROCEDURE ReforceSettings(a: Handle);
*)


PROCEDURE SetNextPneg(a: Handle; n: Node) =
  VAR t: Handle := a;
  BEGIN
    REPEAT
      SetPneg(t, n);
      t := prevE(t);
    UNTIL t = a;
  END SetNextPneg;

TYPE
  FaceEdge <: PublicFaceEdge;
  PublicFaceEdge = Map3D.FaceEdge OBJECT
      mark: BOOLEAN := FALSE;        (* Mark the fe, used by RefineT.*)
      vh: Vertex;                    (* Medial vertex associated to fe.
                                        This two last atributes will be used
                                        for the refinement process.         *)
    METHODS
      init(): FaceEdge;
        (* initializes "self" as an isolated fe of triangulation *)
    END;




PROCEDURE SetNextPpos(a : Handle; n: Node) =
  VAR t : Handle := Flip(a);
  BEGIN
    SetAllPneg(t,n);
    REPEAT
      SetAllPneg(prevF(t),n);
      t := prevE(t); 
    UNTIL t = Flip(a); 
  END SetNextPpos;



PROCEDURE FaceEdgeswso(READONLY a: Handle): REF ARRAY OF Handle =

  CONST IniStackSize = 10000;  (* 1024 *)
  VAR stack := NEW(REF ARRAY OF Handle, IniStackSize);
      nstack: CARDINAL := 0;
      estack : REF ARRAY OF Handle;
      degree : CARDINAL;
    BEGIN
      estack := Map3D.NumberEdgesForDegree(ARRAY OF Handle{a});
      degree := NUMBER(estack^);
      FOR i := 0 TO degree-1 DO
        VAR c: Handle := estack[i];
            ct : Handle;
        BEGIN
          IF Org(a) = Org(c)
            THEN
              ct := c;
              REPEAT
                stack[nstack] := ct;
                INC(nstack);
                ct := prevF(ct);
              UNTIL (ct = c );
            ELSE
              ct := Flip(c);
              REPEAT
                stack[nstack] := ct;
                INC(nstack);
                ct := prevF(ct);
              UNTIL (ct = Flip(c) );
           END;
        END;
      END;
      WITH
        r = NEW(REF ARRAY OF Handle, nstack)
      DO
        r^ := SUBARRAY(stack^, 0, nstack);
      RETURN r
    END;
  END FaceEdgeswso;

PROCEDURE RescueOrigins(a,b: Handle) =
  BEGIN
    WITH
     na = OrgV(a),
     cb = FaceEdgeswso(b)
    DO
      FOR j := 0 TO LAST(cb^) DO
        SetAllOrgs(cb[j],na);
      END
    END;
  END RescueOrigins;

PROCEDURE RescueOnFace(a,b: Handle) =
  BEGIN
    WITH
      a0 = a,           b0 = b, 
      a1 = nextE(a),    b1 = nextE(b),
      a2 = prevE(a),  b2 = prevE(b)
    DO
      RescueOrigins(a0,b0);
      RescueOrigins(a1,b1);
      RescueOrigins(a2,b2);
    END
  END RescueOnFace;

PROCEDURE RescueOnTriangFace(a,b: Handle; n: CARDINAL) =
  VAR
    ta,tb: ARRAY [0..99] OF Handle;
  BEGIN
    (* sanity check *)
    <* ASSERT n >= 1 *>
    ta[0] := a; tb[0] := b;

    IF n > 1 THEN
      FOR i := 1 TO n-1 DO
        ta[i] := Flip(prevE(prevF(nextE(ta[i-1]))));
        tb[i] := Flip(prevE(nextF(nextE(tb[i-1]))));
        <* ASSERT ta[i] # a *>
        <* ASSERT tb[i] # b *>
      END
    END;

    
    FOR j := 0 TO n-1 DO
      RescueOrigins(ta[j],tb[j]);
    END;

    RescueOrigins(nextE  (ta[n-1]),nextE  (tb[n-1]));  
    RescueOrigins(prevE(ta[n-1]),prevE(tb[n-1])); 
    
  END  RescueOnTriangFace;

PROCEDURE DegreeOfFaceEdges(a: Handle) : CARDINAL =
  BEGIN
    WITH r = FaceEdgeswso(a), n = NUMBER(r^) DO
      RETURN n;
    END;
  END DegreeOfFaceEdges;

PROCEDURE Edgeswso(READONLY a: Handle): REF ARRAY OF Handle =
  CONST IniStackSize = 1024;
  VAR stack := NEW(REF ARRAY OF Handle, IniStackSize);
      nstack: CARDINAL := 0;
      estack : REF ARRAY OF Handle;
      degree : CARDINAL;

    BEGIN
      estack := Map3D.NumberEdgesForDegree(ARRAY OF Handle{a});
      degree := NUMBER(estack^);

      (* Sanity Check *)
      FOR i := 0 TO degree-1 DO
        VAR ei := estack[i];
        BEGIN
          FOR j := i+1 TO degree-1 DO
            <* ASSERT ei # estack[j] *>
          END;
        END;
      END;
      
      FOR i := 0 TO degree-1 DO
        VAR c: Handle := estack[i];
            ct : Handle;
        BEGIN
          IF Org(a) = Org(c)
            THEN
              ct := c;
              stack[nstack] := ct;
              INC(nstack);
            ELSE
              ct := Flip(c);
              stack[nstack] := ct;
              INC(nstack);  
          END;
        END;
      END;
      WITH
        r = NEW(REF ARRAY OF Handle, nstack)
      DO
        r^ := SUBARRAY(stack^, 0, nstack);
      RETURN r
    END;
  END Edgeswso;

  PROCEDURE Present(c: Handle) : BOOLEAN =
    (* Return TRUE if org(c) its on the stack, FALSE c.c. *)
    VAR nstack1 : CARDINAL := nstack;
    BEGIN
      WHILE nstack1 > 0 DO
        nstack1 := nstack1 - 1;
        IF ostack[nstack1] = Org(c) THEN RETURN TRUE END;
      END;
      RETURN FALSE;
    END Present;

  PROCEDURE VisitAndMarkIni(c: Handle)=
    (* If org(c) is diferent that the origins handles stacked, stack the handle
       "c".
    *)
    VAR cn : Handle;
    BEGIN
      WITH
        p = Org(c),
        q = NIL
      DO
        IF p # q THEN
          visit(c);
          
          (* stack *)
          IF nstack >= NUMBER(ostack^) THEN DoubleStack() END;
          ostack[nstack] := Org(c);
          nstack := nstack + 1;

	  Stack(c);
          cn := c;
          REPEAT
            Stack(cn);
            cn := Onext(cn);
          UNTIL (cn = c)
        END
      END;
  END VisitAndMarkIni;

PROCEDURE ReforceSettings(a: Handle) = 
  
  PROCEDURE Visit(c: Handle) =
  
    PROCEDURE Visit_(c: Handle) =
      VAR cn,dn: Handle;  
      BEGIN
        WITH v = Org(c) DO
          cn := c;
          REPEAT
            SetOrg(cn,v);
            dn := Flip(prevE(cn));
            REPEAT
              SetOrg(dn, Org(cn));
              dn := nextF(dn);
            UNTIL ( dn = Flip(prevE(cn)) );
            cn := nextF(cn);
          UNTIL ( cn = c );
        END;
      END Visit_;

    VAR cn: Handle := c;
    BEGIN
      REPEAT        
        Visit_(cn);
        cn := Onext(cn);
      UNTIL cn = c;
    END Visit;

  BEGIN
    EnumVertices(a, Visit);  
  END  ReforceSettings;

  
PROCEDURE RescueOrigins(a,b: Handle);
  (* In the gluing procedure the handle "b" will be killed, so we need
     updating the origin relations with the origins of the handle "a".
     A general way, (perhaps not optimized) is done by the procedure
     above.
  *)

PROCEDURE RescueOnFace(a,b: Handle);
  (* Realize the "RescueOrigins" procedure above for each fe on the
     face a.fe.face.
  *)

PROCEDURE RescueOnTriangFace(a,b: Handle; order: CARDINAL);
  (* Realize the "RescueOnFace" procedure above for each fe on the
     triangulated face, when "a.fe.face" is a corner face.
  *)

PROCEDURE FaceEdgeswso(READONLY a: Handle) : REF ARRAY OF Handle ;
  (* Return the handles fes with the same origin of handle "a". *)

PROCEDURE FaceEdgeswso(READONLY a: Handle) : REF ARRAY OF Handle ;
  (* Return the handles fes with the same origin of handle "a". *)