/*
 * TwoWireGraphAnalyzer.cpp
 *
 *  Created on: Apr 24, 2015
 *      Author: Hovanes Egiyan
 */

#include "TwoWireGraphAnalyzer.hh"
map<string, map<string,double> > TwoWireGraphAnalyzer::twgaFitRange = TwoWireGraphAnalyzer::InitRanges();

TwoWireGraphAnalyzer::TwoWireGraphAnalyzer(TGraphErrors* graph, string sName, string label ) :
GraphAnalyzer( graph, sName, label ),
twgaFitFun(0),
twgaSubGraph(0)
{
	cout << "In TwoWireGraphAnalyzer::TwoWireGraphAnalyzer()" << endl;
	CreateSubGraph();
	return;
}

TwoWireGraphAnalyzer::TwoWireGraphAnalyzer( const TwoWireGraphAnalyzer& analyzer ) :
	GraphAnalyzer(analyzer) {
	twgaFitFun		= new TF1(*analyzer.twgaFitFun);
	twgaSubGraph	= new TGraphErrors(*analyzer.twgaSubGraph);
	return;
}

TwoWireGraphAnalyzer::~TwoWireGraphAnalyzer() {
	if( twgaFitFun != 0 ) delete twgaFitFun;
	if( twgaSubGraph != 0 ) delete twgaSubGraph;
	return;
}

void TwoWireGraphAnalyzer::CreateSubGraph() {
	unsigned long nSubPoints = 0;
	unsigned long minIndex = 0;
	unsigned long maxIndex = 0;
	vector<double> subX, subY, subErrX, subErrY;
	for( Int_t iPoint=0; iPoint < fNpoints; iPoint++ ) {
		if( fX[iPoint] >= twgaFitRange[gaLabel]["low"] && nSubPoints == 0  ) {
			minIndex = iPoint;
		}
		if( fX[iPoint] > twgaFitRange[gaLabel]["high"] && iPoint > 0 ) {
			maxIndex = (iPoint - 1);
		}
		if( fX[iPoint] >= twgaFitRange[gaLabel]["low"] &&
				fX[iPoint] <=  twgaFitRange[gaLabel]["high"] ) {
			subX.push_back( fX[iPoint] );
			subY.push_back( fY[iPoint] );
			subErrX.push_back( fEX[iPoint] );
			subErrY.push_back( fEY[iPoint] );
		}
		twgaSubGraph = new TGraphErrors( subX.size(), &subX[0], &subY[0], &subErrX[0], &subErrY[0] );
	}
	return;
}

map<string, map<string,double> > TwoWireGraphAnalyzer::InitRanges() {
	map<string, map<string,double> > localMap;
	{
		map<string,double> tmpMap;
		tmpMap["low"] 	= 30.0;
		tmpMap["high"] 	= 55.0;
		tmpMap["dir"]	= -1.0;
		localMap["X"] = tmpMap;
	}
	{
		map<string,double> tmpMap;
		tmpMap["low"] 	= 65.0;
		tmpMap["high"] 	= 90.0;
		tmpMap["dir"]	= -1.0;
		localMap["Y"] = tmpMap;
	}
	return localMap;
}

void TwoWireGraphAnalyzer::FindInitialParValues() {
	cout << "In TwoWireGraphAnalyzer::FindInitialParValues() for " << fName << ":" << gaLabel << endl;
	double xMin, xMax, yMin, yMax;
	twgaSubGraph->ComputeRange( xMin, yMin, xMax, yMax );
	twgaFitFun->SetParLimits(0, 0.05*yMax, 5.0*yMax );

	double minVal = TMath::Min( twgaFitRange[gaLabel]["dir"] * twgaFitRange[gaLabel]["low"]/sqrt(2),
			twgaFitRange[gaLabel]["dir"] * twgaFitRange[gaLabel]["high"]/sqrt(2) );
	double maxVal = TMath::Max( twgaFitRange[gaLabel]["dir"] * twgaFitRange[gaLabel]["low"]/sqrt(2),
			twgaFitRange[gaLabel]["dir"] * twgaFitRange[gaLabel]["high"]/sqrt(2) );
	twgaFitFun->SetParLimits(1,  minVal, maxVal );

	twgaFitFun->SetParLimits(2, 0.010, 0.5*(twgaFitRange[gaLabel]["high"] - twgaFitRange[gaLabel]["low"]));
	twgaFitFun->SetParLimits(3, 0.0, 2*yMax );
	twgaFitFun->SetParameter(0, yMax );
	twgaFitFun->SetParameter(1, twgaFitRange[gaLabel]["dir"] * twgaSubGraph->GetMean() / sqrt(2.0) ) ;
	twgaFitFun->SetParameter(2, twgaSubGraph->GetRMS() / sqrt(2.0) / 5.0 );
	twgaFitFun->SetParameter(3, 0.001 * yMax );
	cout << "Parameter "<<  twgaFitFun->GetParName(0) << " is " << twgaFitFun->GetParameter(0) << endl;
	cout << "Parameter "<<  twgaFitFun->GetParName(1) << " is " << twgaFitFun->GetParameter(1) << endl;
	cout << "Parameter "<<  twgaFitFun->GetParName(2) << " is " << twgaFitFun->GetParameter(2) << endl;
	cout << "Parameter "<<  twgaFitFun->GetParName(3) << " is " << twgaFitFun->GetParameter(3) << endl;
	cout << "Initial parameters set for " << fName << ":" << gaLabel << endl;
	return;
}


