import java.util.*;
import java.io.*;
import org.biojava.bio.symbol.*;
import org.biojava.bio.seq.*;
import org.biojava.bio.seq.io.*;

/** 
 * Class that manages the terminator search.
 * @author Gard Thomassen
 * @see TerminatorSearchManager
 * @see ProgramManager
 * @see TerminatorCandidate
 */


public class TerminatorSearcher{

    Sequence seq;
    U_residue u_res;
    private Vector tmpCand = new Vector(50,10);
    String ungappedStem;
      

    /** 
     * This method initiates the first part of the terminator
     * search by calling a method that searches
     * the given string for U residues. 
     * @param s the sequence to searchsee TerminatorCandidate
     * @return returns a vector with only the best termiantor candidates
    */

    public Vector startTermSearch(Sequence s){
	seq = s;
	searchUresidues();
	return(clean());
    }
    
    /** 
     * This method searches a string for three consecutive U's
     */
    protected void searchUresidues(){
	String tmp;
		
	if(seq.length()>=14){
	    for(int i = 14; i < (seq.length()-1); i++ ) {
		tmp = seq.subStr(i , (i+2));
		tmp = tmp.toLowerCase();
		
		 if(tmp.charAt(0) == 't' &&
		   tmp.charAt(1) == 't' &&
		   tmp.charAt(2) == 't'){
		     u_res = new U_residue(i , (i + 2));
		     ungappedPermutation(seq);
		     leftGappedPermutation(seq);
		     rightGappedPermutation(seq);
		 }
	    }
	}
    }
        
    /** 
     * This method searches the string upstream from the three consecutive U's 
     * for a possible hairpin structure. One gap on the right stem is allowed in this search.
     */


    protected void rightGappedPermutation(Sequence s){
	Sequence seq = s;
	int u_start = u_res.getStart();
	int u_stopp = u_res.getStop();
	int length = 0;
	String stem1, loop, stem2, tmp, tmp2;
	String u_sequence;
	boolean shortSeq = false;


	if(u_start < 56){
	    shortSeq = true;
	}
	
	//iteration of all possible gapped terminators.
	if(!shortSeq){
	    int j = 0;
	    int k = 0;
	    int l = 0;
	    int m = 0;
	    	    
	    //variating length between U-residue and right stem
	    for(int i=0; i < 2; i++){
		tmp = seq.subStr((u_start-50-i-1), (u_start-i-1));
		int j_decreasement = 0;
		
		//setting the 15 nucleotide U residue
		if(seq.length() >= (u_start+14))
		    u_sequence = seq.subStr(u_start-i, (u_start + 14-i));
		else
		    u_sequence = seq.subStr(u_start-i, seq.length());
				
		//gapping on right stem	
		for(j=21; j>4; j--){
		    for(m=1; m<(20-j_decreasement); m++){
			stem2 = tmp.substring((tmp.length()-j),(tmp.length()-j+m));
			stem2 = stem2.concat(tmp.substring((tmp.length()-j+m+1),(tmp.length())));
			ungappedStem = tmp.substring((tmp.length()-j),tmp.length());
			
			for(k=11; k>3; k--){
			    loop = tmp.substring((tmp.length()-stem2.length()-k),
						 (tmp.length()-stem2.length()-1));
			    stem1 =tmp.substring((tmp.length()-j-loop.length()-j+1),
						 (tmp.length()-j-loop.length())); 
			     
			    scoreTerm(stem1, loop, stem2, 'R', u_sequence, i, u_start, ungappedStem); 
			}
		    }
		    //to keep track of decreasement og stemlength
		    j_decreasement++;
		}
	    }
	}
    }
    /** 
     * This method searches the string upstream from the three consecutive U's 
     * for a possible hairpin structure. One gap on the left stem is allowed in this search.
     */
    protected void leftGappedPermutation(Sequence s){
	Sequence seq = s;
	int u_start = u_res.getStart();
	int u_stopp = u_res.getStop();
	int length = 0;
	String stem1, loop, stem2, tmp, tmp2;
	String u_sequence;
	boolean shortSeq = false;
		
	if(u_start < 56){
	    shortSeq = true;
	}
	
	//iteration of all possible gapped
	if(!shortSeq){
	    int j = 0;
	    int k = 0;
	    int l = 0;
	    int m = 0;
		    
	    //variating length between U-residue and right stem
	    for(int i = 0; i < 2; i++){
		tmp = seq.subStr((u_start-50-i-1), (u_start-i-1));
			
		//setting the 15 nucleotide U residue
		if(seq.length() >= (u_start+14))
		    u_sequence = seq.subStr(u_start-i, (u_start + 14 -i));
		else
		    u_sequence = seq.subStr(u_start-i, seq.length());
				
		//int k_decreasement;
		int j_decreasement = 0;
		//gapping on left stem		
		for(j=20; j>4; j--){
		    stem2 = tmp.substring((tmp.length()-j),(tmp.length()));
		  		    
		    for(k=11; k>3; k--){
			loop = tmp.substring((tmp.length()-stem2.length()-k+1),
					     (tmp.length()-stem2.length()));
			
			for(l=1; l<(20-j_decreasement); l++){
			    stem1 =tmp.substring((tmp.length()-j-k-j),(tmp.length()-j-k-j+l));			    
			    stem1 = stem1.concat(tmp.substring((tmp.length()-j-k-j+l+1),(tmp.length()-j-k+1)));
			    ungappedStem = tmp.substring((tmp.length()-j-k-j),(tmp.length()-j-k+1));
			    
			    scoreTerm(stem1, loop, stem2, 'L', u_sequence, i , u_start, ungappedStem );
			}
		    }
		    j_decreasement++;
		}
	    }
	}
    }// end gapSequence
   
