
	/* go through all the files on the command line (which should be
	   standard-star .pht files) and read in information on airmass,
	   time, filter, etc.  then, from the part of the file after
	   the header, read in all the information on each star into
	   a structure.

	   pick the 'best' aperture for each image, and note the
	   instrumental mag and estimated uncertainty through that aperture.

	   after that, decide on a standard coordinate system and transfer
	   all stars to it, shifting them by appropriate amounts so that
	   they all match up.

	   finally, write out information for all 'acceptable' images, and
	   all 'acceptable' stars, to stdout (or a file, as shown by 
	   a command-line switch), which serves as a 'master list' of all
	   stars in the entire set of images.  */
	


#include <stdio.h>
#include <math.h>
#include <string.h>
#include "pf.h"
#include "multipht.h"

#define SAMEPOS   0.10		/* if two stars have both x- and y-centroids */
                      		/*   that differ by less than this much, */
                     		/*   consider them duplicate entries for the */
                      		/*   same actual star */ 


	/* definitions of static functions in this file */
static int do_listfile( /* char *file, int aper */ );
static int do_pht( /* char *file, int aper */ );
static int slurp_header( /* char *file, struct s_star *star */ );
static int slurp_pht_data( /* char *file, struct s_star *star */ );
static int read_pht_line( /* char *line, s_pht *pht, int *naper */ );
static void print_star( /* struct s_star *s */ );
static void proc_pht( /* int flags, char *template */ );
static void do_print( );
static void sortem( /* int flag */ );
static int comp_name( /* struct s_star *s1, *s2 */ );
static int comp_time( /* struct s_star *s1, *s2 */ );
static int num_in_file( /* char *file */ );
static void check_it( /* */ );
static int identical( /* struct s_image *im, int n */ );


char progname[] = "multipht";
static char usage[] =    "multipht [debug] [print] [auto] [aper= ] [list= ] [outfile=] \n           [template= ] file.pht ... ";

main(argc, argv)
int argc;
char *argv[];
{
	int i, flags, aper, numargfile, *argfile, listflag;
	char file[100], outfile[100], template[100];
	struct s_star *s;

	aper = DEF_APER;
	flags = 0;
	debug_flag = 0;
	listflag = 0;
	num_images = 0;
	numargfile = 0;
	strcpy(outfile, "");
	strcpy(template, "");

	if (argc == 1) {
		fprintf(stderr, "usage: %s\n", usage);
		exit(0);
	}

	/* allocate enough pointers for all the images we'll read */
	if (argc >= MAX_IMAGE) {
		fprintf(stderr, "number of args >= MAX_IMAGE = %d\n", MAX_IMAGE);
		exit(-1);
	}
	if ((image_arr = (struct s_image *) calloc(MAX_IMAGE, sizeof(struct s_image))) == NULL) {
			fprintf(stderr, "can't calloc for image array \n");
			exit(-1);
	}
	if ((argfile = (int *) calloc(argc, sizeof(int))) == NULL) {
		fprintf(stderr, "can't calloc for argfile list?! \n");
		exit(-1);
	}
	
	/* read in data from all the files */
	for (i = 1; i < argc; i++) {
		if (strcmp(argv[i], "help") == 0) {
			fprintf(stderr, "usage: %s\n", usage);
			exit(0);
		}
		if (strcmp(argv[i], "debug") == 0) {
			debug_flag = 1;
			continue;
		}
		if (strcmp(argv[i], "print") == 0) {
			flags |= PRINT;
			continue;
		}
		if (strcmp(argv[i], "auto") == 0) {
			flags |= AUTO;
			continue;
		}
		if (strncmp(argv[i], "aper=", 5) == 0) {
			if (sscanf(argv[i] + 5, "%d", &aper) != 1) {
				fprintf(stderr, "usage: %s\n", usage);
				exit(1);
			}
			aper--;
			if ((aper < 0) || (aper >= (MAX_APER - 1))) {
				fprintf(stderr, "invalid aperture out of range [1-%d]\n",
								MAX_APER);
				exit(-1);
			}
			continue;
		}
		if (strncmp(argv[i], "outfile=", 8) == 0) {
			if (sscanf(argv[i] + 8, "%s", outfile) != 1) {
				fprintf(stderr, "usage: %s\n", usage);
				exit(1);
			}
			continue;
		}
		if (strncmp(argv[i], "template=", 9) == 0) {
			if (sscanf(argv[i] + 9, "%s", template) != 1) {
				fprintf(stderr, "usage: %s\n", usage);
				exit(1);
			}
			continue;
		}
		if (strncmp(argv[i], "list=", 5) == 0) {
			continue;
		}
		/* assume that this argument is a .pht filename */
		argfile[numargfile++] = i;
	}

	/* take care of command-line .pht files */
	for (i = 0; i < numargfile; i++) {
		strcpy(file, argv[argfile[i]]);
		if (strncmp(file + (strlen(file) - 3), "pht", 3) != 0)
			strcat(file, ".pht");
		do_pht(file, aper);
	}
	/* and now take care of listfiles */
	for (i = 0; i < argc; i++) {
		if (strncmp(argv[i], "list=", 5) == 0) {
			if (sscanf(argv[i] + 5, "%s", file) != 1) {
				fprintf(stderr, "usage: %s\n", usage);
				exit(1);
			}
			do_listfile(file, aper);
			continue;
		}
	}

	check_it();

	/* now that we've read in all the data, time to process it */
	proc_pht(flags, template);

	/* and, finally, write the output to the named file (or stdout, by
	   default) */
	write_output(outfile, flags);

	exit(0);
}

	/* go through all the files listed in the given file, taking input
	   from each in turn */

