/*#define DEBUG*/
/*
 * xpboxclass.c
 *
 * classify feature files from x-probed chromosomes
 */

#include <stdio.h>
#include <wstruct.h>
#include <chromanal.h>
#include <otypes.h>
#include "xprobe.h"
#include <vcv.h>

main(argc,argv)
char **argv;
{
	int i,j;
	struct xfeats xf[200],*xfp;
	FILE *in;
	int number,xfpcomp(),vcvcomp(),rank,max2d,xrank;
	char *fname = "stdin";
	
	if (argc > 1) {
		fname = argv[1];
		in = fopen(fname,"r");
		if (in == NULL) {
			fprintf(stderr,"Cannot open file %s\n",fname);
			exit(2);
		}
	}
	else
		in = stdin;
	
	/*
	 * read feature data file
	 */
	xfp = xf;
	max2d = 0;
	while (fscanf(in,"%d",&xfp->number) != EOF) {
		fscanf(in,"%d",&xfp->Cbtype);
		fscanf(in,"%d",&xfp->relarea);
		fscanf(in,"%d",&xfp->cindexa);
		fscanf(in,"%d",&xfp->pindexa);
		fscanf(in,"%d",&xfp->difpos);
		fscanf(in,"%d",&xfp->profpeak);
		fscanf(in,"%d",&xfp->deriv2);
		fscanf(in,"%d",&xfp->objpeak);
		fscanf(in,"%d",&xfp->objmedian);
		xfp->profpeak -= xfp->objmedian;
		xfp->fname = fname;
isitfn(xfp);
		if (xfp->relarea > TOOSMALL && xfp->relarea < TOOBIG) {
			if (xfp->deriv2 > max2d)
				max2d = xfp->deriv2;
			number++;
			xfp++;
		}
	}
fprintf(stderr,"-----\n");

	/*
	 * massage features as necessary
	 */
	xfp = xf;
	for (i=0; i<number; i++,xfp++)
		xfp->deriv2 = 100 * xfp->deriv2 / max2d;

	/*
	 * sort by density profile peak, code rank.
	 */
	qsort(xf,number,sizeof(struct xfeats),xfpcomp);

	/*
	 * Classification:
	 * (1) First MAXRANK only of those which satisfy wide pindexa bounds
	 * (2) relarea limits
	 * (3) cindexa limits, or pindexa limits if cindexa == -1
	 * (4) difpos limit if cindexa != -1
	 */
	xfp = xf;
	for (i=1, rank=1, xrank=1; xrank<=XRANK && i<=number && rank<=MAXRANK; i++,rank++,xfp++) {
		xfp->btype = 2;
		if (xfp->pindexa < PILB || xfp->pindexa > PIUB) {
			rank--;
			continue;
		}
		if (xfp->relarea < AREAMIN || xfp->relarea > AREAMAX ||
		    xfp->deriv2 < D2MIN) {
			isitfn(xfp);
			continue;
		}
		if (xfp->cindexa != -1) {
			if (xfp->cindexa < CIMIN) {
				isitfn(xfp);
				continue;
			}
			if (xfp->cindexa > CIMAX) {
				if (xfp->pindexa > PIMAX || xfp->difpos < 0 ||
				    xfp->difpos > DIFF) {
					isitfn(xfp);
					continue;
				}
			}
			if (xfp->difpos < -DIFF || xfp->difpos > DIFF) {
				isitfn(xfp);
				continue;
			}
		}
		else if (xfp->pindexa < PIMIN || xfp->pindexa > PIMAX) {
			isitfn(xfp);
			continue;
		}
		/*
		 * Potential probed X found
		 */
		xfp->xrank = xrank++;
		xfp->btype = 1;
		xfp->rank = rank;
		writexpboxfeat(stdout,xfp);
			isitfn(xfp);
	}
	for (; i<=number; i++,xfp++)
		isitfn(xfp);
}

int xfpcomp(x1,x2)
register struct xfeats *x1,*x2;
{
	return(x2->profpeak - x1->profpeak);
}

isitfn(x)
struct xfeats *x;
{
	if (x->Cbtype == XPROBE)
		writexpboxfeat(stderr,x);
}