// cc -o reoffset reoffset.c  -lnetpbm -lm
// ./reoffset --xstep 40.0 --ystep 66.2 AW-003.pnm

// BUG: UPDATED PARAMS: IPA-02_1_1.pnm 54.8 14.200000 71 -0.1 0.100000

// This program takes an image which has already been aligned with the
// horizontal axis, using 'deskew', and a given line and character
// spacing in x and y, and determines the x and y starting offset that
// aligns the grid as best as possible between lines and characters.
//
// This is specifically for images of fixed-pitch/fixed character width
// printouts such as old listings of computer printouts.
// Actual segmentation is performed by 'makegrid' (or 'reoffset', in
// the early stages of development)

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>

#ifndef FALSE
#define FALSE (0 != 0)
#endif
#ifndef TRUE
#define TRUE (0 == 0)
#endif

#include <netpbm/pnm.h>
#include <netpbm/shhopt.h>

struct cmdline_info {
  const char *infilename; /* Filespec of input file */
  char *basename;
  float ystep, xstep; // forced
  unsigned int help, verbose, xq, yq;
};

FILE *ifp;
xel *xelrow;
xel **xelrows;
xel bgxel;
int rows, cols, format;
int row;
xelval maxval;

struct cmdline_info cmdline;

static void parse_command_line(int argc, char **argv,
                               struct cmdline_info *cmdlineP) {
  optStruct3 opt;
  unsigned int option_def_index = 0;
  optEntry *option_def = malloc(100 * sizeof(optEntry));

  opt.opt_table = option_def;
  opt.short_allowed = FALSE;
  opt.allowNegNum = TRUE;

  cmdlineP->verbose = FALSE;
  cmdlineP->help = FALSE;

  OPTENT3('h', "help",    OPT_FLAG,   NULL, &cmdlineP->help, 0);
  OPTENT3(  0, "info",    OPT_FLAG,   NULL, &cmdlineP->help, 0);
  OPTENT3('v', "verbose", OPT_FLAG,   NULL, &cmdlineP->verbose, 0);
  OPTENT3('x', "xstep",   OPT_FLOAT,  &cmdlineP->xstep, &cmdlineP->xq, 0);
  OPTENT3('y', "ystep",   OPT_FLOAT,  &cmdlineP->ystep, &cmdlineP->yq, 0);

  pm_optParseOptions3(&argc, argv, opt, sizeof(opt), 0);
  if (argc - 1 < 1) {
    cmdlineP->infilename = "-";
    if (cmdlineP->verbose)
      fprintf(stderr, "Processing /dev/stdin\n");
    cmdlineP->infilename = "temp";
    cmdlineP->basename = strdup("temp");
  } else {
    char *p;
    cmdlineP->infilename = argv[1];
    if ((p = strrchr(argv[1], '/')) != NULL) {
      cmdlineP->basename = strdup(p + 1);
    } else
      cmdlineP->basename = strdup(argv[1]);
    if ((p = strrchr(cmdlineP->basename, '.')) != NULL) {
      *p = '\0';
    }

    if (argc - 1 > 1) {
      pm_error("too many arguments (%d).  "
               "The only positional argument is the filespec.",
               argc - 1);
    }
    if (cmdlineP->verbose)
      fprintf(stderr, "Processing %s\n", cmdlineP->infilename);
  }
}

void read_image(int *maxx, int *maxy) {
  /* Read a bitmap, return its size;
     Numbers are inclusive, ie bounds are 0:maxx-1, 0:maxy-1 */

  ifp = pm_openr(cmdline.infilename);

  pnm_readpnminit(ifp, &cols, &rows, &maxval, &format);
  *maxx = cols;
  *maxy = rows;
  xelrows = malloc(rows * sizeof(xel *));

  for (row = 0; row < rows; row += 1) {
    xelrow = pnm_allocrow(cols);
    if (row == 0)
      bgxel = pnm_backgroundxelrow(xelrow, cols, maxval, format);
    pnm_readpnmrow(ifp, xelrow, cols, maxval, format);
    xelrows[row] = xelrow;
  }
}