static int
do_listfile(file, aper)
char *file;
int aper;
{
	int i, count;
	char name[LINELEN], line[LINELEN];
	FILE *fp;

	if ((fp = fopen(file, "r")) == NULL) {
		fprintf(stderr, "do_listfile: can't open file %s\n", file);
		exit(-1);
	}
	/* first, just count the number of images we'll need */
	count = 0;
	while (fgets(line, LINELEN, fp) != NULL)
		count++;
	if (count >= MAX_IMAGE) {
		fprintf(stderr, "do_listfile: %d images is more than MAX = %d\n",
					count, MAX_IMAGE);
		exit(-1);
	}

	rewind(fp);
	while (fgets(line, LINELEN, fp) != NULL) {
		if (sscanf(line, "%s", name) != 1)
			continue;
		if (strncmp(name + (strlen(name) - 3), "pht", 3) != 0)
			strcat(name, ".pht");
		do_pht(name, aper);
	}
	fclose(fp);

	return(0);
}


	/* do all the necessary things for a .pht file, adding the new star
	   into the array if possible. Use the 'aper'th aperture for all
	   measurements, if the arguement is >= 0. Otherwise, use the 
	   aperture which is indicated by FWHM on a star-by-star basis. 

	   return 0 if all went well, or -1 if there's a problem. */

static int
do_pht(file, aper)
char *file;
int aper;
{
	struct s_image *im;
	int ret;

	im = &(image_arr[num_images]);
	im->index = num_images;		/* this gets reset later */
	im->num_match = 0;
	im->include = 1;

	if ((ret = slurp_header(file, im)) == -1) {
		fprintf(stderr, "%s.do_pht: slurp_header returns -1 for file %s\n",
					progname, file);
		return(-1);
	}

	if (slurp_pht_data(file, im) != 0) {
		fprintf(stderr, "%s.do_pht: slurp_pht_data returns non-zero for file %s\n",
					progname, file);
		return(-1);
	}
	if (slurp_coo(file, im) != 0) {
		fprintf(stderr, "%s.do_pht: slurp_coo_data returns non-zero for file %s\n",
					progname, file);
	}
	pick_aper(im, aper);
	num_images++;

	return(0);
}

	/* read in the appropriate quantities from a pseudo-FITS header.
	   note that we copy only the portion of the star's name before
	   a parenthesis, thus converting the name
	                     SA 101 311 (V)
	   to
	                     SA 101 311
	   we read in filter information later, so we don't need it in the
	   name.

	   return 0 if all went well or -1 otherwise. */

static int
slurp_header(file, image)
char *file;
struct s_image *image;
{
	int i, day, month, year;
	char value[PF_LEN + 1], str[PF_LEN + 1];
	pf_handle pfp;

