MODULE Tridimensional;

(*
   This module contains procedures that computing geometric operations
   onto tridimensional mesh. Also, contain procedures for write and read
   a 3D state.

   Also we include some procedures defined by R. L. W. L. and slightly
   modified in her tools "Animacao Dinamica de Corpos Elasticos".

   See the copyright and authorship futher down.

   Last modification: 19-02-200.

*)

IMPORT Wr, Rd, Lex, Triangulation, Mis, FileRd, FloatMode, FileWr, Fmt,
       Thread, Text, OSError, LR3, Math, FileFmt;

FROM Triangulation IMPORT Topology;
FROM FileFmt IMPORT WriteFooter;
FROM Mis IMPORT WriteCommentsJS, WriteLong, WritePoint3D, WriteInt;

CONST
  (* constants for use with the R.L.W.L. procedures *)
  TE_VERSION = "99-08-25";
  ST_VERSION = "99-08-25";
  EN_VERSION = "99-08-25";
  AlphaChars = Mis.AlphaChars;

PROCEDURE Barycenter3D(
  READONLY top: Topology;
  READONLY c3: Coord3D;
  ): LR3.T =
  VAR
    B: LR3.T := LR3.T{0.0d0, ..};
    N: CARDINAL := 0;
  BEGIN
    FOR i := 0 TO LAST(c3) DO 
      WITH
        v = top.vertex[i]
      DO
        IF v.exists THEN B := LR3.Add(B, c3[i]); INC(N) END;
      END
    END;
    RETURN LR3.Scale(1.0d0/FLOAT(N, LONGREAL), B)
  END Barycenter3D;

PROCEDURE MeanVertexDistance3D(
  READONLY top: Topology;
  READONLY c3: Coord3D;
  ): LONGREAL =
  VAR
    S: LONGREAL := 0.0d0;
    N: CARDINAL := 0;
  BEGIN
    FOR i := 0 TO LAST(c3) DO 
      WITH v = top.vertex[i] DO
        IF v.exists THEN 
          S := S + LR3.NormSqr(c3[i]);
          INC(N)
        END
      END
    END;
    RETURN Math.sqrt(S/FLOAT(N,LONGREAL))
  END MeanVertexDistance3D;

PROCEDURE Displace3D(
  READONLY top: Topology;
  d: LR3.T;
  VAR c3: Coord3D;
  ) =
  BEGIN
    FOR i := 0 TO LAST(c3) DO 
      IF top.vertex[i].exists THEN
        WITH vc = c3[i] DO vc := LR3.Add(vc, d) END
      END
    END
  END Displace3D;

PROCEDURE Scale3D(
  READONLY top: Topology;
  s: LONGREAL;
  VAR c3: Coord3D;
  ) =
  BEGIN
    FOR i := 0 TO LAST(c3) DO 
      IF top.vertex[i].exists THEN
        WITH vc = c3[i] DO vc := LR3.Scale(s,vc) END;
      END
    END
  END Scale3D;

PROCEDURE NormalizeVertexDistance3D(
  READONLY top: Topology;
  VAR c3: Coord3D;
  ) =
  BEGIN
    WITH
      b = Barycenter3D(top, c3)
    DO
      Displace3D(top, LR3.Neg(b), c3)
    END;
    WITH
      s = MeanVertexDistance3D(top, c3)
    DO
      Scale3D(top, 1.0d0/s, c3)
    END;
  END NormalizeVertexDistance3D;

PROCEDURE WriteState3D(
    name: TEXT;
    READONLY top: Topology;
    READONLY c: Coord3D;
    comments: TEXT := " ";
    ) =
  <* FATAL Wr.Failure, Thread.Alerted, OSError.E *>
  BEGIN
    WITH 
      st3 = FileWr.Open(name & ".st3"),
      vWidth = Mis.NumDigits(top.NV - 1)
    DO    
      FileFmt.WriteHeader(st3,"state3D","99-08-25");
      Wr.PutText(st3, "vertices ");
      Wr.PutText(st3, Fmt.Int(top.NV) & "\n");
      IF NOT Text.Empty(comments) THEN
        Mis.WriteCommentsJS(st3, comments & "\n", '|')
      END;
      Mis.WriteCommentsJS(st3, "\nVertex data:\n", '|');
      FOR i := 0 TO top.NV-1 DO
        WITH v = top.vertex[i] DO
          (* state 3D *)
          Wr.PutText(st3, Fmt.Pad(Fmt.Int(v.num), vWidth)); 
          Wr.PutText(st3,"  ");
          WritePoint3D(st3,c[v.num]);
          Wr.PutText(st3,"\n");
        END       
      END;
      FileFmt.WriteFooter(st3,"state3D");
      Wr.Close(st3);
    END;     
  END WriteState3D;