unsigned int samplepixel(int x, int y) {
  xel *el = &xelrows[y][x];
  int rp = el->r & 255;
  int gp = el->g & 255;
  int bp = el->b & 255;
#ifdef PREV_VERSION
  int c = (rp + gp + bp) / 3;
  /* get the value of a pixel.  For now, apply threshholding and return 0/1 */
  if (c >= 64)
    return (1); /* white */
  return (0);   /* black */
  //  return(2); /* our drawn line for display purposes */
#endif
  return rp + gp + bp;
}

//
// The code below is a DDA trivially adapted from a line-drawing algorithm;
// I have faster algorithms than this - even an anti-aliased one - but
// this is rock-solid reliable code that can be trusted during the debugging
// phase (and also it appears to be fast enough already!)
//
// For this program we don't need a general purpose line - we only use
// horizontal and vertical lines.  However I might as well leave it in as I may
// decide at some futute date that I don't want to straighten the image
// externally first.
//

int count_white(int x1, int y1, int x2, int y2) {
  /* calculate a DDA from (xl,yb) to (xr,yt). */
  /* instead of plotting pixels, add them up somehow */
  int temp, dx, dy, x, y, x_sign, y_sign, flag;
  int tot = 0;

  dx = abs(x2 - x1);               // Delta of X
  dy = abs(y2 - y1);               // Delta of Y
  if (((dx >= dy) && (x1 > x2)) || // Make sure that first coordinate
      ((dy > dx) && (y1 > y2)))    // is the one with least value
  {
    temp = x1;
    x1 = x2;
    x2 = temp;
    temp = y1;
    y1 = y2;
    y2 = temp;
  };

  if ((y2 - y1) < 0)
    y_sign = -1; // The direction into which Y-coord shall travel
  else
    y_sign = 1; // Same for X

  if ((x2 - x1) < 0)
    x_sign = -1; // ---- " ----
  else
    x_sign = 1; // ---- " ----

  if (dx >= dy) // Which one of the deltas is the greatest one
  {
    for (x = x1, y = y1, flag = 0; x <= x2; x++, flag += dy) // From x1 to x2
    {                 // Also increase the
      if (flag >= dx) // Increase/decrease     // flag (displacement value)
      {               // y!
        flag -= dx;
        y += y_sign;
      };
      tot += samplepixel(x, y); // Plot the pixel
    };
  } else {
    /* This is the same as above, just with x as y and vice versa.*/
    for (x = x1, y = y1, flag = 0; y <= y2; y++, flag += dx) {
      if (flag >= dy) {
        flag -= dy;
        x += x_sign;
      };
      tot += samplepixel(x, y);
    };
  };
  return (tot);
}

static inline int iround(float num) {
  return (num - floor(num) > 0.5) ? ceil(num) : floor(num);
}

