package org.fda.contiggenerator;

import java.io.File;
import java.util.Random;
import org.apache.commons.math3.distribution.GammaDistribution;
import org.apache.commons.math3.random.RandomGeneratorFactory;
import org.apache.commons.math3.random.RandomGenerator;
import org.fda.data.Utilities;
import org.fda.data.Enums.Distribution;

/**
 *
 * @author Gokhan.Yavas
 */
public abstract class ContigGeneratorWDistribution extends ContigGenerator{
    protected final Random rand;
    protected final GammaDistribution grand;
    protected double mean;
    protected double sd;
    protected final Distribution dist;
    protected boolean setseed=false;

    public double getMean() {
        return mean;
    }

    public double getSd() {
        return sd;
    }

    public void setMean(double mean) {
        this.mean = mean;
    }

    public void setSd(double sd) {
        this.sd = sd;
    }
    
    
    public ContigGeneratorWDistribution(File infile, File outfile, double mean, double sd, Distribution dist){
        super(infile, outfile);
        this.dist = dist;
        this.mean = mean;
        this.sd = sd;
        rand = new Random();
        if(dist==Distribution.GAMMA){
            RandomGenerator r = RandomGeneratorFactory.createRandomGenerator(rand);
            grand = new GammaDistribution(r, this.sd, this.mean);
        }
        else{
            grand=null;
        }
        
    }
    
    public ContigGeneratorWDistribution(File infile, File outfile, double mean, double sd, Distribution dist, boolean createreport){
        super(infile, outfile, createreport);
        this.dist = dist;
        this.mean = mean;
        this.sd = sd;
        rand = new Random();
        if(dist==Distribution.GAMMA){
            RandomGenerator r = RandomGeneratorFactory.createRandomGenerator(rand);
            grand = new GammaDistribution(r, this.sd, this.mean);
        }
        else{
            grand=null;
        }
        
    }
    
    public ContigGeneratorWDistribution(File infile, File outfile, double mean, double sd, Distribution dist, long seed){
        super(infile, outfile);
        this.dist = dist;
        this.mean = mean;
        this.sd = sd;
        rand = new Random(seed);
        if(dist==Distribution.GAMMA){
            RandomGenerator r = RandomGeneratorFactory.createRandomGenerator(rand);
            grand = new GammaDistribution(r, this.sd, this.mean);
        }
        else{
            grand=null;
        }
        
    }
    public ContigGeneratorWDistribution(File infile, File outfile, double mean, double sd, Distribution dist, long seed, boolean createreport){
        super(infile, outfile, createreport);
        this.dist = dist;
        this.mean = mean;
        this.sd = sd;
        rand = new Random(seed);
        if(dist==Distribution.GAMMA){
            RandomGenerator r = RandomGeneratorFactory.createRandomGenerator(rand);
            grand = new GammaDistribution(r, this.sd, this.mean);
        }
        else{
            grand=null;
        }        
    }
    
    public ContigGeneratorWDistribution(File infile, File outfile, double mean, double sd, Distribution dist, boolean setseed, boolean createreport){
        this(infile, outfile, mean, sd, dist, createreport);
        this.setseed=setseed;        
    }
    
    
    private int getGammaDistRandomLength(long l){
        int contlen;
        while(true){
            contlen = (int)grand.sample();
            if(contlen>=Utilities.determineAllowedMinSimContigLength(l))
                break;                                        
        }
        
        return contlen;
    }
    private int getNormalDistRandomLength(long l){                
        int contlen;        
        while(true){
            contlen = (int)Math.round(rand.nextGaussian()*this.sd+this.mean);
            if(contlen>=Utilities.determineAllowedMinSimContigLength(l))
                break;                                        
        }
        return contlen;
    }
    private int getUniformDistRandomLength(long l){                
        int contlen;        
       while(true){
            contlen = rand.nextInt((int)this.mean+1);
            if(contlen>=Utilities.determineAllowedMinSimContigLength(l))
                break;                                        
        }
        return contlen;
    }

    
    private int getGammaDistRandomLength(long l, double target_cov, long current_length){
        int contlen;
        while(true){
            contlen = (int)grand.sample();
            if(contlen>=Utilities.determineAllowedMinSimContigLength(l)){
                double posscov=(double)(contlen+current_length)/l;
                //if(posscov/target_cov<=1.1 && posscov/target_cov>=0.9)
                if(posscov/target_cov<=(1+Utilities.simulation_coverage_ratio_threshold))
                    break;                                                        
            }
        }
        
        return contlen;
    }
    private int getNormalDistRandomLength(long l, double target_cov, long current_length){                
        int contlen;        
        while(true){
            contlen = (int)Math.round(rand.nextGaussian()*this.sd+this.mean);
            if(contlen>=Utilities.determineAllowedMinSimContigLength(l)){
                double posscov=(double)(contlen+current_length)/l;
                //if(posscov/target_cov<=1.1 && posscov/target_cov>=0.9)
                if(posscov/target_cov<=(1+Utilities.simulation_coverage_ratio_threshold))
                    break;                                        
            }
        }
        return contlen;
    }
    private int getUniformDistRandomLength(long l, double target_cov, long current_length){                
        int contlen;        
       while(true){
            contlen = rand.nextInt((int)this.mean+1);            
            if(contlen>=Utilities.determineAllowedMinSimContigLength(l)){
                double posscov=(double)(contlen+current_length)/l;
                //if(posscov/target_cov<=1.1 && posscov/target_cov>=0.9)
                if(posscov/target_cov<=(1+Utilities.simulation_coverage_ratio_threshold))
                    break;                                        
            }
        }
        return contlen;
    }
    
    
    protected int getRandomLength(long ref_cont_length){
        if(null==dist)
            return getUniformDistRandomLength(ref_cont_length);
        else switch (dist) {
            case NORMAL:
                return getNormalDistRandomLength(ref_cont_length);
            case GAMMA:
                return getGammaDistRandomLength(ref_cont_length);
            default:
                return getUniformDistRandomLength(ref_cont_length);
        }
    }                
    protected int getRandomLength(long ref_cont_length, double target_cov, long current_length){
        if(null==dist)
            return getUniformDistRandomLength(ref_cont_length, target_cov, current_length);
        else switch (dist) {
            case NORMAL:
                return getNormalDistRandomLength(ref_cont_length, target_cov, current_length);
            case GAMMA:
                return getGammaDistRandomLength(ref_cont_length, target_cov, current_length);
            default:
                return getUniformDistRandomLength(ref_cont_length, target_cov, current_length);
        }
    }                
    
}
