#include<fstream.h>
#include <stdlib.h>
#include <time.h>
#include <string.h>

#define SIZE 300
#define NSQ  50  // no of sub sequences
void Subsequence(int, int);
void splicer();
void LCS(char *, char*, char *, int, int&); 
void substring(char* ,char* , int , int );
void mergeSub(char* ,char* ,char* , char * , int);
void align(char**, int);

void main()
{	
	
	char newSeq[SIZE], nextSeq[SIZE], overlap[SIZE], mergedSeq[4*SIZE];
	int noSeq;
	int overlapSize,oldOverSize=0;
	int type;

	int i,j;
	ifstream ifp;
	ifp.open("subseq.dat");
	//ifp>>noSeq;  //same as NSQ, for test data

	char **data;
	data= new char*[NSQ];
	for(i=0; i<NSQ; i++)
	{
		data[i]= new char[SIZE+1];
	}
	
/*	// generate test data
	Subsequence(SIZE,NSQ);

	
	for( i=0; i<NSQ; i++)
	{
		ifp.getline(data[i], SIZE, '\n');
	}

	//aling sequences
	align(data, NSQ);
*/
	
	//Hello world! test cases
	LCS("hello lopjyeyvb world", "world is sfzxvbjkpo great",overlap,2, type);
	cout<<"Sequence 1: hello world\n Sequence 2: world is great\n";
	cout<<endl<<"Over lap:"<<overlap<<endl;
	mergeSub("hello world", "world is great",overlap,mergedSeq,type);
	cout<<"\nMerged Sequence is :"<<mergedSeq;

	LCS("great world", " Its a great",overlap,2, type);
	cout<<"Sequence 1: great world\n Sequence 2: Its a great\n";
	cout<<endl<<"Over lap:"<<overlap<<endl;
	mergeSub("great world", "Its a great",overlap,mergedSeq,type);
	cout<<"\nMerged Sequence is :"<<mergedSeq;


}

/**********************************************

Function to chop the sequence into subsequences

**********************************************/
void Subsequence(int seqSize, int noSubSeq)
{
	
	srand(time(NULL));

	char nucleo[4]={'A', 'T', 'C', 'G'};	
	char seq[SIZE];
	int start, endse, temp;
	
	//output data to a file
	ofstream fout;
	fout.open("subseq.dat");

	cout<<"Sequence .........................\n";

	for(int i=0; i<seqSize; i++)
	{
		seq[i]=nucleo[rand()%4];	
		cout<<seq[i];
	}


	//fout<<noSubSeq<<endl;

	for(i=0; i<noSubSeq; i++)
	{
		cout<<"\nSub Sequence "<<i<<" .........................\n";
		start=rand()%300; 
		endse=rand()%300;
		if(start > endse)
		{
			temp=start;
			start=endse;
			endse=temp;
		}
		cout<<"Start "<<start<<" End "<<endse<<endl;

		for(int j=start; j<endse; j++)	
		{ 
			cout<<seq[j];
			fout<<seq[j];
		}
		fout<<endl;
	}
	
	fout.close();

}

/**********************************************

Function to aling sequences, exact alignment 

**********************************************/	
void align(char **data, int nsq)
{
		
	char newSeq[SIZE], nextSeq[SIZE], overlap[SIZE], mergedSeq[4*SIZE];
	int noSeq;
	int overlapSize,oldOverSize=0;
	int type;

	int i,j;
	i=0;

	while(i<nsq)
	{
		 //next one onwards
		for( j=i+1; j<nsq; j++)
		{
			//find LCS
			LCS(data[i],data[j],overlap,2, type);
			
			overlapSize=strlen(overlap);
			
			if(overlapSize>0 && type>0)	//has a match, creates an empty mergedSubst if not checked, to be used later to define size of overlap  
			{
				mergeSub(data[i], data[j],overlap, mergedSeq, type);
			
				//show only new sequences
				if((strlen(mergedSeq) > strlen(data[i])) && (strlen(mergedSeq) > strlen(data[j])))
				{
					cout<<"\n..........................................................."<<endl;
					cout<<"\n\nSequence 1  Length:"<<strlen(data[i])<<"\n"<<data[i];
					cout<<"\n\nSequence 2 Length:"<<strlen(data[j])<<"\n"<<data[j];
					cout<<endl<<i<<":"<<j<<" Over lap "<<overlap<<endl;
					cout<<"\n\nMerged Sequence Length :"<<strlen(mergedSeq)<<"\n"<<mergedSeq;				
				}
	
			}
		
		
		}
		i++;
	}
}

