/*
 * dcdemo.c		16-07-87		Jim Piper
 *
 * (C) Medical Research Council 1988
 *
 * dcdemo1 [-p] [-q <QUANTUM>] [-P] [-B] [-M] [-R] [-C n] [-X n]
 *    [-D dist] [-d dist] [-E dist ] [-m PC1 PC2] [-F n f1 .. fn] [file]
 *
 * FLAGS:
 * -----
 * -p		turn on printing
 * -q quantum	quantum (factor) for minimum detection in profile
 * -P		extract profile generated candidates
 * -B		extract boundary generated candidates
 * -M		refine position of boundary candidate crossing profile
 *		if there is an nearby profile local minimum.
 *
 * -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 10, 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 6.
 * -R		Incorporate a final review stage of the positives, to reduce the
 *		false positives.
 * -m PC1 PC2	setbox classifier percentages (see write-up).
 * -F n f1 ... fn select n features (numbered from one, as in dicen.h)
 * -X n		box classifier method (default 3)
 */

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


int QUANTUM, PROF_CANDS, BOUND_CANDS;
int DISPLAY,PRINT,BPOLMOD;
int POSPRINT = 0;
int CANDSEL = 1;
int OTHNEARDIST = 10;
int SELNEARDIST = 8;
int ENDDIST = 6;
int BOXMETHOD = 3;
int REVIEWP = 0;
int Cbtype[MAXOBJS];
static struct pframe f;
int PC1 = 80, PC2 = 140;
/* which features (numbering from one, as in dicen.h) */
int nfeat = 4;
int selfeat[MAXFEATURES] = { 9,13,18,19};
int selllim[MAXFEATURES] = { 1, 1, 1, 1};

