//
// This file is part of MorphoGraphX - http://www.MorphoGraphX.org
// Copyright (C) 2012-2015 Richard S. Smith and collaborators.
//
// If you use MorphoGraphX in your work, please cite:
//   http://dx.doi.org/10.7554/eLife.05864
//
// MorphoGraphX is free software, and is licensed under under the terms of the 
// GNU General (GPL) Public License version 2.0, http://www.gnu.org/licenses.
//
#include <HistogramMaker.hpp>
#include <GraphUtils.hpp>
#include <MeshProcessHeatMap.hpp>
#include <math.h> 
#include <iostream>
#include <vector>
#include <fstream>
#include <limits>
#include <Progress.hpp>
#include <QFileDialog>
#include <QFileDialog>


namespace mgx
{

// return the attr map for a specified measure 
void HistogramMaker::getAttrMap(AttrMap<int, double>& data, QString measure)
{
  Mesh* m1 = currentMesh();

  if(measure == "Current Heat Map"){
  	data.clear();
  	forall(auto p, m1->labelHeat()){
  	  data[p.first] = p.second;
  	}
  	return;
  }

  QString prefix = "Measure Label Double ";
  QString proName = "Mesh/Heat Map/Measures";
  QString measureAttr = measure;
  measureAttr.replace( " ", "" );
  data = m1->attributes().attrMap<int, double>(prefix+measureAttr);

  if(data.size() == 0){ // calculate it
    qDebug() << "Attr Map \t" << prefix+measure << " doesnt exist\n";
    QStringList parms;
    Process *pro = makeProcess(proName+measure);
    if(pro){
    qDebug() << "Process call\t" << proName+measure << "\n";
    getLastParms(proName+measure, parms);
    pro->run();

    data = m1->attributes().attrMap<int, double>(prefix+measureAttr);
    }
  }
 // clMap.addCellFeature(measure, data, selectedCells);
}


 // fill the combo boxes with all available measures and/or attr maps
void HistogramMaker::fillMeasures(){

    ui.xDimCmb->clear();
    ui.yDimCmb->clear();

    QStringList defaultEntries;
    defaultEntries << "none";

    int idx = 0;
    int selIdxX = 0;
    int selIdxY = 0;

    // find measure processes
    QStringList procs = mgx::listProcesses();

    QStringList addProcs;
    addProcs << "Current Heat Map";
    std::set<QString> procsSet;
    QStringList addAttr;
    forall(const QString& name, procs) {
      mgx::ProcessDefinition* def = mgx::getProcessDefinition(name);
      QStringList list = def->name.split("/");

      QString subfolder = "Measures";

      if(list[0] == "Mesh" and list[1] == "Heat Map" and list[2] == subfolder){
        QString newProc = "/" + list[3] + "/" + list[4];
        addProcs << newProc;
        procsSet.insert(newProc);
        if(newProc == selectedX) selIdxX = idx;
        if(newProc == selectedY) selIdxY = idx;
        idx++;
      }
      
    }

    // now the attr maps
    Mesh* m = currentMesh();
    Attributes *attributes;
    attributes = &m->attributes();
    QStringList attr = attributes->getAttrList(); 
    QStringList measureAttr;

    forall(const QString &text, attr) {
      QStringList list = text.split(" ");

      if(list.size() < 4) continue;
      if(list[0] != "Measure") continue;
      QString name;
      for(int i = 3; i<list.size(); i++){
        name = name + list[i] + " ";
      }
      name = name.trimmed();
      if(procsSet.find(name) == procsSet.end()){
        addProcs << name;
          // if(name == selectedX) selIdxX = idx;
          // if(name == selectedY) selIdxY = idx;
          // idx++;
      }
    }

    addProcs.sort(Qt::CaseSensitive);

    forall(const QString &name, attr){
      if(name == selectedX) selIdxX = idx;
      if(name == selectedY) selIdxY = idx;
      idx++;
    }

    ui.xDimCmb->addItems(addProcs);
    ui.yDimCmb->addItems(addProcs);

    ui.xDimCmb->setCurrentIndex(selIdxX+1);
    ui.yDimCmb->setCurrentIndex(selIdxY+1);

  }


void HistogramMaker::changeValue(){

  bin_num =  ui.spinBoxBins->value();
  x_max =  ui.spinBoxXMax->value();
  x_min =  ui.spinBoxXMin->value();
  y_max =  ui.spinBoxYMax->value();
  y_min =  ui.spinBoxYMin->value();

  selectedX = ui.xDimCmb->currentText();
  selectedY = ui.yDimCmb->currentText();



/*
  clMap.customMinMax = true;

  clMap.xMinMax.x() = ui.spinBoxSigmaXMin->value();
  clMap.xMinMax.y() = ui.spinBoxSigmaXMax->value();
  clMap.yMinMax.x() = ui.spinBoxSigmaYMin->value();
  clMap.yMinMax.y() = ui.spinBoxSigmaYMax->value();

  clMap.heatMax = ui.spinBoxScaleHeat->value();

  changeHeatmap();
*/
}


bool HistogramMaker::initialize(QWidget* parent)
  {

    dlg = new QDialog(parent);
    ui.setupUi(dlg);
    this->dlg = dlg;

    fillMeasures();

    // check which is selected
    connect(ui.yDimCmb, SIGNAL(currentIndexChanged(QString)), this, SLOT(changeValue()));
    connect(ui.xDimCmb, SIGNAL(currentIndexChanged(QString)), this, SLOT(changeValue()));


    // calculate measure & update heatmap
    ui.spinBoxBins->setValue(bin_num);
    ui.spinBoxXMax->setValue(x_max);
    ui.spinBoxXMin->setValue(x_min);
    ui.spinBoxYMax->setValue(y_max);
    ui.spinBoxYMin->setValue(y_min);

    connect(ui.spinBoxBins, SIGNAL(valueChanged(double)), this, SLOT(changeValue()));
    connect(ui.spinBoxXMax, SIGNAL(valueChanged(double)), this, SLOT(changeValue()));
    connect(ui.spinBoxXMin, SIGNAL(valueChanged(double)), this, SLOT(changeValue()));
    connect(ui.spinBoxYMax, SIGNAL(valueChanged(double)), this, SLOT(changeValue()));
    connect(ui.spinBoxYMin, SIGNAL(valueChanged(double)), this, SLOT(changeValue()));
    

    if(dlg->exec() == QDialog::Accepted && parent){
      return OutputHistogram(stringToBool(parm("Separate by Parent")));
    }
    else
      return false;

  }

