/*
 * Write a fits image to a file descriptor; useful for talking to DS9
 *
 * This version knows about psLib data structures
 */
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <assert.h>
#include <ctype.h>
#include <unistd.h>
#include <fcntl.h>
#include "pslib.h"
#include "poisErrorCodes.h"

#include "simpleFits.h"
   
#define FITS_SIZE 2880

typedef struct {
    char keyword[9];
    enum {
	UNKNOWN_CARD = -1,
	LOGICAL_CARD,
	STRING_CARD,
	INTEGER_CARD,
	DOUBLE_CARD,
	FLOAT_CARD,
    } type;
    union {
	char *l;
	char *s;
	int i;
	double d;
	float f;      
    } val;
    char *commnt;
} CARD;
      

/*
 * Write a string value
 */
static int
write_card_s(int fd,
	     int ncard,
	     char *record,
	     const char *keyword,
	     const char *val,
	     const char *commnt)
{
    char *card = &record[80*ncard];
    char value[20];			/* blank-padded val, if needed */

    if(strlen(val) < 8) {		/* FITS requires at least 8 chars */
	sprintf(value, "%-8s", val);
	val = value;
    }

    if(*keyword == '\0' ||
       !strcmp(keyword,"COMMENT") || !strcmp(keyword,"END") || !strcmp(keyword,"HISTORY")) {
	sprintf(card,"%-8.8s%-72s",keyword,val);
    } else {
	sprintf(card,"%-8.8s= '%s' %c%-*s",keyword,val,
		(commnt[0] == '\0' ? ' ' : '/'),
		(int)(80 - 14 - strlen(val)), commnt);
    }
/*
 * Write record if full
 */
    if(++ncard == 36) {
	if(write(fd, record, FITS_SIZE) != FITS_SIZE) {
	    psError(POIS_ERR_FITS, true,"Cannot write header record\n");
	}
	ncard = 0;
    }
   
    return ncard;
}   

/*****************************************************************************/
/*
 * Write an integer (%d)
 */
static int
write_card_i(int fd,
	     int ncard,
	     char *record,
	     const char *keyword,
	     const int val,
	     const char *commnt)
{
    char *card = &record[80*ncard];

    sprintf(card,"%-8.8s= %20d %c%-48s",keyword,val,
	    (commnt[0] == '\0' ? ' ' : '/'),commnt);
/*
 * Write record if full
 */
    if(++ncard == 36) {
	if(write(fd, record, FITS_SIZE) != FITS_SIZE) {
	    psError(POIS_ERR_FITS, true,"Cannot write header record\n");
	}
	ncard = 0;
    }
   
    return ncard;
}   

/*****************************************************************************/
/*
 * Write a logical value
 */
static int
write_card_l(int fd,
	     int ncard,
	     char *record,
	     const char *keyword,
	     const char *val,
	     const char *commnt)
{
    char *card = &record[80*ncard];

    if(strcmp(val,"T") != 0 && strcmp(val,"F") != 0) {
	psError(POIS_ERR_FITS, true,"Invalid logical %s for keyword %s\n",val,keyword);
	val = "?";
    }
    sprintf(card,"%-8.8s= %20s %c%-48s",keyword,val,
	    (commnt[0] == '\0' ? ' ' : '/'),commnt);
/*
 * Write record if full
 */
    if(++ncard == 36) {
	if(write(fd, record, FITS_SIZE) != FITS_SIZE) {
	    psError(POIS_ERR_FITS, true,"Cannot write header record\n");
	}
	ncard = 0;
    }
   
    return ncard;
}   

/*****************************************************************************/
/*
 * Write a double floating value
 */
static int
write_card_d(int fd,
	     int ncard,
	     char *record,
	     const char *keyword,
	     const double val,
	     const char *commnt)
{
    char *card = &record[80*ncard];

    sprintf(card,"%-8.8s= %20g %c%-48s",keyword,val,
	    (commnt[0] == '\0' ? ' ' : '/'),commnt);
/*
 * Write record if full
 */
    if(++ncard == 36) {
	if(write(fd, record, FITS_SIZE) != FITS_SIZE) {
	    psError(POIS_ERR_FITS, true,"Cannot write header record\n");
	}
	ncard = 0;
    }
   
    return ncard;
}   

