MODULE SceneShade;

IMPORT Wr, Fmt, Thread;
IMPORT Scene, H3, R4;

FROM Scene IMPORT 
  FaceList, MaxLights,
  CountFaces, CountEdges, CountVertices;

(* DEBUG: *)
FROM Stdio IMPORT stderr;
CONST Debug = FALSE;

PROCEDURE ClipToWindow(
    VAR sc: FaceList;   (* Scene to clip. *) 
    READONLY w2i: H3.PMap;  (* World to Image transform matrix. *)
    xmin, xmax: REAL;       (* Image X range. *)
    ymin, ymax: REAL;       (* Image Y range. *)
    zmax: REAL;             (* Near clipping plane. *)
  ) =
  VAR plane: ARRAY [0..5] OF H3.Plane;
      side: ARRAY [0..5] OF H3.Sign;
      ins, bds, ots: Scene.T;
  BEGIN
    (* Build image volume: *)
    plane[0] := H3.Plane{f := R4.T{-xmin, +1.0, 00.0, 00.0}};
    plane[1] := H3.Plane{f := R4.T{+xmax, -1.0, 00.0, 00.0}};
    plane[2] := H3.Plane{f := R4.T{-ymin, 00.0, +1.0, 00.0}};
    plane[3] := H3.Plane{f := R4.T{+ymax, 00.0, -1.0, 00.0}};
    plane[4] := H3.Plane{f := R4.T{+zmax, 00.0, 00.0, -1.0}};
    FOR i := 0 TO 4 DO
      plane[i] := H3.InvMapPlane(plane[i], w2i);
      side[i] := +1;
    END;
    plane[5] := H3.Plane{f := R4.T{1.0, 00.0, 00.0, 00.0}};
    side[5] := +1;
    
    (* Clip against image volume: *)
    ins.t := sc; ins.e := NIL; ins.v := NIL;
    Scene.CloseScene(ins);
    bds := Scene.Empty; 
    ots := Scene.Empty;
    Scene.ClipSceneToRegion(
      ins, bds, ots, plane, side, TRUE, TRUE, FALSE, debug := TRUE
    );
    
    (* Merge the "bds" and "ins" lists: *)
    <* ASSERT ots = Scene.Empty *>
    sc := ins.t;
    Scene.MoveFaceList(bds.t, sc);
  END ClipToWindow;
  
PROCEDURE Visible(
    VAR tri: FaceList;  (* Faces to paint *)
    opq: FaceList;      (* Opaque triangles *)
    READONLY obs: H3.Point;
  ) =
  VAR vis, inv: FaceList;
  BEGIN
    vis := NIL; inv := NIL;
    AuxVisible(tri, opq, obs, FALSE, vis, inv);
    <* ASSERT inv = NIL *>
    tri := vis;
  END Visible;
  
PROCEDURE Shadow(
    VAR tri: FaceList;       (* Faces to illuminate *)
    opq: FaceList;           (* Opaque triangles *)
    READONLY lightPos: H3.Point;
    lightNo: [0..MaxLights-1];
  ) =
  VAR vis, inv, tmp: FaceList;
  BEGIN
    vis := NIL; inv := NIL;
    AuxVisible(tri, opq, lightPos, TRUE, vis, inv);
    tri := NIL;
    WHILE vis # NIL DO
      vis.t.lights[lightNo] := TRUE; 
      tmp := vis.rest; vis.rest := tri; tri := vis; 
      vis := tmp
    END;
    WHILE inv # NIL DO
      inv.t.lights[lightNo] := FALSE; 
      tmp := inv.rest; inv.rest := tri; tri := inv; 
      inv := tmp
    END;
  END Shadow;
  
TYPE
  Tetrahedron = RECORD
      plane: ARRAY [0..3] OF H3.Plane;
      point: ARRAY [0..3] OF H3.Point;
      side:  ARRAY [0..3] OF H3.Sign;
    END;

