[u/mrichter/AliRoot.git] / TPC / CalibMacros / CalibTimeGain.C

/// \file CalibTimeGain.C
///
/// Macro to generate and update OCDB entries for a given run:
/// this is a TObjArray which has at 0 the MIP-Spline and at 1 the Fermi-Plateau-Spline ...
///
/// \author marian.ivanov@cern.ch
/// \author A.Kalweit@gsi.de
///
/// How to use it locally:
///
/// ~~~{.cpp}
/// // Load libraries
/// gSystem->Load("libANALYSIS");
/// gSystem->Load("libSTAT");
/// gSystem->Load("libTPCcalib");
/// gSystem->AddIncludePath("-I$ALICE_ROOT/TPC");
/// .L $ALICE_ROOT/TPC/CalibMacros/CalibTimeGain.C+
///
/// //Make calibration
/// CalibTimeGain("CalibObjectsTrain1.root",0,120000);
/// TBrowser b;
/// b.Add(gainArray);
/// ~~~

#if !defined(__CINT__) || defined(__MAKECINT__)
#include "TObjArray.h"
#include "TGraphErrors.h"
#include "THnSparse.h"
#include "TCanvas.h"
#include "TFile.h"

#include "AliTPCcalibTimeGain.h"
#include "AliSplineFit.h"
#include "AliCDBMetaData.h"
#include "AliCDBId.h"
#include "AliCDBManager.h"
#include "AliCDBStorage.h"
#endif


TGraphErrors * graphMIP    = 0;       // graph time dependence of MIP
TGraphErrors * graphCosmic = 0;       // graph time dependence at Plateu
AliSplineFit * fitMIP = 0;            // fit of dependence - MIP
AliSplineFit * fitCosmic = 0;         // fit of dependence - Plateu
TObjArray    * gainArray = new TObjArray(4); // array to be stored in the OCDB
AliTPCcalibTimeGain * gainMIP =0;     // calibration component for MIP
AliTPCcalibTimeGain * gainCosmic =0;  // calibration component for cosmic


void UpdateOCDBGain(Int_t  startRunNumber, Int_t endRunNumber, const char* storagePath);
void ReadGainGlobal(Char_t* fileName="CalibObjectsTrain1.root");
void MakeQAPlot(Float_t  FPtoMIPratio);
Bool_t AnalyzeGain(Int_t startRunNumber, Int_t endRunNumber, Int_t minEntriesGaussFit = 500, Float_t FPtoMIPratio = 1.43); 


void CalibTimeGain(Char_t* fileName="CalibObjectsTrain1.root", Int_t startRunNumber=0, Int_t endRunNumber=AliCDBRunRange::Infinity(),  TString  ocdbStorage=""){
  /// Update OCDB gain

  ReadGainGlobal(fileName);
  AnalyzeGain(startRunNumber,endRunNumber, 1000,1.43);
  MakeQAPlot(1.43);  
  if (ocdbStorage.Length()==0) ocdbStorage+="local://"+gSystem->GetFromPipe("pwd")+"/OCDB";
  UpdateOCDBGain( startRunNumber, endRunNumber, ocdbStorage.Data());
}