/*****************************************************************************/
/*
 * Write a floating value
 */
static int
write_card_f(int fd,
	     int ncard,
	     char *record,
	     const char *keyword,
	     const float val,
	     const char *commnt)
{
    char *card = &record[80*ncard];

    sprintf(card,"%-8.8s= %20g %c%-48s",keyword,val,
	    (commnt[0] == '\0' ? ' ' : '/'),commnt);
/*
 * Write record if full
 */
    if(++ncard == 36) {
	if(write(fd, record, FITS_SIZE) != FITS_SIZE) {
	    psError(POIS_ERR_FITS, true,"Cannot write header record\n");
	}
	ncard = 0;
    }
   
    return ncard;
}   

/*****************************************************************************/
/*
 * Write a CARD
 */
static int
write_card(int fd,
	   int ncard,
	   char *record,
	   const CARD *card)
{
    switch (card->type) {
      case UNKNOWN_CARD:		// we should already have complained
	return ncard;
      case LOGICAL_CARD:
	return write_card_l(fd, ncard, record, card->keyword, card->val.l, card->commnt);
      case STRING_CARD:
	return write_card_s(fd, ncard, record, card->keyword, card->val.s, card->commnt);
      case INTEGER_CARD:
	return write_card_i(fd, ncard, record, card->keyword, card->val.i, card->commnt);
      case DOUBLE_CARD:
	return write_card_d(fd, ncard, record, card->keyword, card->val.d, card->commnt);
      case FLOAT_CARD:
	return write_card_f(fd, ncard, record, card->keyword, card->val.f, card->commnt);
    }
    printf("Impossible type for card %s, %d\n", card->keyword, card->type);
    abort();
    return ncard;			/* NOTREACHED */
}

/*****************************************************************************/
/*
 * Byte swap ABABAB -> BABABAB in place
 */
static void swap_2(char *arr,		// array to swap
		   int n)		// number of bytes
{
    char *end,
	t;

    if(n%2 != 0) {
	psError(POIS_ERR_FITS, true,"Attempt to byte swap odd number of bytes: %d\n",n);
	n = 2*(int)(n/2);
    }

    for(end = arr + n;arr < end;arr += 2) {
	t = arr[0];
	arr[0] = arr[1];
	arr[1] = t;
    }
}

/*
 * Byte swap ABCDABCD -> DCBADCBA in place (e.g. sun <--> vax)
 */
static void swap_4(char *arr,		// array to swap
		   int n)		// number of bytes
{
    char *end,
	t;

    if(n%4 != 0) {
	psError(POIS_ERR_FITS, true,"Attempt to byte swap non-multiple of 4 bytes: %d\n",n);
	n = 4*(int)(n/4);
    }

    for(end = arr + n;arr < end;arr += 4) {
	t = arr[0];
	arr[0] = arr[3];
	arr[3] = t;
	t = arr[1];
	arr[1] = arr[2];
	arr[2] = t;
    }
}

/*
 * Byte swap ABCDEFGH -> HGFEDCBA in place (e.g. sun <--> vax)
 */
static void swap_8(char *arr,		// array to swap
		   int n)		// number of bytes
{
    char *end,
	t;

    if(n%8 != 0) {
	psError(POIS_ERR_FITS, true,"Attempt to byte swap non-multiple of 8 bytes: %d\n",n);
	n = 8*(int)(n/8);
    }

    for(end = arr + n;arr < end;arr += 8) {
	t = arr[0];
	arr[0] = arr[7];
	arr[7] = t;
	t = arr[1];
	arr[1] = arr[6];
	arr[6] = t;
	t = arr[2];
	arr[2] = arr[5];
	arr[5] = t;
	t = arr[3];
	arr[3] = arr[4];
	arr[4] = t;
    }
}

/*****************************************************************************/