PROCEDURE AuxVisible(
    VAR tri: FaceList;   (* Faces to analyze *)
    opq: FaceList;       (* Opaque triangles *)
    READONLY obs: H3.Point;
    invF: BOOLEAN;
    VAR vis, inv: FaceList;
  ) = 
  (* See SceneShadeF *)
  VAR tvis, tinv: FaceList;
      blk: FaceList;
      antiobs: H3.Point := H3.Point{c := R4.Neg(obs.c)};
      bnum, tnum: CARDINAL;
  <* FATAL Wr.Failure, Thread.Alerted *>
  BEGIN
    AuxComputeEdgePlanes(tri, opq, obs);
    tnum := 0;
    WHILE tri # NIL DO 
      (* Remove first triangle from list "tri", put in "tvis" *)
      tvis := tri; tri := tvis.rest; tvis.rest := NIL;
      tinv := NIL;
      WITH 
        t = tvis.t,
        tvol = BuildTetrahedron(t^, obs)
      DO
        IF Debug THEN 
          Wr.PutText(stderr, "t[");
          Wr.PutText(stderr, Fmt.Pad(Fmt.Int(tnum), 4, '0'));
          Wr.PutText(stderr, "] = ");
          Scene.PrintFace(stderr, t);
          Wr.PutText(stderr, "\n");
        END;
        IF NOT Degenerate(tvol) THEN
          (* look for opaque triangles that may cover "t": *)
          blk := opq;
          bnum := 0;
          WHILE blk # NIL AND tvis # NIL DO
            IF blk.t.org # t.org THEN
              WITH 
                b = blk.t,
                bvol = BuildTetrahedron(b^, antiobs)
              DO
                IF TetraIntersect(bvol, tvol) THEN
                  AuxAuxVisible(tvis, bvol, invF, tvis, tinv);
                  IF Debug THEN 
                    Wr.PutText(stderr, "b[");
                    Wr.PutText(stderr, Fmt.Pad(Fmt.Int(bnum), 4, '0'));
                    Wr.PutText(stderr, "] = ");
                    Scene.PrintFace(stderr, t);
                    Wr.PutText(stderr, "  tvis = ");
                    Wr.PutText(stderr, Fmt.Pad(Fmt.Int(Scene.CountFaces(tvis)), 4));
                    Wr.PutText(stderr, "  tinv = ");
                    Wr.PutText(stderr, Fmt.Pad(Fmt.Int(Scene.CountFaces(tinv)), 4));
                    Wr.PutText(stderr, "\n");
                  END
                END
              END;
            END;
            blk := blk.rest; INC(bnum);
          END;
        END;
        Scene.MoveFaceList(tvis, vis);
        IF invF THEN 
          Scene.MoveFaceList(tinv, inv) 
        ELSE 
          tinv := NIL
        END
      END;
      INC(tnum);
    END
  END AuxVisible;
  
PROCEDURE AuxComputeEdgePlanes(
    tri, opq: FaceList;
    READONLY obs: H3.Point;
  ) =
  VAR sc: Scene.T := Scene.Empty;
  BEGIN
    Scene.CopyFaceList(tri, sc.t);
    Scene.CopyFaceList(opq, sc.t);
    Scene.CloseScene(sc);
    ComputeEdgePlanes(sc.e, obs)
  END AuxComputeEdgePlanes;
  
PROCEDURE ComputeEdgePlanes(
    els: Scene.EdgeList;
    READONLY obs: H3.Point
  ) =
  (*
    Computes the "opl" field of all edges in "els", for the
    given observer "obs". *)
  BEGIN
    (* Unmark all lines: *)
    VAR ep := els; 
    BEGIN
      WHILE ep # NIL DO ep.e.line.mark := FALSE; ep := ep.rest END
    END;
    
    (* Compute planes: *)
    VAR ep := els; 
    BEGIN
      WHILE ep # NIL DO 
        WITH
          e = ep.e,
          ln = e.line^
        DO
          IF NOT ln.mark THEN 
            ln.opl := H3.PlaneFromLineAndPoint(ln.ln, obs);
            ln.mark := TRUE
          END
        END;
        ep := ep.rest
      END
    END
  END ComputeEdgePlanes;