void ReadGainGlobal(Char_t* fileName){
  /// read calibration entries from file

  TFile fcalib(fileName);
  TObjArray * array = (TObjArray*)fcalib.Get("TPCCalib");
  if (array){
    gainMIP    = ( AliTPCcalibTimeGain *)array->FindObject("calibTimeGain");
    gainCosmic = ( AliTPCcalibTimeGain *)array->FindObject("calibTimeGainCosmic");
  }else{
    gainMIP    = ( AliTPCcalibTimeGain *)fcalib.Get("calibTimeGain");
    gainCosmic = ( AliTPCcalibTimeGain *)fcalib.Get("calibTimeGainCosmic");
  }
  TH1 * hisT=0;
  Int_t firstBinA =0, lastBinA=0;

  if (gainCosmic){ 
    hisT= gainCosmic->GetHistGainTime()->Projection(1);
    firstBinA = hisT->FindFirstBinAbove(2);
    lastBinA  = hisT->FindLastBinAbove(2);    
    gainCosmic->GetHistGainTime()->GetAxis(1)->SetRange(firstBinA,lastBinA);
    delete hisT;
  }

  if (gainMIP){ 
    hisT= gainCosmic->GetHistGainTime()->Projection(1);
    firstBinA = hisT->FindFirstBinAbove(2);
    lastBinA  = hisT->FindLastBinAbove(2);    
    gainMIP->GetHistGainTime()->GetAxis(1)->SetRange(firstBinA,lastBinA);
    delete hisT;
  }

}


Bool_t AnalyzeGain(Int_t startRunNumber, Int_t endRunNumber, Int_t minEntriesGaussFit,  Float_t FPtoMIPratio){
  ///

  gainMIP->GetHistGainTime()->GetAxis(5)->SetRangeUser(startRunNumber, endRunNumber);
  // 1.) try to create MIP spline
  gainMIP->GetHistGainTime()->GetAxis(2)->SetRangeUser(1.51,2.49); // only beam data
  gainMIP->GetHistGainTime()->GetAxis(4)->SetRangeUser(0.39,0.51); // only MIP pions
  //
  graphMIP = AliTPCcalibBase::FitSlices(gainMIP->GetHistGainTime(),0,1,minEntriesGaussFit,10,0.1,0.7);
  if (graphMIP->GetN()==0) graphMIP = 0x0;
  if (graphMIP) fitMIP = AliTPCcalibTimeGain::MakeSplineFit(graphMIP);
  if (graphMIP) graphMIP->SetName("TGRAPHERRORS_MEAN_GAIN_BEAM_ALL");// set proper names according to naming convention
  gainArray->AddAt(fitMIP,0);
  

  // 2.) try to create Cosmic spline
  gainCosmic->GetHistGainTime()->GetAxis(2)->SetRangeUser(0.51,1.49); // only cosmics
  gainCosmic->GetHistGainTime()->GetAxis(4)->SetRangeUser(20,100);    // only Fermi-Plateau muons
  //
  graphCosmic = AliTPCcalibBase::FitSlices(gainCosmic->GetHistGainTime(),0,1,minEntriesGaussFit,10);
  if (graphCosmic->GetN()==0) graphCosmic = 0x0;
  //
  if (graphCosmic) {
    for(Int_t i=0; i < graphCosmic->GetN(); i++) {
      graphCosmic->GetY()[i] /= FPtoMIPratio;
      graphCosmic->GetEY()[i] /= FPtoMIPratio;
    }
  }
  //
  if (graphCosmic) fitCosmic = AliTPCcalibTimeGain::MakeSplineFit(graphCosmic);
  if (graphCosmic) graphCosmic->SetName("TGRAPHERRORS_MEAN_GAIN_COSMIC_ALL"); // set proper names according to naming convention
  gainArray->AddAt(fitCosmic,1);
  // with naming convention and backward compatibility
  gainArray->AddAt(graphMIP,2);
  gainArray->AddAt(graphCosmic,3);
  cout << "graphCosmic: " << graphCosmic << " graphMIP " << graphMIP << endl;
  return kTRUE;

}


void UpdateOCDBGain(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
  /// Update OCDB entry

  AliCDBMetaData *metaData= new AliCDBMetaData();
  metaData->SetObjectClassName("TObjArray");
  metaData->SetResponsible("Alexander Kalweit");
  metaData->SetBeamPeriod(1);
  metaData->SetAliRootVersion("05-24-00"); //root version
  metaData->SetComment("Calibration of the time dependence of the gain due to pressure and temperature changes.");
  AliCDBId id1("TPC/Calib/TimeGain", startRunNumber, endRunNumber);
  AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
  gStorage->Put(gainArray, id1, metaData);    
}