static int
write_fits_hdr(int fd,
	       int bitpix,
	       int naxis,
	       int *naxes,
	       CARD **cards,		/* extra header cards */
	       const psMetadata *metaData, /* metadata to write; or NULL */
	       int primary)		/* is this the primary HDU? */
{
    int i;
    int ncard;				/* number of cards written */
    char record[FITS_SIZE + 1];		/* write buffer */
   
    ncard = 0;
    if(primary) {
	ncard = write_card_l(fd,ncard,record,"SIMPLE","T","");
    } else {
	ncard = write_card_s(fd,ncard,record,"XTENSION","IMAGE","");
    }
    ncard = write_card_i(fd,ncard,record,"BITPIX",bitpix,"");
    ncard = write_card_i(fd,ncard,record,"NAXIS",naxis,"");
    for(i = 0; i < naxis; i++) {
	char key[] = "NAXIS.";
	sprintf(key, "NAXIS%d", i + 1);
	ncard = write_card_i(fd,ncard,record,key,naxes[i],"");
    }
    if(primary) {
	ncard = write_card_l(fd,ncard,record,
			     "EXTEND","T","There may be extensions");
    }
/*
 * Write extra header cards
 */
    if(cards != NULL) {
	while(*cards != NULL) {
	    ncard = write_card(fd,ncard,record,*cards);
	    cards++;
	}
    }
    /*
     * Is there metadata to write?
     */
    if (metaData != NULL) {
	CARD card;			    // used to write metadata to
	card.keyword[8] = '\0';		    // we\'ll be using strncpy, so need to guarantee NUL termination
	psMetadataItem *meta = NULL;	    // used to traverse a set of metadata
	int new_error = 1;		    // is this a new error?


	psListIterator *iter = psListIteratorAlloc(metaData->list, PS_LIST_HEAD, false); 
	while ((meta = psListGetAndIncrement(iter)) != NULL) {
	    strncpy(card.keyword, meta->name, 8);
	    card.commnt = meta->comment;
	    switch (meta->type) {
	      case PS_META_F64:
		card.type = DOUBLE_CARD;
		card.val.d = meta->data.F64;
		break;
	      case PS_META_F32:
		card.type = FLOAT_CARD;
		card.val.f = meta->data.F32;
		break;
	      case PS_META_S32:
		card.type = INTEGER_CARD;
		card.val.i = meta->data.S32;
		break;
	      case PS_META_STR:
		card.type = STRING_CARD;
		card.val.s = meta->data.V;
		break;
	      default:
		card.type = UNKNOWN_CARD;
		psError(POIS_ERR_FITS, new_error,
			"I don't know how to write metadata of type %d", meta->type);
		new_error = 0;
		break;
	    }
	    if (card.type != UNKNOWN_CARD) {
		ncard = write_card(fd, ncard, record, &card);
	    }
	}
	psFree(iter);
    }

    ncard = write_card_s(fd,ncard,record,"END","","");
    while(ncard != 0) {
	ncard = write_card_s(fd,ncard,record,"","","");
    }

    return 0;
}

/*
 * Pad out to a FITS record boundary
 */
static void pad_to_fits_record(int fd,		// output file descriptor
			       int npixel,	// number of pixels already written to HDU
			       int bitpix)	// bitpix for this datatype
{
    const int bytes_per_pixel = (bitpix > 0 ? bitpix : -bitpix)/8;
    int nbyte = npixel*bytes_per_pixel;
    
    if(nbyte%FITS_SIZE != 0) {
	char record[FITS_SIZE + 1];	/* write buffer */
	
	nbyte = FITS_SIZE - nbyte%FITS_SIZE;
	memset(record, ' ', nbyte);
	if(write(fd, record, nbyte) != nbyte) {
	    psError(POIS_ERR_FITS, true, "error padding file to multiple of fits block size\n");
	}	       
    }
}