/**********************************************

Function to put back the subsequences, using LCS -Dynamic programming

Notes:Conventions used -1 for up, -2 for left and -3 for up and left in the table

  x and y are the sub sequences to be compared, y is the overlap, style is the type of overlap
  contiguous for 2 and 1 for non-contiguous.  

**********************************************/
void LCS(char *x, char *y, char *z, int style, int &type)
{

//variables	
int i,j;

int x_l, y_l;

x_l=strlen(x);   //lenght of string x and y
y_l=strlen(y);


int **c, **b;		//allocate space
c= new int* [x_l+1];
b= new int* [x_l+1];

for(i=0; i<=x_l; i++)
{
	b[i]=new int[y_l+1];
	c[i]=new int[y_l+1];

	for(j=0; j<=y_l; j++)
		b[i][j]=0;
}


for(i=0; i<=x_l; i++)   //initialize the tables b and c
{
		c[i][0]=0;
		b[i][0]=-1;   //-1 indicates up
}


for(i=0; i<=y_l; i++)   //initialize the tables b and c
{
		c[0][i]=0;
		b[0][i]=-2;  //-2 indicates left
}


//dynamic programming, constructs the table that can be used to reconstruct the LCS

for(i=1; i<=x_l; i++)
	for(j=1; j<=y_l; j++)
	{
		if(x[i-1]==y[j-1])  //two characters in the strings match
		{
			c[i][j]=c[i-1][j-1]+1;
		
				b[i][j]=-3;           //indicates up and left
		
		}
		else if(c[i-1][j]>=c[i][j-1])
		{
			c[i][j]=c[i-1][j];
			b[i][j]=-1;
		}
		else
		{
			c[i][j]=c[i][j-1];
			b[i][j]=-2;
		}
	
	}

	//the largest no is the length of the common string
	int len=c[x_l][y_l];
	int start, stop=0, max, oldstart, oldstop, ctr, len2=0;

	int ystart, yend, jcopy;

	
	z[len]='\0';

	i=x_l; j=y_l;

//back track
	oldstop=0; stop=0; ctr=0;

while(i>0 || j>0)
{
	if(b[i][j]==-3)
	{
		i--; j--;
		len--;
		z[len]=x[i];

		
		if(ctr==0)
		{
			start=i;
			jcopy=j;
			stop=0;
			ctr=1;
		}
		else
		{

			stop++;

			if(i==0)
			{
				if(stop>oldstop)
				{
					oldstop=stop;
					oldstart=start;
					yend=jcopy;
				}

			}
		}
	}
	else if(b[i][j]==-1)
	{
		i--;
	
		if(stop>oldstop)
		{
			oldstop=stop;
			oldstart=start;
			yend=jcopy;

		}
		ctr=0;
		
	}
	else if(b[i][j]==-2)
	{
		j--;

		if(stop>oldstop)
		{
			oldstop=stop;
			oldstart=start;
			yend=jcopy;
		}	
		ctr=0;
		
	}

}

//////////////////////////Testing
cout<<"new value "<<yend<<oldstart;
int i1,j1, count, oldcount=0;
for(i=0; i<x_l; i++)
for(j=0; j<y_l; j++)
{
	count=0;
	if(b[i][j]==-3)
	{
		i1=i+1; j1=j+1;
		count=1;
		while((i1<x_l) &&(j1<y_l) &&(b[i1][j1]==-3))
		{
			if(b[i1][j1]==-3)
			{
				cout<<"=="<<endl;
				i1++; j1++;
				count++;
			}
		}
	}//if
	if(count >= oldcount)
	{
		yend=j1;  
		oldstart=i1; //xend
		oldstop=count;
		oldstop--;
	}
}
cout<<"new value "<<yend<<oldstart;
//////////////////////////Testing


//cout<<"\nLongest common Sub sequence: "<<z<<endl;

if(oldstop ==0)
{
	z[0]='\0';
	type=0;   // no match
}
else if(oldstop>0 && style==2) 
{	
	int xstart=oldstart-oldstop;
	int xend=oldstart;

	ystart=yend-oldstop;
	
	substring(z, x, xstart, oldstart); //we only want it if its more than 1 char
	//cout<<"\nLongest common Substring "<<z<<endl;

	
	//Determine if the common substring is at the ends
	if((xend ==x_l-1) &&  (ystart==0)) //suffix of x prefix of y
	{
		type=1;
	}
	else if((xstart==0) && (yend==y_l-1))  //suffix of y prefix of x
	{
		type=2;
	}
	else
		type =0; // not a suffix or prefix
}


for(i=0; i<=x_l; i++)
{
	cout<<i<<" :";
	for(j=0; j<=y_l; j++)
	{
		cout<<c[i][j]<<"*"<<b[i][j]<<" ";
	}
cout<<endl;
}


//free space

for(i=0; i<=x_l; i++)
{  
	delete[] b[i];
	delete[] c[i];
}
delete[] b;
delete[] c;



  }

/**********************************************

Function to merge two subsequences based on 
an exact match at the edges.

 
**********************************************/
void mergeSub(char* x,char* y,char* overlap, char * z, int type)
{
	int i,j;
	int len_m=strlen(overlap);
	int len_x=strlen(x);
	int len_y=strlen(y);

	if(type==1) // suffix of x prefix of y
	{
		strcpy(z,x); //copy x to new string
		for(i=len_x, j=len_m; j<len_y; i++,j++)
		{
			z[i]=y[j];
		}
		z[i]='\0';

	}
	else if (type ==2)  //suffix of y prefix of x
	{
		strcpy(z,y); //copy x to new string
		for(i=len_y, j=len_m; j<len_x; i++,j++)
		{
			z[i]=x[j];
		}
		z[i]='\0';
	}
	else if (type ==0)
	{
		z[0]='\0';
	}

}

/**********************************************

Function to put back the subsequences-Not in Use

**********************************************/
void splicer(int type)
{

	int i, j, nosub;
	char reconst[SIZE];
	char temp[SIZE];

	ifstream fin;
	fin.open("subseq.dat");
	fin>>nosub;

	switch(type)
	{
	
	case 1:					//similarity based exact matching
							// this would be the simplest and redundant/slow way of matching
		int overlap_size=24;

		for(i=0; i<nosub; i++)
		{
			for(j=0; j<nosub; j++)
			{
			//	for()
			//	{

		//		}
			}
		}
		break;

	}

}
/************************************************************/
void substring(char* y,char* x, int start, int len)
{
  int i,j;
	for(i=start, j=0; i<=len; i++,j++)
	{
		y[j]=x[i]; 
	}
	y[j]='\0';
}