    /** 
     * This method searches the string upstream from the three consecutive U's 
     * for a possible hairpin structure. No gaps are allowed in this search.
     */

    protected void ungappedPermutation(Sequence s){
	Sequence seq = s;
	int u_start = u_res.getStart();
	int u_stopp = u_res.getStop();
	String stem1, loop, stem2, tmp;
	String u_sequence;
	boolean shortSeq = false;

	//if sequence is shorter than max terminator length
	if(u_start < 55 ){
	    //tmp = seq.subStr(0, (u_start-1));
	    shortSeq = true;
	}
	
	//ungapped iteration of all permutations of terminator sequence
	if(!shortSeq){
	    int j = 0;
	    int k = 0;
	    ungappedStem = null;
	    
	    for(int i = 0; i < 2; i ++){
		tmp = seq.subStr((u_start-50-i), (u_start-i-1));
		
		if(seq.length() >= (u_start+14))
		    u_sequence = seq.subStr(u_start-i, (u_start + 14-i));
		
		else{
		    u_sequence = seq.subStr(u_start-i, seq.length());
		}
		
		for(j=21; j>4; j--){
		    stem2 = tmp.substring((tmp.length()-j+1),(tmp.length()));
		    for(k=11; k>3; k--){
			loop = tmp.substring((tmp.length()-stem2.length()-k+1),
					     (tmp.length()-stem2.length()));
			
			stem1 =tmp.substring((tmp.length()-stem2.length()-loop.length()-j+1),
					     (tmp.length()-stem2.length()-loop.length())); 
			
			scoreTerm(stem1, loop, stem2, 'N', u_sequence, i , u_start, ungappedStem);
			
		    }
		}
	    }
	}// end if(!shortSeq)
	
    }

    /** 
     * This method scores a terminator candidate suggested by different mutation methods. Scoring is done
     * according to the algorithm form Ermaoleava.
     * @param s1 the left stem in the hairpin
     * @param l the loop of the hairpin
     * @param s2 the right stem of the hairpin
     * @param gap a character telling whether there is a gap or not.
     * @param u the 15 nucleotides following the hairpin
     * @param u_st the position of the start of the U_residue
     * @param k the distance between the hairpin end and the u_residue start
     */
    
