/*
 * License: FreeBSD (Berkeley Software Distribution)
 * Copyright (c) 2013, Sara Sheehan, Kelley Harris, Yun Song
 * 
 * based on:
 * Copyright (c) 2011, Joshua Paul, Matthias SteinrŸcken, Yun Song
 */

package edu.berkeley.utility;

import java.io.BufferedReader;
import java.io.FileReader;
import java.io.IOException;

import Jama.EigenvalueDecomposition;
import Jama.Matrix;

// the parameters for a finite-site, finite-alleles coalescent with recombination
public class ParamSet {
	
	// mutation parameters
	private final int numAlleles;
	private final double mutRate;
	private final Matrix mutMatrix;
	private final Matrix mutDiagonal;
	private final Matrix mutU;
	private final Matrix mutUinverse;
	private final double[] logStationaryProb;
		
	// recombination parameters
	private final double recRate;
	
	// read in the parameters from parameter file
	public ParamSet(String paramFile) throws IOException {
		
		BufferedReader paramReader = new BufferedReader(new FileReader(paramFile));

		String[] thetaStrings = Utility.readLine(paramReader).split(";");
		assert (thetaStrings.length == 1);
		this.mutRate = Double.parseDouble(thetaStrings[0]);
			
		String[] matrixStrings = Utility.readLine(paramReader).split(";");
		assert (matrixStrings.length == 1);
		this.mutMatrix = Utility.readMatrix(matrixStrings[0]);
		this.numAlleles = this.mutMatrix.getRowDimension();
			
		String[] rhoStrings = Utility.readLine(paramReader).split(";");
		assert (rhoStrings.length == 1);
		this.recRate = Double.parseDouble(rhoStrings[0]);
		
		// diagonalize mutMatrix (do this only once)
		EigenvalueDecomposition eigenD = new EigenvalueDecomposition(this.mutMatrix);
		this.mutDiagonal = eigenD.getD();
		this.mutU = eigenD.getV();
		this.mutUinverse = this.mutU.inverse();
		
		// compute the stationary distribution (do this only once)
		this.logStationaryProb = StationaryProbability.getLogStationaryDistribution(this.mutMatrix.getArray());
	}
		
	// mutation getters
	public int numAlleles() { return this.numAlleles; }
	public double getMutationRate() { return this.mutRate; }
	public Matrix getMutationMatrix() { return this.mutMatrix; }
	public Matrix getMutationDiagonal() { return this.mutDiagonal; }
	public Matrix getMutationU() { return this.mutU; }
	public Matrix getMutationUinverse() { return this.mutUinverse; }
	public double[] getLogStationaryProb() {return this.logStationaryProb; }
		
	// recombination getter
	public double getRecombinationRate() { return this.recRate; }
}