void TwoWireGraphAnalyzer::Draw( Option_t* opt) {
	TGraphErrors::Draw(opt);
	GetXaxis()->SetRangeUser( twgaFitRange[gaLabel]["low"], twgaFitRange[gaLabel]["high"]);
	GetYaxis()->SetTitleOffset( 1.2 );
	TGraphErrors::Draw(opt);
	return;
}

TFitResultPtr TwoWireGraphAnalyzer::FitGraph() {
	if( twgaFitFun == 0 ) {
		string funName = string(static_cast<const char*>(fName) ) + string("_fun");
		CreateFitFun( funName, twgaFitRange[gaLabel]["low"], twgaFitRange[gaLabel]["high"]);
	}
	this->FindInitialParValues();
//	twgaFitFun->ReleaseParameter(0);
//	return Fit( this->twgaFitFun, "", "", twgaFitRange[gaLabel]["low"], twgaFitRange[gaLabel]["high"] );
//	twgaFitFun->ReleaseParameter(1);
//	return Fit( this->twgaFitFun, "", "", twgaFitRange[gaLabel]["low"], twgaFitRange[gaLabel]["high"] );
//	twgaFitFun->ReleaseParameter(2);
//	return Fit( this->twgaFitFun, "EX0", "", twgaFitRange[gaLabel]["low"], twgaFitRange[gaLabel]["high"] );
//	twgaFitFun->ReleaseParameter(3);
	return Fit( this->twgaFitFun, "EX0", "", twgaFitRange[gaLabel]["low"], twgaFitRange[gaLabel]["high"] );
}

void TwoWireGraphAnalyzer::CreateFitFun( string funName, double xMin, double xMax ) {
//	void* funPtr = reinterpret_cast<void*>( TwoWireGraphAnalyzer::TwoWireScanFitFun );
	twgaFitFun = new TF1( funName.c_str(), TwoWireScanFitFun, xMin, xMax, 5);
	cout << "Fit function array is " << hex << showbase << twgaFitFun << endl;
	twgaFitFun->SetParName( 0, "Amplitude");
	twgaFitFun->SetParName( 1, (string("Mean ") + gaLabel).c_str());
	twgaFitFun->SetParName( 2, (string("Sigma ") + gaLabel).c_str());
	twgaFitFun->SetParName( 3, "Background");
	twgaFitFun->SetParName( 4, "Direction" );

	twgaFitFun->FixParameter(4, twgaFitRange[gaLabel]["dir"] );
}


double TwoWireGraphAnalyzer::TwoWireScanFitFun( double* xArray, double* param ) {
//	cout << "In TwoWireGraphAnalyzer::TwoWireScanFitFun" << endl;
	// Define a fit function with a flat background and a gaussian
	// the parameters are defined for the X or Y directions, not the
	// direction of the motion of the motor.
	double ampl = param[0];
	double mean = param[1];
	double sigma = param[2];

	double bkg0 = param[3];

	double dirSign = param[4];

	mean *= dirSign;

	// Scale by square root of too since the motion of the harp is at 45 degrees.
	double x = xArray[0] / sqrt( 2.0 );
//	double x 		= xArray[0];

	// Calculate non-distorted value
	double value = ampl * exp( -(x - mean) * (x - mean) / (2.0 * sigma * sigma) ) + bkg0;

	return value;
}