  // HistogramMaker GUI process for HeatMap clustering
  REGISTER_PROCESS(HistogramMaker);


  bool HistogramMaker::OutputHistogram(bool separateParents){
 
       AttrMap<int, double> dataX, dataY;
        getAttrMap(dataX,selectedX); 
        getAttrMap(dataY,selectedY);
	std::cout<<"Creating histogram with: "<<selectedX.toUtf8().constData()<<" "<<selectedY.toUtf8().constData()<<std::endl;
	std::cout<<"And binning parameters: "<<x_min<<"-"<<x_max<<", "<<y_min<<"-"<<y_max<<" in "<<bin_num<<" bins."<<std::endl;

    std::set<int> parentLabels;
    Mesh* m1 = currentMesh();

    if(separateParents){
	    forall(auto p, m1->parents()){
	      if(p.second < 1) continue;
	      parentLabels.insert(p.second);
	    }
    } else {
    	parentLabels.insert(0);
    }

 

	std::vector<double> bins;
//	std::vector<std::vector<int>> bin_labels;
	std::vector<double> bin_bounds;
	std::vector<double> bin_variance;
	std::vector<double> counts;
	double under_bin=0;
	double under_count=0;
	double over_bin=0;
	double over_count=0;
	double under_y=0;
	double under_y_count = 0;
	double over_y=0;
	double over_y_count =0;
	int in_x_not_y =0;
	bins.resize(bin_num);
	bin_variance.resize(bin_num);
	bin_bounds.resize(bin_num+1);
	counts.resize(bin_num);
	if(x_min >= x_max){
		x_min = 1e20;
		x_max = -1e20;
		forall(IntDouble x,dataX){
			if(x.second>x_max)
				x_max = x.second;
			if(x.second<x_min)
				x_min = x.second;
		}
		//this will miss max and min value sometimes
		x_min -= 0.01*(x_min);
		x_max += 0.01*(x_max);
	}
	if(y_min >= y_max){
		y_min = -1e20;
		y_max = 1e20;
	}

    QString filename;
//	    if(/*filename.isEmpty() and*/ parent)
      filename = QFileDialog::getSaveFileName(0, "Choose spreadsheet file to save", QDir::currentPath(), "CSV files (*.csv)");
    if(filename.isEmpty())
      return false;
    if(!filename.endsWith(".csv", Qt::CaseInsensitive))
      filename += ".csv";
//	    parms[0] = filename;

    QFile file(filename);
    if(!file.open(QIODevice::WriteOnly)) {
      setErrorMessage(QString("File '%1' cannot be opened for writing").arg(filename));
      return false;
    }
    QTextStream out(&file);

	//NEED TO COUNT PRESENT IN Y BUT NOT X
	//OPTION FOR SUM OR AVERAGE
	//OPTION FOR WHAT TO DO WITH BIN_BOUNDS FOR FIRST AND LAST SAMPLE
	//NEED TO GET CSV FILE NAME!
//	for(int i=0;i<NUM_X;i++){
	forall(int parentLabel, parentLabels){


		for(int i=0;i<bin_num;i++){
		  float cur_bound = ((float)i/(float)bin_num)*(x_max-x_min)+x_min;
		  bin_bounds[i]=cur_bound;
		  counts[i]=0;	
		  bins[i] = 0;
		  bin_variance[i] = 0;
		}
		bin_bounds[bin_num] = x_max;


		forall(IntDouble x,dataX){
			if(separateParents and m1->parents()[x.first] != parentLabel) continue;
			
			float x_val = x.second;
			if(dataY.find(x.first) == dataY.end() or std::isnan(dataY[x.first]) or std::isinf(dataY[x.first])){
				in_x_not_y++;
				continue;
			}
			float y_val =dataY[x.first] ;
			if(y_val < y_min){
				under_y+=y_val;
				under_y_count++;
				continue;
			}
			if(y_val > y_max){
				over_y+=y_val;
				over_y_count++;
				continue;
			}


			if(x_val<bin_bounds[0]){
				under_bin+=y_val;
				under_count++;
			}
			else if(x_val>=bin_bounds[bin_num]){
				over_bin+=y_val;
				over_count++;
			}
			else{
				for(int j=0;j<bin_num;j++){
					if(x_val>=bin_bounds[j] and x_val<bin_bounds[j+1]){
						bins[j]+=y_val;
						//Construct first term of variance
						bin_variance[j]+=y_val*y_val;
						counts[j]++;
						//bin_labels[j].push_back(x.first);
					}
				}
			}

		}

		//Compute the variance of each bin
		for(int i=0;i<bin_num;i++){			
			if(counts[i]<=0)
				continue;
	    else if(counts[i] == 1){
	      bin_variance[i] = 0;
	      continue;
	    }
			double average = bins[i]/counts[i];
			bin_variance[i]/=counts[i];
			bin_variance[i]-=average*average;

		}

		//Output CSV file here
		out<<"Values of "<<selectedY.toUtf8().constData()<<" binned by "<<selectedX.toUtf8().constData()<<endl;
		if(separateParents) out<<"Parent Label "<<parentLabel<<endl;
		out<<"Bin number, Interval start, Interval end, Count, Accumulated value, Average, Std-dev"<<endl;
		for(int i=0;i<bin_num;i++)
			out<<i<<","<<bin_bounds[i]<<","<<bin_bounds[i+1]<<","<<counts[i]<<","<<bins[i]<<","<<(counts[i]>0? bins[i]/counts[i] : 0)<<","<<(counts[i]>0? sqrt(bin_variance[i]) : 0)<<endl;

	    out<<endl;
		out<<"Below min bin"<<","<<"-inf"<<","<<bin_bounds[0]<<","<<under_count<<","<<under_bin<<","<<(under_count>0? under_bin/under_count : 0)<<",NA"<<endl;

	    out<<endl;
		out<<"Over max bin"<<","<<bin_bounds[bin_num]<<","<<"+inf"<<","<<over_count<<","<<over_bin<<","<<(over_count>0? over_bin/over_count : 0)<<",NA"<<endl;

	    out<<endl;
		out<<"Below min thresh"<<","<<"-inf"<<","<<y_min<<","<<under_y_count<<","<<under_y<<","<<(under_y_count>0? under_y/under_y_count : 0)<<",NA"<<endl;

	    out<<endl;
		out<<"Over max thresh"<<","<<y_max<<","<<"+inf"<<","<<over_y_count<<","<<over_y<<","<<(over_y_count>0? over_y/over_y_count : 0)<<",NA"<<endl;

		out<<"Total Num of cells in not binned,"<<in_x_not_y<<endl;
		if(separateParents) out<<endl;
    }
	return true;
  }

}
