/*
 * Jalview - A Sequence Alignment Editor and Viewer ($$Version-Rel$$)
 * Copyright (C) $$Year-Rel$$ The Jalview Authors
 * 
 * This file is part of Jalview.
 * 
 * Jalview is free software: you can redistribute it and/or
 * modify it under the terms of the GNU General Public License 
 * as published by the Free Software Foundation, either version 3
 * of the License, or (at your option) any later version.
 *  
 * Jalview is distributed in the hope that it will be useful, but 
 * WITHOUT ANY WARRANTY; without even the implied warranty 
 * of MERCHANTABILITY or FITNESS FOR A PARTICULAR 
 * PURPOSE.  See the GNU General Public License for more details.
 * 
 * You should have received a copy of the GNU General Public License
 * along with Jalview.  If not, see <http://www.gnu.org/licenses/>.
 * The Jalview Authors are detailed in the 'AUTHORS' file.
 */
package jalview.analysis;

import jalview.api.analysis.ScoreModelI;
import jalview.api.analysis.SimilarityParamsI;
import jalview.datamodel.AlignmentView;
import jalview.math.MatrixI;

import java.io.PrintStream;

/**
 * Performs Principal Component Analysis on given sequences
 */
public class PCA implements Runnable
{
  /*
   * inputs
   */
  final private AlignmentView seqs;

  final private ScoreModelI scoreModel;

  final private SimilarityParamsI similarityParams;

  /*
   * outputs
   */
  private MatrixI symm;

  private MatrixI eigenvector;

  private String details;

  /**
   * Constructor given the sequences to compute for, the similarity model to
   * use, and a set of parameters for sequence comparison
   * 
   * @param sequences
   * @param sm
   * @param options
   */
  public PCA(AlignmentView sequences, ScoreModelI sm, SimilarityParamsI options)
  {
    this.seqs = sequences;
    this.scoreModel = sm;
    this.similarityParams = options;
  }

  /**
   * Returns Eigenvalue
   * 
   * @param i
   *          Index of diagonal within matrix
   * 
   * @return Returns value of diagonal from matrix
   */
  public double getEigenvalue(int i)
  {
    return eigenvector.getD()[i];
  }

  /**
   * DOCUMENT ME!
   * 
   * @param l
   *          DOCUMENT ME!
   * @param n
   *          DOCUMENT ME!
   * @param mm
   *          DOCUMENT ME!
   * @param factor
   *          DOCUMENT ME!
   * 
   * @return DOCUMENT ME!
   */
  public float[][] getComponents(int l, int n, int mm, float factor)
  {
    float[][] out = new float[getHeight()][3];

    for (int i = 0; i < getHeight(); i++)
    {
      out[i][0] = (float) component(i, l) * factor;
      out[i][1] = (float) component(i, n) * factor;
      out[i][2] = (float) component(i, mm) * factor;
    }

    return out;
  }

  /**
   * DOCUMENT ME!
   * 
   * @param n
   *          DOCUMENT ME!
   * 
   * @return DOCUMENT ME!
   */
  public double[] component(int n)
  {
    // n = index of eigenvector
    double[] out = new double[getHeight()];

    for (int i = 0; i < out.length; i++)
    {
      out[i] = component(i, n);
    }

    return out;
  }

  /**
   * DOCUMENT ME!
   * 
   * @param row
   *          DOCUMENT ME!
   * @param n
   *          DOCUMENT ME!
   * 
   * @return DOCUMENT ME!
   */
  double component(int row, int n)
  {
    double out = 0.0;

    for (int i = 0; i < symm.width(); i++)
    {
      out += (symm.getValue(row, i) * eigenvector.getValue(i, n));
    }

    return out / eigenvector.getD()[n];
  }

  /**
   * Answers a formatted text report of the PCA calculation results (matrices
   * and eigenvalues) suitable for display
   * 
   * @return
   */
  public String getDetails()
  {
    return details;
  }

  /**
   * Performs the PCA calculation
   */
  @Override
  public void run()
  {
    /*
     * print details to a string buffer as they are computed
     */
    StringBuilder sb = new StringBuilder(1024);
    sb.append("PCA calculation using ").append(scoreModel.getName())
            .append(" sequence similarity matrix\n========\n\n");
    PrintStream ps = wrapOutputBuffer(sb);

    try
    {
      eigenvector = scoreModel.findSimilarities(seqs, similarityParams);

      sb.append(" --- OrigT * Orig ---- \n");
      eigenvector.print(ps, "%8.2f");

      symm = eigenvector.copy();

      eigenvector.tred();

      sb.append(" ---Tridiag transform matrix ---\n");
      sb.append(" --- D vector ---\n");
      eigenvector.printD(ps, "%15.4e");
      ps.println();
      sb.append("--- E vector ---\n");
      eigenvector.printE(ps, "%15.4e");
      ps.println();

      // Now produce the diagonalization matrix
      eigenvector.tqli();
    } catch (Exception q)
    {
      q.printStackTrace();
      sb.append("\n*** Unexpected exception when performing PCA ***\n"
              + q.getLocalizedMessage());
      sb.append(
              "*** Matrices below may not be fully diagonalised. ***\n");
    }

    sb.append(" --- New diagonalization matrix ---\n");
    eigenvector.print(ps, "%8.2f");
    sb.append(" --- Eigenvalues ---\n");
    eigenvector.printD(ps, "%15.4e");
    ps.println();

    details = sb.toString();
  }

  /**
   * Returns a PrintStream that wraps (appends its output to) the given
   * StringBuilder
   * 
   * @param sb
   * @return
   */
  protected PrintStream wrapOutputBuffer(StringBuilder sb)
  {
    PrintStream ps = new PrintStream(System.out)
    {
      @Override
      public void print(String x)
      {
        sb.append(x);
      }

      @Override
      public void println()
      {
        sb.append("\n");
      }
    };
    return ps;
  }

  /**
   * Answers the N dimensions of the NxN PCA matrix. This is the number of
   * sequences involved in the pairwise score calculation.
   * 
   * @return
   */
  public int getHeight()
  {
    // TODO can any of seqs[] be null?
    return seqs.getSequences().length;
  }
}