	if ((pfp = pf_open(file, 1, PF_EXIST)) != 0) {
		fprintf(stderr, "%s.slurp_header: can't pf_open file %s\n", 	
					progname, file);
		return(-1);
	}
	if (pf_getval(pfp, "DATID", value) == NULL) {
		if (pf_getval(pfp, "REQID", value) == NULL) {
			fprintf(stderr, "%s.slurp_header: can't read DATID or REQID from file %s\n",
					progname, file);
			pf_close(pfp);
			return(-1);
		}
	}
	if (pf_valtostr(value, image->name) != 0) {
		fprintf(stderr, "%s.slurp_header: can't valtostr DATID/REQID ..%s.. from file %s\n",
					progname, value, file);
		pf_close(pfp);
		return(-1);
	}

	if (pf_getval(pfp, "OBJECT", value) == NULL) {
		fprintf(stderr, "%s.slurp_header: can't read OBJECT from file %s\n",
					progname, file);
		pf_close(pfp);
		return(-1);
	}
	if (pf_valtostr(value, image->object) != 0) {
		fprintf(stderr, "%s.slurp_header: can't valtostr OBJECT ..%s.. from file %s\n",
					progname, value, file);
		pf_close(pfp);
		return(-1);
	}

	if (pf_getval(pfp, "FILTER", value) == NULL) {
		fprintf(stderr, "%s.slurp_header: can't read FILTER from file %s\n",
					progname, file);
		pf_close(pfp);
		return(-1);
	}
	if (pf_valtostr(value, image->filter) != 0) {
		fprintf(stderr, "%s.slurp_header: can't valtostr FILTER ..%s.. from file %s\n",
					progname, value, file);
		pf_close(pfp);
		return(-1);
	}

	if (pf_getval(pfp, "EXPTIME", value) == NULL) {
		fprintf(stderr, "%s.slurp_header: can't read EXPTIME from file %s\n",
					progname, file);
		pf_close(pfp);
		return(-1);
	}
	if (pf_valtodouble(value, &(image->exptime)) != 0) {
		fprintf(stderr, "%s.slurp_header: can't valtodouble EXPTIME ..%s.. from file %s\n",
					progname, value, file);
		pf_close(pfp);
		return(-1);
	}

	if (pf_getval(pfp, "AIRMASS", value) == NULL) {
		fprintf(stderr, "%s.slurp_header: can't read AIRMASS from file %s\n",
					progname, file);
		pf_close(pfp);
		return(-1);
	}
	if (pf_valtodouble(value, &(image->airmass)) != 0) {
		fprintf(stderr, "%s.slurp_header: can't valtodouble AIRMASS ..%s.. from file %s\n",
					progname, value, file);
		pf_close(pfp);
		return(-1);
	}

	if (pf_getval(pfp, "UT", value) == NULL) {
		fprintf(stderr, "%s.slurp_header: can't read UT from file %s\n",
					progname, file);
		pf_close(pfp);
		return(-1);
	}
	if (pf_valtostr(value, str) != 0) {
		fprintf(stderr, "%s.slurp_header: can't valtostr UT ..%s.. from file %s\n",
					progname, value, file);
		pf_close(pfp);
		return(-1);
	}
	if (read_timestr(str, &(image->ut)) != 0) {
		fprintf(stderr, "%s.slurp_header: can't read_timestr UT ..%s.. from file %s\n",
					progname, str, file);
		pf_close(pfp);
		return(-1);
	}

	if (pf_getval(pfp, "DATE-OBS", value) == NULL) {
		fprintf(stderr, "%s.slurp_header: can't read DATE-OBS from file %s\n",
					progname, file);
		pf_close(pfp);
		return(-1);
	}
	if (read_date(value, &day, &month, &year) != 0) {
		fprintf(stderr, "%s.slurp_header: can't read_date DATE-OBS ..%s.. from file %s\n",
					progname, value, file);
		pf_close(pfp);
		return(-1);
	}
	get_jd(day, month, year, image->ut, &(image->mjd));
	image->mjd = MJD(image->mjd);

	set_apertures(pfp);

	pf_close(pfp);

