/* See dg_pack.h */
/* Last edited on 2016-04-01 01:11:28 by stolfilocal */

#define _GNU_SOURCE
#include <affirm.h>
#include <bool.h>
#include <ix.h>
#include <stdlib.h>
#include <assert.h>
#include <dg_grid.h>
  
#include <dg_pack.h>

/* CELL PACKS */

uint64_t dg_pack_slot_count(dg_dim_t d, dg_grid_size_t psz[])
  { uint64_t n = 1;
    int i;
    for (i = 0; i < d; i++) { n *= psz[i]; }
    return n;
  }

void dg_pack_cells
  ( dg_dim_t d, 
    dg_cell_id_t idb,
    dg_grid_size_t psz[], 
    dg_cell_id_t id[]
  )
  { dg_grid_pos_t rpos[d]; /* Relative position of cell in pack. */
    /* Indexing steps: */
    ix_step_t pst[d];
    pst[0] = 1;
    { dg_axis_t i; for (i=1; i <d; i++) { pst[i] = pst[i-1]*psz[i-1]; } }
    ix_fill(d, rpos, 0);
    int nC = 0;
    do
      { id[nC] = dg_cell_shift(d, idb, rpos);  nC++; }
    while (! ix_next(d, rpos, psz, ix_order_F, pst, NULL, NULL, NULL, NULL, NULL));
  }

dg_pack_slot_index_t dg_pack_find_cell
  ( dg_dim_t d, 
    dg_cell_id_t id, 
    dg_cell_id_t idb,
    dg_grid_size_t psz[]
  )
  {
    dg_rank_t idr = dg_cell_rank(id);
    dg_rank_t br = dg_cell_rank(idb);
    if (idr != br) { return dg_NOT_SLOT; }

    dg_grid_pos_t idp[dg_dim_MAX];
    dg_cell_position(d, id, idp);
    
    dg_grid_pos_t bp[dg_dim_MAX];
    dg_cell_position(d, idb, bp);
 
    int i;
    dg_pack_slot_index_t pkx = 1;
    dg_pack_slot_index_t ix = 0;
    for (i = 0; i < d; i++)
      { demand(psz[i] > 0, "bad pack size");
        dg_grid_pos_t gszi = dg_max_grid_pos(d, br, (dg_axis_t)i);
        /* Get relative position of {id} in pack along axis {i}: */
        dg_grid_pos_t ei = (idp[i] - bp[i] + gszi) % gszi;
        if (ei >= psz[i]) { return dg_NOT_SLOT; }
        ix += ei * pkx;
        pkx *= psz[i];
      }
    return pkx;
  }

bool_t dg_pack_intersects_cell
  ( dg_dim_t d, 
    dg_cell_id_t id, 
    dg_cell_id_t idb,
    dg_grid_size_t psz[]
  )
  {
    dg_rank_t br = dg_cell_rank(idb);
    dg_grid_pos_t bp[dg_dim_MAX];
    dg_cell_position(d, idb, bp);
    
    /* If {id} is deeper than {idb}, get its ancestor of same rank: */
    dg_rank_t idr = dg_cell_rank(id);
    if (idr > br) { id = id >> (idr - br); idr = br; }
    /* Get posns of min and max descendants of {id} on level {Er}: */
    dg_grid_pos_t idp_lo[dg_dim_MAX];
    dg_cell_position(d, id << (br - idr), idp_lo);
    dg_grid_pos_t idp_hi[dg_dim_MAX];
    dg_cell_position(d, ((id + 1) << (br - idr)) - 1, idp_hi);

    /* Get grid sizes on level {br}: */
    dg_grid_size_t bgsz[dg_dim_MAX];
    dg_grid_size(d, br, bgsz);
    
    int i;
    for (i = 0; i < d; i++)
      { demand(psz[i] > 0, "bad pack size");
        dg_grid_size_t bgszi = bgsz[i];
        /* Get rel position range {plop..phip} of pack in level {br}: */ 
        dg_grid_pos_t plop = bp[i];
        dg_grid_pos_t phip = (plop + psz[i] - 1) % bgszi;
        /* Make sure {phip >= plop}: */
        if (phip < plop) { phip += bgszi; }
        /* Get rel position range {idlop..idhip} of {id} in the same level: */
        dg_grid_pos_t idlop = idp_lo[i];
        dg_grid_pos_t idhip = idp_hi[i];
        assert(idlop <= idhip); /* A single cell can't wrap around. */
        /* Check whether the ranges intersect (with or without wrap-around): */
        if (idlop < plop)
          { if ((idhip < plop) && (phip < idlop + bgszi)) { return FALSE; } }
        else if (plop < idlop)
          { if (idlop > phip) { return FALSE; }}
      }
    return TRUE;
  }

void dg_pack_box_root_relative
  ( dg_dim_t d, 
    dg_cell_id_t idb,
    dg_grid_size_t psz[], 
    interval_t B[]
  )
  { dg_cell_box_root_relative(d, idb, B);
    int i;
    for (i = 0; i < d; i++)
      { demand(psz[i] > 0, "bad pack size");
        HI(B[i]) = HI(B[i]) + ((double)(psz[i]-1))*(HI(B[i]) - LO(B[i]));
      }
  }