int main(int argc, char **argv) {
  int x, y, maxx, maxy;
  int *xscore, *yscore, *start, *end, *weight, *midpoint, *step;
  int xstep, ystep;
  int highest, lowest;

  pnm_init(&argc, argv);
  parse_command_line(argc, argv, &cmdline);

  if (cmdline.help) {
    fprintf(stderr, "syntax: reoffset [options] filename.png\n\n");
    exit(EXIT_SUCCESS);
  }
  // &cmdlineP->xstep, &cmdlineP->xq
  if (!cmdline.xq || !cmdline.yq) {
    fprintf(stderr, "syntax: reoffset --xstep nn.n --ystep nn.n filename.png\n\n");
    exit(EXIT_FAILURE);
  }

  xstep = iround(cmdline.xstep * 10.0);  
  ystep = iround(cmdline.ystep * 10.0);
  
  read_image(&maxx, &maxy); /* get width and height of image */
  pm_close(ifp);

  if (cmdline.verbose) {
    fprintf(stderr, "I will be verbose\n");
  }

  // add --lines and --cols for estimating rough starting point by dividing page
  // height by lines etc.
  if (cmdline.verbose) {
    fprintf(stderr, "%s %dx%d grid %d.%0dx%d.%0d\n",
            cmdline.infilename,
            maxx, maxy,
            xstep/10,xstep%10, ystep/10,ystep%10);
  }

  // I have deliberately duplicated the code below so that X and Y are handled
  // separately. They *could* be merged into one by adding another dimension to
  // most arrays for [x] and [y] indicies, but I thought it more readable to
  // keep them separate, and also I'm expecting that some tweaking might be
  // necessary for handling rows that is not applicable to handling colums, in
  // the future.

  // ============================================================================================================
  // ROUGHLY LOCATE Y SEPARATION

  xscore = malloc((maxx + 1) * sizeof(int));
  yscore = malloc((maxy + 1) * sizeof(int));
  start = malloc((maxy + 1) * sizeof(int));
  end = malloc((maxy + 1) * sizeof(int));
  weight = malloc((maxy + 1) * sizeof(int));
  midpoint = malloc((maxy + 1) * sizeof(int));
  step = malloc((maxy + 1) * sizeof(int));

  highest = -1;
  lowest = 0x7FFFFFFF;

  for (y = 0; y < maxy; y++) {
    /* draw a line from 0,0 to maxx,y for all y */
    yscore[y] = count_white(0 /*maxx/10*/, y, maxx /*- maxx/10*/ - 1,
                            y); // skip left and right margins
    // if (cmdline.verbose) fprintf(stderr, "yscore[%d] = %d\n", y, yscore[y]);
    if (yscore[y] > highest)
      highest = yscore[y];
    if (yscore[y] < lowest)
      lowest = yscore[y];
  }
  if (cmdline.verbose)
    fprintf(stderr, "lowest = %d\n", lowest);

  // This hack or something like it is needed in order not to find a better
  // result from sampling only every second line by accident :-/
  for (y = 0; y < maxy; y++) {
    yscore[y] -= lowest;
    //yscore[y] = (yscore[y] * 80) / (highest - lowest);
    //yscore[y] -= 70;
    if (yscore[y] < 0)
      yscore[y] = 0;
  }

  // A page of maxy pixels in height is likely to have something between 30 and
  // 200 lines of text so we'll be looking at roughly maxy/200 to maxy/30 pixels
  // per line. Then to allow .1 pixel resolution, we'll use 10 times as many
  // sample separations.
  {
    int gap, lines;
    int mingap, maxgap;
    //int best_count[maxgap - mingap + 1];
    int best_offset[1];
    int best_lines[1];
    int all_best_gap = -1, all_best_count = -1;
    mingap = maxgap = gap = ystep;
    {
      int offset;
      int this_best_count = 0, this_best_offset = -1, this_best_lines = 0;
      for (offset = 0; offset < gap; offset++) {
        int count = 0, y = offset;
        lines = 0;
        for (;;) {
          int truncy = y / 10;
          // if ((truncy < maxy/10) || (truncy > (maxy-maxy/10))) { // skip
          // headers and footers
          count += yscore[truncy];
          //}
          lines += 1;
          y += gap;
          if (y / 10 >= maxy)
            break;
        }
        if (count > this_best_count) {
          fprintf(stderr, "Offset %d: Updating count (y) %d\n", offset, count);
          this_best_count = count;
          this_best_offset = offset;
          this_best_lines = lines;
        }
      }
      this_best_count = (this_best_count) / lines;
      //best_count[gap - mingap] = this_best_count;
      // fprintf(stderr, "%0.1f: %d\n", gap/10.0, this_best_count);
      best_offset[gap - mingap] = this_best_offset;
      best_lines[gap - mingap] = this_best_lines;
      if (this_best_count > all_best_count) {
        fprintf(stderr, "Offset %d: Updating all best count (y) %d\n", this_best_offset, this_best_count);
        all_best_count = this_best_count;
        all_best_gap = gap;
      }
    }

    fprintf(stdout,
            "%s %0.1f %f %d ",
            cmdline.infilename, all_best_gap / 10.0,
            best_offset[all_best_gap - mingap] / 10.0,
            best_lines[all_best_gap - mingap]);
  }

  // ------------------------------------------------------------------------------------------------------------
  // NOW DO THE SAME FOR X!

  // This goes wrong in roughly 1 in 40 pages, where it usually finds a column
  // separation that is 50% or 33% of what it should be. This can be corrected
  // on a second pass after generated the majority of the pages - we can go back
  // and regenerate any pages whose separations were far from the mean.

  if (start) {
    free(start);
    start = NULL;
  }
  if (end) {
    free(end);
    end = NULL;
  }
  if (weight) {
    free(weight);
    weight = NULL;
  }
  if (midpoint) {
    free(midpoint);
    midpoint = NULL;
  }
  if (step) {
    free(step);
    step = NULL;
  }

  start = malloc((maxx + 1) * sizeof(int));
  end = malloc((maxx + 1) * sizeof(int));
  weight = malloc((maxx + 1) * sizeof(int));
  midpoint = malloc((maxx + 1) * sizeof(int));
  step = malloc((maxx + 1) * sizeof(int));

  highest = -1;
  lowest = 0x7FFFFFFF;

  for (x = 0; x < maxx; x++) {
    /* draw a line from x,0 to x,maxy for all x */
    xscore[x] = count_white(x, 0, x, maxy - 1);
    // if (cmdline.verbose) fprintf(stderr, "xscore[%d] = %d\n", x, xscore[x]);
    if (xscore[x] > highest)
      highest = xscore[x];
    if (xscore[x] < lowest)
      lowest = xscore[x];
  }
  for (x = 0; x < maxx; x++) {
    xscore[x] -= lowest;
    //xscore[x] = (xscore[x] * 80) / (highest - lowest);
    //xscore[x] -= 70;
    if (xscore[x] < 0)
      xscore[x] = 0;
  }
  {
    // Pages were probably 80 or 132 chars per line, but we'll look between 50
    // and 200 so we'll be looking at roughly maxy/200 to maxy/50 pixels per
    // char. Then to allow .1 pixel resolution, we'll use 10 times as many
    // sample separations.
    int gap, lines;
    int mingap, maxgap;
    
    //int best_count[maxgap - mingap + 1];
    int best_offset[1];
    int best_lines[1];
    int all_best_gap = -1, all_best_count = -1;
    mingap = maxgap = gap = xstep;
    {
      int offset;
      int this_best_count = 0, this_best_offset = -1, this_best_lines = 0;
      for (offset = 0; offset < gap; offset++) {
        int count = 0, x = offset;
        lines = 0;
        for (;;) {
          int truncx = x / 10;
          // if ((truncy < maxy/10) || (truncy > (maxy-maxy/10))) { // skip
          // headers and footers
          count += xscore[truncx];
          //}
          lines += 1;
          x += gap;
          if (x / 10 >= maxx)
            break;
        }
        if (count > this_best_count) {
          fprintf(stderr, "Offset %d: Updating count (x) %d\n", offset, count);
          this_best_count = count;
          this_best_offset = offset;
          this_best_lines = lines;
        }
      }
      this_best_count = (this_best_count) / lines;
      //best_count[gap - mingap] = this_best_count;
      // fprintf(stderr, "%0.1f: %d\n", gap/10.0, this_best_count);
      best_offset[gap - mingap] = this_best_offset;
      best_lines[gap - mingap] = this_best_lines;
      if (this_best_count > all_best_count) {
        fprintf(stderr, "Offset %d: Updating all best count (x) %d\n", this_best_offset, this_best_count);
        all_best_count = this_best_count;
        all_best_gap = gap;
      }
    }

    fprintf(stdout,
            "%0.1f %f %d\n",
            all_best_gap / 10.0,
            best_offset[all_best_gap - mingap] / 10.0,
            //best_count[all_best_gap - mingap],
            best_lines[all_best_gap - mingap]);
  }

  exit(EXIT_SUCCESS);
  return EXIT_FAILURE;
}
