/* See {Map.h} */
#include <Map.h>

/* 3D map operations pecialized for triangulations -- maps with tetrahedral cells.*/
/* Last edited on 2024-11-14 05:17:43 by stolfi */

#define Map_C_COPYRIGHT \
  "Copyright © 2000 by the State University of Campinas (UNICAMP)"

#define _GNU_SOURCE

#include <assert.h>
#include <stdint.h>
#include <stdlib.h>
#include <stdio.h>

#include <fget.h>
#include <filefmt.h>
#include <jsfile.h>
#include <jsmath.h>

#include <OctfEnum.h>
#include <OctfData.h>
#include <Octf.h>

ElemTableRec_t MakeElemTable(Place_vec_t root)
  {
    ElemTableRec_t top;
    
    /* Step 1: Renumber all wedges in the submap and make a vector of them: */
    top.wedge = EnumWedges(root, NULL);
    uint32_t NW = top.wedge.ne;
    
    /* Step 2: Collect one place for each existing element record: */
    
    auto void SaveElem(Place_vec_t *tbl, uint32_t num, uint32_t *NP, Place_t p, uint32_t dim);
      /* Assumes that {p} is a place on some element with dimension
        {dim} and serial number {num} among such elements. Stores {p}
        into {tbl->e[num]}, the table of {dim}-dimensional elements,
        updates {*NP} to be the max {num} seen plus one,
        and performs some consistency checks. */
    
    auto void CloseTable(Place_vec_t *tbl, uint32_t N, uint32_t dim);
      /* Trims {. */
    
    /* Collect one place on each distinct element of the submap: */
    uint32_t NV = 0; top.node = Place_vec_new(0);
    uint32_t NE = 0; top.edge = Place_vec_new(0);
    uint32_t NF = 0; top.wall = Place_vec_new(0);
    uint32_t NP = 0; top.cell = Place_vec_new(0);
    int i;
    for (i = 0; i < NW; i++)
      { Wedge_t w = top.wedge.e[i];
        Place_t p = WedgeBase(w);
        assert(w->num == i);
        /* Collect all element slots of the wedge {w}: */
        if(OrgV(p)  != NULL) { SaveElem(&(top.node), OrgV(p) ->num, &NV, p,        0); }
        if(DesV(p)  != NULL) { SaveElem(&(top.node), DesV(p) ->num, &NV, Clock(p), 0); }
        if(PEdge(p) != NULL) { SaveElem(&(top.edge), PEdge(p)->num, &NE, p,        1); }
        if(PWall(p) != NULL) { SaveElem(&(top.wall), PWall(p)->num, &NF, p,        2); }
        if(PnegP(p) != NULL) { SaveElem(&(top.cell), PnegP(p)->num, &NP, p,        3); }        
        if(PposP(p) != NULL) { SaveElem(&(top.cell), PposP(p)->num, &NP, Clock(p), 3); }
      }
    /* Check whether serial numbers were consecutive, for each dimension: */ 
    CloseTable(&(top.node), NV, 0);
    CloseTable(&(top.edge), NE, 1);
    CloseTable(&(top.wall), NF, 2);
    CloseTable(&(top.cell), NP, 3);
    
    /* Step 3: provide new data records for the NULL slots: */
    
    ProvideElemRecords(root, 0, &(top.node));
    ProvideElemRecords(root, 1, &(top.edge));
    ProvideElemRecords(root, 2, &(top.wall));
    ProvideElemRecords(root, 3, &(top.cell));
    
    return top;
    
    void SaveElem(Place_vec_t *tbl, uint32_t num, uint32_t *NP, Place_t p, uint32_t dim)
      { assert(p != NULL);
        Elem_t pe = PElem(p,dim);
        /* Each wedge may have at most two {dim}-elements, so: */
        demand(num < 2*top.wedge.ne, "numbers are not consecutive");
        /* Expand table if needed: */
        Place_vec_expand(tbl, num);
        /* Update {*NP} and clear the new table entries: */
        while ((*NP) <= num) { tbl->e[*NP] = NULL; (*NP)++; }
        /* Check for duplicate numbers: */
        { Place_t q = tbl->e[num];
          if (q != NULL) 
            { Elem_t qe = PElem(q,dim);
              demand(qe == pe, "found two element records with same number");
            }
        }
        tbl->e[num] = p;
        /* Check {p.flp[i].data[dim] == p.data[dim]} for all {i != dim}: */
        int i;
        for (i = 0; i < 4; i++)
          { if (i != dim)
              { Place_t q = Flip(p, i);
                Elem_t qe = PElem(q,dim);
                demand(qe == pe, "found two adjacent places with distinct elem records");
              }
          }
      }
    
    void CloseTable(Place_vec_t *tbl, uint32_t N, uint32_t dim)
      { int i;
        for (i = 0; i < N; i++)
          { demand(tbl->e[i] != NULL, "element numbers are not consecutive"); }
        Place_vec_trim(tbl, N);
      }
  }