    protected void scoreTerm(String s1, String l, String s2, char gap, String u, int k, int u_st, String ungappedStem ){
	
	double termScore, tailScore;
	int GC=0,AU=0,GU=0,MM=0,GP=0,LP=0;
	String leftLeg = s1.toLowerCase();
	String rightLeg = reverseString(s2.toLowerCase());
	String realRightLeg = s2;
	String loop = l;
	String u_seq = u.toLowerCase();
	int u_start = u_st;
	int distance = k; 
		
	LP = l.length();
	
	if (gap == 'R' || gap == 'L') 
	    GP = 1;
	else 
	    GP = 0;
	
	for(int i=0; i<leftLeg.length(); i++) {
	    
	    if(leftLeg.charAt(i) == 'g' &&
	       rightLeg.charAt(i) == 'c' 
	       || leftLeg.charAt(i) == 'c' &&
	       rightLeg.charAt(i) == 'g') 
		GC++;
	    
	    else if ((leftLeg.charAt(i) == 'a' &&
		      rightLeg.charAt(i) == 't') 
		      ||(leftLeg.charAt(i) == 't' &&
			  rightLeg.charAt(i) == 'a'))
		AU++;
	    
	    else if((leftLeg.charAt(i) == 'g' &&
		     rightLeg.charAt(i) == 't') 
		    ||(leftLeg.charAt(i) == 't' &&
		    rightLeg.charAt(i) == 'g'))
		GU++;
	    
	    else 
		MM++;
	}

	//Score with parameters from Ermoleava
	termScore =-(GC * 2.3)-(AU * 0.9)+(GU * 1.3)+(MM * 3.5)+(GP * 6.0)+(1.0 * LP)-5.7;
	
	if(termScore < -6 ){
	    tailScore = tail_Scorer(u_seq);
	
	    if(tailScore < -2){
		TerminatorCandidate terminatorCandidate = new TerminatorCandidate (seq.getName(), termScore, tailScore, realRightLeg, leftLeg, loop, u_seq, distance, u_res.getStart()-2, gap , ungappedStem);
		
		tmpCand.addElement(terminatorCandidate);
	    }
	}
    }
    
    /**
     * This method scores a given tail sequence (max 15 nucleotides) 
     * to check if it scores above the threshold of the tail-score.
     * @param s the sequence to score
     * @return the tail score
     */
    
    protected double tail_Scorer(String s){
	
	double uScore = 0;
	double last_u_score = 1;
	
	for(int z=0; z<s.length(); z++){
		if(s.charAt(z) =='t'){
		    uScore = uScore + (last_u_score * 0.9);
		    last_u_score = last_u_score * 0.9;
		}
		else{
		    uScore = uScore + (last_u_score * 0.6);
		    last_u_score = last_u_score * 0.6;
		}
	}
	return (0 - uScore);
    }
    
    /**
     * This method is used to reverse one stem of the hairpin so that
     * the computation of matches is correct.
     * @param s the string to reverse
     * @return the reversed string
     */
    
    public String reverseString(String s){
	
	char[] tmp = new char [s.length()];
    	int n, m = 0;
		
	for(n = s.length()-1; n >= 0; n--){
	    tmp[m] = s.charAt(n);
	    m++;
	}
	
	return(new String(tmp));
    }	
    
    /**
     * This method is used to find almost similar termiantors and to keep only the best of these
     * Two termiantors are similar if their start and stop points are no more than +/- 10
     * positions away from each other.
     */
    

    private Vector clean(){
       	
	String name;
	int start = 0, stop = 0 ;
	double candScore = 0;
	TerminatorCandidate cand,tmp;
	int k;
		
	for(int i = 0; i < tmpCand.size(); i ++){
	    cand = (TerminatorCandidate)tmpCand.elementAt(i);
	    name = cand.seq_name;
	    start = cand.relativeStart;
	    stop = cand.relativeEnd;
	    k = i+1;
	    candScore = cand.hairpinEnergy + cand.u_Score ;
	    	    
	    if(k < tmpCand.size()){
		tmp = (TerminatorCandidate)tmpCand.elementAt(k);
				
		while(k < tmpCand.size() && tmp.seq_name.equals(name)){
		    tmp = (TerminatorCandidate)tmpCand.elementAt(k);
		    if(((((start-tmp.relativeStart) <= 10) && (start-tmp.relativeStart) >= 0 ) ||
			(((tmp.relativeStart - start) <= 10) && ((tmp.relativeStart - start) >= 0)))
		       && ((((stop-tmp.relativeEnd) <= 10) && (stop-tmp.relativeEnd) >= 0) ||
			   (((tmp.relativeEnd - stop) <= 10) && ((tmp.relativeEnd - stop) >= 0))))
			{
			    if((tmp.hairpinEnergy + tmp.u_Score) < candScore){
				tmp = (TerminatorCandidate)tmpCand.remove(k);
				k--;
			    }
			    else{
				tmp = (TerminatorCandidate)tmpCand.remove(i);
				k =  tmpCand.size();
				i--;
			    }
			}
		    k++;
		}
	    }
	}
	return tmpCand;
    }
}