	if (debug_flag) {
		printf("image %4d: %10s %20s %10s %6.1lf %5.2lf %5.2lf %8.3lf \n",
			image->index, image->name, image->object, image->filter, 
			image->exptime, image->airmass, image->ut, image->mjd);
	}
	return(0);
}
	

		
	/* read in the data from a .pht file, sticking the instrumental
	   magnitudes (together with their estimated uncertainties) into
	   the given star structure. 

	   return 0 if all OK, or -1 if there's a problem. */

static int 
slurp_pht_data(file, im)
char *file;
struct s_image *im;
{
	int i, j, n, found, naper;
	char line[LINELEN];
	FILE *fp;

	if ((im->num_star = num_in_file(file)) < 1) {
		im->stars = (struct s_star *) NULL;
		return(0);
	}
	else if ((im->stars = (struct s_star *) calloc(im->num_star, 
						sizeof(struct s_star))) == NULL) {
		fprintf(stderr, "%s.slurp_pht_data: can't calloc for stars in file %s\n",
					progname, file);
		exit(-1);
	}

	/* first, skip past all lines in the file until we come to the
	   one that consists of 'END' - leave us just past that line */
	if ((fp = fopen(file, "r")) == NULL) {
		fprintf(stderr, "%s.slurp_pht_data: can't open file %s\n",
					progname, file);
		return(-1);
	}
	found = 0;
	while (fgets(line, LINELEN, fp) != NULL) {
		if (strcmp(line, "END\n") == 0) {
			found = 1;
			break;
		}
	}
	if (found == 0) {
		fclose(fp);
		fprintf(stderr, "%s.slurp_pht_data: file %s has no END\n",
					progname, file);
		return(-1);
	}

	/* now that we're here, read each line in */
	for (i = 0; i < im->num_star; i++) {
		im->stars[i].image = im;
		im->stars[i].image_id = -1;
		im->stars[i].id = -1;
		im->stars[i].xc = -1.0;
		im->stars[i].yc = -1.0;
		im->stars[i].sky = -1.0;
		im->stars[i].skysig = -1.0;
		im->stars[i].n_aper = 0;
		for (j = 0; j < MAX_APER; j++) {
			im->stars[i].mag[j] = NODATA;
			im->stars[i].err[j] = NODATA;
		}
		im->stars[i].w3 = 1;
		im->stars[i].w4 = 0.0;
		im->stars[i].saturated = 0;
	}
	n = 0;
	while (fgets(line, LINELEN, fp) != NULL) {
		if (read_pht_line(line, &(im->stars[n]), &naper) != 0) {
			fprintf(stderr, "%s.slurp_pht_data: bad data line %d\n", 
						progname, n);
		}
		if (identical(im, n)) {
			im->num_star--;
		}
		else {
			n++;
		}
	}
	fclose(fp);

	return(0);
}
		

	/* read in one line's worth of .pht data. place the	number of 
	   apertures found into the 'naper' parameter, and print a warning
	   message to stdout if there are more than MAX_APER. 

	   return 0 if all is well, or -1 if there's a problem.  */
	
static int 
read_pht_line(line, s)
char *line;
struct s_star *s;
{
	int i, pos;
	double foo, bar;

	if (sscanf(line, "%d %lf %lf %lf %lf ", &(s->image_id), &(s->xc),
			&(s->yc), &(s->sky), &(s->skysig)) != 5) {
		fprintf(stderr, "%s.read_pht_line: bad line ..%s..\n",
				progname, line);
		return(-1);
	}
	
	pos = PHT_POS;
	for (i = 0; i < MAX_APER; i++) {
		if (sscanf(line + pos, "%lf %lf ", &(s->mag[i]), &(s->err[i])) != 2) {
			if (i == 0) {
				fprintf(stderr, "%s.read_pht_line: no mag data in line ..%s..\n",
						progname, line);
				return(-1);
			}
			else {
				break;
			}
		}
		if (s->err[i] < 0) {
			s->err[i] = 0.0 - s->err[i];
			s->saturated = 1;
		}
		pos += PHT_LEN;
	}
	s->n_aper = i;

	/* check to see if there's more data - issue a warning message if there
	   is so that the user knows we're skipping some of it */
	if (sscanf(line + pos, "%lf %lf ", &foo, &bar) == 2) {
		fprintf(stderr, "%s.read_pht_line: mag/err pairs after number %d ignored\n",
				progname, MAX_APER);
	}
	return(0);
}

	/* print out the info in a star structure */

