(* Oct-edge data structure. *)
(* Last edited by Rober Marcone Rosi on Mon Oct 18 1993 *)

MODULE Oct;

IMPORT Lex, Convert, Word, Wr, Rd, Fmt, Thread;
(* IMPORT Stdio; *)

(* The arc Onext(a) is given by arc(a.edge.enext[rbit], a.edge.bnext[rbit]),
   provided that the "flip' bit of a is 0. Otherwise Onext(a) is computed by 
   the formula: Onext(a) = Flip(Rot(Onext(Rot(Flip(a)))))
*)
VAR NoEdge: CARDINAL := 0;
TYPE 
  MarkBit = (* BITS 1 FOR *) BOOLEAN;

PROCEDURE SenseBit(a: Arc): SBit =
  BEGIN RETURN Word.RightShift(a.bits, 2) END SenseBit;

(* If a = flip^f(rot^r(e)), where e is the reference arc of edge,
   r in [0..3], and f in [0..1], then a.bits = (2r + f).
   Therefore,
   
     a.flip.bits = XOR(a.bits, 1)
     
     If fbit is 0
       a.rot.bits = (a.bits + 2) MOD 8 
     otherwise
       a.rot.bits = (a.bits + 6) MOD 8
*)

PROCEDURE FlipBit (a: Arc): FBit =
  BEGIN RETURN Word.And(a.bits, 1) END FlipBit;
  
PROCEDURE DualBit (a: Arc): DBit =
  BEGIN RETURN Word.And(Word.RightShift(a.bits, 1), 1) END DualBit;

PROCEDURE RotBits (a: Arc): RBits =
  BEGIN RETURN Word.And( Word.RightShift(a.bits, 1), 3) END RotBits;
  
PROCEDURE GetArcNo (a: Arc): ArcNo =
  BEGIN RETURN 8 * a.edge.no + a.bits END GetArcNo;

PROCEDURE GetEdgeNo(e: Edge): EdgeNo =
  BEGIN RETURN e.no END GetEdgeNo;
  
PROCEDURE SetEdgeNo(e: Edge; no: EdgeNo) =
  BEGIN e.no := no END SetEdgeNo;
  
PROCEDURE GetOrgNo (a: Arc): VertexNo =
  BEGIN RETURN a.edge.orgno[RotBits(a)] END GetOrgNo;
  
PROCEDURE SetOrgNo (a: Arc; no: VertexNo) =
  BEGIN a.edge.orgno[RotBits(a)] := no END SetOrgNo;

PROCEDURE Rot(a: Arc): Arc =
  BEGIN
    RETURN Arc{edge := a.edge, 
               bits := Word.And( a.bits + 2 + Word.LeftShift( 
               Word.And(a.bits, 1), 2), 7) };
  END Rot;

PROCEDURE Flip(a: Arc): Arc =
  BEGIN
    RETURN Arc{edge := a.edge, bits := Word.Xor( a.bits, 1 )};
  END Flip;
  
PROCEDURE Sym(a: Arc): Arc =
  BEGIN
    RETURN Arc{edge := a.edge, bits := Word.And( (a.bits + 4), 7 )};
  END Sym;
  
PROCEDURE Tor(a: Arc): Arc =
  BEGIN
    RETURN Arc{edge := a.edge,
               bits := Word.And( a.bits + 6 + Word.LeftShift( 
                       Word.And(a.bits, 1), 2), 7) };    
  END Tor;

PROCEDURE Onext(a: Arc): Arc =
  BEGIN
    WITH s = Word.And( Word.RightShift( (a.bits + 1), 1 ), 3 ) DO     
      WITH 
      r = a.edge.bnext[s] 
      DO
        a.edge := a.edge.enext[s];
        IF ( Word.And( a.bits, 1 ) = 1) THEN 
          a.bits := Word.Xor( Word.And( r + 2 + Word.LeftShift( 
                    Word.And(r, 1), 2), 7), 1) ;
        ELSE
          a.bits := r;
        END;
      END
    END;
    RETURN a;
  END Onext;

PROCEDURE Oprev(a: Arc): Arc =
  BEGIN RETURN Rot(Onext(Rot(a))) END Oprev;
  
  
PROCEDURE Dnext(a: Arc): Arc =
  BEGIN RETURN Sym(Onext(Sym(a)))  END Dnext;

