# include "gastro2.h"

void gheader (char *file, CmpCatalog *Target) {

  double dR;
  Header header;
  FILE *f, *g;
  int i, oldsize, nbytes, status;
  char line[1024];

  if (!fits_read_header (file, &header)) {
    fprintf (stderr, "ERROR: can't find image file %s (3)\n", file);
    exit(0);
  }
  oldsize = header.size;

  if (Target[0].answer.N < 6) {
    fits_modify (&header, "NASTRO", "%d", 1, 0);
    fits_modify (&header, "NASTRO", "%C", 1, "number of stars used for astrometry");
    goto skipstuff;
  }
  
  fits_modify (&header, "NASTRO", "%d", 1, Target[0].answer.N);
  fits_modify (&header, "NASTRO", "%C", 1, "number of stars used for astrometry");

  if (Target[0].coords.Npolyterms > 1) {
    fits_modify (&header, "CTYPE1",   "%s",  1, "RA---PLY");
    fits_modify (&header, "CTYPE2",   "%s",  1, "DEC--PLY");
  } else {
    fits_modify (&header, "CTYPE1",   "%s",  1, "RA---TAN");
    fits_modify (&header, "CTYPE2",   "%s",  1, "DEC--TAN");
  }    
  fits_modify (&header, "CDELT1",   "%le", 1, Target[0].coords.cdelt1); 
  fits_modify (&header, "CDELT2",   "%le", 1, Target[0].coords.cdelt2);
  fits_modify (&header, "CRVAL1",   "%lf", 1, Target[0].coords.crval1);
  fits_modify (&header, "CRVAL2",   "%lf", 1, Target[0].coords.crval2);  
  fits_modify (&header, "CRPIX1",   "%lf", 1, Target[0].coords.crpix1);
  fits_modify (&header, "CRPIX2",   "%lf", 1, Target[0].coords.crpix2);
  fits_modify (&header, "PC001001", "%le", 1, Target[0].coords.pc1_1);
  fits_modify (&header, "PC001002", "%le", 1, Target[0].coords.pc1_2);
  fits_modify (&header, "PC002001", "%le", 1, Target[0].coords.pc2_1);
  fits_modify (&header, "PC002002", "%le", 1, Target[0].coords.pc2_2);
  fits_modify (&header, "NPLYTERM", "%d", 1, Target[0].coords.Npolyterms);
  if (Target[0].coords.Npolyterms > 1) {
    /* RA Terms */
    fits_modify (&header, "PCA1X2Y0", "%le", 1, Target[0].coords.polyterms[0][0]);   /* polyterms[0]); */
    fits_modify (&header, "PCA1X1Y1", "%le", 1, Target[0].coords.polyterms[1][0]);   /* polyterms[1]); */
    fits_modify (&header, "PCA1X0Y2", "%le", 1, Target[0].coords.polyterms[2][0]);   /* polyterms[2]); */

    if (Target[0].coords.Npolyterms > 2) {
      fits_modify (&header, "PCA1X3Y0", "%le", 1, Target[0].coords.polyterms[3][0]);   /* polyterms[3]); */
      fits_modify (&header, "PCA1X2Y1", "%le", 1, Target[0].coords.polyterms[4][0]);   /* polyterms[4]); */
      fits_modify (&header, "PCA1X1Y2", "%le", 1, Target[0].coords.polyterms[5][0]);   /* polyterms[5]); */
      fits_modify (&header, "PCA1X0Y3", "%le", 1, Target[0].coords.polyterms[6][0]);   /* polyterms[6]); */
    }
    /* Dec Terms */
    fits_modify (&header, "PCA2X2Y0", "%le", 1, Target[0].coords.polyterms[0][1]);   /* polyterms[7]); */
    fits_modify (&header, "PCA2X1Y1", "%le", 1, Target[0].coords.polyterms[1][1]);   /* polyterms[8]); */
    fits_modify (&header, "PCA2X0Y2", "%le", 1, Target[0].coords.polyterms[2][1]);   /* polyterms[9]); */

    if (Target[0].coords.Npolyterms > 2) {
      fits_modify (&header, "PCA2X3Y0", "%le", 1, Target[0].coords.polyterms[3][1]);   /* polyterms[10]); */
      fits_modify (&header, "PCA2X2Y1", "%le", 1, Target[0].coords.polyterms[4][1]);   /* polyterms[11]); */
      fits_modify (&header, "PCA2X1Y2", "%le", 1, Target[0].coords.polyterms[5][1]);   /* polyterms[12]); */
      fits_modify (&header, "PCA2X0Y3", "%le", 1, Target[0].coords.polyterms[6][1]);   /* polyterms[13]); */
    }
  }

  dR = fabs (Target[0].answer.dR*Target[0].coords.cdelt1*3600.0);
  fits_modify (&header, "CERROR", "%lf", 1, dR);
  fits_modify (&header, "CERROR", "%C", 1, "scatter in astrometry soln (arcsec)");
  fits_modify (&header, "CPRECISE", "%lf", 1, dR / sqrt(1.0*Target[0].answer.N));
  fits_modify (&header, "CPRECISE", "%C", 1, "precision of astrometry soln (arcsec)");
  fits_modify (&header, "EQUINOX", "%lf", 1, 2000.0);
  /* we force equinox to be 2000.0 for all images */

skipstuff:
  if (header.size > oldsize) {
    if (VERBOSE) fprintf (stderr, "header expanded, creating new copy\n");
    sprintf (line, "mv %s %s~", file, file);
    status = system (line);
    if (status) {
      fprintf (stderr, "ERROR: unable to create %s~, exiting\n", file);
      exit (0);
    }
    sprintf (line, "%s~", file);
    f = fopen (line, "r");
    g = fopen (file, "w");
    if (f == NULL) {
      fprintf (stderr, "ERROR: can't find image file %s (4)\n", line);
      exit(0);
    }
    if (g == NULL) {
      fprintf (stderr, "ERROR: can't open output image file %s (4)\n", file);
      exit(0);
    }
    nbytes = fwrite (header.buffer, 1, header.size, g);
    fseek (f, oldsize, SEEK_SET);
    for (i = 0; (nbytes = fread (header.buffer, 1, header.size, f)) > 0; i++) {
      if (nbytes != fwrite (header.buffer, 1, nbytes, g)) {
	fprintf (stderr, "ERROR: failure writing output data file\n");
	exit (0);
      }
    }
    fclose (f);
    fclose (g);
  } else {
    f = fopen (file, "r+");
    if (f == NULL) {
      fprintf (stderr, "ERROR: can't find image file %s (4)\n", file);
      exit(0);
    }
    
    fseek (f, 0, SEEK_SET);
    nbytes = fwrite (header.buffer, 1, header.size, f);
    
    fclose (f);
  }
  free (header.buffer);

}