static void
print_star(s)
struct s_star *s;
{
	int i;
	struct s_image *im;

	if ((im = s->image) != NULL) {
		printf("%-10s %-3s %6.2lf %4.2lf %5.2lf ", im->name, im->filter,
				im->exptime, im->airmass, im->ut);
	}
	printf("%5.1lf %5.1lf  %1d  %2d ", s->xc, s->yc, s->saturated, s->n_aper);
	for (i = 0; i < s->n_aper; i++) {
		printf("%5.2lf %5.3lf  ", s->mag[i], s->err[i]);
	}

	printf("\n");

}

	/* process all the .pht data; use the given file as the template */

static void 
proc_pht(flags, template)
int flags;
char *template;
{
	int i, max, max_image;

	if (flags & PRINT)
		do_print();

	do_match(template);
	assign_ids();

	do_weights(flags);

}

	/* print out the data on all the stars in all images, using the 
	   NODATA value if there's no good data */

static void
do_print()
{
	int i, j;
	struct s_image *image;

	for (i = 0; i < num_images; i++) {
		printf("           image %3d: %s\n", i, image_arr[i].name);
		for (j = 0; j < image_arr[i].num_star; j++) {
			print_star(&(image_arr[i].stars[j]));
		}
	}
}

	
	/* sort stars by some key (or keys) */

static void
sortem(flag)
int flag;
{
	if (flag == BY_TIME) {
	}
	else if (flag == BY_NAME) {
	}
	else {
		fprintf(stderr, "%s.sortem: unknown sort flag %d\n", progname, flag);
	}
}
		

	/* compare two stars by name and secondarily, by filter, and 
	   tertially, by time. */

static int
comp_name(s1, s2)
struct s_star *s1, *s2;
{
	int ret;

}

	/* compare two stars by UT time */

static int
comp_time(s1, s2)
struct s_star *s1, *s2;
{
}

	/* return the number of stars in an image */

static int
num_in_file(file)
char *file;
{
	int i, found, naper;
	char line[LINELEN];
	FILE *fp;

	/* first, skip past all lines in the file until we come to the
	   one that consists of 'END' - leave us just past that line */
	if ((fp = fopen(file, "r")) == NULL) {
		fprintf(stderr, "%s.num_in_file: can't open file %s\n",
					progname, file);
		return(-1);
	}
	found = 0;
	while (fgets(line, LINELEN, fp) != NULL) {
		if (strcmp(line, "END\n") == 0) {
			found = 1;
			break;
		}
	}
	if (found == 0) {
		fclose(fp);
		fprintf(stderr, "%s.slurp_pht_data: file %s has no END\n",
					progname, file);
		return(-1);
	}
	i = 0;
	while (fgets(line, LINELEN, fp) != NULL) {
		i++;
	}
	fclose(fp);
	return(i);
}

	/* print warning messages if an image's object name and/or filter does 
	   not match the first image's */

static void
check_it()
{
	int i, j;

	for (i = 1; i < num_images; i++) {
		if (strcmp(image_arr[i].object, image_arr[0].object) != 0) {
			fprintf(stderr, "check: image %3d ..%s.. doesn't match image %3d ..%s..\n", 
					image_arr[i].index, image_arr[i].object, image_arr[0].index,
					image_arr[0].object);
			fprintf(stderr, "");
		}
		if (strcmp(image_arr[i].filter, image_arr[0].filter) != 0) {
			fprintf(stderr, "check: image %3d filter ..%s.. doesn't match image %3d ..%s..\n", 
					image_arr[i].index, image_arr[i].filter, image_arr[0].index,
					image_arr[0].filter);
			fprintf(stderr, "");
		}
	}
}

	/* return 1 if the n'th star in the given image has exactly the same
	   position as the previous star - that is, if it is an erroneous
	   second entry for the same star.  return 0 if not */

static int
identical(im, n)
struct s_image *im;
int n;
{
	int i;
	double xc, yc;

	if (n < 1)
		return(0);
	xc = fabs(im->stars[n].xc - im->stars[n - 1].xc);
	yc = fabs(im->stars[n].yc - im->stars[n - 1].yc);
	if ((xc < SAMEPOS) && (yc < SAMEPOS))
		return(1);
	return(0);
}
