/*
 * dcrule.c		Simon Towers		05-04-88
 *
 * perform candidate selection by rules
 *
 * dcclass [-C n] [-D dist] [-d dist] [-E dist ] [-p] [file]
 *
 * -C n		selects whether classification on type 1 candidates (n=1),
 *		type 2 candidates (n=2), or some other candidate selection
 *		method (to be determined).  Default 1.
 * -D dist	If flag -D is set, then
 *		candidates of the "other" type are also included if their
 *		nearest candidate of the selected type is more than
 *		"dist" away (feature neardist).  Default 8, i.e. will
 *		reject "other" candidates within 1 micron of selected candidate.
 * -d dist	If two candidates of the SAME selected type are closer than
 *		"dist", reject one.  Default 8.
 * -E dist	reject any candidate less than "dist" to medial axis end.
 *		Default 8.
 * -p		turn on printing
 */

#include <stdio.h>
#include <math.h>
#include <wstruct.h>
#include <chromanal.h>
#include <dicen.h>
#include <otypes.h>


main(argc,argv)
char **argv;
{
	struct chromcands *chcp,chc[MAXOBJS];
	struct cand *cndp,cnd[MAXCANDS];
	int i, j, number, numberused;
	char *cell,s[256];
	FILE *infile, *pf;
	struct dcresults r;
	int PRINT = 0;
	int CANDSEL = 1;
	int OTHNEARDIST = 8;
	int SELNEARDIST = 8;
	int ENDDIST = 8;

	r.TP = r.FP = r.FN = r.UTP = r.UFN = 0;
	r.compoFP = r.twistedFP = r.acroFN = 0;
	r.onecandFN = r.nchrom = r.nwithcands = 0;
	r.areaFN = r.ovFP = r.ovFN = 0;

	while (argc > 1 && argv[1][0] == '-') {
		switch (argv[1][1]) {
		case 'p':
			PRINT++;
			break;
		case 'C':
			CANDSEL = atoi(argv[2]);
			argc--;
			argv++;
			if (CANDSEL != 1 && CANDSEL != 2) {
				fprintf(stderr,"Illegal -C value\n");
				exit(1);
			}
			break;
		case 'D':
			OTHNEARDIST = atoi(argv[2]);
			argc--;
			argv++;
			break;
		case 'E':
			ENDDIST = atoi(argv[2]);
			argc--;
			argv++;
			break;
		case 'd':
			SELNEARDIST = atoi(argv[2]);
			argc--;
			argv++;
			break;
		default:
			fprintf(stderr,"unrecognised switch -%c\n",argv[1][1]);
			break;
		}
		argc--;
		argv++;
	}
	if (argc > 1) {
		cell = argv[1];
		infile = fopen(cell,"r");
		if (infile == NULL) {
			fprintf(stderr,"Can't open %s\n",cell);
			exit(1);
		}
	} else {
		cell = "<stdin";
		infile = stdin;
	}

	/*
	 * read data
	 */
	numberused = dcread(chc,cnd,infile,CANDSEL,OTHNEARDIST,SELNEARDIST,ENDDIST,&r);

	/*
	 * apply the rules for selection
	 */
	dcrule(CANDSEL, chc, numberused);

	/*
	 * display the found true and false positives
	 * and the false negatives
	 */
	if (PRINT)
		dcresprint(cell, chc, numberused);
	else
		dcresaccum(chc, numberused, stdout, &r, 0);
}


dcrule(candtype, chc, nobjs)
struct chromcands *chc;
{
	register int i;
	int l, cmerecount, compare(), sum, prf, min, mediana, areas[MAXOBJECTS];
	int SET, dw, dp, near, def, mod;
	struct chromcands *chcp;
	register struct cand *cnd;
	
	/*
	 * find median area
	 */
	chcp = chc;
	for (i=0 ; i<nobjs ; i++,chcp++)
		areas[i] = chcp->area;
	qsort(areas, nobjs, sizeof(int), compare);
	mediana = areas[nobjs/2];
	
