calTransition47.cxx

#include <iostream>
#include <vector>       // std::vector
#include <algorithm>

using namespace std;

void calTransition47()
{

	//input files
	TString inHistFile = "../parameters/histOutput47.root";
	TString inHalfVoltPosFile = "../parameters/halfVoltPositions680.par";
//	TString inCalParFile = "../parameters/cal1e_HV1000_nonLin.par";
	TString inCalParFile = "../parameters/cal1e_HV950_nonLin.par";

	//output files
	TString outPosFile = "parameters/threeVoltPositions47.par";
//	TString outCalFile = "parameters/cal1e_47_HV1000_nonLin.par";
	TString outCalFile = "parameters/cal1e_47_HV950_nonLin.par";

	TFile *fr = new TFile(inHistFile.Data());
	if (fr->IsZombie()) {
		Error("calTransition47.cxx", "File %s was not open correctly.", inHistFile.Data());
		return;
	}

	ifstream inCalFile;
	inCalFile.open(inCalParFile.Data());
	if (!inCalFile.is_open()) {
		Error("calTransition47.cxx", "File with input parameters %s was not open.", inCalParFile.Data());
		return;
	}
	Info("calTransition47.cxx", "File with 680 pF calpars %s will be loaded.", inCalParFile.Data());

	ifstream inHalfVoltFile;
	inHalfVoltFile.open(inHalfVoltPosFile.Data());
	if (!inHalfVoltFile.is_open()) {
		Error("cal1ePeaks.cxx", "File with input parameters %s was not open.", inHalfVoltPosFile.Data());
		return;
	}
	Info("calTransition47.cxx", "File with 0.5 V positions at 680 pF %s will be loaded.", inHalfVoltPosFile.Data());


	//output files
	ofstream calPosFile;
	calPosFile.open(outPosFile.Data());
	if (calPosFile.is_open()) {
		Info("calTransition47.cxx", "Positions of 3.0V peak at 4.7 pF will be saved in %s.", outPosFile.Data());
	} else {
		Error("calTransition47.cxx", "File %s was not open.", outPosFile.Data());
		return;
	}

	ofstream calFile;
	calFile.open(outCalFile.Data());
	if (calFile.is_open()) {
		Info("calTransition47.cxx", "Calibration parameters for 4.7 nF will be saved in %s.", outCalFile.Data());
	} else {
		Error("calTransition47.cxx", "File %s was not open.", outCalFile.Data());
		return;
	}


	TH1D *hch[64];
	TGraph *gr[64];
	TString hName;

	for (Int_t i = 0; i <64; i++) {
		hName.Form("ch%d", i);
		hch[i] = (TH1D*)fr->Get(hName.Data());
	}
	


	TString cName;
	TString cTitle;
	
	TCanvas *cWork[11];
	for (Int_t i = 0; i < 11; i++) {
		cWork[i] = new TCanvas();
		cName.Form("can%d", i);
		cTitle.Form("canvas %d", i);
		cWork[i]->SetName(cName.Data());
		cWork[i]->SetTitle(cTitle.Data());
		cWork[i]->Divide(2, 3);
		cWork[i]->cd(1);
	}



	Double_t x, y;
	Double_t relHight = 0.5;

	string line;
	Int_t ch = 0;
	Double_t calPar = 0.;
	Double_t halfVoltPos = 0.;

	for (Int_t i = 0; i < 64; i++) {
		hch[i]->GetXaxis()->SetRangeUser(100, 800);
		Int_t cNumber = i/6;
		cWork[cNumber]->cd(i-6*cNumber+1);
		hch[i]->Draw();

		hch[i]->ShowPeaks(2., "", relHight);

		TList *functions = hch[i]->GetListOfFunctions();
		TPolyMarker *pm = (TPolyMarker*)functions->FindObject("TPolyMarker");

		const Int_t noPeaks = pm->GetN();
		if (noPeaks != 1) {
			Error("calTransition47.cxx", "%d peaks found in hist %s", noPeaks, hch[i]->GetName());
			continue;
		}

//		Info("calTransition47.cxx", "%d peak found in hist %s", noPeaks, hch[i]->GetName());

		Double_t *energy = pm->GetX();
		calPosFile << i << "\t" << energy[0] << endl;

		std::getline(inCalFile, line);
		std::istringstream(line,ios_base::in) >> ch >> calPar;

		std::getline(inHalfVoltFile, line);
		std::istringstream(line,ios_base::in) >> ch >> halfVoltPos;

//		cout << i << "\t" << energy[0] << "\t" << calPar << "\t\t" << halfVoltPos << "\t" << energy[0]*calPar / (6*halfVoltPos) << endl;
//		cout << energy[0]*calPar / (6*halfVoltPos) << endl;
		calFile << i << "\t" << energy[0]*calPar / (6*halfVoltPos) << endl;

	}


	calFile.close();
	inCalFile.close();
	inHalfVoltFile.close();

//	fr->Close();

	return;
}