/*
 Program From Progressive Alignment of Amino Acid
 Sequences and Construction of Phylogenetic Trees
 from Them, by D. F. Feng and R. F. Doolittle,
 Methods of Enzymology, v266, p368 (1996)
				2/7/1996
 Read "readme.doc" before using.
*/

#include <stdio.h>
int hist[100],igroup;

main(ac,av)
int ac;
char *av[];
{
	FILE *fd1, *fd2, *fd3, *fd4, *fd5, *fd6, *fd7, *fd8, *fd9, *fopen();
	char dummyf[500],dummyb[500];
	int i,j,k,len,newlen,left_pos,right_pos;
	
	if(ac<3)
		{
		fprintf(stderr,"\nPurpose: Analyze Clusters (nodes) from");
		fprintf(stderr," BTREE for Agreement with Initial Tree.\n");
		fprintf(stderr,"Usage: ");
		fprintf(stderr,"%s rf1 rf2 [>wf]\n",av[0]);
		fprintf(stderr,"\trf1: Output file (from SETREE.)\n");
		fprintf(stderr,"\trf2: rboot (from BTREE.)\n");
		fprintf(stderr,"\twf: Results\n\n");
		exit(0);
		}

	fd1 = fopen(av[1],"r");
	fgets(dummyf, sizeof dummyf, fd1);
	fclose(fd1);
	fd2 = fopen("xyxy","w");
	fclose(fd2);
	fd2 = fopen("xyxy","a");
/*	printf("%s\n",dummyf);
	printf("Original Branching Order:\n"); */
	len = strlen(dummyf);

	for(i=0; i<100; ++i)
		hist[i] = 0;
	j=0;
	for(i=len; i>-1;)
		{
		if(dummyf[i]==')' || dummyf[i]=='(' || dummyf[i]==',')
			dummyb[j++] = dummyf[i--];
		else if(dummyf[i]==':'&&dummyf[i-1]!=')')
			{
			i--;
			for(k=0; k<4; ++k)
				dummyb[j++] = dummyf[i--];
			}
		else
			i--;
		}
	dummyb[j]='\0';
	len = strlen(dummyb);
	j=0;
	for(i=len-1; i>-1; i--)
		dummyf[j++] = dummyb[i];
	dummyf[j]='\0';
/*	printf("%s\n",dummyf); */

	for(;;)
		{
		len = strlen(dummyf);
	
		left_pos = find_left(len,dummyf);
		if(left_pos==0)
			break;
		right_pos = extract(fd2,left_pos,dummyf);
		j = 0;
		for(i=0; i<len;++i)
			{
			if(i!=left_pos && i!=right_pos)
				dummyb[j++] = dummyf[i];	
			}
		dummyb[j] = '\0';
		newlen = strlen(dummyb);
		strcpy(dummyf,dummyb);
		}
	fclose(fd2);
	fd3 = fopen("xyxy","r");
	fd4 = fopen("original","w");
	while((fscanf(fd3,"%s",dummyf))!=EOF)
		{
		len=strlen(dummyf);
		for(i=0; i<len; ++i)
			{
			if(dummyf[i]==',')
				fprintf(fd4," ");
			else
				fprintf(fd4,"%c",dummyf[i]);
			}
		fprintf(fd4,"\n");
		}
	fclose(fd3);
	fclose(fd4);
	
	fd5 = fopen("rboot","r");
	while(fgets(dummyf, sizeof dummyf, fd5) != NULL)
		{
		fd6 = fopen("fromboot","w");
	/*	printf("\nBoot Branching Order:\n"); */
		for(;;)
			{
			len = strlen(dummyf);
			left_pos = find_left(len,dummyf);
			if(left_pos==0)
				break;
			right_pos = extract(fd6,left_pos,dummyf);
			j = 0;
			for(i=0; i<len;++i)
				{
				if(i!=left_pos && i!=right_pos)
					dummyb[j++] = dummyf[i];	
				}
			dummyb[j] = '\0';
			newlen = strlen(dummyb);
			strcpy(dummyf,dummyb);
			}
		fclose(fd6);
	
		if((fd7 = fopen("fromboot","r")) == NULL)
			exit(fprintf(stderr,"can't open fromboot\n"));
		if((fd8 = fopen("boot","w")) == NULL)
			exit(fprintf(stderr,"can't open boot\n"));
		while((fscanf(fd7,"%s",dummyf))!=EOF)
			{
			len = strlen(dummyf);
			for(i=0; i<len; )
				{
				fprintf(fd8,"%c",dummyf[i++]);
				if(i%4==0 && i!=len)
					fprintf(fd8," ");
				}
			fprintf(fd8,"\n");
			}
		fclose(fd7);
		fclose(fd8);
		statistics();
		}
	fclose(fd5);
	fd9=fopen("original","r");
	for(i=0; i<igroup; ++i)
		{
		fgets(dummyf, sizeof dummyf, fd9);
		printf("%3d,  %s",hist[i],dummyf);
		}
	fclose(fd9);
	remove("xyxy");
	remove("original");
	remove("fromboot");
	remove("boot");
}

find_left(length,string)
int length;
char *string;
{
	int max,left,i,pos;

	max = -10;
	left = 0;
	for(i=0; i<length; ++i)
		{
		if(string[i] == '(')
			left++;
		else if(string[i] == ')')
			left--;
		if(left>max)
			{
			max = left;
			pos = i;
			}
		}
	return(pos);
}

extract(fd,left_pos,string)
FILE *fd;
int left_pos;
char *string;
{
	int i,right_pos;

	i=left_pos+1;
	while(string[i]!=')')
		fprintf(fd,"%c",string[i++]);
	fprintf(fd,"\n");
	right_pos = i;
	return(right_pos);
}

statistics()
{
	FILE *ref, *test, *fopen();
	int i,j,k,c,match,ri,ti;
	char ref_entry[52][5],test_entry[52][5];

	ref = fopen("original","r");
	test = fopen("boot","r");
	igroup = 0;

	for(i=0; i<52; ++i)
		{
		c=fgetc(ref);
		if(c==EOF)
			{
			fclose(ref);
			fclose(test);
			return;
			}
		else if(c=='\n')
			{
			igroup++;
		/*	printf("# ref_entry=%d\n",i);
			printf("\t"); */
			for(j=0; j<52; ++j)
				{
				c=fgetc(test);
				if(c==EOF)
					{
				/*	printf("igroup=%d\n",igroup); */
					fclose(test);
					test=fopen("boot","r");
					i = -1;
					goto skip;
					}
				else if(c=='\n')
					{
				/*	printf("# test_entry=%d\n\t",j); */
					match=0;
					if(i==j)
						{
						for(ri=0; ri<i; ++ri)
							for(ti=0; ti<j; ++ti)
								if(strcmp(ref_entry[ri],test_entry[ti])==0)
									match++;
						}
					if(match==i)
						hist[igroup-1]++;
					j=-1;
					goto skip1;
					}
				else
					ungetc(c,test);
				fscanf(test,"%s",test_entry[j]);
			/*	printf("%s ",test_entry[j]); */
				skip1:
					;
				}
			}
		else
			ungetc(c,ref);
		fscanf(ref,"%s",ref_entry[i]);
	/*	printf("%s ",ref_entry[i]); */
		skip:
			;
		}
	return;
}