/****************************************/
/*	mstat		valdes 8 31 82	*/
/*					*/
/*	print statistics of 		*/
/*	matched catalog			*/
/****************************************/
#include	"focas1.h"
#include	"match1.h"


/* MSTAT -- Print statistics of matched catalog */

mstat (catalog)
char	*catalog;
{
    int     i, j, objs;
    long    entries, objects;
    long    counts[MAXMATCH][MAXMATCH], pairs[MAXMATCH][MAXMATCH];
    float   avgwt[MAXMATCH][MAXMATCH], avgwt2[MAXMATCH][MAXMATCH];
    float   pavgwt[MAXMATCH][MAXMATCH], pavgwt2[MAXMATCH][MAXMATCH];
    float   awt, awt2, wt;
    struct sgparm   msp[MAXMATCH];

    cfd = mcatopen (catalog, 0);
    for (i = 0; i < ncats; i++) {
	rdcathdr (cfd, 1);
	msp[i] = sp;
    }

    entries = 0;
    objects = 0;
    for (i = 0; i < ncats; i++)
	for (j = 0; j < ncats; j++) {
	    counts[i][j] = 0;
	    avgwt[i][j] = 0.;
	    avgwt2[i][j] = 0.;
	    pairs[i][j] = 0;
	    pavgwt[i][j] = 0.;
	    pavgwt2[i][j] = 0.;
	}

    for (;;) {
	if (rdmcatob (cfd, &mentry, mob))
	    break;
	objs = charbits (mentry.entries);
	entries += 1;
	objects += objs;
	wt = mentry.wt & 0377;
	for (i = 0; i < ncats; i++) {
	    for (j = 0; j < ncats; j++) {
		if ((mob[i].entnum < 0) || (mob[j].entnum < 0))
		    continue;
		pairs[i][j]++;
		pavgwt[i][j] += wt;
		pavgwt2[i][j] += wt * wt;
	    }
	}
	for (i = 0; i < ncats; i++) {
	    if (mob[i].entnum < 0)
		continue;
	    counts[i][objs - 1]++;
	    avgwt[i][objs - 1] += wt;
	    avgwt2[i][objs - 1] += wt * wt;
	}
    }
    fclose (cfd);

    printf ("\n\n	Statistics of matched file\n\n");
    for (i = 0; i < ncats; i++)
	printf ("Catalog %d: %s   %s", i + 1, msp[i].pltnm, msp[i].dtestr);
    printf ("\nEntries = %d  Objects = %d\n\n", entries, objects);
    for (j = 0; j < ncats; j++) {
	printf ("\n		%d Matches\n\n", j);
	for (i = 0; i < ncats; i++) {
	    if (counts[i][j] == 0)
		continue;
	    awt = avgwt[i][j] / counts[i][j];
	    if (counts[i][j] > 1) {
		awt2 = (avgwt2[i][j]-counts[i][j]*awt*awt) / (counts[i][j] - 1);
		awt2 = awt2 > 0. ? sqrt (awt2) : 0.;
	    } else
		awt2 = 0.;
	    printf (
	    "Catalog %1d:  Objects = %5d  Avg wt = %5.1f  RMS wt = %5.1f\n",
		i + 1, counts[i][j], awt, awt2);
	}
    }
    for (i = 0; i < ncats; i++) {
	printf ("\n          Catalog %d: %5d objects\n\n", i + 1, pairs[i][i]);
	for (j = 0; j < ncats; j++) {
	    if (j == i)
		continue;
	    if (pairs[i][j] == 0)
		continue;
	    awt = pavgwt[i][j] / pairs[i][j];
	    if (pairs[i][j] > 1) {
		awt2 = (pavgwt2[i][j] - pairs[i][j] * awt * awt) /
		    (pairs[i][j] - 1);
		awt2 = awt2 > 0. ? sqrt (awt2) : 0.;
	    } else
		awt2 = 0.;
	    printf (
"Matches with catalog %1d:  Matches = %5d  Avg wt = %5.1f  RMS wt = %5.1f\n",
		j + 1, pairs[i][j], awt, awt2);
	}
    }
}