static int write_fits_data(int fd,
			   int bitpix,
			   int npix,
			   char *data,
			   int pad_record) // Pad data out to a FITS record boundary?
{
    char *buff = NULL;			/* I/O buffer */
    const int bytes_per_pixel = (bitpix > 0 ? bitpix : -bitpix)/8;
    int nbyte = npix*bytes_per_pixel;
    int swap_bytes = 0;			/* the default */
    static int warned = 0;		/* Did we warn about BZERO/BSCALE? */

#if PS_BIGENDIAN == 0			/* we'll need to byte swap FITS */
    if(bytes_per_pixel > 1) {
	swap_bytes = 1;
    }
#endif

    if(swap_bytes) {
	if((buff = psAlloc(FITS_SIZE*bytes_per_pixel)) == NULL) {
	    psError(POIS_ERR_FITS, true,"Can't allocate storage for buff\n");
	    return -1;
	}
    }
    
    if(bytes_per_pixel == 2 && !warned) {
	warned = 1;
	fprintf(stderr,"Worry about BZERO/BSCALE\n");
    }
   
    while(nbyte > 0) {
	int nwrite = (nbyte > FITS_SIZE) ? FITS_SIZE : nbyte;
      
	if(swap_bytes) {
	    memcpy(buff, data, nwrite);
	    if(bytes_per_pixel == 2) {
		swap_2((char *)buff, nwrite);
	    } else if(bytes_per_pixel == 4) {
		swap_4((char *)buff, nwrite);
	    } else if(bytes_per_pixel == 8) {
	        swap_8((char *)buff, nwrite);
	    } else {
		fprintf(stderr,"You cannot get here\n");
		abort();
	    }
	} else {
	    buff = data;
	}
      
	if(write(fd, buff, nwrite) != nwrite) {
	    perror("Error writing image: ");
	    break;
	}

	nbyte -= nwrite;
	data += nwrite;
    }
    if (nbyte == 0) {			/* no write errors */
	if (pad_record) {
	    pad_to_fits_record(fd, npix, bitpix);
	}
    }
   
    if(swap_bytes) {
	psFree(buff);
    }

    return (nbyte == 0 ? 0 : -1);
}