main(argc,argv)
char **argv;
{
	int i,objnum,number,nonum;
	int scale,cellw,celll;
	int numberused, cellkol1,cellline1,celllastkl,celllastln;
	struct chromosome *objlist[MAXOBJS], *obj, *readobj();
	struct chromcands *chc, *chcc;
	struct cand *cnd;
	char s[512], *cell;
	FILE *infile, *fin;
	int FIP=1;

#ifdef MIZAR
	setmemalloc(48*1024);
	fprintf(stderr,"Demonstration dicentric finding program\n");
#endif MIZAR
	chc = (struct chromcands *)Calloc(sizeof(struct chromcands),MAXOBJS);
	chcc = (struct chromcands *)Calloc(sizeof(struct chromcands),MAXOBJS);
	/*
	 * defaults
	 */
	DISPLAY = PRINT = BPOLMOD = PROF_CANDS = BOUND_CANDS = 0;
	QUANTUM = 12;
	scale = -1;
	nonum = 1;
#ifndef MIZAR
	setout("stosun -p270,0 -d1 -r870,870",'p');
#endif MIZAR

	while (argc > 1 && argv[1][0] == '-') {
		switch (argv[1][1]) {
		case 'p':
			PRINT++;
			break;
		case 'q':
			QUANTUM = atoi(argv[2]);
			argc--;
			argv++;
			break;
		case 'B':
			BOUND_CANDS++;
			break;
		case 'P':
			PROF_CANDS++;
			break;
		case 'M':
			BPOLMOD++;
			break;
		case 'R':
			REVIEWP++;
			break;
		case 'm':
			PC1 = atoi(argv[2]);
			PC2 = atoi(argv[3]);
			argc -= 2;
			argv += 2;
			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;
		case 'F':
			nfeat = atoi(argv[2]);
			argc--;
			argv++;
			for (i=0; i<nfeat; i++) {
				selfeat[i] = atoi(argv[2]);
				argc--;
				argv++;
			}
			break;
		case 'X':
			BOXMETHOD = atoi(argv[2]);
			argc--;
			argv++;
			if (BOXMETHOD < 0 || BOXMETHOD > 7) {
				fprintf(stderr,"Unknown box method\n");
				exit(1);
			}
			break;
		default:
			fprintf(stderr,"unrecogised 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;
	}
	if (PROF_CANDS == 0 && BOUND_CANDS == 0)
		PROF_CANDS = BOUND_CANDS = 1;

	number = 0;
	while((obj = readobj(infile)) != NULL) {
		if (obj->type == 1 && obj->plist == NULL) {
			fprintf(stderr,"NO PROPERTY LISTS\n");
			exit(1);
		}
#ifndef MIZAR
		if (obj->type == 1)
			swabplist(obj->plist);
#endif MIZAR
		if (obj->type == 1 && obj->plist->otype <= 1) {
			objlist[number++] = obj;
			if (number == 1) {
				cellline1 = obj->idom->line1;
				cellkol1 = obj->idom->kol1;
				celllastln = obj->idom->lastln;
				celllastkl = obj->idom->lastkl;
			} else {
				cellline1 = min(cellline1,obj->idom->line1);
				cellkol1 = min(cellkol1,obj->idom->kol1);
				celllastln = max(celllastln,obj->idom->lastln);
				celllastkl = max(celllastkl,obj->idom->lastkl);
			}
		}
	}
#ifndef MIZAR
	normaliseset(objlist,number);
#endif MIZAR

	f.type = 60;
	f.dx = 3000;
	f.dy = 2200;
	f.ix = f.iy = 1;
	f.ox = (cellkol1 + celllastkl)/2;
	f.oy = (cellline1 + celllastln)/2;

	if (scale < 0) {
		celll = celllastln - cellline1 + 1;
		if (FIP)
			celll = 3*celll/4;
		cellw = celllastkl - cellkol1 + 1;
		f.scale = 3840/max(celll,cellw);
		if (f.scale > 8)
		f.scale = 8;
	}
	if (FIP) {
		f.ix = f.scale;
		f.iy = (3*f.scale+2)/4;
		f.scale = 1;
	}

	if (number > 0) {
#ifdef MIZAR
		ddgsopen(4,0);		/* VIRTUOSO-INVERT */
/*		ddgsopen(2,0);		/* VME512 */
#endif MIZAR
		dinit();
		clear();
		lwidth(2);
		colourlut();
		colour(248);
		for (objnum=1; objnum<=number; objnum++)
			picframe(objlist[objnum-1],&f);
		tsize(50,60);
		
		/*
		 * measure features
		 */
		colour(2);
		dispstring("Profile candidate metacentric centromeres",32,4400);
		colour(4);
		dispstring("Boundary candidate metacentric centromeres",32,4300);
		numberused = dccfmain(objlist,number,chc);
		update();
#ifdef SUN
		fin = fopen("/dev/tty","r");
		fprintf(stderr,"Type return to continue\n");
		getc(fin);
#endif SUN
		/*
		 * classify
		 */
		if (BOUND_CANDS != 1 || PROF_CANDS != 1) {
			if (PROF_CANDS)
				CANDSEL = 1;
			else if (BOUND_CANDS)
				CANDSEL = 2;
		}
		colour(1);
		dispstring("Selected centromeres",32,4200);
		dispstring("Selected centromeres",40,4200);
		dcclassmain(chcc,chc,numberused);
#ifdef SUN
		fprintf(stderr,"Type return to exit\n");
		getc(fin);
		setrgb(0,255,255,180);	/* restore LUT */
		update();
#endif SUN
		finish();
	}
}


max(i,j)
{
	return ( i>j? i: j);
}


min(i,j)
{
	return ( i<j? i: j);
}


#define SCALE	32
#define FIP 1			/* FIP digitisation assumed */


dccfmain(objl,totalnumber,chc)
struct chromcands *chc;
struct chromosome **objl;
{
	struct chromcands *chcp = chc;
	int i, number, numberused, sernum;
	FILE *pf;

	pf = stdout;

	numberused = number = sernum = 0;
	while (sernum < totalnumber)
		if (dcfm(objl[sernum++],chcp,pf,&f,BOUND_CANDS,PROF_CANDS,BPOLMOD,QUANTUM,FIP,PRINT)) {
			chcp++;
			numberused++;
		}

	/*
	 * print summary about candidates
	 */
	dccandsummary(chc,numberused,pf);
	if (PRINT) {
		chcp = chc;
		for (i=0; i<numberused; i++)
			dcprint(chcp++,pf);
	}

	/*
	 * write to disc and exit
	 */
	return(numberused);
}


dcclassmain(chc,chco,numberused)
struct chromcands *chc;
{
	int i;
	FILE *pf, *infile;
	struct dcresults r;
	struct cand *cnd = (struct cand *)Malloc(300*sizeof(struct cand));
	struct chromcands *chcp;

	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;


	pf = stdout;
	if (PRINT >= 1) {
		fprintf(pf,"FEATURES USED BY BOX CLASSIFIER: ");
		for (i=0; i<nfeat; i++)
			fprintf(pf,"%2d ",selfeat[i]);
		fprintf(pf,"\n\n");
	}

	/*
	 * select required data for classification
	 */
	numberused = dcselcand(chc,cnd,chco,numberused,CANDSEL,OTHNEARDIST,SELNEARDIST,ENDDIST,&r);
	if (PRINT) {
		fprintf(stderr,"%d objects\n",numberused);
		chcp = chc;
		for (i=0; i<numberused; i++)
			dcprint(chcp++,pf);
	}

	/*
	 * normalise features
	 */
	dcnormalise(chc,numberused);

	/*
	 * do the classification
	 */
	dcboxclass(CANDSEL,BOXMETHOD,chc,numberused,nfeat,selfeat,selllim,PC1,PC2,PRINT,pf);

	/*
	 * review positives
	 */
	if (REVIEWP)
		dcreview(chc,numberused);

	/*
	 * display the found true and false positives
	 * and the false negatives
	 */
	dcselcanddisp(chc,numberused,&f);
	dcresaccum(chc,numberused,pf,&r,POSPRINT);
}