#define TP_FILE_TYPE "topology"
#define TP_FILE_VERSION "99-08-25"
  /* File type and version for {WriteTopology} files (Lozada'a 1999 format) */

void WriteTopology(char *name, char *tag, ElemTableRec_t *top, char *cmt /* DF " " */)
    char *fileName = jsprintf("%s%s.tp", name, tag);
    FILE *wr = open_write(fileName, TRUE);
    FWriteTopology(wr, top, cmt);
    fclose(wr);
  }

TopReadData_t ReadTopology(char *name, char *tag)
    char *fileName = jsprintf("%s%s.tp", name, tag);
    FILE *rd = open_read(fileName, TRUE);
    TopReadData_t tc = FReadTopology(rd);
    fclose(rd);
    return tc;
  }

void FWriteTopology(FILE *wr, ElemTableRec_t *top, char *cmt /* DF " " */)
  { /* [!!! Should write the wedges only !!!] */
    /* Get some numbers: */
    uint32_t NE = top->edge.ne;
    uint32_t NF = top->wall.ne;
    uint32_t NV = top->node.ne;
    uint32_t NP = top->cell.ne;
    uint32_t NW = top->wedge.ne;
    
    /* Compute max elem index widths, for nicely aligned output: */
    uint32_t NVP = (uint32_t)umax(NP, NV);
    uint32_t NEF = (uint32_t)umax(NE, NF);
    uint32_t nVEFP = (uint32_t)umax(NVP, NEF);
    
    /* Write the file header line and the client comments: */
    filefmt_write_header(wr, TP_FILE_TYPE, TP_FILE_VERSION);
    filefmt_write_comment(wr, cmt, '|');

    /* Write the element counts and other global fields: */
    uint32_t nWidth = (nVEFP == 0 ? 1 : digits(nVEFP));   /* Max digits in element counts. */
    fprintf(wr, "vertices   %*d\n", nWidth, NV);  
    fprintf(wr, "edges      %*d\n", nWidth, NE);
    fprintf(wr, "faces      %*d\n", nWidth, NF);
    fprintf(wr, "polyhedra  %*d\n", nWidth, NP);
    fprintf(wr, "facetedges %d\n", NW);
    fprintf(wr, "der        %d\n", 0 /* {top->der} */);
    fprintf(wr, "bdr        %d\n", 0 /* {top->bdr} */);
    
    int i;
    
    /* Write the edge and wall places: */
    uint32_t wWidth = (NW <= 1 ? 1 : digits(NW-1)); /* Max digits in wedge indices. */
    uint32_t efWidth = (NEF <= 1 ? 1 : digits(NEF-1)); /* Max digits in edge/wall indices. */
    filefmt_write_comment(wr, "\nEdge data:\n", '|');
    for (i = 0; i < NE; i++)
      { Place_t e = top->edge.e[i];
        fprintf(wr, "%*d ", efWidth, i); 
        PrintPlace(wr, e, wWidth, TRUE);  
      }
    filefmt_write_comment(wr, "\nFace data:\n", '|');
    for (i = 0; i < NF; i++)
      { Place_t f = top->wall.e[i];
        fprintf(wr, "%*d ", efWidth, i); 
        PrintPlace(wr, f, wWidth, TRUE);  
      }

    uint32_t vpWidth = (NVP <= 1 ? 1 : digits(NVP-1)); /* Max digits in node/cell indices. */
    
    auto void PrintOrgV(Place_t p);
    auto void PrintPnegP(Place_t p);
    
    void PrintOrgV(Place_t p)
      { Node_t v = OrgV(p);
        if (v == NULL)
          { fprintf(wr, " %*s ", vpWidth, "-"); }
        else
          { fprintf(wr, " %*dv", vpWidth, v->num); }
      }
          
    void PrintPnegP(Place_t p)
      { Cell_t c = PnegP(p);
        if (c == NULL)
          { fprintf(wr, " %*s ", vpWidth, "-"); }
        else
          { fprintf(wr, " %*dp", vpWidth, c->num); }
      }
    
    /* Write the wedge data: */
    filefmt_write_comment(wr, "\nFacetEdge data:\n", '|');
    for (i = 0; i < NW; i++)
      { Wedge_t w = top->wedge.e[i];
        /* Print the wedge's index (and serial number): */
        assert(w->num == i);
        fprintf(wr, "%*d ", wWidth, i);
        PrintWedge(wr, w, wWidth);
        /* Print the node and cell slots: */
        fprintf(wr, " ");
        Place_t p = WedgeBase(w);
        Place_t q = Clock(p);
        PrintOrgV(p); PrintPnegP(p);
        PrintOrgV(q); PrintPnegP(q);
        /* Print the wall and edge slots: */
        fprintf(wr, " ");
        Wall_t f = w->wall;
        fprintf(wr, " %*df", efWidth, f->num);
        Edge_t e = w->edge;
        fprintf(wr, " %*de", efWidth, e->num);
        fprintf(wr, "\n");
      }

    /* Finish off: */
    filefmt_write_footer(wr, TP_FILE_TYPE);
    fflush(wr);
  }
  
