/* See {salamic_slicer_Minetto.h} */
/* Last edited on 2016-04-21 18:47:23 by stolfilocal */

#define _GNU_SOURCE
#include <stdio.h>
#include <stdint.h>
#include <string.h>
#include <assert.h>
#include <stdlib.h>

#include <bool.h>
#include <vec.h>
#include <affirm.h>
#include <binsearch_int32.h>
#include <group_sort_uint32.h>

#include <stmesh.h>
/* #include <stmesh_utils.h> */

#include <salamic_slicer.h>
#include <salamic_closer.h>
#include <salamic_stats.h>

#include <salamic_slicer_Minetto.h>

/* INTERNAL PROTOTYPES */

void salamic_slicer_Minetto_group_faces
  ( stmesh_t mesh,
    bool_t preSorted,
    int_vec_t *planeZ, 
    bool_t uniform,
    stmesh_face_t f[],   /* (OUT) Faces, sorted by group. */
    uint32_t gstart[]    /* (OUT) beg and end of each group. */
  );   
  /* Fills {f[0..nf-1]} with all the faces of {mesh},
    sorted by increasing slicing group.

    There are {ng = np+1} slicing groups, numbered {0..ng-1} or {0..np}.
    For {ig} in {1..np-1}, group {ig} has every face {t} that start
    between planes {ig-1} and {ig}, that is, such that
    {planeZ.e[ig-1] < t.minZ < planeZ.e[ig]}. Group 0 has all faces
    that start below the first plane ({t.minZ < planeZ.e[0]}) and
    group {np} has all faces that start after the last plane
    ({planeZ.e[np-1] < t.minZ}).
    
    Also sets {gstart[0..ng-1]} to indicate the beginning and end of
    each group in {f}. Namely, for any group {ig} in {0..ng-1}, the
    faces in group {ig} are {f[k0..k1-1]}, where
    {k0=gstart[ig]} and {k1=gstart[ig+1]}. Note that {gstart} must
    have {ng+1 = np+2} elements.
    
    Assumes (requires) that the {Z} coordinates of planes are distinct from 
    the {Z}-coordinates of all mesh vertices; so that 
    no vertex lies on a slicing plane. 
    
    The flags {preSorted} and {uniform} can be set to save some
    processing time in special cases. If {preSorted} is true, the
    procedure assumes (requires) that the faces of {mesh} are already
    sorted by increasing {.minZ}. If {uniform} is true, it assumes that
    the slicing planes (if there is more than one) are uniformly spaced. */

void salamic_slicer_Minetto_group_faces_presorted
  ( stmesh_t mesh,
    int_vec_t *planeZ, 
    stmesh_face_t f[],   /* (OUT) Face indices, sorted by group. */
    uint32_t gstart[]    /* (OUT) Beg and end of each group. */
  );   
  /* Same as {salamic_slicer_Minetto_group_faces}, for the special case that faces in {mesh} are already
    sorted by increasing {.minZ}.  The planes may be arbitrarily spaced. */

void salamic_slicer_Minetto_group_faces_unsorted
  ( stmesh_t mesh,
    int_vec_t *planeZ, 
    bool_t uniform,
    stmesh_face_t f[],   /* (OUT) Face indices, sorted by group. */
    uint32_t gstart[]    /* (OUT) Beg and end of each group. */
  );   
  /* Same as {salamic_slicer_Minetto_group_faces}, for the general case 
    of unsorted faces.  If {uniform} is true, assumes (and requires) that 
    the slicing planes be uniformly spaced. */

void salamic_slicer_Minetto_face_group_sort
  ( stmesh_t mesh,
    uint32_t ng,        /* Number of groups. */
    uint32_t fgroup[],  /* Group assigned to each face. */
    stmesh_face_t f[],  /* (OUT) Faces sorted by group. */
    uint32_t gstart[]   /* (OUT) Beg and end of each group. */
  );
  /* Sorts the {nf} faces of {mesh} into {ng} groups. 
    
    Assumes that the groups are identified by /group indices/ in {0..ng-1},
    and that each face of {mesh} with index {uxf} belongs to the group with index {fgroup[uxf]}.
    It stores the faces into {f[0..nf-1]}, sorted by increasing group index. 
    
    The procedure also sets {gstart[0..ng]} to indicate the beginning and end of
    each group in {f}. Namely, for any group index {ig} in {0..ng-1}, the
    faces in group {ig} are {f[k0..k1-1]}, where
    {k0=gstart[ig]} and {k1=gstart[ig+1]}. Note that {gstart} must
    have {ng+1} elements, not {ng}. */
    

/* IMPLEMENTATIONS */

