/* This program creates .fits images from rdo input files.   R. Stiening September 17, 2005
USAGE:  ./raw2fits infile.rdo frame#
ACTION: 3 fits files (r1, r2, and r2-r1) are created from frame# frame of infile.rdo
COMPILE SCRIPT: gcc raw2fits.c -lcfitsio -lm -oraw2fits
*/

#include <stdio.h>
#include <fitsio.h>
#include <errno.h>
#include <string.h>

#define pixels   65536

unsigned short int   r1_buffer[pixels];
unsigned short int   r2_buffer[pixels];
int  		    r21_buffer[pixels];
int *array[256];
unsigned short int frames=0;
int i,j;
int size_read;
int size_writ;
int status = 0;
int true = 1;
int bitpix=LONG_IMG;
unsigned short int selected_frame;
long fpixel,nelements,exposure;
long naxis = 2;
long naxes[2] = {256,256};
void printerror(int status);
char fname1[25];
char fname2[25];
char fname3[25];

int main(int argc, char *argv[])
{
if(argc != 3) 
	{
	printf("Command line does not have two arguments.\n");
	printf("Usage raw2fits filename.rdo frame \n");
	return 1;
	}
	selected_frame=atoi(argv[2]);
	FILE *fin;
	fin=fopen(argv[1],"r");
	if(fin == (FILE *)NULL)
		{
		printf("%s: %s\n",argv[1],strerror(errno));
		exit(-1);
		}
       for(i=0;i<pixels;i++) 
                {
                r1_buffer[i] = 0;
                r2_buffer[i] = 0;
		r21_buffer[i]= 0;
                }
	while(frames < selected_frame)
        {
	size_read=fread(r1_buffer,1,sizeof(r1_buffer),fin);
	if(size_read != sizeof(r1_buffer)) break;
	size_read=fread(r2_buffer,1,sizeof(r2_buffer),fin);
	if(size_read != sizeof(r2_buffer)) break;	
	frames++;
	}
	printf("Frames read= %i\n",frames);
	fclose(fin);
for(i=0;i<8;i++)
	{
	*(fname1+i) = *(argv[1]+i);
	*(fname2+i) = *(argv[1]+i);
	*(fname3+i) = *(argv[1]+i);
	}
	printf("%s\n",argv[1]);
	printf("%s\n",fname1);
	sprintf(fname1,"%s_r1_%i.fits",fname1,selected_frame);
	sprintf(fname2,"%s_r2_%i.fits",fname2,selected_frame);
	sprintf(fname3,"%s_r2-r1_%i.fits",fname3,selected_frame);
	printf("%s  %s  %s\n",fname1,fname2,fname3);
//Selectd frame is now in the array
//Compute r2-r1 array using signed integers
	for(i=0;i<pixels;i++)
		{
		r21_buffer[i]=(int) r2_buffer[i]- (int) r1_buffer[i];
		}
//Allocate memory for fits array
	array[0] = (int *)malloc(naxes[0]*naxes[1]*sizeof(int));
//Initialize pointers to start of each row of the image
	for(i=1;i<naxes[1];i++)
		array[i]=array[i-1] + naxes[0];
// Load fits array and unscramble as needed
	for (j=0;j<naxes[1];j++)
		{
		for(i=0;i < naxes[0]; i++)
			{
				array[j][i]=(int) r1_buffer[i+256*j];
			}	
		}
// Here we write the r2-r1 image to disk
	fitsfile *fptr;
	status = 0;
	remove(fname1);
	remove(fname2);
	remove(fname3);
	if(fits_create_file(&fptr,fname1,&status)) printerror(status);
	if(fits_create_img (fptr, bitpix,naxis,naxes,&status)) printerror(status);
	fpixel=1;
	nelements = naxes[0]*naxes[1];
	if(fits_write_img(fptr,TINT,fpixel,nelements,array[0],&status)) printerror(status);
	if(fits_close_file(fptr,&status)) printerror(status);
	for (j=0;j<naxes[1];j++)
		{
		for(i=0;i < naxes[0]; i++)
			{
				array[j][i]=(int) r2_buffer[i+256*j];
			}	
		}
	if(fits_create_file(&fptr,fname2,&status)) printerror(status);
	if(fits_create_img (fptr, bitpix,naxis,naxes,&status)) printerror(status);
	fpixel=1;
	nelements = naxes[0]*naxes[1];
	if(fits_write_img(fptr,TINT,fpixel,nelements,array[0],&status)) printerror(status);
	if(fits_close_file(fptr,&status)) printerror(status);
	for (j=0;j<naxes[1];j++)
		{
		for(i=0;i < naxes[0]; i++)
			{
				array[j][i]=r2_buffer[i+256*j]-r1_buffer[i+256*j];
			}	
		}
	if(fits_create_file(&fptr,fname3,&status)) printerror(status);
	if(fits_create_img (fptr, bitpix,naxis,naxes,&status)) printerror(status);
	fpixel=1;
	nelements = naxes[0]*naxes[1];
	if(fits_write_img(fptr,TINT,fpixel,nelements,array[0],&status)) printerror(status);
	if(fits_close_file(fptr,&status)) printerror(status);
	free(array[0]);
	return 0;
}
void printerror( int status)
{
    /*****************************************************/
    /* Print out cfitsio error messages and exit program */
    /*****************************************************/
    if (status)
    {
       fits_report_error(stderr, status); /* print error report */

       exit( status );    /* terminate the program, returning error status */
    }
    return;
}