/*
 * xmatch.c
 *
 * match auto- and manually-split chromosomes
 */
#include <stdio.h>
#include <wstruct.h>
#include <chromanal.h>

main(argc,argv)
char **argv;
{
	struct chromosome *c1[1000],*c2[1000],*o, *readobj();
	FILE *f1, *f2, *f3;
	int i,j,k,l;
	int n1, n2, nm, nj;
	int nmatch,tmatch;

	f1 = fopen(argv[1],"r");
	f2 = fopen(argv[2],"r");
	if (f1 == NULL || f2 == NULL)
		exit(1);
	if (argc > 3)
		f3 = fopen(argv[3],"w");
	else
		f3 = NULL;

printf("Matching %s with %s",argv[1],argv[2]);
if (argc > 3)
	printf(", output to %s",argv[3]);
printf("\n");
	n1 = n2 = 0;

	while ((o = readobj(f1)) != NULL) {
		swabplist(o->plist);
		if (o->plist->Cotype <= CHROMOSOME)
			c1[n1++] = o;
	}
	while ((o = readobj(f2)) != NULL) {
		swabplist(o->plist);
		if (o->plist->Cotype <= CHROMOSOME)
			c2[n2++] = o;
	}

printf("%d, %d chromosomes read\n",n1,n2);
	for (i=0; i<n1; i++) {
		c1[i]->plist->area = area(c1[i]);
		c1[i]->plist->mass = mass(c1[i]);
	}
	for (i=0; i<n2; i++) {
		c2[i]->plist->area = area(c2[i]);
		c2[i]->plist->mass = mass(c2[i]);
	}
	tmatch = 0;
	for (i=0; i<n1; i++) {
		nmatch = 0;
		for (j=0; j<n2; j++) {
			nm = match(c1[i],c2[j],i+1,j+1);
			nmatch += nm;
			if (nm == 1)
				nj = j;
		}
		if (nmatch > 1)
			printf("Chr %d - %d matches\n",i+1,nmatch);
		else if (nmatch == 1) {
			tmatch++;
			c1[i]->plist->Cbtype = c2[nj]->plist->Cbtype;
		}
	}
printf("---------%d matches----------------------\n",tmatch);
	if (f3 != NULL)
		for (i=0; i<n1; i++)
			writeobj(f3,c1[i]);
}

match(o1,o2,n1,n2)
struct chromosome *o1, *o2;
{
	struct chromplist *p1, *p2;
	struct intervaldomain *i1, *i2;
	int a1,a2,d;

	i1 = o1->idom;
	i2 = o2->idom;
	p1 = o1->plist;
	p2 = o2->plist;

	d = i1->line1 - i2->line1;
	if (d < -10 || d > 10)
		return(0);
	d = i1->lastln - i2->lastln;
	if (d < -10 || d > 10)
		return(0);
	d = i1->kol1 - i2->kol1;
	if (d < -10 || d > 10)
		return(0);
	d = i1->kol1 - i2->kol1;
	if (d < -10 || d > 10)
		return(0);

	a1 = p1->area;
	a2 = p2->area;
	if (a1 < 50*a2/100 || a2 < 50*a1/100)
		return(0);
/*
	a1 = mass(o1);
	a2 = mass(o2);
	if (a1 < 50*a2/100 || a2 < 50*a1/100)
		return(0);
*/

/*	printf("********** %d matches %d\n",n1,n2);	*/
	return(1);
}