package org.umdnj.JDNA;

import java.util.*;
import java.lang.Math;
/*  Implementation of the local alignment algorithm developed by Smith & Waterman
 *
 *  Details taken from J. Mol. Biol. (1981), Volume 147, pages 195-197
 *
 *  Author: Ryan Golhar (golharam@umdnj.edu)
 */

public class SmithWaterman {

    private static double s(char a, char b) {
        if (a==b)
            return 1;
        else
            return (-0.25);
    }

    private static double W(int k) {
        return 1.0 + (0.25)*k;
    }

    public static void align(String A, String B, Alignment al) {
        int n = A.length();
        int m = B.length();
        double[][] H = new double[1+n][1+m];
        int k, l;
        int i, j;
        double c1, c2, c3, c4; // case 1, case 2, case 3, case 4

        // start algorithm
        for (k = 0; k <= n; k++)
            H[k][0] = 0;

        for (l = 0; l <= m; l++)
            H[0][l] = 0;

        for (i = 1; i <= n; i++) {
            for (j = 1; j <= m; j++) {
                c1 = c2 = c3 = c4 = 0;

                // case 1
                c1 = H[i-1][j-1] + s(A.charAt(i-1), B.charAt(j-1));

                // case 2 (max of H[i-k][j] - W[k]), k >= 1

                // v1:
                c2 = Math.max(c2, H[i-1][j] - W(1));

                // v1:
                c3 = Math.max(c3, H[i][j-1] - W(1));

                // v2:
/*                k = i;
                c2 = 0;
                while (k >= 1) {
                    c2 = Math.max(c2, H[i-k][j] - W(k));
                    k--;
                }

                // case 3 (max of H[i][j-l] - W(l)), l >= 1

                // v2:
                l = j;
                c3 = 0;
                while (l >= 1) {
                    c3 = Math.max(c3, H[i][j-l] - W(l));
                    l--;
                }
*/

                // case 4
                //c4 = 0;

                H[i][j] = Math.max(Math.max(c1, c2), Math.max(c3, c4));
            }
        }

        // output
        makealignment(A, B, H, al);//, aTraceback, bTraceback);
    }

    private static void makealignment(String A, String B, double[][] H, Alignment al) { //, int[][] aTraceback, int[][] bTraceback) {
        int maxrow, maxcol;
        int row, col;
        int rowlen = A.length();
        int collen = B.length();
        StringBuffer outA, outB, midline;
        boolean skipA, skipB;
        double diag, left, top;

//        System.out.println("SIDE: " + A);
//        System.out.println("TOP: " + B);
         maxrow = maxcol = 0;

        // print matrix
        // find maximum value cell
/*        System.out.print("\t");
        for (col = 1; col <= collen; col++) {
                System.out.print("\t " + B.charAt(col-1));
        }
        System.out.println();
*/

        for (row = 0; row <= rowlen; row++) {
 //           if (row > 0)
 //               System.out.print(A.charAt(row-1));

            for (col = 0; col <= collen; col++) {
 //               System.out.print("\t" + Double.toString(H[row][col]));
                if (H[row][col] > H[maxrow][maxcol]) {
                    maxrow = row;
                    maxcol = col;
                }
            }
//            System.out.println();
        }

        // print alignments
        // maxrow, maxcol gives us the end of the path,
        // so we need to trace back to the beginning of the path
        // as we do that, add the base to the beginning of the string as
        // we are "tracing back".  Once done, concat the beginning and end.
        outA = new StringBuffer();
        outB = new StringBuffer();
        skipA = skipB = false;

        row = maxrow;
        col = maxcol;
        while (H[row][col] > 0) {
            diag = H[row-1][col-1];
            top = H[row-1][col];
            left = H[row][col-1];

            if ((!skipA) && (!skipB)) {
                outA.insert(0, A.charAt(row-1));
                outB.insert(0, B.charAt(col-1));
            } else if (skipA) {
                outA.insert(0, '-');
                outB.insert(0, B.charAt(col-1));
            } else {
                outA.insert(0, A.charAt(row-1));
                outB.insert(0, '-');
            }

            if ((diag >= top) && (diag >= left)) {
                row--;
                col--;
                skipA = skipB = false;
            } else
            if ((top >= diag) && (top >= left)) {
                row--;
                skipA = false;
                skipB = true;
            } else {
                col--;
                skipA = true;
                skipB = false;
            }
        }

        // insert the beginning of the sequences
/*        if (row >= 0) {
            outA.insert(0, A.substring(0, row+1));
        }

        if (col >= 0) {
            outB.insert(0, B.substring(0, col+1));
        }

        if (maxrow < rowlen) {
            outA.append(A.substring(maxrow+1));
        }

        if (maxcol < collen) {
            outB.append(B.substring(maxcol+1));
        }
    */
        int i;
	midline = new StringBuffer();

	for (i = 0; i < outA.length(); i++) {
            if (outA.charAt(i) == outB.charAt(i))
                midline.append("|");
	    else
	        midline.append(" ");
	}

        //return new Alignment(outA.toString(), outB.toString(), midline.toString());
        al.setAlignment(outA.toString(), outB.toString(), midline.toString());
//      System.out.println("Alignments:");
//	System.out.println(outA);
//	System.out.println(midline);
//	System.out.println(outB);
//	System.out.println();
    }

/*
    public static void main(String args[]) {
        SmithWaterman smwm = new SmithWaterman();
        // sw (side, top)
        //smwm.sw("AAUGCCAUUGACGG", "CAGCCUCGCUUAG");
        //smwm.sw("ACCCCAGGCTTTACACAT", "TTGACACCCTCCCAATTGTA");
        smwm.sw(args[0], args[1]);
    }
*/
}
