//
// pgmiso.cxx
//
// Generates a set of isolines from a portable grey map (or PGM) file.
// Limitations: PGM must be binary, not ascii, and maxval should be 255,
// which means an 8 bit greyscale. You can force a maxval with the 
// pnmdepth utility. For more info on PGM, see the NetPBM package.
// 
// There are no dependencies, simply build with:
// c++ -o pgmiso pgmiso.cxx
//
// Copyright by Bram Stolk, bram at sara.nl
// Licensed under the Gnu Public License, or GPL.
// See gnu.org website for more information on this license.
//

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



static void drawline(unsigned char *dstimg, int w, int h, unsigned char col[3], int x0, int y0, int x1, int y1)
{
  int dx=0;
  int dy=0;
  if (x1<x0) dx=-1;
  if (x1>x0) dx=1;
  if (y1<y0) dy=-1;
  if (y1>y0) dy=1;
  do
  {
    memcpy(dstimg+3*x0+3*y0*w, col, 3);
    x0 += dx;
    y0 += dy;
  } while (x0 != x1 || y0 != y1);
  memcpy(dstimg+3*x0+3*y0*w, col, 3);
}


int main(int argc, char *argv[])
{
  if (argc!=4)
  {
    fprintf(stderr,"Usage: %s nr-of-isolines input.pgm output.ppm\n", argv[0]);
    exit(1);
  }
  int isocnt = atoi(argv[1]);
  FILE *f = fopen(argv[2],"rb");
  assert(f);

  int w,h,maxval;
  char ws;
  int retval=fscanf(f,"P5 %d %d %d%c",&w,&h,&maxval,&ws);
  assert(retval==4);
  assert(maxval==255);
  unsigned char *srcimg = new unsigned char [w*h];

  retval=fread(srcimg, w*h, 1, f);
  assert(retval==1);
  fclose(f);

  int step=4;
  int hstep=step/2;
  unsigned char *dstimg = new unsigned char [w*h*3];

  // Copy the source image to the destination image.
  // You may want to skip this if you want the isolated isolines only,
  // without the original greyscale image in the background.
  unsigned char *reader = srcimg;
  unsigned char *writer = dstimg;
  for (int i=0; i<w*h; i++)
  {
    *writer++ = *reader;
    *writer++ = *reader;
    *writer++ = *reader++;
  }

  for (int isonr = 0; isonr < isocnt; isonr++)
  {
    float parm = isonr / (isocnt-1.0);
    unsigned char col[3] = {(unsigned char) ((1.0-parm)*255), (unsigned char)(parm*255), 0};
    float isodelta = 255 / isocnt;
    int isoval = (int) ((isonr+0.5) * isodelta);

    fprintf(stderr,"parameter %5.3f: isoval %3d, colour %02x %02x %02x\n", parm, isoval, col[0], col[1], col[2]);
    for (int y=0; y<h-step; y+=step)
    {
      for (int x=0; x<w-step; x+=step)
      {
        int xx = x+step;
        int yy = y+step;
        unsigned char b=0;
        if (srcimg[x+y*w]   > isoval) b |= 1;
        if (srcimg[xx+y*w]  > isoval) b |= 2;
        if (srcimg[x+yy*w]  > isoval) b |= 4;
        if (srcimg[xx+yy*w] > isoval) b |= 8;
        if (b==0 || b==15)
        {
        }
        if (b==1 || b==14)
        {
          drawline(dstimg, w, h, col, x+hstep, y, x, y+hstep);
        }
        if (b==2 || b==13)
        {
          drawline(dstimg, w, h, col, x+hstep, y, x+step, y+hstep);
        }
        if (b==3 || b==12)
        {
          drawline(dstimg, w, h, col, x, y+hstep, x+step, y+hstep);
        }
        if (b==4 || b==11)
        {
          drawline(dstimg, w, h, col, x, y+hstep, x+hstep, y+step);
        }
        if (b==5 || b==10)
        {
          drawline(dstimg, w, h, col, x+hstep, y, x+hstep, y+step);
        }
        if (b==6)
        {
          drawline(dstimg, w, h, col, x+hstep, y, x, y+hstep);
          drawline(dstimg, w, h, col, x+hstep, y+step, x+step, y+hstep);
        }
        if (b==7 || b==8)
        {
          drawline(dstimg, w, h, col, x+hstep, y+step, x+step, y+hstep);
        }
        if (b==9)
        {
          drawline(dstimg, w, h, col, x+hstep, y, x+step, y+hstep);
          drawline(dstimg, w, h, col, x, y+hstep, x+hstep, y+step);
        }
      }
    }
  }
  
  f=fopen(argv[3],"wb");
  assert(f);
  fprintf(f,"P6 %d %d %d\n", w,h,maxval);
  retval=fwrite(dstimg, 3*w*h, 1, f);
  assert(retval==1);
  fclose(f);
  fprintf(stderr, "image with %d isolines saved as %s\n", isocnt, argv[3]);
}