/* See pz_image.h */
/* Last edited on 2017-06-20 21:04:19 by stolfilocal */

#include <pz_image.h>

#include <pz_types.h>
#include <jspnm.h>
#include <uint16_image.h>
#include <jsstring.h>
#include <affirm.h>

#include <stdio.h>

/* INTERNAL PROTOTYPES */

/* For these procedures, a pixel {(x,y)} is ``just inside the object'' if
    {img[y,x] >maxbg}, and {img[y,x-1] <= maxbg}. */

void pz_image_find_pixel_inside
  ( pz_image_t *img, 
    int *xp, int *yp, 
    uint16_t maxbg
  );
  /*
    Finds a pixel {(x,y)} of {img} that lies just inside the object.
    If there is none, returns {x=img.cols, y=img.rows}. */

r3_vec_t pz_image_follow_boundary(
    pz_image_t *img, 
    int x,
    int y,
    uint16_t maxbg,
    pz_side_t side 
  );
  /*
    Given an image {img}, and a pixel {(x,y)} that is just inside the
    object, this procedure walks counterclockwise along the object's
    boundary. The points of the boundary.  The points of the selected
    curve (inner, outer, or interpolated) are returned as a {r3_vec_t}. */

/* IMPLEMENTATIONS */

pz_image_t *pz_image_read_pgm (char *fName)
  { char *fn = txtcat(fName, ".pgm");
    pz_image_t *img = uint16_image_read_pnm_named(fn, TRUE);
    free(fn);
    return img;
  }

void pz_image_write_pgm(char *fName, pz_image_t *img)
  { char *fn = txtcat(fName, ".pgm");
    uint16_image_write_pnm_named(fn, img, FALSE, TRUE);
    free(fn);
  }

r3_vec_t pz_image_extract_boundary
  ( pz_image_t *img, 
    uint16_t maxbg,
    pz_side_t side 
  )
  { int x = -1, y = -1;
    pz_image_find_pixel_inside(img, &x, &y, maxbg);
    fprintf(stderr, "found pixel inside = (%d,%d)\n", x, y);
    if ((x >= img->cols) || (y >= img->rows))
      { return r3_vec_new(0); }
    else
      { return pz_image_follow_boundary(img, x, y, maxbg, side); }
  }

void pz_image_find_pixel_inside
  ( pz_image_t *img, 
    int *xp, int *yp, 
    uint16_t maxbg
  )
  { 
    int x, y = 0; 
    for (y = 0; y < img->rows; y++)
      { for (x = 0; x < img->cols; x++)
          { if (img->smp[y][x] > maxbg) { *xp = x; *yp = y; return; } }
      }
    *xp = img->cols; *yp = img->rows;
  }

r3_vec_t pz_image_follow_boundary(
    pz_image_t *img, 
    int x,
    int y,
    uint16_t maxbg,
    pz_side_t  side 
  )
  {
    int xi = x,   yi = y;  /* The inner foot. */
    int xo = x-1, yo = y;  /* The outer foot. */
    
    r3_vec_t f = r3_vec_new(2000);
    int n = 0;
    
#   define Inside(X,Y) (img->smp[Y][X] > maxbg)
#   define Outside(X,Y) (img->smp[Y][X] <= maxbg)

    do 
      { 
        /* Store the current outline point: */
        double xm = 0, ym = 0;
        switch(side)
          { 
            case pz_side_interpolated:
              { /* Find {(xm,ym}) between {(xo,yo)} and {(xi,yi)} */
                /* where {img[xm,ym] = grey}, by linear interpolation: */
                double vo = (double)(maxbg + 0.5 - img->smp[yo][xo]);
                double dv = (double)(img->smp[yi][xi] - img->smp[yo][xo]);
                double alpha = vo/dv;
                /* Just to be sure: */
                if (alpha > 1) { alpha = 1; }
                if (alpha < 0) { alpha = 0; }
                xm = (1 - alpha)*xo + alpha*xi;
                ym = (1 - alpha)*yo + alpha*yi;
              }
              break;
            case pz_side_inner:
              xm = (double)xi; ym = (double)yi; 
            break;
            case pz_side_outer:
              xm = (double)xo; ym = (double)yo; 
            break;
          }
        fprintf(stderr, " point = (%10.3f,%10.3f)\n", xm, ym);
        r3_vec_expand(&f, n);
        f.el[n] = (r3_t){{xm, ym, 0}};
        n++;   
        
        /* Take one step counterclockwise along outline: */
        if (yo < yi)
          { if (Inside(xi+1, yi))
              { if (Outside(xo+1, yo))
                  { xi = xi+1; xo = xo+1; }
                else
                  { xi = xi+1; yi = yi-1; }
              }
            else
              { xo = xo+1; yo = yo+1; }
          }
        else if (yo > yi)
          { if (Inside(xi-1, yi))
              { if (Outside(xo-1, yo))
                  { xi = xi-1; xo = xo-1; }
                else
                  { xi = xi-1; yi = yi+1; }
              }
            else
              { xo = xo-1; yo = yo-1; }                
          }
        else if (xo >xi)
          { if (Inside(xi, yi+1))
              { if (Outside(xo, yo+1))
                  { yi = yi+1; yo = yo+1; }
                else
                  { xi = xi+1; yi= yi+1; }
              }
            else
              { xo = xo-1; yo = yo+1; } 
          }
        else if (xo < xi)
          { if (Inside(xi, yi-1))
              { if (Outside(xo, yo-1))
                  { yi = yi-1; yo = yo-1; }
                else
                  { xi = xi-1; yi = yi-1; }
              }
            else
              { xo = xo+1; yo = yo-1; }
          }
        else
         { affirm(FALSE, "crossed legs"); }
      }
    while ((xi != x) || (yi != y) || (xo != x-1) || (yo != y));
  
    r3_vec_trim(&f, n);
    return f;
  }

/*
  Copyright © 2001 Universidade Estadual de Campinas (UNICAMP).
  Authors: Helena C. G. Leitão and Jorge Stolfi.
  
  This file can be freely distributed, used, and modified, provided
  that this copyright and authorship notice is preserved, and that any
  modified versions are clearly marked as such.
  
  This software has NO WARRANTY of correctness or applicability for
  any purpose. Neither the authors nor their employers chall be held
  responsible for any losses or damages that may result from its use.
*/