int
rhlWriteFits(int fd,			/* file descriptor to write to */
	     const psImage *data,	/* The data to write */
	     const psMetadata *meta,	/* associated metadata */
	     const char *WCS)		/* name of WCS for pixel coordinates */
{
    CARD **cards = NULL;		/* extra header info */
    int CRVAL1 = 0, CRVAL2 = 0;		/* card indices for keywords */
    int i;
    const int naxis = 2;		// == NAXIS
    int naxes[2];			/* values of NAXIS1 etc */
    int ncard = 12;			/* number of extra header cards */

    psErrorClear();			// clear error stack
    /*
     * Does the data have a WCS? If so, call our pixel-based WCS "A"
     */
    if (WCS != NULL && *WCS == '\0') {
	WCS = NULL;
    }
    if (WCS == NULL &&
	meta != NULL && psMetadataLookup((psMetadata *)meta, "CTYPE1") != NULL) {
	WCS = "A";
    } else {
	psErrorClear();			// set by psMetadataLookup XXX
	WCS = "";
    }
    /*
     * What sort if image is it?
     */
    int bitpix;				// == BITPIX
    char **rows = NULL;
    switch (data->type.type) {
      case PS_TYPE_F64:
	bitpix = -64;
	rows = (char **)data->data.F64;
	break;
      case PS_TYPE_F32:
	bitpix = -32;
	rows = (char **)data->data.F32;
	break;
      case PS_TYPE_U16:
	bitpix = -16;
	rows = (char **)data->data.U16;
	break;
      case PS_TYPE_S8:
	bitpix = 8;
	rows = (char **)data->data.S8;
	break;
      case PS_TYPE_U8:
	bitpix = 8;
	rows = (char **)data->data.U8;
	break;
      default:
	bitpix = 0;			// make gcc happy
	psAbort(__func__, "Unsupported image data type: %d", data->type.type);
	break;				// NOTREACHED
    }
    assert (rows != NULL);
    
    /*
     * Allocate cards for FITS headers
     */
    if (WCS != NULL) {
       cards = psAlloc((ncard + 1)*sizeof(CARD *));
       assert(cards != NULL);
       
       cards[0] = psAlloc((ncard + 1)*sizeof(CARD));
       assert(cards[0] != NULL);
       
       for(i = 0; i < ncard; i++) {
	  cards[i] = cards[0] + i;
       }
       cards[i] = NULL;
       /*
	* Generate cards for WCS, so that pixel (0,0) is correctly labelled
	*/
       i = 0;
       sprintf(cards[i]->keyword, "CRVAL1%s", WCS); CRVAL1 = i;
       cards[i]->type = INTEGER_CARD;
       cards[i]->val.i = 0;
       cards[i]->commnt = "(output) Column pixel of Reference Pixel";
       
       i++;
       sprintf(cards[i]->keyword, "CRVAL2%s", WCS); CRVAL2 = i;
       cards[i]->type = INTEGER_CARD;
       cards[i]->val.i = 0;
       cards[i]->commnt = "(output) Row pixel of Reference Pixel";
       
       i++;
       sprintf(cards[i]->keyword, "CRPIX1%s", WCS); CRVAL1 = i;
       cards[i]->type = INTEGER_CARD;
       cards[i]->val.i = 1;
       cards[i]->commnt = "Column Pixel Coordinate of Reference";
       
       i++;
       sprintf(cards[i]->keyword, "CRPIX2%s", WCS); CRVAL1 = i;
       cards[i]->type = INTEGER_CARD;
       cards[i]->val.i = 1;
       cards[i]->commnt = "Row Pixel Coordinate of Reference";
       
       i++;
       sprintf(cards[i]->keyword, "CTYPE1%s", WCS); CRVAL1 = i;
       cards[i]->type = STRING_CARD;
       cards[i]->val.s = "LINEAR";
       cards[i]->commnt = "Type of projection";
       
       i++;
       sprintf(cards[i]->keyword, "CTYPE1%s", WCS); CRVAL1 = i;
       cards[i]->type = STRING_CARD;
       cards[i]->val.s = "LINEAR";
       cards[i]->commnt = "Type of projection";
       
       i++;
       sprintf(cards[i]->keyword, "CUNIT1%s", WCS); CRVAL1 = i;
       cards[i]->type = STRING_CARD;
       cards[i]->val.s = "PIXEL";
       cards[i]->commnt = "Column unit";
       
       i++;
       sprintf(cards[i]->keyword, "CUNIT2%s", WCS); CRVAL1 = i;
       cards[i]->type = STRING_CARD;
       cards[i]->val.s = "PIXEL";
       cards[i]->commnt = "Row unit";
      
       cards[++i] = NULL;
       assert(i <= ncard);
    }

    naxes[0] = data->numCols;
    naxes[1] = data->numRows;
    
    write_fits_hdr(fd, bitpix, naxis, naxes, cards, meta, 1);
    for (int r = 0; r < data->numRows; r++) {
	if(write_fits_data(fd, bitpix, data->numCols, rows[r], ((r==data->numRows - 1) ? 1 : 0)) < 0) {
	    psError(POIS_ERR_UNKNOWN, 1, "Error writing data for row %d", r);
	    break;
	}
    }

    if (cards != NULL) {
       psFree(cards[0]);
       psFree(cards);
    }

    psErr *err = psErrorLast();
    int code = err->code == PS_ERR_NONE ? 0 : -1;
    psFree(err);
    
    return code;
}   

/************************************************************************************************************/

int
rhlWriteFitsFile(const char *filename,	/* file to write to (or "| cmd"); NULL => "|xpaset ds9 fits"*/
		 const psImage *data,	/* The data to write */
		 const psMetadata *meta, /* associated metadata */
		 const char *WCS)	/* name of WCS for pixel coordinates */
{
    if (filename == NULL) {
	filename =  "| xpaset ds9 fits"; // i.e. display on DS9
    }

    int fd;

    if (filename[0] == '|') {		// a command
	const char *cmd = filename + 1;
	while (isspace(*cmd)) {
	    cmd++;
	}

	fd = fileno(popen(cmd, "w"));
    } else {
	fd = creat(filename, 777);
    }

    if (fd < 0) {
	psError(POIS_ERR_FITS, 1, "Cannot open \"%s\"", filename);
	return -1;
    }

    int ret = rhlWriteFits(fd, data, meta, WCS);

    (void)close(fd);

    return ret;
}