void salamic_slicer_Minetto_slice
  ( stmesh_t mesh,
    bool_t preSorted,
    int_vec_t *planeZ, 
    bool_t uniform,
    char *closer, 
    char *outPrefix, 
    salamic_stats_t *st
  )
  {
    uint32_t nf = stmesh_face_count(mesh);   /* Number of triangles. */
    uint32_t np = planeZ->ne; /* Number of slicing planes. */
    uint32_t ng = np + 1;     /* Number of face groups. */
    float eps = stmesh_get_eps(mesh);
    bool_t verbose = TRUE;
    
    /* Sort the triangles by slicing group: */
    stmesh_face_t *f = notnull(malloc(nf*sizeof(stmesh_face_t)), "no mem"); /* All faces, group-sorted */
    uint32_t *gstart = notnull(malloc((ng+1)*sizeof(uint32_t)), "no mem");  /* Start of each group in {f}. */
    salamic_slicer_Minetto_group_faces(mesh, preSorted, planeZ, uniform, f, gstart);    

    /* Sweep over the mesh with a horizontal plane: */
    stmesh_face_t *a = notnull(malloc(nf*sizeof(stmesh_face_t)), "no mem");  /* Active faces, group-sorted */
    uint32_t na = 0; /* Number of active faces. */
    uint32_t ixp;
    for (ixp = 0; ixp < np; ixp++)
      { /* Get the quantized Z-coordinate of next slicing plane: */
        int32_t pZ = planeZ->e[ixp];
        /* Get the corresponding group index {ig}: */
        uint32_t ig = ixp;
        /* Append to the active list all faces in group {ig}: */
        uint32_t na_plus = gstart[ig+1] - gstart[ig]; /* Number of faces in group. */
        uint32_t kf;
        for (kf = gstart[ig]; kf < gstart[ig+1]; kf++)
          { assert(na < nf);
            a[na] = f[kf];
            na++; 
          }
        /* Remove from the active list all triangles that do not cross the plane: */
        uint32_t na_new = 0; /* Faces that remain active: */
        uint32_t ja;
        for (ja = 0; ja < na; ja++)
          { /* Get the next active face {f}: */
            stmesh_face_t aj = a[ja];
            int32_t minZ, maxZ;
            stmesh_face_get_zrange(aj, &minZ, &maxZ);
            assert(minZ < pZ);
            if (maxZ > pZ)
              { /* Still active, keep it: */
                a[na_new] = aj;
                na_new++; 
              }
            else
              { /* Became inactive, skip it: */
                assert(maxZ < pZ); /* No vertices on slicing planes. */
              }
          }
        uint32_t na_minus = na_new - na; /* Number of inactivated faces. */
        na = na_new;
        if (verbose) 
          { fprintf(stderr, "  plane %4u at Z = %+11d (%+.8f mm): %u faces", ixp, pZ, pZ*(double)eps, na_new);
            fprintf(stderr, " (added %d removed %d)\n", na_plus, na_minus);
          }
        if (strcmp(closer, "NONE") != 0)
          { salamic_closer_close(mesh, pZ, na, a, closer, outPrefix, st); }
      }
      
    /* Free storage:*/
    free(f);
    free(gstart);
    free(a);
  }

void salamic_slicer_Minetto_group_faces
  ( stmesh_t mesh,
    bool_t preSorted,
    int_vec_t *planeZ, 
    bool_t uniform,
    stmesh_face_t f[],  /* (OUT) Face indices, sorted by group. */
    uint32_t gstart[]   /* (OUT) beg and end of each group. */
  )   
  {
    if (preSorted)
      { salamic_slicer_Minetto_group_faces_presorted(mesh, planeZ, f, gstart); }
    else 
      { salamic_slicer_Minetto_group_faces_unsorted(mesh, planeZ, uniform, f, gstart); }
  }
    
void salamic_slicer_Minetto_group_faces_presorted
  ( stmesh_t mesh,
    int_vec_t *planeZ, 
    stmesh_face_t f[],  /* (OUT) Face indices, sorted by group. */
    uint32_t gstart[]   /* (OUT) Beg and end of each group. */
  )   
  {
    uint32_t nf = stmesh_face_count(mesh);   /* Number of triangles. */
    uint32_t np = planeZ->ne; /* Number of planes. */
    uint32_t ng = np + 1;     /* Number of groups. */
    
    /* Scan the faces in index order, noting the group breaks. */
    assert(ng < INT32_MAX);   /* Paranoia - so that we can use {int32_t} for {ig}. */
    int32_t ig = -1;          /* Index of current group. */
    int32_t pZ = INT32_MIN;   /* Quantized {Z} of plane {ig} */
    int32_t minZ_next = INT32_MIN; /* Mininmum {.minZ} for next face. */
    uint32_t uxf;
    for (uxf = 0; uxf < nf; uxf++)
      {  /* Save the face the {f} list: */
        f[uxf] = stmesh_get_face(mesh, uxf);
        /* Get the face record and its min{Z} coordinate {fZ}: */
        int32_t fminZ, fmaxZ;
        stmesh_face_get_zrange(f[uxf], &fminZ, &fmaxZ);
        demand(fminZ >= minZ_next, "mesh faces are not pre-sorted by {.minZ}");
        demand(fminZ < INT32_MAX, "excessive vertex {Z} coord");

        /* Update {ig} so that {fminZ} is less than group's top coord {pZ}; set {gstart}: */
        while (fminZ > pZ) 
          { /* Get to next group: */
            ig++;
            assert ((ig >= 0) && (ig < ng));
            /* Group {ig} starts here: */
            gstart[ig] = uxf;
            /* Get the group's max Z coordinate in {pZ}: */
            if (ig < np)
              { int32_t pZ_new = planeZ->e[ig];
                demand(pZ_new >= pZ + 2, "plane {Z} coords not properly spaced");
                demand(pZ_new < INT32_MAX, "excessive plane {Z} coord");
                pZ = pZ_new;
              }
            else
              { pZ = INT32_MAX; }
          }
        demand(fminZ < pZ, "bug: vertex lies on slicing plane");
        minZ_next = fminZ;
      }
    /* Complete {gstart} if needed: */
    while (ig < ng) { ig++; gstart[ig] = nf; } 
  }