void MakeQAPlot(Float_t  FPtoMIPratio) {
  /// Make QA plot to visualize results

  if (graphCosmic) {
    TCanvas * canvasCosmic = new TCanvas("gain Cosmic", "time dependent gain QA histogram cosmic");
    canvasCosmic->cd();
    TH2D * gainHistoCosmic = gainCosmic->GetHistGainTime()->Projection(0,1);
    gainHistoCosmic->SetDirectory(0);
    gainHistoCosmic->SetName("GainHistoCosmic");
    gainHistoCosmic->GetXaxis()->SetTimeDisplay(kTRUE);
    gainHistoCosmic->GetXaxis()->SetTimeFormat("#splitline{%d/%m}{%H:%M}");
    gainHistoCosmic->Draw("colz");
    graphCosmic->SetMarkerStyle(25);
    graphCosmic->Draw("lp");
    graphCosmic->SetMarkerStyle(25);
    TGraph * grfitCosmic = fitCosmic->MakeGraph(graphCosmic->GetX()[0],graphCosmic->GetX()[graphCosmic->GetN()-1],50000,0);
    if (grfitCosmic) {
      for(Int_t i=0; i < grfitCosmic->GetN(); i++) {
 	grfitCosmic->GetY()[i] *= FPtoMIPratio;	
      }
      for(Int_t i=0; i < graphCosmic->GetN(); i++) {
 	graphCosmic->GetY()[i] *= FPtoMIPratio;	
      }
    }
    graphCosmic->Draw("lp");
    grfitCosmic->SetLineColor(2);
    grfitCosmic->Draw("lu");
    gainArray->AddLast(gainHistoCosmic);
    gainArray->AddLast(canvasCosmic->Clone());
    delete canvasCosmic;    
  }
  if (fitMIP) {
    TCanvas * canvasMIP = new TCanvas("gain MIP", "time dependent gain QA histogram MIP");
    canvasMIP->cd();
    TH2D * gainHistoMIP    = gainMIP->GetHistGainTime()->Projection(0,1);
    gainHistoMIP->SetName("GainHistoCosmic");
    gainHistoMIP->SetDirectory(0);
    gainHistoMIP->GetXaxis()->SetTimeDisplay(kTRUE);
    gainHistoMIP->GetXaxis()->SetTimeFormat("#splitline{%d/%m}{%H:%M}");
    gainHistoMIP->Draw("colz");
    graphMIP->SetMarkerStyle(25);
    graphMIP->Draw("lp");
    TGraph * grfitMIP = fitMIP->MakeGraph(graphMIP->GetX()[0],graphMIP->GetX()[graphMIP->GetN()-1],50000,0);
    grfitMIP->SetLineColor(2);
    grfitMIP->Draw("lu");    
    gainArray->AddLast(gainHistoMIP);
    gainArray->AddLast(canvasMIP->Clone());
    delete canvasMIP;    
  }  
}
Commit	Line	Data
35d81915	1	/// \file CalibTimeGain.C
	2	///
	3	/// Macro to generate and update OCDB entries for a given run:
	4	/// this is a TObjArray which has at 0 the MIP-Spline and at 1 the Fermi-Plateau-Spline ...
	5	///
	6	/// \author marian.ivanov@cern.ch
	7	/// \author A.Kalweit@gsi.de
	8	///
	9	/// How to use it locally:
	10	///
	11	/// ~~~{.cpp}
	12	/// // Load libraries
	13	/// gSystem->Load("libANALYSIS");
	14	/// gSystem->Load("libSTAT");
	15	/// gSystem->Load("libTPCcalib");
	16	/// gSystem->AddIncludePath("-I$ALICE_ROOT/TPC");
	17	/// .L $ALICE_ROOT/TPC/CalibMacros/CalibTimeGain.C+
	18	///
	19	/// //Make calibration
	20	/// CalibTimeGain("CalibObjectsTrain1.root",0,120000);
	21	/// TBrowser b;
	22	/// b.Add(gainArray);
	23	/// ~~~
	24