PROCEDURE Dprev(a: Arc): Arc =
  BEGIN RETURN Tor(Onext(Tor(a)))  END Dprev;
  
PROCEDURE Lnext(a: Arc): Arc =
  BEGIN RETURN Rot(Onext(Tor(a)))  END Lnext;
  
PROCEDURE Lprev(a: Arc): Arc =
  BEGIN RETURN Sym(Onext(a))  END Lprev;
  
PROCEDURE Rnext(a: Arc): Arc =
  BEGIN RETURN Tor(Onext(Rot(a)))  END Rnext;

PROCEDURE Rprev(a: Arc): Arc =
  BEGIN RETURN Onext(Sym(a))  END Rprev;
  
PROCEDURE Make (no: EdgeNo): Arc =
  BEGIN
    WITH 
      e = NEW(Edge).init(no := no)
    DO
      RETURN Arc{edge := e, bits := 0}
    END;
  END Make;

<*UNUSED*>
PROCEDURE Debug() =
  BEGIN 
  END Debug;

PROCEDURE Splice (a, b: Arc) =
<* FATAL Wr.Failure, Thread.Alerted *>

  BEGIN
  <* ASSERT b # Flip(Onext(a)) *>
  IF a # b THEN 
    WITH 
      ta = Onext(a), tb = Onext(b),
      c = Rot(ta), d = Rot(tb), 
      tc = Onext(c), td = Onext(d)
    DO
      IF ( Word.And(a.bits, 1) = 0 ) THEN 
        WITH
          rba = RotBits(a)
        DO
          a.edge.enext[rba] := tb.edge;
          a.edge.bnext[rba] := tb.bits;
        END
      ELSE
        WITH
          ra = Rot(Flip(a)), fd = Flip(d), rbra = RotBits(ra)
        DO
          ra.edge.enext[rbra] := fd.edge;
          ra.edge.bnext[rbra] := fd.bits;
        END
      END;
        
      IF ( Word.And(b.bits, 1) = 0 ) THEN 
        WITH
          rbb = RotBits(b)
        DO
          b.edge.enext[rbb] := ta.edge;
          b.edge.bnext[rbb] := ta.bits;          
        END
      ELSE
        WITH
          rb = Rot(Flip(b)), fc = Flip(c), rbrb = RotBits(rb)            
        DO
          rb.edge.enext[rbrb] := fc.edge;
          rb.edge.bnext[rbrb] := fc.bits;
        END
      END;
      

      IF ( Word.And(c.bits, 1) = 0 ) THEN 
        WITH
          rbc = RotBits(c) 
        DO
          c.edge.enext[rbc] := td.edge;
          c.edge.bnext[rbc] := td.bits;
        END
      ELSE
        WITH
          rc = Rot(Flip(c)), fb = Flip(b), rbrc = RotBits(rc)
        DO
          rc.edge.enext[rbrc] := fb.edge;
          rc.edge.bnext[rbrc] := fb.bits;
        END          
      END;

      IF ( Word.And(d.bits, 1) = 0 ) THEN
      WITH
        rbd = RotBits(d)
      DO
        d.edge.enext[rbd] := tc.edge;
        d.edge.bnext[rbd] := tc.bits;          
      END
      ELSE
        WITH
          rd = Rot(Flip(d)), fa = Flip(a), rbrd = RotBits(rd) 
        DO
          rd.edge.enext[rbrd] := fa.edge;
          rd.edge.bnext[rbrd] := fa.bits;
        END
      END;
    END;
  END;  
  END Splice;

PROCEDURE SetOnext(a, b: Arc) =
  BEGIN
    IF Onext(a) # b THEN Splice(a, Oprev(b)) END
  END SetOnext;

PROCEDURE Enum(a: Arc; visit: VisitProc; edges: BOOLEAN := FALSE)=
  CONST stacksize = 16*1024; (* Guess: max edges in a connected component *)
  VAR stack: ARRAY [0..stacksize-1] OF Edge;
      bitstack: ARRAY [0..stacksize-1] OF (* BITS 3 FOR *) [0..7];
      top: [0..stacksize];
      
  PROCEDURE VisitandMark (c: Arc)=
    (* If edge(c) is unmarked: visit, mark, and stack it. *)
    BEGIN
      IF NOT c.edge.marks[FlipBit(c)] THEN 