	/*
	 * For each candidate check whether the rules are satisfied.
	 * Also check if area exceeds the cell median area, otherwise
	 * cannot be dicentric.
	 */
	chcp = chc;
	for (i=0 ; i<nobjs ; i++,chcp++) {
		cmerecount = 0;
		near = 0;
		cnd = chcp->cands;
		for (l=0 ; l<chcp->ncands ; l++,cnd++) {
			sum = cnd->u.sy.curvsum;
			prf = cnd->u.sy.nprofval;
			min = cnd->u.sy.curvmin;
			if ((sum<151 && prf<81) || (sum<156 && prf<71)
				|| (sum<161 && prf<51)) cnd->ctype = 1;
			else
				cnd->ctype = 0;
			if (cnd->u.sy.modality>1 && cnd->u.sy.ncvhdeficiency>20)
				cnd->ctype = 0;
			if (cnd->u.sy.ncvhdeficiency > 30) cnd->ctype = 0;
			if (min > (205-prf/2)) cnd->ctype = 0;
			if (min == 180) cnd->ctype = 0;
			if (cnd->ctype)	{
				cmerecount++;
				if (cnd->neardist <= 4) near++;
			}
		}
		
		/* if 1 doubly selected centromere then relax conditions on others */
		if (cmerecount==1 && near) {
			cnd = chcp->cands;
			for (l=0 ; l<chcp->ncands ; l++,cnd++) {
				if (cnd->ctype) continue;
				def = cnd->u.sy.ncvhdeficiency;
				mod = cnd->u.sy.modality;
				if (cnd->u.sy.curvsum<161 &&
					cnd->u.sy.nprofval<81 &&
					(def<21 || (def<31 && mod==1))) {
					cnd->ctype++;
					cmerecount++;
					break;
				}
			}
		}

		/* do some comparisons between centromere candidates */
		if (cmerecount >= 2) {
			SET = 0;
			cnd = chcp->cands;
			for (l=0 ; l<chcp->ncands ; l++,cnd++) {
				if (!cnd->ctype) continue;
				if (!SET) {
					SET++;
					dw = cnd->u.sy.width;
					dp = cnd->u.sy.nprofval;
				}
				else {
					if ((abs(dw-cnd->u.sy.width)<5) &&
						(abs(dp-cnd->u.sy.nprofval)<21))
						break;
					else
						cmerecount--;
				}
			}
		}

		switch (cmerecount) {
			case 0:
				chcp->btype = ACENTRIC;
				break;
			case 1:
				chcp->btype = MONOCENTRIC;
				break;
			case 2:
			default:
				if (chcp->area >= mediana)
					chcp->btype = DICENTRIC;
				else
					chcp->btype = DICENTRIC | UNSHURE;
				break;
		}
	}
}


static compare(a, b)
register int *a, *b;
{
	return(*a - *b);
}


dcresprint(cell, chc, number)
char *cell;
struct chromcands *chc;
{
	int i;
	int PRINTED = 0;

	for (i=0; i<number; i++,chc++) {
		if (chc->porigin & 022) continue;
		if (chc->btype & UNSHURE) continue;
		if (chc->Cbtype & UNSHURE) continue;
		if (chc->Cbtype & ACROCENTRIC) continue;
		if (chc->btype&DICENTRIC) {
			if (!PRINTED) {
				printf("DATA FILE %s \n", cell);
				PRINTED++;
			}
			dcprint(chc, stdout);
			if (chc->Cbtype&DICENTRIC)
				printf("TRUE DICENTRIC: %d\n",chc->number);
			else if (chc->Cbtype&CROSSARMS)
				printf("FALSE DICENTRIC (crossed arms): %d\n",chc->number);
			else if (chc->Cotype == 2 || chc->Cotype == 3)
				printf("FALSE DICENTRIC (true composite missed): %d\n",chc->number);
			else
				printf("FALSE DICENTRIC (no excuse): %d\n",chc->number);
		}
		else if (chc->Cbtype&DICENTRIC) {
			if (!PRINTED) {
				printf("DATA FILE %s \n", cell);
				PRINTED++;
			}
			dcprint(chc, stdout);
			printf("MISSED DICENTRIC (%d candidates): %d\n",chc->ncands,chc->number);
		}
	}
	if (PRINTED) printf("------------------------------------------------\n\n");
}