12068303	25	#if !defined(__CINT__) \|\| defined(__MAKECINT__)
b6d506d5	26	#include "TObjArray.h"
	27	#include "TGraphErrors.h"
	28	#include "THnSparse.h"
	29	#include "TCanvas.h"
12068303	30	#include "TFile.h"
b6d506d5	31
	32	#include "AliTPCcalibTimeGain.h"
	33	#include "AliSplineFit.h"
	34	#include "AliCDBMetaData.h"
	35	#include "AliCDBId.h"
	36	#include "AliCDBManager.h"
	37	#include "AliCDBStorage.h"
12068303	38	#endif
	39
	40
	41
	42
	43	TGraphErrors * graphMIP = 0; // graph time dependence of MIP
	44	TGraphErrors * graphCosmic = 0; // graph time dependence at Plateu
	45	AliSplineFit * fitMIP = 0; // fit of dependence - MIP
	46	AliSplineFit * fitCosmic = 0; // fit of dependence - Plateu
	47	TObjArray * gainArray = new TObjArray(4); // array to be stored in the OCDB
	48	AliTPCcalibTimeGain * gainMIP =0; // calibration component for MIP
	49	AliTPCcalibTimeGain * gainCosmic =0; // calibration component for cosmic
	50
	51
	52	void UpdateOCDBGain(Int_t startRunNumber, Int_t endRunNumber, const char* storagePath);
	53	void ReadGainGlobal(Char_t* fileName="CalibObjectsTrain1.root");
	54	void MakeQAPlot(Float_t FPtoMIPratio);
	55	Bool_t AnalyzeGain(Int_t startRunNumber, Int_t endRunNumber, Int_t minEntriesGaussFit = 500, Float_t FPtoMIPratio = 1.43);
	56
	57
	58
bade7d8d	59	void CalibTimeGain(Char_t* fileName="CalibObjectsTrain1.root", Int_t startRunNumber=0, Int_t endRunNumber=AliCDBRunRange::Infinity(), TString ocdbStorage=""){
35d81915	60	/// Update OCDB gain
35d81915	61
12068303	62	ReadGainGlobal(fileName);
	63	AnalyzeGain(startRunNumber,endRunNumber, 1000,1.43);
	64	MakeQAPlot(1.43);
	65	if (ocdbStorage.Length()==0) ocdbStorage+="local://"+gSystem->GetFromPipe("pwd")+"/OCDB";
	66	UpdateOCDBGain( startRunNumber, endRunNumber, ocdbStorage.Data());
	67	}