PROCEDURE ReadState3D(name: TEXT): REF Coord3D =
  <* FATAL Rd.Failure,Thread.Alerted,FloatMode.Trap,Lex.Error,OSError.E *>
  VAR
    c: REF Coord3D;
    comments: TEXT;
    nv: CARDINAL;
  BEGIN
    WITH 
      rs = FileRd.Open(name & ".st3")
    DO
      (* Read Headers File Formats*)
      FileFmt.ReadHeader(rs,"state3D","99-08-25");

      (* Element counts: *)
      Lex.Skip(rs, cs := AlphaChars);
      nv := Lex.Int(rs);
      Lex.Skip(rs);
      comments := Mis.ReadCommentsJS(rs, '|');
      
      c := NEW(REF Coord3D, nv);
      (* Read vertex data state: *)  
      FOR j := 0 TO nv-1 DO 
        Lex.Skip(rs);
        WITH 
          nv = Lex.Int(rs) 
        DO 
          WITH cv = c[nv] DO
            cv[0] := Lex.LongReal(rs);
            Lex.Skip(rs);
            cv[1] := Lex.LongReal(rs);
            Lex.Skip(rs);
            cv[2] := Lex.LongReal(rs);
          END
        END 
      END;

      FileFmt.ReadFooter(rs,"state3D");
      Rd.Close(rs);        
      RETURN c;
    END;  
  END ReadState3D;