PROCEDURE BuildTetrahedron(
    READONLY base: Scene.Face;
    READONLY apex: H3.Point;
  ): Tetrahedron =
  (*
    Collects the planes and vertices of the 
    tetrahedron with given "base" and "apex".
    
    Vertex "point[i]" is opposite to face "plane[i]"
    and lies on the side "side[i]" of that plane.
    
    Assumes that the "line.opl" fields of all edges
    have been computed for the given "apex" as observer.
  *)
  VAR tet: Tetrahedron;
  BEGIN
    FOR i := 0 TO 3 DO
      WITH 
        pl = tet.plane[i], 
        pt = tet.point[i], 
        s = tet.side[i] 
      DO
        IF i = 0 THEN
          pl := base.plane.pl;
          pt := apex
        ELSE
          WITH
            j = i -1,
            k = (j + 1) MOD 3,
            ej = base.edge[j],
            ek = base.edge[k],
            dk = base.edir[k]
          DO
            pl := ej.line.opl;
            pt := ek.vert[1-dk].wcs
          END
        END;
        s := H3.Side(pt, pl);
      END
    END;
    RETURN tet
  END BuildTetrahedron;

PROCEDURE Degenerate(READONLY vol: Tetrahedron): BOOLEAN =
  BEGIN
    FOR i := 0 TO 3 DO 
      IF vol.side[i] = 0 THEN RETURN TRUE END 
    END;
    RETURN FALSE
  END Degenerate;
  
PROCEDURE ComparePlanes(READONLY a, b: H3.Plane): H3.Sign =
  CONST rtol = 1.0d-6;
  BEGIN
    IF R4.Sin(a.f, b.f) > rtol THEN 
      RETURN 0
    ELSIF R4.Cos(a.f, b.f) > 0.0d0 THEN
      RETURN +1
    ELSE
      RETURN -1
    END
  END ComparePlanes;
  
PROCEDURE TetraIntersect(READONLY va, vb: Tetrahedron): BOOLEAN =
  VAR sep: BOOLEAN;
  <* FATAL Wr.Failure, Thread.Alerted *>
  BEGIN
    (* Check for common faces: *)
    FOR i := 0 TO 3 DO 
      WITH
        apl = va.plane[i], ain = va.side[i],
        bpl = vb.plane[i], bin = vb.side[i],
        rel = ComparePlanes(apl, bpl)
      DO
        IF rel # 0 THEN
          IF Debug THEN
            Wr.PutText(stderr, "    apl = ");
            PrintPlane(stderr, apl);
            Wr.PutText(stderr, "\n");
            Wr.PutText(stderr, "    bpl = ");
            PrintPlane(stderr, bpl);
            Wr.PutText(stderr, "\n");
            Wr.PutText(stderr, "    rel = ");
            Wr.PutText(stderr, Fmt.Int(rel));
            Wr.PutText(stderr, "\n");
          END;
          IF ain * bin # rel THEN RETURN FALSE END
        END
      END
    END;
    
    (* Look for separating plane: *)
    FOR i := 0 TO 3 DO 
      WITH
        apl = va.plane[i], ain = va.side[i],
        bpl = vb.plane[i], bin = vb.side[i]
      DO
        (* Try separating with "apl": *)
        sep := TRUE;
        FOR j := 0 TO 3 DO 
          IF H3.Side(vb.point[j], apl) = ain THEN sep := FALSE END;
        END;
        IF sep THEN RETURN FALSE END;
        (* Try separating with "bpl": *)
        sep := TRUE;
        FOR j := 0 TO 3 DO 
          IF H3.Side(va.point[j], bpl) = bin THEN sep := FALSE END
        END;
        IF sep THEN RETURN FALSE END;
      END
    END;
    RETURN TRUE
  END TetraIntersect;
  