TopReadData_t FReadTopology(FILE *rd)
  { 
    ElemTableRec_t top; 
    char *cmt;
    /* [!!! Must free all discarded {cmt}s !!!] */
     
    /* Read the file header and comments: */
    filefmt_read_header(rd, TP_FILE_TYPE, TP_FILE_VERSION);
    cmt = filefmt_read_comment(rd, '|');

    auto uint32_t ReadCount(char *name);
    
    uint32_t ReadCount(char *name)
     { fget_skip_spaces(rd); 
       fget_match(rd, name);
       uint32_t N = fget_uint32(rd, 10);
       fget_eol(rd);
       return N;
     }
    
    /* Element counts: */
    uint32_t NV = ReadCount("vertices");
    uint32_t NE = ReadCount("edges");
    uint32_t NF = ReadCount("faces");
    uint32_t NP = ReadCount("polyhedra");
    uint32_t NW = ReadCount("facetedges");
    (void) ReadCount("der"); /* Used to be {top->der} */
    (void) ReadCount("bdr"); /* Used to be {top->bdr} */

    int i, j, k;
    
    /* Create the wedge table and the wedge records: */
    top.wedge = Wedge_vec_new(NW);
    for (i = 0; i < NW; i++)
      { Wedge_t w = MakeWedge();
        top.wedge.e[i] = w;
        w->num = i;
      }

    /* Create the element tables: */
    top.node = Place_vec_new(NV);
    top.edge = Place_vec_new(NE);
    top.wall = Place_vec_new(NF);
    top.cell = Place_vec_new(NP);
    
    /* Read edge data: */
    (void)filefmt_read_comment(rd, '|');
    for (j = 0; j < NE; j++)
      { uint32_t en = fget_uint32(rd, 10); /* Edge index. */ 
        demand(en == j, "edge index mismatch");
        Place_t p = ReadPlace(rd, &(top.wedge)); /* A place on the edge. */
        top.edge.e[en] = p; 
        Edge_t e = MakeEdge(en);
        PWedge(p)->edge = e; e->pa = p;
        fget_eol(rd);
      }

    /* Read wall data: */
    (void)filefmt_read_comment(rd, '|');
    for (j = 0; j < NF; j++)
      { uint32_t fn = fget_uint32(rd, 10); /* Wall index. */ 
        demand(fn == j, "wall index mismatch");
        Place_t p = ReadPlace(rd, &(top.wedge)); /* A place on the wall. */
        top.wall.e[fn] = p; 
        Wall_t f = MakeWall(fn);
        PWedge(p)->wall = f; f->pa = p; 
        fget_eol(rd);
      }

    /* Read wedge data: */
    (void)filefmt_read_comment(rd, '|');
    for (j = 0; j < NW; j++)
      { /* Read the wedge index: */
        uint32_t wn = fget_uint32(rd, 10); /* Wedge index. */
        Wedge_t w = top.wedge.e[wn];
        demand(wn == j, "wedge index mismatch");
        /* Read the wedge-to-wedge pointers {w->next[0..3]}: */
        ReadWedge(rd, w, &(top.wedge));
        /* Read the node/cell pointers {w->org[0..3]}: */
        for (k = 0; k < 4; k++)
          { fget_skip_spaces(rd);
            int ch = fgetc(rd);
            if (ch == '-')
              { w->org[k] = NULL; }
            else
              { ungetc(ch, rd);
                uint32_t vpn = fget_uint32(rd, 10); /* Node or cell index. */
                ch = fgetc(rd);
                if (ch == 'v')
                  { Node_t v = MakeNode(vpn);
                    w->org[k] = (NodeOrCell_t)v;
                    top.node.e[vpn] = PlaceInWedge(w, (k << 1));
                  }
                else if (ch == 'p')
                  { Cell_t c = MakeCell(vpn);
                    w->org[k] = (NodeOrCell_t)c;
                    top.cell.e[vpn] = Tors(PlaceInWedge(w, (k << 1)));
                  }
                else
                  { demand(FALSE, "invalid suffix of node/cell number"); }
              }
          }
        /* Read the wall and edge pointers {w->wall}, {w->edge}: */
        uint32_t fn = fget_uint32(rd, 10); /* Wall index. */
        w->wall = PWedge(top.wall.e[fn])->wall;
        demand(fgetc(rd) == 'f', "invalid suffix of wall number");
        uint32_t en = fget_uint32(rd, 10); /* Edge index. */
        w->edge = PWedge(top.edge.e[en])->edge;
        demand(fgetc(rd) == 'e', "invalid suffix of edge number");
        /* Skip the end-of-line character: */
        fget_eol(rd);
      }

    filefmt_read_footer(rd, TP_FILE_TYPE);
    return (TopReadData_t){top, cmt};
  }

void SetExWallsOfEdge(Place_t p, bool_t exists)
  { Place_t q = p;
    do { PWall(q)->exists = exists; q = NextF(q); } while (q != p);
  }

void SetExNode(Place_t p, bool_t exists)
  { OrgV(p)->exists = exists; }

void SetExEdge(Place_t p, bool_t exists)
  { PEdge(p)->exists = exists; }

void SetExWall(Place_t p, bool_t exists)
  { PWall(p)->exists = exists; }

void SetExCell(Place_t p, bool_t exists)
  { PnegP(p)->exists = exists; }

bool_t WallOnBoundary(Place_t p)
  { return (PposP(p) == NULL) || (PnegP(p) == NULL); }

#define Map_C_author \
  "Created 2000 by L. A. P. Lozada, based on 1996 procedures by\n" \
  "  R. M. Rosi and J. Stolfi.\n" \
  "Revisions:\n" \
  "  2000-10-27 : Modified {SetCubeProperties} procedure.\n" \
  "  2007-01-26 : Converted to C by J.Stolfi.\n" \
  "  2007-01-31 : Moved from {Triangulation.c} to here the\n" \
  "               procedure {MakeElemTable} and related ones.\n" \
  "  2007-02-03 Added {ProvideAllElemRecords}."