/*
 * ThreeWireGraphAnalyzer.cpp
 *
 *  Created on: May 11, 2022
 *      Author: Hovanes Egiyan
 */

#include "ThreeWireGraphAnalyzer.hh"

using namespace std;

map<string, map<string,double> > ThreeWireGraphAnalyzer::twgaFitRange = ThreeWireGraphAnalyzer::InitRanges();

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

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

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

void ThreeWireGraphAnalyzer::CreateSubGraph() {
	vector<double> subX, subY, subErrX, subErrY;
	for( Int_t iPoint=0; iPoint < fNpoints; iPoint++ ) {
		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] );
		twgaSubGraph->SetNameTitle( GetName(), GetTitle() );
		twgaSubGraph->GetXaxis()->SetTitle(this->GetXaxis()->GetTitle());
		twgaSubGraph->GetYaxis()->SetTitle(this->GetYaxis()->GetTitle());
	}
	return;
}

map<string, map<string,double> > ThreeWireGraphAnalyzer::InitRanges() {
	map<string, map<string,double> > localMap;
	{
		map<string,double> tmpMap;
		tmpMap["low"] 	= 30.0;
		tmpMap["high"] 	= 58.0;
		tmpMap["dir"]	= -1.0;
		tmpMap["scale"] = sqrt(2) ;
		localMap["X"] = tmpMap;
	}
	{
		map<string,double> tmpMap;
		tmpMap["low"] 	= 55.0;
		tmpMap["high"] 	= 74.0;
		tmpMap["dir"]	= -1.0;
		tmpMap["scale"] = 1.0  ;
		localMap["U"] = tmpMap;
	}
	{
		map<string,double> tmpMap;
		tmpMap["low"] 	= 73.0;
		tmpMap["high"] 	= 90.0;
		tmpMap["dir"]	= -1.0;
		tmpMap["scale"] = sqrt(2) ;
		localMap["Y"] = tmpMap;
	}
	return localMap;
}

void ThreeWireGraphAnalyzer::FindInitialParValues() {
	cout << "In ThreeWireGraphAnalyzer::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"] / twgaFitRange[gaLabel]["scale"],
			twgaFitRange[gaLabel]["dir"] * twgaFitRange[gaLabel]["high"] / twgaFitRange[gaLabel]["scale"] );
	double maxVal = TMath::Max( twgaFitRange[gaLabel]["dir"] * twgaFitRange[gaLabel]["low"] / twgaFitRange[gaLabel]["scale"],
			twgaFitRange[gaLabel]["dir"] * twgaFitRange[gaLabel]["high"] / twgaFitRange[gaLabel]["scale"] );
	twgaFitFun->SetParLimits(1,  minVal, maxVal );

	twgaFitFun->SetParLimits(2, 0.010, 0.5*(twgaFitRange[gaLabel]["high"] - twgaFitRange[gaLabel]["low"]) / twgaFitRange[gaLabel]["scale"] );
	twgaFitFun->SetParLimits(3, 0.0, 2*yMax );
	twgaFitFun->SetParameter(0, yMax );
	twgaFitFun->SetParameter(1, twgaFitRange[gaLabel]["dir"] * twgaSubGraph->GetMean() / twgaFitRange[gaLabel]["scale"] ) ;
	twgaFitFun->SetParameter(2, twgaSubGraph->GetRMS() / twgaFitRange[gaLabel]["scale"] / 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 ThreeWireGraphAnalyzer::Draw( Option_t* opt) {
//	if( this->GetMinimum() < 0.1 ) this->SetMinimum(0.1);
//	TGraphErrors::Draw(opt);
//	GetXaxis()->SetRangeUser( twgaFitRange[gaLabel]["low"], twgaFitRange[gaLabel]["high"]);
//	GetYaxis()->SetTitleOffset( 1.2 );
//	TGraphErrors::Draw(opt);
//	return;
//}

void ThreeWireGraphAnalyzer::Draw( Option_t* opt) {
	if(  TMath::MinElement( twgaSubGraph->GetN(),twgaSubGraph->GetY() ) < 0.1 ) twgaSubGraph->SetMinimum(0.1);
	twgaSubGraph->SetMaximum( 5.0*TMath::MaxElement( twgaSubGraph->GetN(), twgaSubGraph->GetY() ) );
	twgaSubGraph->Draw(opt);
	twgaSubGraph->GetXaxis()->SetRangeUser( twgaFitRange[gaLabel]["low"], twgaFitRange[gaLabel]["high"]);
//	twgaSubGraph->GetYaxis()->SetRangeUser( twgaSubGraph->GetMinimum(),  twgaSubGraph->GetMaximum() );
	twgaSubGraph->GetYaxis()->SetTitleOffset( 1.2 );
	twgaSubGraph->Draw(opt);
	return;
}


TFitResultPtr ThreeWireGraphAnalyzer::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 twgaSubGraph->Fit( this->twgaFitFun, "EX0", "", twgaFitRange[gaLabel]["low"], twgaFitRange[gaLabel]["high"] );
}

void ThreeWireGraphAnalyzer::CreateFitFun( string funName, double xMin, double xMax ) {
//	void* funPtr = reinterpret_cast<void*>( ThreeWireGraphAnalyzer::ThreeWireScanFitFun );
	twgaFitFun = new TF1( funName.c_str(), ThreeWireScanFitFun, xMin, xMax, 6);
	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->SetParName( 5, "Scale" );

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


double ThreeWireGraphAnalyzer::ThreeWireScanFitFun( double* xArray, double* param ) {
//	cout << "In ThreeWireGraphAnalyzer::ThreeWireScanFitFun" << 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];
	double scale = param[5];

	mean *= dirSign;

	// Scale by square root of too since the motion of the harp is at 45 degrees for x and y and 90 degrees for u.
	double x = xArray[0] / scale;
//	double x 		= xArray[0];

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

	return value;
}