void salamic_slicer_Minetto_group_faces_unsorted
  ( stmesh_t mesh,
    int_vec_t *planeZ, 
    bool_t uniform,
    stmesh_face_t f[], /* OUT: Face indices, sorted by group. */
    uint32_t gstart[]  /* OUT: Beg and end of each group. */
  )
  {
    uint32_t nf = stmesh_face_count(mesh);   /* Number of triangles. */
    uint32_t np = planeZ->ne; /* Number of planes. */
    uint32_t ng = np + 1;     /* Number of groups. */
    bool_t verbose = TRUE;
    
    /* Get quantized {Z} coordinates {bZ,eZ} of first and last planes: */
    int32_t bZ = planeZ->e[0];
    int32_t eZ = planeZ->e[np-1];
    if (verbose) { fprintf(stderr, "grouping unsorted faces for %u slicing planes in %+d..%+d\n", np, bZ, eZ); }
    demand(eZ >= bZ, "invalid plane {Z} coordinates");

    /* Get quantized step {dZ} of {Z} coordinate between planes, or 0 if not uniform: */
    int32_t dZ;
    if ((uniform) && (np >= 2))
      { dZ = planeZ->e[1] - bZ;
        demand((((eZ - bZ) % dZ) == 0), "plane spacing is not uniform");
        demand((((eZ - bZ)/dZ) == np - 1), "inconsistent plane spacing");
      }
    else
      { dZ = 0; }
    
    /* Assigns group index {fgroup[uxf]} to each face with index {uxf}: */
    uint32_t *fgroup = notnull(malloc(nf*sizeof(uint32_t)), "no mem"); /* Group of each face. */
    uint32_t uxf;
    for (uxf = 0; uxf < nf; uxf++)
      { /* Get the face record and its {.minZ} coordinate {fZ}: */
        stmesh_face_t fk = stmesh_get_face(mesh, uxf);
        int32_t fminZ, fmaxZ;
        stmesh_face_get_zrange(fk, &fminZ, &fmaxZ);
        demand((fminZ & 1) == 0, "mesh vertices are not even");

        /* Compute the group index {ig}: */
        uint32_t ig;
        if (fminZ < bZ)
          { ig = 0; }
        else if (fminZ > eZ)
          { ig = np; }
        else
          { if (uniform)
              { ig = (fminZ - bZ)/dZ + 1; }
            else
              { ig = (uint32_t)binsearch_int32(fminZ, np, planeZ->e); }
            demand((ig >= 1) && (fminZ > planeZ->e[ig-1]), "bug: vertex on lower slicing plane");
            demand((ig < np) && (fminZ < planeZ->e[ig]), "bug: vertex on upper slicing plane");
          }
        assert(ig < ng);
        
        /* Save and count it: */
        fgroup[uxf] = ig;
      }
      
    /* Put in {f[0..nf-1]} the face indices, sorted by group, and set {gstart[0..ng]}: */
    salamic_slicer_Minetto_face_group_sort(mesh, ng, fgroup, f, gstart);

    free(fgroup);
  }
  
void salamic_slicer_Minetto_face_group_sort
  ( stmesh_t mesh,
    uint32_t ng,        /* Number of groups. */
    uint32_t fgroup[],  /* Group assigned to each face. */
    stmesh_face_t f[],  /* (OUT) Faces sorted by group. */
    uint32_t gstart[]   /* (OUT) Beg and end of each group. */
  )
  {
    uint32_t nf = stmesh_face_count(mesh);   /* Number of triangles. */
    
    /* Sort the face indices by group: */
    stmesh_face_unx_t *uxf = notnull(malloc(nf*sizeof(stmesh_face_unx_t)), "no mem");
    group_sort_uint32(nf, ng, fgroup, uxf, gstart);
    
    /* Pick the faces in that order: */
    uint32_t k;
    for (k = 0; k < nf; k++) { f[k] = stmesh_get_face(mesh, uxf[k]); }

    /* Paranoia check: */
    assert(gstart[ng] == nf);
    
    free(uxf);
  }