#include <stdio.h>
#include <stdlib.h>
#include <string.h>
//#include <libshp/shapefil.h>
#include <shapefil.h>


char buf[1000];
int dbrecnr = 0;

static
void usage(char *myname, char *err) {
  if (err) fprintf(stderr, "%s\n", err);
  fprintf(stderr, "usage: %s [-noclose] [-<deg>] infilename outbasename\n", myname);
  exit(1);
}

static
void readonepoly(FILE *infile, SHPHandle shpfile, DBFHandle dbfhandle, int dbfieldnr, double shift, int do_close) {
  /* upon entering, the first data line is already in buf */
  double *xa, *ya;
  int npairs = 0;
  int mallocsize;

  fprintf(stderr, "starting polygon nr %d ", dbrecnr);
  xa = malloc(100*sizeof(double));
  ya = malloc(100*sizeof(double));
  mallocsize = 100;
  do {
    double x, y;
    int err;
    if ('>' == buf[0]) break;
    err = sscanf(buf, "%lg %lg", &x, &y);
    if (err != 2) {
      fprintf(stderr, "sscanf did not match 2 doubles\n");
      break;
    }
    //if (x < -180) fprintf(stderr, "\nfound small X %g\n", x);
    // avoid wrap over date line
    if (x < 0) x += 360;
    x += shift;
    if (x > 360) x -= 360;
    if (x < 0) x += 360;
    if (npairs == mallocsize) {
      xa = realloc(xa, (mallocsize+100)*sizeof(double));
      ya = realloc(ya, (mallocsize+100)*sizeof(double));
      mallocsize += 100;
    }
    xa[npairs] = x;
    ya[npairs] = y;
    npairs++;
    fgets(buf, sizeof(buf), infile);
  } while (!feof(infile));

  /* sigh. check if the poly line is closed. if not, close it */
  if (npairs > 0) {
    SHPObject *shpobj;
    if (xa[0] != xa[npairs-1] && ya[0] != ya[npairs-1]) {
      if (do_close) {
        fprintf(stderr, " closing ");
        if (npairs == mallocsize) {
          xa = realloc(xa, (mallocsize+1)*sizeof(double));
          ya = realloc(ya, (mallocsize+1)*sizeof(double));
          mallocsize += 1;
        }
        xa[npairs] = xa[0];
        ya[npairs] = ya[0];
        npairs++;
      } else {
        fprintf(stderr, "unclosed ");
      }
    } else {
      fprintf(stderr, "closed ");
    }
    fprintf(stderr, "got %d pairs\n", npairs);
    shpobj = SHPCreateSimpleObject(SHPT_POLYGON, npairs, xa, ya, NULL);
    SHPWriteObject(shpfile, -1, shpobj);
    SHPDestroyObject(shpobj);
    free(xa);
    free(ya);
    if (!DBFWriteDoubleAttribute(dbfhandle, dbrecnr, dbfieldnr, (double)dbrecnr)) printf("DBFWriteNULLAttribute failed\n");
    dbrecnr++;
  }
}

int main(int argc, char **argv) {
  FILE *infile;
  SHPHandle shpfile;
  DBFHandle dbfhandle;
  int dbfieldnr;
  char *p;
  double shift = 0;
  int do_close = 1;
  char *myname = argv[0];

  if (argc < 3) { usage(myname, "need input and output file names"); }
  if (!strcmp("-noclose", argv[1])) {
    do_close = 0;
    argv++;
    argc--;
  }
  if (argv[1][0] == '-') {
    shift = strtod(argv[1]+1, NULL);
    fprintf(stderr, "shifting longitudes by additional %g degrees\n", shift);
    argv++;
    argc--;
  }
  if (argc != 3) { usage(myname, "need input and output file names"); }
  infile = fopen(argv[1], "r");
  if (!infile) { perror(argv[1]); exit(1); }
  shpfile = SHPCreate(argv[2], SHPT_POLYGON);
  if (shpfile < 0) { perror(argv[2]); exit(1); }
  dbfhandle = DBFCreate(argv[2]);
  if (dbfhandle < 0) { perror(argv[2]); exit(1); }
  dbfieldnr = DBFAddField(dbfhandle, "dummy", FTDouble, 10, 4);

  while (!feof(infile)) {
    p = fgets(buf, sizeof(buf), infile);
    if (!p) break;
    if ('>' == *p) {
      //fprintf(stderr, "skipping %s\n", buf);
      continue;
    }
    readonepoly(infile, shpfile, dbfhandle, dbfieldnr, shift, do_close);
  }
  SHPClose(shpfile);
  DBFClose(dbfhandle);
  fprintf(stderr, "found %d polygons\n", dbrecnr);
  return 0;
}