b6d506d5	68
b6d506d5	69
b6d506d5	70
b6d506d5	71
12068303	72	void ReadGainGlobal(Char_t* fileName){
35d81915	73	/// read calibration entries from file
35d81915	74
12068303	75	TFile fcalib(fileName);
	76	TObjArray * array = (TObjArray*)fcalib.Get("TPCCalib");
	77	if (array){
	78	gainMIP = ( AliTPCcalibTimeGain *)array->FindObject("calibTimeGain");
	79	gainCosmic = ( AliTPCcalibTimeGain *)array->FindObject("calibTimeGainCosmic");
	80	}else{
	81	gainMIP = ( AliTPCcalibTimeGain *)fcalib.Get("calibTimeGain");
	82	gainCosmic = ( AliTPCcalibTimeGain *)fcalib.Get("calibTimeGainCosmic");
	83	}
	84	TH1 * hisT=0;
	85	Int_t firstBinA =0, lastBinA=0;
	86
	87	if (gainCosmic){
	88	hisT= gainCosmic->GetHistGainTime()->Projection(1);
	89	firstBinA = hisT->FindFirstBinAbove(2);
	90	lastBinA = hisT->FindLastBinAbove(2);
	91	gainCosmic->GetHistGainTime()->GetAxis(1)->SetRange(firstBinA,lastBinA);
	92	delete hisT;
	93	}
	94
	95	if (gainMIP){
	96	hisT= gainCosmic->GetHistGainTime()->Projection(1);
	97	firstBinA = hisT->FindFirstBinAbove(2);
	98	lastBinA = hisT->FindLastBinAbove(2);
	99	gainMIP->GetHistGainTime()->GetAxis(1)->SetRange(firstBinA,lastBinA);
	100	delete hisT;
	101	}
	102
	103	}
	104
	105
	106
	107	Bool_t AnalyzeGain(Int_t startRunNumber, Int_t endRunNumber, Int_t minEntriesGaussFit, Float_t FPtoMIPratio){
35d81915	108	///
35d81915	109
12068303	110	gainMIP->GetHistGainTime()->GetAxis(5)->SetRangeUser(startRunNumber, endRunNumber);
b6d506d5	111	// 1.) try to create MIP spline
12068303	112	gainMIP->GetHistGainTime()->GetAxis(2)->SetRangeUser(1.51,2.49); // only beam data
12068303	113	gainMIP->GetHistGainTime()->GetAxis(4)->SetRangeUser(0.39,0.51); // only MIP pions
b6d506d5	114	//
12068303	115	graphMIP = AliTPCcalibBase::FitSlices(gainMIP->GetHistGainTime(),0,1,minEntriesGaussFit,10,0.1,0.7);
b6d506d5	116	if (graphMIP->GetN()==0) graphMIP = 0x0;
b6d506d5	117	if (graphMIP) fitMIP = AliTPCcalibTimeGain::MakeSplineFit(graphMIP);
b6d506d5	118	if (graphMIP) graphMIP->SetName("TGRAPHERRORS_MEAN_GAIN_BEAM_ALL");// set proper names according to naming convention
12068303	119	gainArray->AddAt(fitMIP,0);
b6d506d5	120
12068303	121
b6d506d5	122	// 2.) try to create Cosmic spline
12068303	123	gainCosmic->GetHistGainTime()->GetAxis(2)->SetRangeUser(0.51,1.49); // only cosmics
12068303	124	gainCosmic->GetHistGainTime()->GetAxis(4)->SetRangeUser(20,100); // only Fermi-Plateau muons
b6d506d5	125	//
12068303	126	graphCosmic = AliTPCcalibBase::FitSlices(gainCosmic->GetHistGainTime(),0,1,minEntriesGaussFit,10);
b6d506d5	127	if (graphCosmic->GetN()==0) graphCosmic = 0x0;
12068303	128	//
	129	if (graphCosmic) {
	130	for(Int_t i=0; i < graphCosmic->GetN(); i++) {
	131	graphCosmic->GetY()[i] /= FPtoMIPratio;
	132	graphCosmic->GetEY()[i] /= FPtoMIPratio;
	133	}
	134	}
	135	//
b6d506d5	136	if (graphCosmic) fitCosmic = AliTPCcalibTimeGain::MakeSplineFit(graphCosmic);
b6d506d5	137	if (graphCosmic) graphCosmic->SetName("TGRAPHERRORS_MEAN_GAIN_COSMIC_ALL"); // set proper names according to naming convention
12068303	138	gainArray->AddAt(fitCosmic,1);
b6d506d5	139	// with naming convention and backward compatibility
12068303	140	gainArray->AddAt(graphMIP,2);
12068303	141	gainArray->AddAt(graphCosmic,3);
b6d506d5	142	cout << "graphCosmic: " << graphCosmic << " graphMIP " << graphMIP << endl;
12068303	143	return kTRUE;
	144
	145	}
	146
	147
	148
	149	void UpdateOCDBGain(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
35d81915	150	/// Update OCDB entry
35d81915	151
b6d506d5	152	AliCDBMetaData *metaData= new AliCDBMetaData();
	153	metaData->SetObjectClassName("TObjArray");
	154	metaData->SetResponsible("Alexander Kalweit");
	155	metaData->SetBeamPeriod(1);
	156	metaData->SetAliRootVersion("05-24-00"); //root version
	157	metaData->SetComment("Calibration of the time dependence of the gain due to pressure and temperature changes.");
	158	AliCDBId id1("TPC/Calib/TimeGain", startRunNumber, endRunNumber);
	159	AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
12068303	160	gStorage->Put(gainArray, id1, metaData);
12068303	161	}
b6d506d5	162
12068303	163	void MakeQAPlot(Float_t FPtoMIPratio) {
35d81915	164	/// Make QA plot to visualize results
35d81915	165
12068303	166	if (graphCosmic) {
	167	TCanvas * canvasCosmic = new TCanvas("gain Cosmic", "time dependent gain QA histogram cosmic");
	168	canvasCosmic->cd();
	169	TH2D * gainHistoCosmic = gainCosmic->GetHistGainTime()->Projection(0,1);
	170	gainHistoCosmic->SetDirectory(0);
	171	gainHistoCosmic->SetName("GainHistoCosmic");
	172	gainHistoCosmic->GetXaxis()->SetTimeDisplay(kTRUE);
	173	gainHistoCosmic->GetXaxis()->SetTimeFormat("#splitline{%d/%m}{%H:%M}");
	174	gainHistoCosmic->Draw("colz");
	175	graphCosmic->SetMarkerStyle(25);
	176	graphCosmic->Draw("lp");
	177	graphCosmic->SetMarkerStyle(25);
	178	TGraph * grfitCosmic = fitCosmic->MakeGraph(graphCosmic->GetX()[0],graphCosmic->GetX()[graphCosmic->GetN()-1],50000,0);
	179	if (grfitCosmic) {
	180	for(Int_t i=0; i < grfitCosmic->GetN(); i++) {
	181	grfitCosmic->GetY()[i] *= FPtoMIPratio;
	182	}
	183	for(Int_t i=0; i < graphCosmic->GetN(); i++) {
	184	graphCosmic->GetY()[i] *= FPtoMIPratio;
	185	}
b6d506d5	186	}
12068303	187	graphCosmic->Draw("lp");
	188	grfitCosmic->SetLineColor(2);
	189	grfitCosmic->Draw("lu");
	190	gainArray->AddLast(gainHistoCosmic);
	191	gainArray->AddLast(canvasCosmic->Clone());
	192	delete canvasCosmic;
b6d506d5	193	}
12068303	194	if (fitMIP) {
	195	TCanvas * canvasMIP = new TCanvas("gain MIP", "time dependent gain QA histogram MIP");
	196	canvasMIP->cd();
	197	TH2D * gainHistoMIP = gainMIP->GetHistGainTime()->Projection(0,1);
	198	gainHistoMIP->SetName("GainHistoCosmic");
	199	gainHistoMIP->SetDirectory(0);
	200	gainHistoMIP->GetXaxis()->SetTimeDisplay(kTRUE);
	201	gainHistoMIP->GetXaxis()->SetTimeFormat("#splitline{%d/%m}{%H:%M}");
	202	gainHistoMIP->Draw("colz");
	203	graphMIP->SetMarkerStyle(25);
	204	graphMIP->Draw("lp");
	205	TGraph * grfitMIP = fitMIP->MakeGraph(graphMIP->GetX()[0],graphMIP->GetX()[graphMIP->GetN()-1],50000,0);
	206	grfitMIP->SetLineColor(2);
	207	grfitMIP->Draw("lu");
	208	gainArray->AddLast(gainHistoMIP);
	209	gainArray->AddLast(canvasMIP->Clone());
	210	delete canvasMIP;
	211	}
b6d506d5	212	}
	213
	214