(* procedures associated to the R.L.W.L.'s tools  *)

PROCEDURE WriteSpace(wr: Wr.T; n: CARDINAL) =
  <* FATAL Wr.Failure, Thread.Alerted *>
  BEGIN
    FOR i := 1 TO n DO
      Wr.PutChar(wr,' ')
    END
  END WriteSpace;

PROCEDURE WriteEOL(wr: Wr.T) =
  <* FATAL Wr.Failure, Thread.Alerted *>
  BEGIN
    Wr.PutChar(wr,'\n')
  END WriteEOL;

PROCEDURE WriteStates(
    wr: Wr.T;
    time: LONGREAL;
    READONLY vertex: ARRAY OF Vertex;
    prec: CARDINAL := 4;
    comments: TEXT := "";
   ) =
  
  PROCEDURE WriteCoord3(x: LONGREAL) =
    <* FATAL Wr.Failure, Thread.Alerted *>    
    BEGIN
      Wr.PutText(wr, Fmt.LongReal(x, Fmt.Style.Sci, prec := prec))
    END WriteCoord3;

  PROCEDURE WritePoint3(READONLY p: LR3.T) =   
    <* FATAL Wr.Failure, Thread.Alerted *>    
    BEGIN
      WriteCoord3(p[0]);
      Wr.PutText(wr, " ");
      WriteCoord3(p[1]);
      Wr.PutText(wr, " ");
      WriteCoord3(p[2]);
      Wr.PutText(wr, " ");
    END WritePoint3;

  <* FATAL Wr.Failure, Thread.Alerted *>    
  BEGIN
    WITH
      NV = NUMBER(vertex),
      vWidth = Mis.NumDigits(NV-1)
    DO
      FileFmt.WriteHeader(wr, "state3D", "99-08-25");
      Wr.PutText(wr, "vertices " & Fmt.Int(NUMBER(vertex)) & "\n");
      WriteCommentsJS(wr, comments & "\n", '|');
      Mis.WriteCommentsJS(wr, "\nVertex data:\n", '|');
      Wr.PutText(wr,"t " &
                 Fmt.Pad(Fmt.LongReal(time,Fmt.Style.Fix,prec:=2),4)&"\n");
      FOR i := 0 TO LAST(vertex) DO
        Wr.PutText(wr, Fmt.Pad(Fmt.Int(i), vWidth)); 
        Wr.PutText(wr, " ");
        WritePoint3(vertex[i]); WriteEOL(wr);
      END;
      WriteFooter(wr, "state3D");
      Wr.Close(wr)
    END
  END WriteStates;


(* Notice: actually WriteStates = WriteState3D

*)


PROCEDURE WriteTetrahedra(
    name: TEXT;
    READONLY top: Topology;   
    READONLY cell: ARRAY OF Tetrahedron;
    comments: TEXT := " ";
   ) =
  <* FATAL Wr.Failure, Thread.Alerted, OSError.E *>
  BEGIN
    WITH
      wr = FileWr.Open(name & ".te"),
      pWidth = Mis.NumDigits(top.NP - 1)
    DO
      FileFmt.WriteHeader(wr, "tetrahedra","99-08-25");
      Wr.PutText(wr, "tetrahedra " & Fmt.Int(NUMBER(cell)) & "\n");
      WriteCommentsJS(wr, comments & "\n", '|');
      FOR i := 0 TO LAST(cell) DO 
        Wr.PutText(wr, Fmt.Pad(Fmt.Int(i), pWidth)); 
        Wr.PutText(wr," ");
        WriteTetrahedron(wr, cell[i], 0);
      END;
      WriteFooter(wr, "tetrahedra");
      Wr.Close(wr)
    END
  END WriteTetrahedra; 

PROCEDURE WriteTetrahedron
   (wr: Wr.T;
    READONLY t: Tetrahedron; 
    nbase: CARDINAL;
   ) =
  BEGIN
    WriteSpace(wr,0);
    WriteLong(wr, t.A[0,0]); WriteSpace(wr,0);
    WriteLong(wr, t.A[0,1]); WriteSpace(wr,0);
    WriteLong(wr, t.A[0,2]); WriteSpace(wr,1);
    WriteLong(wr, t.A[1,0]); WriteSpace(wr,0);
    WriteLong(wr, t.A[1,1]); WriteSpace(wr,0);
    WriteLong(wr, t.A[1,2]); WriteSpace(wr,1);
    WriteLong(wr, t.A[2,0]); WriteSpace(wr,0);
    WriteLong(wr, t.A[2,1]); WriteSpace(wr,0);
    WriteLong(wr, t.A[2,2]); WriteEOL(wr);
    
    WriteSpace(wr, 3);
    WriteInt(wr, nbase + t.p0); WriteSpace(wr,0);
    WriteInt(wr, nbase + t.p1); WriteSpace(wr,0);
    WriteInt(wr, nbase + t.p2); WriteSpace(wr,0);
    WriteInt(wr, nbase + t.p3); WriteEOL(wr);
    
    WriteSpace(wr, 2);
    WriteLong(wr, t.density); WriteSpace(wr,1);
    WriteLong(wr, t.alpha);   WriteSpace(wr,1);
    WriteLong(wr, t.beta);    WriteSpace(wr,1);
    WriteEOL(wr);
  END WriteTetrahedron;

PROCEDURE ReadTetrahedron(rd: Rd.T; VAR t: Tetrahedron) =
  <* FATAL Rd.Failure, Thread.Alerted, FloatMode.Trap, Lex.Error *>
  BEGIN
    t.A[0,0] := Lex.LongReal(rd); Lex.Skip(rd);
    t.A[1,0] := Lex.LongReal(rd); Lex.Skip(rd);
    t.A[2,0] := Lex.LongReal(rd); Lex.Skip(rd);
    t.A[0,1] := Lex.LongReal(rd); Lex.Skip(rd);
    t.A[1,1] := Lex.LongReal(rd); Lex.Skip(rd);
    t.A[2,1] := Lex.LongReal(rd); Lex.Skip(rd);
    t.A[0,2] := Lex.LongReal(rd); Lex.Skip(rd);
    t.A[1,2] := Lex.LongReal(rd); Lex.Skip(rd);
    t.A[2,2] := Lex.LongReal(rd); Lex.Skip(rd);
  
    t.p0 := Lex.Int(rd); Lex.Skip(rd);
    t.p1 := Lex.Int(rd); Lex.Skip(rd);
    t.p2 := Lex.Int(rd); Lex.Skip(rd);
    t.p3 := Lex.Int(rd); Lex.Skip(rd);

    t.density := Lex.LongReal(rd); Lex.Skip(rd);
    t.alpha   := Lex.LongReal(rd); Lex.Skip(rd);
    t.beta    := Lex.LongReal(rd); Lex.Skip(rd);
  END ReadTetrahedron;

PROCEDURE ReadVectors(rd: Rd.T; n: CARDINAL; pos: Vectors3D) =
 <* FATAL Rd.EndOfFile, Rd.Failure,Thread.Alerted, FloatMode.Trap, Lex.Error *>
  BEGIN
    FOR i := 0 TO n-1 DO
      EVAL Lex.Int(rd); EVAL Rd.GetChar(rd);
      pos[i][0] := Lex.LongReal(rd); Lex.Skip(rd);
      pos[i][1] := Lex.LongReal(rd); Lex.Skip(rd);
      pos[i][2] := Lex.LongReal(rd); Lex.Skip(rd);
    END
  END ReadVectors;

PROCEDURE ReadHeader(rd: Rd.T; type, param: TEXT): CARDINAL
  RAISES {Rd.EndOfFile} =
  <* FATAL Rd.Failure, Thread.Alerted, FloatMode.Trap, Lex.Error *>
  VAR n := 0;
  BEGIN
    TRY Lex.Match(rd,"begin " & type & " (format of " & GetVersion(type) & ")")
    EXCEPT ELSE RAISE Rd.EndOfFile END;
    Lex.Skip(rd);
    IF param # NIL THEN
      Lex.Match(rd, param);
      n := Lex.Int(rd); Lex.Skip(rd);
    END;
    EVAL Mis.ReadCommentsJS(rd,'|');
    RETURN n;
  END ReadHeader;

PROCEDURE ReadFooter(rd: Rd.T; type: TEXT) =
  <* FATAL Rd.Failure, Thread.Alerted, Lex.Error *>
  BEGIN
    Lex.Skip(rd);
    Lex.Match(rd, "end " & type)
  END ReadFooter;

PROCEDURE GetVersion(type: TEXT): TEXT =
  BEGIN
    IF  Text.Equal(type, "tetrahedra")   THEN RETURN TE_VERSION
    ELSIF Text.Equal(type, "state3D")    THEN RETURN ST_VERSION
    ELSIF Text.Equal(type, "energies")   THEN RETURN EN_VERSION
    ELSE <* ASSERT FALSE *>
    END
  END GetVersion;

PROCEDURE WriteHeader(wr: Wr.T; type, param: TEXT; n: CARDINAL) =
  <* FATAL Wr.Failure, Thread.Alerted *>
  BEGIN
    Wr.PutText(wr,"begin " & type & " (format of " & GetVersion(type) & ")\n");
    IF param # NIL THEN
      Wr.PutText(wr, param);
      WriteInt(wr, n);
      WriteEOL(wr);
    END
  END WriteHeader;

BEGIN
END Tridimensional.

(**************************************************************************)
(*                                                                        *)
(* Copyright (C) 2000     Universidade Estadual de Campinas (UNICAMP)     *)
(*                                                                        *)
(* Authors:                                                               *)
(*   L. P. Lozada  & J. Stolfi   - UNICAMP <lozada@dcc.unicamp.br>        *)
(*                                                                        *)
(* This file can be freely used, distributed, and modified, provided that *)
(* this copyright and authorship notice is included in every copy or      *)
(* derived version.                                                       *)
(*                                                                        *)
(* DISCLAIMER: This software is offered ``as is'', without any guarantee  *)
(*   as to fitness for any particular purpose.  Neither the copyright     *)
(*   holder nor the authors or their employers can be held responsible    *)
(*   for any damages that may result from its use.                        *)
(*                                                                        *)
(**************************************************************************)