PROCEDURE AuxAuxVisible(
    VAR tri: FaceList;        (* Faces to analyze *)
    READONLY tetra: Tetrahedron;  (* Shadow volume *)
    invF: BOOLEAN;
    VAR vis, inv: FaceList;
   ) =
  (*
    Given a list "tri" of triangles, prepends onto "vis"
    and "inv" the parts of those triangles that are outside or
    inside the region "tetra", respectively.
    
    The invisible parts are returned only if "invF" is TRUE.

    The procedure will reuse the list nodes of "tri".
  *)
  VAR ins, bds, ots: Scene.T := Scene.Empty;
  BEGIN
    ins.t := tri; tri := NIL;
    IF Degenerate(tetra) THEN 
      Scene.MoveFaceList(ins.t, vis)
    ELSE
      Scene.CloseScene(ins);
      Scene.ClipSceneToRegion(
        ins, bds, ots, 
        tetra.plane, tetra.side, 
        invF, invF, TRUE
      );
      IF invF THEN
        Scene.MoveFaceList(bds.t, inv);
        Scene.MoveFaceList(ins.t, inv)
      END;
      Scene.MoveFaceList(ots.t, vis)
    END;
  END AuxAuxVisible;

<*UNUSED*>
PROCEDURE DebugPlotScene( 
    prefix: TEXT;
    call: CARDINAL;
    sc: Scene.T;
    READONLY obs: H3.Point;
  ) =
  BEGIN
    WITH
      file = "debug-" & prefix & "-" & Fmt.Pad(Fmt.Int(call), 5, '0'),
      foc = H3.Point{c := R4.T{1.0, 0.0, 0.0, 0.0}},
      upp = H3.Point{c := R4.T{0.0, 1.0e-8, 1.0e-4, 1.0}},
      w2i = H3.PerspMap(obs, foc, upp)
    DO
      Scene.PlotScene(file, sc, -2.0d0, 2.0d0, -2.0d0, 2.0d0, w2i)
    END
  END DebugPlotScene;

PROCEDURE DebugTetra(READONLY vol: Tetrahedron) =
  <* FATAL Wr.Failure, Thread.Alerted *>
  BEGIN
    FOR i := 0 TO 3 DO
      Wr.PutText(stderr, "  | ");
      Wr.PutText(stderr, "pl = ");
      PrintPlane(stderr, vol.plane[i]);
      Wr.PutText(stderr, "  ");
      Wr.PutText(stderr, "sign = ");
      Wr.PutText(stderr, Fmt.Pad(Fmt.Int(vol.side[i]), 2, '0'));
      Wr.PutText(stderr, "\n");
    END;
    Wr.PutText(stderr, "\n")
  END DebugTetra;

<*UNUSED*>
PROCEDURE PrintPoint(wr: Wr.T; READONLY p: H3.Point) =
  <* FATAL Wr.Failure, Thread.Alerted *>
  BEGIN
    FOR i := 0 TO 3 DO
      IF i = 0 THEN Wr.PutText(wr, "[") ELSE Wr.PutText(wr, ", ") END;
      PrintReal(wr, p.c[i])
    END;
    Wr.PutText(wr, "]");
  END PrintPoint;
  
PROCEDURE PrintPlane(wr: Wr.T; READONLY p: H3.Plane) =
  <* FATAL Wr.Failure, Thread.Alerted *>
  BEGIN
    FOR i := 0 TO 3 DO
      IF i = 0 THEN Wr.PutText(wr, "<") ELSE Wr.PutText(wr, ", ") END;
      PrintReal(wr, p.f[i])
    END;
    Wr.PutText(wr, ">");
  END PrintPlane;
  
PROCEDURE PrintReal(wr: Wr.T; x: REAL) =
  CONST SigDigits = 6;
  <* FATAL Wr.Failure, Thread.Alerted *>
  BEGIN
    Wr.PutText(wr, Fmt.Real(x, prec := SigDigits, style := Fmt.Style.Auto))
  END PrintReal;
  
<*UNUSED*>
PROCEDURE PrintCounts(wr: Wr.T; sc: Scene.T) =
  <* FATAL Wr.Failure, Thread.Alerted *>
  BEGIN
    Wr.PutText(wr, "  { v = ");
    Wr.PutText(wr, Fmt.Int(CountVertices(sc.v)));
    Wr.PutText(wr, "  e = ");
    Wr.PutText(wr, Fmt.Int(CountEdges(sc.e)));
    Wr.PutText(wr, "  t = ");
    Wr.PutText(wr, Fmt.Int(CountFaces(sc.t)));
    Wr.PutText(wr, " }");
  END PrintCounts;

BEGIN 
END SceneShade.