/**************************************************************************
 *                                                                        *
 * Copyright (c) 1996 Michael Richmond and Richard Treffers               *
 *                                                                        *
 *                    This software may be copied and distributed for     *
 *                    educational, research and not for profit services   *
 *                    provided that this copyright and statement are      *
 *                    included in all such copies.                        *
 *                                                                        *
 **************************************************************************/


/*
	Bins data files (i.e. turns 400x400 -> 200x200 image)

		BIN file bin=  [norm]

	BSCALE, CDELT1 and CDELT2 are adjusted in the header of the 
	binned file.

	4/19/92 - changed data type for sun 
			- put exit(0) and added outfile=  -rrt
*/

#include "pcvista.h"
#include "fits.h"
#include <stdio.h>
#include <math.h>
#include <unistd.h>
#include <stdlib.h>
#include <string.h>

#ifdef PROTO
int main(int, char **);
static void symbols(FITS_HANDLE, FITS_HANDLE, int, int);
#else
static void symbols();
#endif

#define USAGE "usage: bin file [bin=] [norm] [outfile=] [help]"

int
main(argc,argv)
int argc;
char *argv[];
{
	int i, j, k, l, m, binsize = 2, row;
	int nrow, ncol, nr, nc, scale = 1;
	FITS_HANDLE fin, fout;
	static int16 data[NMAX];
	static int norm_flag,new_file;
	static long sum[NMAX];
	char *gotit, buf[NBUF],outfile[NBUF];

	if (argc == 1) {
		error(-1, USAGE);
	}
	if (strcmp(argv[1], "help") == 0) {
		printf("%s\n", USAGE);
		exit(0);
	}

	fin = fits_open(argv[1], "r", &nrow, &ncol);
	strcpy(outfile,argv[1]);
	
	for (i = 2; i < argc; i++) {
		if (keyword("help", argv[i])) {
			printf("%s\n", USAGE);
			exit(0);
		}
		if ((gotit = find("bin", argv[i])) != NULL) {
			binsize = (int)evaluate(gotit);
			continue;
		}
		if ((gotit = find("outfile", argv[i])) != NULL) {
			strcpy(outfile,gotit);
			new_file=1;
			continue;
		}
		if (keyword("norm", argv[i])) {
			norm_flag = 1;
			continue;
		}
		sprintf(buf, "Unknown option '%s'", argv[i]);
		error(-1, buf);
	}
	if (binsize <= 1)
		error(-1, "binsize must be greater than 1");

	nr = nrow/binsize;
	nc = ncol/binsize;	

	fout = fits_open(outfile, new_file?"w":"x", &nr, &nc);

	if (norm_flag)
		scale = binsize*binsize;
			
	for (m = 0, row = 0; m < nr; m++) {
		for (j = 0; j < binsize; j++) {
			fits_get_data(fin, row++, 0, data, ncol);
			for (i = 0, k = 0; i < nc; i++ ) {
				for (l = 0; l < binsize; l++, k++)
					sum[i] += data[k];		
			}
		}
		for (i = 0; i < nc; i++) {
			data[i] = (int16)(sum[i]/scale);
			sum[i] = 0;
		}
		fits_put_data(fout, m, 0, data, nc);
	}
	if(!new_file){
		symbols(fin, fout, binsize, norm_flag);
		fits_cut(fout);		/* truncate file to present length */
	}
	fits_close(fout);
	exit (0);
}


static void 
symbols(fin, fout, binsize, norm_flag)		/* modify FITS header as needed */
FITS_HANDLE fin, fout;
int binsize, norm_flag;
{
	char buf[NBUF];
	double delta1, delta2, scale;

	if (fits_get_symbol(fin, "CDELT1", buf) == FITS_PASS) {
		delta1 = atof(buf);
		if (fits_get_symbol(fin, "CDELT2", buf) == FITS_FAIL) {
			error(0, "CDELT2 is missing,  CDELT1 substituted");
			delta2 = delta1;
		}
		else
			delta2 = atof(buf);
		sprintf(buf, "%20.6E", binsize*delta1);	
		fits_put_symbol(fout, "CDELT1", buf);
		sprintf(buf, "%20.6E", binsize*delta2);	
		fits_put_symbol(fout, "CDELT2", buf);
	}
	if (norm_flag && (fits_get_symbol(fin, "BSCALE", buf) == FITS_PASS)) {
		scale = atof(buf);
		sprintf(buf, "%20.6E", binsize*binsize*scale);
		fits_put_symbol(fout, "BSCALE", buf);
	}
}