(*      IF NOT c.edge.marks[0] THEN        *)
        visit(c);
        IF NOT edges THEN visit(Sym(c)) END;
        <* ASSERT top < stacksize *>
        c.edge.marks[FlipBit(c)] := TRUE; 
        stack[top] := c.edge; bitstack[top] := c.bits;
        top := top + 1;
      END;
    END VisitandMark;
  
  VAR seen: [0..stacksize]; (* # of quads whose childeren were looked at *)
      atop: Arc;
  BEGIN
    <* ASSERT NOT a.edge.marks[FlipBit(a)] *>
    top := 0;
    seen := 0;
    TRY
      VisitandMark (a);
      WHILE seen < top DO 
        atop.edge := stack[seen]; atop.bits := bitstack[seen];
        VisitandMark(Onext(atop));
        VisitandMark(Onext(Sym(atop)));
        seen := seen + 1
      END;
    FINALLY
      (* Erase all marks *)
      WHILE top > 0 DO
        top := top - 1;
        atop.edge := stack[top]; atop.bits := bitstack[top];
        atop.edge.marks[FlipBit(atop)] := FALSE; 
      END;

(*      WHILE top > 0 DO
        top := top - 1;
        stack[top].marks[FlipBit(stack[top]) := FALSE
      END;
   *)
    END;
  END Enum;

PROCEDURE EnumVertices (a: Arc; visit: VisitProc)=
  CONST stacksize = 16*1024; (* Guess: max quads in a connected component *)
  VAR stack: ARRAY [0..stacksize-1] OF Edge;
      bitstack: ARRAY [0..stacksize-1] OF (* BITS 3 FOR *) [0..7];
      top: [0..stacksize]; (* # of visited(=marked=stacked) vertices *)
  
  PROCEDURE VisitandMark (c: Arc)=
    (* If org(c) is unmarked: visit, mark, and stack it *)
    VAR cn: Arc;
    BEGIN
      IF NOT c.edge.marks[SenseBit(c)] THEN
        visit(c);
        <* ASSERT top < stacksize *>
        cn := c;
        REPEAT
          cn.edge.marks[SenseBit(cn)] := TRUE;
          cn := Onext(cn);
        UNTIL (cn = c);
        stack[top] := c.edge; bitstack[top] := c.bits;
        top := top + 1;
      END;
    END VisitandMark;
  
  VAR seen: [0..stacksize]; (* # of vertices whose neighbs were looked at *)
      atop, atn: Arc;
  BEGIN
    <* ASSERT NOT a.edge.marks[0] *>
    top := 0;
    seen := 0;
    TRY
      VisitandMark (a);
      WHILE seen < top DO
        atop.edge := stack[seen]; atop.bits := bitstack[seen];
        atn := atop;
        REPEAT
          VisitandMark(Sym(atn));
          atn := Onext(atn);
        UNTIL (atn = atop);
        seen := seen + 1;
      END; 
    FINALLY
      (* Erase all marks *)
      WHILE top > 0 DO
        top := top - 1;
        atop.edge := stack[top]; atop.bits := bitstack[top];
        atn := atop;
        REPEAT
          atn.edge.marks[SenseBit(atn)] := FALSE;
          atn:= Onext(atn);
        UNTIL (atop = atn);
      END;
    END;
  END EnumVertices;
  
PROCEDURE NumberArcs(READONLY a: ARRAY OF Arc; primal: BOOLEAN): REF ARRAY OF Arc =
  VAR 
    sz: CARDINAL := 1024;     (* Number of edges allocated on stack *)
    ne: CARDINAL := 0;        (* Number of edges seen *)
    estack := NEW(REF ARRAY OF Edge, sz);  (* Edge stack *)
    
    na: CARDINAL := 0;    (* Number of arcs seen *)
    astack := NEW(REF ARRAY OF Arc, 8 * sz);  (* Arc stack *)
      
  (* An arc "a" is "marked" if astack[8 * a.edge.no + a.bits] = a. *)

  PROCEDURE VisitAndMark (t: Arc) =
    (* If t is unmarked: visit, mark, and stack it. *)
    BEGIN
      WITH teno = t.edge.no, tno = 8*teno + t.bits DO
        IF teno < ne AND astack[tno] = t THEN
          (* Arc is already marked *)
        ELSE
          IF teno < ne AND estack[teno] = t.edge THEN
            (* New arc of old edge *)
          ELSE
            (* New edge *)
            IF nstack = sz THEN
              WITH newstack = NEW(REF ARRAY OF Arc, 2 * sz) DO
                SUBARRAY(newstack^, 0, sz) := stack^;
                stack := newstack
              END;
            END;
            t.edge.no := nstack;
            stack[nstack] := t;
            INC(nstack);
          t.edge.no := nstack;
          stack[nstack] := t;
          INC(nstack);
        END
      END
    END VisitAndMark;
  
  VAR 
    seen: CARDINAL := 0; (* # of edges whose childeren were looked at *)
    s: Arc;
  BEGIN
    nstack := 0;
    seen := 0;
    FOR i := 0 TO LAST(a) DO
      VisitAndMark (a[i])
    END;
    WHILE seen < nstack DO
      s := stack[seen];
      VisitAndMark(Onext(s));
      VisitAndMark(Onext(Sym(s)));
      seen := seen + 1
    END;
    WITH
      r = NEW(REF ARRAY OF Arc, nstack)
    DO
      r^ := SUBARRAY(stack^, 0, nstack);
      RETURN r
    END;
  END NumberEdges;

PROCEDURE NumberVertices(READONLY e: ARRAY OF Arc): REF ARRAY OF Arc =
  VAR 
    sz: CARDINAL := 1024;
    stack := NEW(REF ARRAY OF Arc, sz);
    
    nstack: CARDINAL := 0;
      
  (* An arc "a" is "marked" if stack[Get].edge = e. *)

  PROCEDURE VisitAndMark (t: Arc) =
    (* If t is unmarked: visit, mark, and stack it. *)
    BEGIN
      IF t.edge.no < nstack AND stack[t.edge.no].edge = t.edge THEN
        (* Edge is already marked *)
      ELSE
        IF nstack = sz THEN
          WITH newstack = NEW(REF ARRAY OF Arc, 2 * sz) DO
            SUBARRAY(newstack^, 0, sz) := stack^;
            stack := newstack
          END;
        END;
        t.edge.no := nstack;
        stack[nstack] := t;
        INC(nstack);
      END
    END VisitAndMark;
  
  VAR 
    seen: CARDINAL := 0; (* # of edges whose childeren were looked at *)
    s: Arc;
  BEGIN
    nstack := 0;
    seen := 0;
    FOR i := 0 TO LAST(a) DO
      VisitAndMark (a[i])
    END;
    WHILE seen < nstack DO
      s := stack[seen];
      VisitAndMark(Onext(s));
      VisitAndMark(Onext(Sym(s)));
      seen := seen + 1
    END;
    WITH
      r = NEW(REF ARRAY OF Arc, nstack)
    DO
      r^ := SUBARRAY(stack^, 0, nstack);
      RETURN r
    END;
  END NumberVertices;

PROCEDURE NumberVerticesAndFaces(
    READONLY e: ARRAY OF Arc
  ): RECORD v, f: REF ARRAY OF Arc END =
  (* 
    Assumes "e[i]" is one of the arcs of edge number "i".
    Numbers the origin and destination vertices of all those arcs,
    and also the left and right faces.
    Returns a vector "v" with one arc out of each vertex,
    and a vector "f" with one arc on the ccw boundary of each face. *)

  VAR nv: CARDINAL := 0;
      nf: CARDINAL := 0;
      
      v: REF ARRAY OF Arc;
      f: REF ARRAY OF Arc;

  BEGIN
    WITH
      NE = NUMBER(e),
      rep = NEW(REF ARRAY OF ArcNo, 8*NE)^
        (* 
          Define #a = GetArcNo(a), a[i] = arc such that #a = i.
          Invariant: rep[#a] = #a if arc "a" is the least arc out of its origin;
          otherwise rep[#a] = #b where b is some lesser arc with Org(a)=Org(b). *)
    DO
      PROCEDURE Identify(a: Arc) =
        (* Records the fact that "a" and "Onext(a)" have the same origin *)
        VAR ia := GetArcNo(a);
            ib := GetArcNo(Onext(a));
        BEGIN
          WHILE rep[ia] # ia DO ia := rep[ia] END;
          WHILE rep[ib] # ib DO ib := rep[ib] END;
          IF ia < ib THEN
            rep[ib] := ia
          ELSIF ib < ia THEN
            rep[ia] := ib
          ELSE
            (* Ok *)
          END
        END Identify;
        
      BEGIN
        (* Initialize "rep" *)
        FOR ia := 0 TO LAST(rep) DO rep[ia] := ia END;
        
        (* Identify arcs with same origin: *)
        FOR ie := 0 TO LAST(e) DO
          VAR ee := e[ie];
          BEGIN
            FOR r := 0 TO 3 DO
              Identify(ee); Identify(Flip(ee));
              ee := Rot(ee)
            END;
          END
        END;
        
        (* Count and collect distinct vertices and faces: *)
        
        PROCEDURE Count(a: Arc; VAR n: CARDINAL) =
          BEGIN
            WITH ia = GetArcNo(a) DO
              IF rep[ia] = ia THEN INC(n) END;
            END
          END Count;
          
        PROCEDURE Collect(a: Arc; VAR n: CARDINAL; VAR r: ARRAY OF Arc) =
          BEGIN
            WITH ia = GetArcNo(a) DO
              IF rep[ia] = ia THEN r[n] := a; INC(n) END;
            END
          END Collect;
          
        BEGIN
          nv := 0; nf := 0;
          FOR ie := 0 TO LAST(e) DO
            WITH ei = e[i] DO
              Count(ei, nv);
              Count(Rot(ei), nf);
              Count(Sym(ei), nv);
              Count(Tor(ei), nf)
            END;
          END;
          v := NEW(REF ARRAY OF Arc, nv);
          f := NEW(REF ARRAY OF Arc, nf);
          nv := 0; nf := 0;
          FOR ie := 0 TO LAST(e) DO
            WITH ei = e[i] DO
            
              Collect(ei, nv, v^);
              Collect(Rot(ei), nf, f^);
              Collect(Sym(ei), nv, v^);
              Collect(Tor(ei), nf, f^)
            END;
          END;
        END
      END
    END;
    RETURN 
  END NumberVerices;

PROCEDURE PrintArc (wr: Wr.T; a: Arc) =
  <* FATAL Wr.Failure, Thread.Alerted *>
  BEGIN
    Wr.PutText (wr, 
      Fmt.Int(a.edge.no) & ":" &
      Fmt.Int(RotBits(a)) & ":" &
      Fmt.Int(FlipBit(a))
    )
  END PrintArc;
  
CONST
  Digits = SET OF CHAR{'0'..'9'};

PROCEDURE ReadArc(rd: Rd.T; READONLY map: ARRAY OF Edge): Arc =
VAR f, n, r: CARDINAL;
  <* FATAL Rd.Failure, Rd.EndOfFile, Thread.Alerted, Convert.Failed, Lex.Error *>
  BEGIN
    Lex.Skip(rd);
    n := Lex.Int(rd, cs := Digits);
    <* ASSERT Rd.GetChar(rd) = ':' *>
    r := Lex.Int(rd, cs := Digits);
    <* ASSERT Rd.GetChar(rd) = ':' *>
    f := Lex.Int(rd, cs := Digits);
    <* ASSERT n < NUMBER(map) *>
    <* ASSERT r < 4 AND f < 2 *>
    RETURN Arc{edge := map[n], bits := 2 * r + f}
  END ReadArc;

PROCEDURE PrintEdge(wr: Wr.T; e: Edge) =
  <* FATAL Wr.Failure, Thread.Alerted *>
  BEGIN
    FOR i := 0 TO 3 DO 
      IF i > 0 THEN Wr.PutChar(wr, ' ') END;
      PrintArc(wr, Arc{edge := e.enext[i], bits := e.bnext[i]});
    END
  END PrintEdge;
  
PROCEDURE ReadEdge(rd: Rd.T; e: Edge; READONLY map: ARRAY OF Edge) =
  BEGIN
    FOR i := 0 TO 3 DO 
      SetOnext(Arc{edge := e, bits := 2 * i}, ReadArc(rd, map))
    END
  END ReadEdge;

BEGIN
END Oct.