* provided "as is" without express or implied warranty. *
**************************************************************************/
-/* $Log$ */
+/* $Log$co: warning: `/* $Log' is obsolescent; use ` * $Log'.
+
+ * Revision 1.1 2007/08/08 15:58:01 hristov
+ * */
+
//_________________________________________________________________________
// Calibration coefficients
#include "AliEMCALGeometry.h"
#include "TH1F.h"
+#include <TString.h>
-TableClassImpl(AliEMCALCalibCoefs,calibCoef)
+ClassImp(calibCoef)
-// Get initial Calib Data from DB
-AliEMCALCalibCoefs* AliEMCALCalibCoefs::GetCalibTableFromDb(const char *tn)
+calibCoef::calibCoef() : absId(-1), cc(-1), eCc(-1)
+{ }
+calibCoef::calibCoef(const Int_t id, const Double_t c, const Double_t ec) :
+absId(id), cc(c), eCc(ec)
+{}
+// ----------------------------------------------------------------------
+
+ClassImp(AliEMCALCalibCoefs)
+
+AliEMCALCalibCoefs::AliEMCALCalibCoefs() : TNamed("",""), fTable(0), fCurrentInd(0), fCalibMethod(0)
+{
+}
+
+AliEMCALCalibCoefs::AliEMCALCalibCoefs(const char* name, const Int_t nrow) : TNamed(name,"table of cell information") , fTable(0), fCurrentInd(0), fCalibMethod(0)
+{
+ fTable = new TObjArray(nrow);
+}
+
+void AliEMCALCalibCoefs::AddAt(calibCoef* r)
{
- // Initial cc with decalibration
- char* dbString = "local:///data/r22b/ALICE/PROD/CALIBRATION_May_2007/PI0/PDSF/10GEV/DECALIB/DeCalibDB";
+ (*fTable)[fCurrentInd] = new calibCoef(*r);
+ fCurrentInd++;
+}
+
+AliEMCALCalibCoefs::~AliEMCALCalibCoefs()
+{
+ if(fTable) {
+ fTable->Delete();
+ delete fTable;
+ }
+}
+
+calibCoef* AliEMCALCalibCoefs::GetTable(Int_t i) const
+{
+ return (calibCoef*)fTable->At(i);
+}
+
+
+// Get initial Calib Data from DB
+AliEMCALCalibCoefs* AliEMCALCalibCoefs::GetCalibTableFromDb(const char *tn, AliEMCALCalibData **calData)
+{
+ //
+ // See ~/macros/ALICE/sim.C for choice of CDB
+ // Get calib. table which was used than calculated rec.points
+ //
+ static char *calibType = "EMCAL/Calib/*";
+ static char *calibTypeData = "EMCAL/Calib/Data";
+ // Initial cc
+ calData[0] = 0;
+
+ AliCDBManager* man = AliCDBManager::Instance();
+ AliCDBStorage* specificStorage = man->GetSpecificStorage(calibType);
+
AliEMCALCalibData* caldata = (AliEMCALCalibData*)
- (AliCDBManager::Instance()->GetStorage(dbString)
- ->Get("EMCAL/Calib/Data",gAlice->GetRunNumber())->GetObject());
+ specificStorage->Get(calibTypeData, gAlice->GetRunNumber())->GetObject();
if(caldata == 0) return 0;
- AliEMCALCalibCoefs *tab = new AliEMCALCalibCoefs(tn);
+ AliEMCALGeometry *g = AliEMCALGeometry::GetInstance();
+
+ AliEMCALCalibCoefs *tab = new AliEMCALCalibCoefs(tn,g->GetNCells());
tab->SetCalibMethod(AliEMCALCalibCoefs::kMC);
- AliEMCALGeometry *g = AliEMCALGeometry::GetInstance();
for(Int_t id=0; id<g->GetNCells(); id++){
Int_t nSupMod=0, nModule=0, nIphi=0, nIeta=0, iphiCell=0, ietaCell=0;
g->GetCellIndex(id, nSupMod, nModule, nIphi, nIeta);
calibCoef r;
r.absId = id;
- r.cc = caldata->GetADCchannel(nSupMod, ietaCell, iphiCell);
+ r.cc = caldata->GetADCchannel(nSupMod, ietaCell, iphiCell); // column(ietaCell) : row(iphiCell)
r.eCc = 0.0;
tab->AddAt(&r);
}
- tab->Purge();
+ // tab->Purge();
tab->SetCalibMethod(AliEMCALCalibCoefs::kPI0);
+ calData[0] = caldata;
+
+ printf("\n <I> AliEMCALCalibCoefs::GetCalibTableFromDb \n name |%s| title |%s| | storage |%s|\n",
+ caldata->GetName(), caldata->GetTitle(), specificStorage->GetURI().Data());
return tab;
}
TH1F* AliEMCALCalibCoefs::GetHistOfCalibTableFromDb(const char *tn)
{
// First SM only
- AliEMCALGeometry *g = AliEMCALGeometry::GetInstance("");
+ AliEMCALCalibData *calData[1]; // unused here
- AliEMCALCalibCoefs* tab = GetCalibTableFromDb(tn);
+ AliEMCALCalibCoefs* tab = GetCalibTableFromDb(tn, calData);
if(tab==0) return 0;
TH1F *h = new TH1F("hCCfirst", " cc first (in MeV)", 70, 12., 19.);
calibCoef *r;
- for(Int_t i=0; i<tab->GetNRows(); i++){
+ for(Int_t i=0; i<tab->GetSize(); i++){
r = tab->GetTable(i);
if(i>=1152) break;
- h->Fill(r->cc*1000.);
+ h->Fill(r->cc*1000.); // GEV ->MeV
}
delete tab;
return h;
}
+AliEMCALCalibData* AliEMCALCalibCoefs::GetCalibTableForDb(const AliEMCALCalibCoefs *tab, const char* dbLocation,
+const char* coment)
+{
+ printf("<I> AliEMCALCalibCoefs::GetCalibTableForDb started \n");
+ // See EMCAL/macros/CalibrationDB/AliEMCALSetCDB.C
+ // Define and save to CDB
+ AliEMCALCalibData* caldata=0;
+ if(tab==0) return caldata;
+
+ Int_t firstRun = 0;
+ Int_t lastRun = 10;
+ Int_t beamPeriod = 1;
+ char* objFormat = "";
+ caldata = new AliEMCALCalibData("EMCAL");
+ caldata->SetTitle(coment);
+
+ AliEMCALGeometry *g = AliEMCALGeometry::GetInstance();
+
+ for(int id=0; id<tab->GetSize(); id++){
+ Float_t ped=0.;
+
+ Int_t nSupMod=0, nModule=0, nIphi=0, nIeta=0, iphiCell=0, ietaCell=0;
+ g->GetCellIndex(id, nSupMod, nModule, nIphi, nIeta);
+ g->GetCellPhiEtaIndexInSModule(nSupMod, nModule, nIphi, nIeta, iphiCell, ietaCell);
+
+ calibCoef *r = tab->GetTable(id);
+ // ietaCell - column; iphiCell - row
+ caldata->SetADCchannel (nSupMod, ietaCell, iphiCell, r->cc);
+ caldata->SetADCpedestal(nSupMod, ietaCell, iphiCell, ped);
+ }
+ printf("<I> Fill AliEMCALCalibData table \n");
+ //Store calibration data into database
+
+ AliCDBMetaData md;
+ md.SetComment(objFormat);
+ md.SetBeamPeriod(beamPeriod);
+ md.SetResponsible("Aleksei Pavlinov");
+
+ AliCDBManager* man = AliCDBManager::Instance();
+ if(man == 0) {
+ printf("<E> AliEMCALCalibCoefs::GetCalibTableForDb : define AliCDBManager, NO saving !! \n");
+ } else {
+ printf("<I> AliCDBManager %p \n", man);
+ AliCDBId id("EMCAL/Calib/Data",firstRun,lastRun); // create in EMCAL/Calib/Data DBFolder
+ TString DBFolder(dbLocation);
+ AliCDBStorage* loc = man->GetStorage(DBFolder.Data());
+ loc->Put(caldata, id, &md);
+ }
+
+ return caldata;
+}
+
calibCoef *AliEMCALCalibCoefs::GetRow(const int absId)
{
calibCoef *r=0;
- for(int id=0; id<GetNRows(); id++){
+ for(int id=0; id<fTable->GetSize(); id++){
r = GetTable(id);
if(r->absId == absId) return r;
}
void AliEMCALCalibCoefs::PrintTable()
{
- printf(" Table : %s : nrows %i \n", GetName(), int(GetNRows()));
- for(int i=0; i<GetNRows(); i++) PrintTable(i);
+ printf(" Table : %s : nrows %i \n", GetName(), int(fTable->GetSize()));
+ for(int i=0; i<fTable->GetSize(); i++) PrintTable(i);
}
void AliEMCALCalibCoefs::PrintTable(const Int_t i)
{
- if(i>=GetNRows()) return;
+ if(i>=fTable->GetSize()) return;
printf("row %i \n", i);
PrintRec(GetTable(i));
}
/* $Id$ */
//_________________________________________________________________________
-// Table of Calibration coefficients
+// Table of Calibration coefficients
+// Should be extended.
//
//*-- Author: Aleksei Pavlinov (WSU, Detroit, USA)
// --- ROOT system ---
-#include <TTable.h>
+#include <TNamed.h>
+#include <TObjArray.h>
// unit is GeV
-struct calibCoef {
+class calibCoef : public TObject {
+ public:
+ virtual const char* GetName() const {return Form("CC%5.5i",absId);}
+ calibCoef();
+ calibCoef(const Int_t id, const Double_t c, const Double_t ec);
+ virtual ~calibCoef() {};
+
Int_t absId; // absolute id of cell
Double_t cc; // Calib. coef
Double_t eCc; // Calib. coef. error
+ //
+ ClassDef(calibCoef,1) // Cell calibration information
};
class TH1F;
+class AliEMCALCalibData;
-class AliEMCALCalibCoefs : public TTable {
+class AliEMCALCalibCoefs : public TNamed {
public:
enum EEmcalCalibType {kMC, kEQUALIZATION, kMIP, kPI0}; // type of EMCAL calibrations
+ AliEMCALCalibCoefs(); // default constractor
+ AliEMCALCalibCoefs(const char* name, const Int_t nrow);
+ virtual ~AliEMCALCalibCoefs();
+
+ AliEMCALCalibCoefs & operator = (const AliEMCALCalibCoefs & /*rvalue*/) {
+ // assignement operator requested by coding convention but not needed
+ Fatal("operator =", "not implemented");
+ return *this;
+ };
+ void AddAt(calibCoef* r);
+ calibCoef* GetTable(Int_t i) const;
+ Int_t GetSize() const {return fTable->GetSize();}
+ Int_t GetNRows() const {return fCurrentInd;}
+ void Purge() {/* nothing */};
+
void SetCalibMethod(Int_t var) {fCalibMethod=var;}
Int_t GetCalibMethod() {return fCalibMethod;}
calibCoef* GetRow(const int absId);
// Get initial Calib Data from DB
- static AliEMCALCalibCoefs *GetCalibTableFromDb(const char *tn="CCIN");
- static TH1F *GetHistOfCalibTableFromDb(const char *tn="CCTMP");
+ static AliEMCALCalibCoefs *GetCalibTableFromDb(const char *tn="CCIN", AliEMCALCalibData **calData=0);
+ // const char* dbLocation="local:///data/r22b/ALICE/PROD/CALIBRATION_May_2007/PI0/PDSF/10GEV/DECALIB/DeCalibDB");
+ static TH1F *GetHistOfCalibTableFromDb(const char *tn="CCIN");
+ //const char* dbLocation="local:///data/r22b/ALICE/PROD/CALIBRATION_May_2007/PI0/PDSF/10GEV/DECALIB/DeCalibDB");
+ // Create DB calib table
+ static AliEMCALCalibData* GetCalibTableForDb(const AliEMCALCalibCoefs *tab,
+ const char* dbLocation="local://.", // current directory
+ const char* coment="pi0 calibration, MC, Jun 2007");
// Menu
void PrintTable(); // *MENU*
void PrintTable(const Int_t i); // *MENU*
void PrintRec(calibCoef *r);
+ //
protected:
+ TObjArray *fTable;
+ Int_t fCurrentInd;
Int_t fCalibMethod; // method of calibration - EEmcalCalibType
- ClassDefTable(AliEMCALCalibCoefs , calibCoef)
- ClassDef(AliEMCALCalibCoefs,1) // Table of Calibration coefficients
+ ClassDef(AliEMCALCalibCoefs,2) // Table of Calibration coefficients
};
#endif // ALIEMCALCalibCoefs_H
* provided "as is" without express or implied warranty. *
**************************************************************************/
-/* $Id:*/
+/*
+$Log$
+*/
//_________________________________________________________________________
-// Top EMCAL folder which will keep all information about EMCAL itself,
-// super Modules (SM), modules, towers, set of hists and so on.
+// Cell folder which will keep all information about cell(tower) itself
+// Initial version was created with TDataSet staf
+// TObjectSet -> TFolder; Sep 6, 2007
//
//*-- Author: Aleksei Pavlinov (WSU, Detroit, USA)
#include "AliEMCALFolder.h"
#include "AliEMCALSuperModule.h"
#include "AliEMCALCalibData.h"
+#include "AliEMCALRecPointsQaESDSelector.h"
#include "AliEMCALCalibCoefs.h"
#include <TH1.h>
#include <TF1.h>
#include <TNtuple.h>
-#include <TObjectSet.h>
typedef AliEMCALHistoUtilities u;
Double_t MPI02 = MPI0*MPI0; // mass**2
AliEMCALCell::AliEMCALCell() :
-TObjectSet(),
+TFolder(),
+ fParent(0),fLh(0),
fAbsId(0),fSupMod(0),fModule(0),fPhi(0),fEta(0),fPhiCell(0),fEtaCell(0),fCcIn(0),fCcOut(0),
fFun(0)
{
}
AliEMCALCell::AliEMCALCell(const Int_t absId, const char* title) :
-TObjectSet(Form("Cell%4.4i",absId)),
+ TFolder(Form("Cell%4.4i",absId),title),
+ fParent(0),fLh(0),
fAbsId(absId),fSupMod(0),fModule(0),fPhi(0),fEta(0),fPhiCell(0),fEtaCell(0),fCcIn(0),fCcOut(0),
fFun(0)
{
- SetTitle(title);
AliEMCALGeometry *g = AliEMCALGeometry::GetInstance();
g->GetCellIndex(fAbsId, fSupMod, fModule, fPhi, fEta);
{
// dtor
}
+//-------------------------------------------------------------------------------------
void AliEMCALCell::SetCCfromDB(AliEMCALCalibData *ccDb)
{
void AliEMCALCell::SetCCfromCCTable(AliEMCALCalibCoefs *t)
{
- if(t == 0) return;
- this->AddObject((TObject*)BookHists(), kTRUE);
+ if(t == 0) {
+ // Dump();
+ return;
+ }
+ if(fLh == 0) {
+ fLh = BookHists();
+ Add(fLh);
+ }
calibCoef *r = t->GetTable(fAbsId);
if(r && r->absId == fAbsId) {
TString optFit(""), OPT(opt);
OPT.ToUpper();
if(h==0) return;
- printf("<I> AliEMCALCell::FitHist : h %p |%s| is started : opt %s\n", h, h->GetName(), opt);
+ printf("<I> AliEMCALCell::FitHist : |%s| is started : opt %s\n", h->GetName(), opt);
TString tit(h->GetTitle());
TF1 *GausPol2 = 0, *g=0, *bg=0;
- if(h->GetListOfFunctions()->GetSize() == 0) {
+ if(h->GetListOfFunctions()->GetSize() == 0 || 1) {
g = u::Gausi(name, 0.0, 0.4, h); // gaus estimation
g->SetParLimits(0, h->Integral()/20., h->Integral());
} else {
GausPol2 = (TF1*)h->GetListOfFunctions()->At(0);
optFit = "IME+";
+ printf("<I> Function is defined alredy : %s optFit %s \n", GausPol2->GetTitle(), optFit.Data());
}
// optFit = "IME+";
h->Fit(GausPol2, optFit.Data(),"", 0.01, 0.28);
gStyle->SetOptFit(111);
u::DrawHist(h,2);
}
- printf("<I> AliEMCALCell::FitHist : h %p |%s| is ended \n\n", h, h->GetName());
+ printf("<I> AliEMCALCell::FitHist : |%s| is ended \n\n", h->GetName());
}
void AliEMCALCell::FitEffMassHist(const char* opt)
{
- AliEMCALFolder* EMCAL = (AliEMCALFolder*)(GetParent()->GetParent()->GetParent());
+ AliEMCALFolder* EMCAL = AliEMCALRecPointsQaESDSelector::GetEmcalFolder();
Int_t it = EMCAL->GetIterationNumber();
TH1* h = (TH1*)GetHists()->At(0);
if(fFun) {
Double_t mpi = fFun->GetParameter(1), mpi2 = mpi*mpi;
Double_t ccTmp = fCcIn * MPI02 / mpi2;
- if(it<=6) {
+ if(it<=1) { // Jul 16, 2007
fCcOut = ccTmp;
} else {
fCcOut = (ccTmp + fCcIn)/2.;
printf(" %s | fCcIn %6.5f -> % 6.5f <- fCcOut \n", GetTitle(), fCcIn , fCcOut);
}
-void AliEMCALCell::Print()
+void AliEMCALCell::PrintInfo()
{
+ printf(" %s %s \n", GetName(), GetTitle());
+ if(fLh == 0 ) return;
TH1* h = (TH1*)GetHists()->At(0);
TF1 *f = (TF1*)h->GetListOfFunctions()->At(0);
- printf(" %s %s \n", GetName(), GetTitle());
if(fFun) printf(" fFun : %s | %s \n", fFun->GetName(), fFun->GetTitle());
else fFun = f;
- if(f) f->Dump();
+ // if(f) f->Dump();
}
TList* AliEMCALCell::BookHists()
gROOT->cd();
TH1::AddDirectory(1);
- AliEMCALFolder* EMCAL = (AliEMCALFolder*)(GetParent()->GetParent()->GetParent());
+ AliEMCALFolder* EMCAL = AliEMCALRecPointsQaESDSelector::GetEmcalFolder();
Int_t it = EMCAL->GetIterationNumber();
new TH1F("01_EffMass", "effective mass of #gamma,#gamma(m_{#pi^{0}}=134.98 MeV) ", 60,0.0,0.3);
//_________________________________________________________________________
// EMCAL cell - keep everyrhing for calibration task
+// Initial version was created with TDataSet staf
+// TObjectSet -> TFolder; Sep 6, 2007
//
//*-- Author: Aleksei Pavlinov (WSU, Detroit, USA)
-#include <TObjectSet.h>
+#include <TFolder.h>
class TList;
class TH1;
class AliEMCALCalibData;
class AliEMCALCalibCoefs;
-class AliEMCALCell : public TObjectSet {
+class AliEMCALCell : public TFolder {
public:
void SetCCfromDB(AliEMCALCalibData *ccDb); // obsolete
void SetCCfromCCTable(AliEMCALCalibCoefs *t);
- TList* GetHists() {return (TList*)fObj;}
+ TList* GetHists() {return fLh;}
+ TObject* GetParent() {return fParent;}
+ void SetParent(TObject *parent) {fParent=parent;}
Int_t GetAbsId() const {return fAbsId;}
Int_t GetSupMod() const {return fSupMod;}
Int_t GetModule() const {return fModule;}
static void FitHist(TH1* h, const char* name="",const char* opt="");
// Menu
void FitEffMassHist(const char* opt=""); //*MENU*
- void Print(); //*MENU*
+ void PrintInfo(); //*MENU*
protected:
+ TObject* fParent; // parent
+ TList* fLh;
+ //
Int_t fAbsId; // abs cell id
Int_t fSupMod; // super module number
Int_t fModule; // module number inside SM
Int_t fEtaCell; // eta number of cell SM
// CC staf
Double_t fCcIn; // input cc in GeV (from Db or table
- Double_t fCcOut; // output cc in GeV (from fir now)
+ Double_t fCcOut; // output cc in GeV (from fit now)
TF1* fFun; //! fitting function - gaus + pol2
//
#include <TH1.h>
#include <TF1.h>
#include <TFile.h>
-#include <TFileIter.h>
#include <TCanvas.h>
#include <TClonesArray.h>
#include <TKey.h>
#include <TNtuple.h>
#include <TLegend.h>
-#include <TArrayS.h>
+#include <TLegendEntry.h>
+#include <TLine.h>
const TString AliEMCALFolder::fgkBaseFolderName("EMCAL");
const TString AliEMCALFolder::fgkCCFirstName("CCFIRST");
const TString AliEMCALFolder::fgkCCinName("CCIN");
const TString AliEMCALFolder::fgkCCoutName("CCOUT");
+ const TString AliEMCALFolder::fgkDirOfRootFiles("$HOME/ALICE/SHISHKEBAB/RF/CALIB/JUL16/");
typedef AliEMCALHistoUtilities u;
TList *lObj = 0;
AliEMCALFolder::AliEMCALFolder() :
- TObjectSet(),
+ TFolder(),
fCounter(0), fGeometry(0), fNumOfCell(0), fLhists(0), fLofCells(0),fPi0SelPar(0),fCalibData(0),
fCellNtuple(0)
{
}
AliEMCALFolder::AliEMCALFolder(const char* name, const char* title, Bool_t putToBrowser) :
- TObjectSet(name),
+ TFolder(name,title),
fCounter(-1), fGeometry(0), fNumOfCell(0), fLhists(0), fLofCells(0),fPi0SelPar(0),fCalibData(0),
fCellNtuple(0)
{
- SetTitle(title);
Init(putToBrowser);
}
AliEMCALFolder::AliEMCALFolder(const Int_t it, const char* title, Bool_t putToBrowser) :
- TObjectSet(Form("%s_%2.2i", AliEMCALFolder::fgkBaseFolderName.Data(),it)),
+ TFolder(Form("%s_%2.2i", AliEMCALFolder::fgkBaseFolderName.Data(),it),title),
fCounter(it), fGeometry(0), fNumOfCell(0), fLhists(0), fLofCells(0),fPi0SelPar(0),fCalibData(0),
fCellNtuple(0)
{
- SetTitle(title);
Init(putToBrowser);
}
{
lObj = new TList;
lObj->SetName("Objects"); // name is good ?
- this->AddObject((TObject*)lObj, kTRUE);
+ this->Add((TObject*)lObj);
+ //this->AddObject((TObject*)lObj, kTRUE);
// Get default geometry - "SHISH_77_TRD1_2X2_FINAL_110DEG"; May 29, 2007
fGeometry = AliEMCALGeometry::GetInstance(); // should be define before
lObj->Add(fGeometry);
// Initial cc with decalibration
- char* dbString = "local:///data/r22b/ALICE/PROD/CALIBRATION_May_2007/PI0/PDSF/10GEV/DECALIB/DeCalibDB"; // hard coded now
- fCalibData = (AliEMCALCalibData*)
- (AliCDBManager::Instance()->GetStorage(dbString)
- ->Get("EMCAL/Calib/Data",gAlice->GetRunNumber())->GetObject());
+ // Jun 13, 2007;
+ // Jul 13 - See ~/macros/ALICE/sim.C for choice of CDB
+ AliEMCALCalibData *calData[1];
+ // Get from defined arrea;
+ // First table is table which used in rec.points finder.
+ Add(AliEMCALCalibCoefs::GetCalibTableFromDb(fgkCCFirstName.Data(),calData));
+ fCalibData = calData[0];
lObj->Add(fCalibData);
- // Jun 13, 2007
- Add(AliEMCALCalibCoefs::GetCalibTableFromDb(fgkCCFirstName.Data()));
- if(GetIterationNumber()==1) {
- Add(AliEMCALCalibCoefs::GetCalibTableFromDb()); // first iteration only
+ if(GetIterationNumber()<=1) {
+ Add(AliEMCALCalibCoefs::GetCalibTableFromDb(fgkCCinName.Data(), calData));
}
+
// Selection Parameter
fPi0SelPar = AliEMCALPi0SelectionParam::Set1();
this->Add(fPi0SelPar);
AliEMCALSuperModule* AliEMCALFolder::GetSuperModule(const Int_t nm)
{
- TDataSet *set = 0;
AliEMCALSuperModule *sm = 0;
- set = FindByName(Form("SM%2.2i",nm));
+ TObject *set = FindObject(Form("SM%2.2i",nm));
if(set) sm = (AliEMCALSuperModule*)set;
return sm;
if(cl1->E() < cl2->E()) cl = cl2; // Get cluster with highest energy
const Int_t nDigits = cl->GetNumberOfDigits();
- // const UShort_t* adc = cl->GetDigitAmplitude(); // for future
- const Short_t* absId = cl->GetDigitIndex()->GetArray();
+ const Short_t* absId = cl->GetDigitIndex()->GetArray();
+
+ AliEMCALCalibCoefs *tFirst = GetCCFirst();
+ calibCoef *rFirst=0;
- int indMax = 0;
- double emax = cl->GetDigitAmplitude()->At(0);
+ int indMax = 0, id=0;
+ id = Int_t(absId[0]);
+ rFirst = tFirst->GetTable(id);
+ double emax = cl->GetTrueDigitEnergy(indMax, rFirst->cc);
if(nDigits > 1) {
for(int i=1; i<nDigits; i++) {
- if(emax < cl->GetDigitAmplitude()->At(i)) {
+ id = Int_t(absId[i]);
+ rFirst = tFirst->GetTable(id);
+ if(emax < cl->GetTrueDigitEnergy(i, rFirst->cc)) {
indMax = i;
- emax = cl->GetDigitAmplitude()->At(i);
+ emax = cl->GetTrueDigitEnergy(i, rFirst->cc);
}
}
}
//
static Int_t nSupModMax, nModuleMax, nIphiMax, nIetaMax, iphiCellMax, ietaCellMax;
static AliEMCALCell* cell;
- static TDataSet *set;
+ static TFolder *set;
static AliEMCALSuperModule *sm;
fGeometry->GetCellIndex(absIdMax, nSupModMax, nModuleMax, nIphiMax, nIetaMax);
cell = new AliEMCALCell(absIdMax, Form("sm%2.2i:phi%2.2i:eta%2.2i(%4.4i)",nSupModMax,iphiCellMax,ietaCellMax,absIdMax));
SetCell(cell, absIdMax);
// For browser
- set = FindByName(Form("SM%2.2i",nSupModMax));
+ set = dynamic_cast<TFolder*>(FindObject(Form("SM%2.2i",nSupModMax)));
if(set==0) {
sm = new AliEMCALSuperModule(nSupModMax);
Add(sm);
+ sm->SetParent(this);
sm->Init();
} else {
- sm = (AliEMCALSuperModule*)set;
+ sm = dynamic_cast<AliEMCALSuperModule*>(set);
+ }
+ if(sm) {
+ sm->AddCellToEtaRow(cell, ietaCellMax);
+ cell->SetParent(sm);
}
- sm->AddCellToEtaRow(cell, ietaCellMax);
//
cell->SetCCfromCCTable(GetCCIn());
} else {
cell = GetCell(absIdMax);
- set = FindByName(Form("SM%2.2i",nm));
+ set = dynamic_cast<TFolder*>(FindObject(Form("SM%2.2i",nm)));
if(set) sm = (AliEMCALSuperModule*)set;
}
if(sm == 0) assert(0);
AliEMCALCalibCoefs* AliEMCALFolder::GetCCTable(const char* name)
{
- TDataSet *ds = FindByName(name);
- if(ds) return (AliEMCALCalibCoefs*)ds;
- else return 0;
+ TObject *obj = FindObject(name);
+ if(obj) return (AliEMCALCalibCoefs*)obj;
+ else return 0;
}
Int_t AliEMCALFolder::GetSMNumber(AliESDCaloCluster* cl)
// Recalibration staf - Jun 18,2007
AliEMCALRecPoint* AliEMCALFolder::GetRecPoint(AliESDCaloCluster *cl, AliEMCALCalibCoefs *tOld,AliEMCALCalibCoefs *tNew,
-TList *l)
+TList *l, Double_t deff, Double_t w0, Double_t phiSlope)
{
//
// Static function;
// Get recalibrated rec.point from ESD cluster
+ // If tNew == 0 -> get ideal calibration with ADCCHANNELEC
//
- // printf(" AliEMCALFolder::GetRecPoint() : RECALIBRATION \n");
static Double_t ADCCHANNELEC = 0.0153; // Update 24 Apr 2007: 250./16/1024 - width of one ADC channel in GeV
- static Float_t ECAW0 = 4.5; // hard coded now - see AliEMCALClusterizerv1::InitParameters()
- Int_t ampDigi=0;
- Double_t eDigiNew=0.0;
+ //static Float_t ECAW0 = 4.5; // hard coded now - see AliEMCALClusterizerv1::InitParameters()
+ static Double_t ECAW0 = 5.5; // Beter case for simulation
+ Int_t ampDigi=0, indMax=-1;
+ Double_t eDigiNew=0.0, eDigiMax=0.0;
+ static TArrayD ed;
+ if(w0 > 0.5) ECAW0 = w0;
AliEMCALRecPoint *rp=0;
calibCoef *rOld=0, *rNew=0;
- if(cl && tOld && tNew){
+ // printf(" AliEMCALFolder::GetRecPoint() : RECALIBRATION : w0 %f ECAW0 %f \n", w0, ECAW0);
+ if(cl && tOld){
// cl->PrintClusterInfo(1);
const Int_t ndg = cl->GetNumberOfDigits();
- //const Int_t nprim = cl->GetNumberOfPrimaries();
+ // const Int_t nprim = cl->GetNumberOfPrimaries();
- const Short_t* prim = cl->GetLabels()->GetArray();
+ const Short_t prim = cl->GetLabel();
const Short_t* dgAbsId = cl->GetDigitIndex()->GetArray();
- // const UShort_t* dgAmp = cl->GetDigitAmplitude()(); // This is energy - bad definition
+ // const UShort_t* dgAmp = cl->GetDigitAmplitude(); // This is energy - bad definition
rp = new AliEMCALRecPoint(""); // opt=""
rp->SetClusterType(AliESDCaloCluster::kClusterv1);
AliEMCALDigit* dg=0;
TClonesArray digits("AliEMCALDigit", ndg);
Int_t absId = 0;
+ ed.Set(ndg); // resize array
for(Int_t i=0; i<ndg; i++){
// Save just abs id and amplitude of digits which will be used for recalculation of
// cluster energy and position
- absId = Int_t(dgAbsId[i]);
- rOld = tOld->GetTable(absId);
- rNew = tNew->GetTable(absId);
- ampDigi = Int_t(cl->GetDigitAmplitude()->At(i) / rOld->cc);
+ absId = Int_t(dgAbsId[i]);
+ rOld = tOld->GetTable(absId);
+ ampDigi = cl->GetTrueDigitAmplitude(i, rOld->cc); // True amplitude
- new(digits[i]) AliEMCALDigit(Int_t(prim[i]),0,absId, ampDigi, 0.0, i, 0.0);
+ new(digits[i]) AliEMCALDigit(Int_t(prim),0, absId, ampDigi, 0.0, i, 0.0);
dg = (AliEMCALDigit*)digits[i];
- eDigiNew = ampDigi * rNew->cc;
- //eDigiNew = ampDigi * ADCCHANNELEC;
- //eDigiNew = Double_t(cl->GetDigitAmplitude()->At(i)); // Copy from ESD for checking
+
+ if(tNew) rNew = tNew->GetTable(absId);
+ if(rNew) {
+ rNew = tNew->GetTable(absId);
+ eDigiNew = Double_t(ampDigi) * rNew->cc; // Recalibrate with new cc
+ } else {
+ eDigiNew = Double_t(ampDigi) * ADCCHANNELEC; // Ideal calibration
+ }
+ //eDigiNew = Double_t(cl->GetTrueDigitEnergy(i)); // Copy from ESD for checking
rp->AddDigit(*dg, Float_t(eDigiNew));
+ ed[i] = eDigiNew;
+ if(eDigiMax<eDigiNew) {
+ eDigiMax = eDigiNew;
+ indMax = i;
+ }
u::FillH1(l, 14, eDigiNew);
u::FillH1(l, 15, Double_t(absId));
//printf("<I> digit %i amp %i rOld->cc %6.5f GeV rNew->cc %6.5f GeV\n", i, ampDigi, rOld->cc, rNew->cc);
}
//printf("<I> recalibration of digits was done ! \n");
// rp->EvalAll(ECAW0, &digits);
- rp->EvalLocalPosition(ECAW0, &digits); // I need just position
+ if(indMax>=0) rp->SetIndMaxDigit(indMax);
+ if(deff > 0.0) { // for fit
+ rp->EvalLocalPositionFit(deff, ECAW0, phiSlope, &digits); // I need just position
+ } else { // get w0 and deff from parametrisation - Sep 4, 2007
+ rp->EvalLocalPosition2(&digits, ed);
+ }
digits.Delete();
}
//rp->Print("print");
{
//
// Jun 5, 2007; See TFileIter and StFMC.cxx
- //
- TString FN(fn);
+ // Jul 16 - added fgkDirOfRootFiles
+ // Sep 7, 2007 - should be changed without TFileIter
+ /*
+ TString FN = fgkDirOfRootFiles;
+ FN += fn;
if(FN.Contains(".root")==0) FN += ".root";
TFileIter f(FN.Data(),opt,"EMCAL object");
UInt_t eventNum = 0; // just one object
- UInt_t runNumber = 0; // 0 now, - may statistics on selectorr
+ UInt_t runNumber = 0; // 0 now, - may statistics on selector
f.NextEventPut(this, eventNum, runNumber);
printf(" Save %s to file %s\n", GetName(), FN.Data());
+ */
}
AliEMCALFolder* AliEMCALFolder::Read(const char *fn, const char *opt)
{
//
// Jun 27, 2007
+ // Jul 16 - added fgkDirOfRootFiles
//
+ printf("<I> AliEMCALFolder::Read(%s,%s)\n",fn,opt);
AliEMCALFolder* EMCAL = 0;
TH1::AddDirectory(0); // this is obligatory
- TFile f(fn,opt);
- TList *l = f.GetListOfKeys();
- printf(" The total number of the objects: %d in file %s\n", l->GetSize(), fn);
+ TString FN = fgkDirOfRootFiles;
+ FN += fn;
+ if(FN.Contains(".root")==0) FN += ".root";
- TKey *key = (TKey*)l->At(0);
- EMCAL = (AliEMCALFolder*)key->ReadObj();
- f.Close();
- if(EMCAL) EMCAL->InitAfterRead();
+ TFile f(FN.Data(),opt);
+ if(f.IsOpen()) {
+ TList *l = f.GetListOfKeys();
+ printf("<I> The total number of the objects: %d \n File %s\n", l->GetSize(), FN.Data());
+ TKey *key = (TKey*)l->At(0);
+ EMCAL = (AliEMCALFolder*)key->ReadObj();
+ f.Close();
+ if(EMCAL) EMCAL->InitAfterRead();
+ }
return EMCAL;
}
void AliEMCALFolder::InitAfterRead()
-{
- lObj = (TList*)fObj;
- fLhists = (TList*)lObj->FindObject("HistsOfEmcal");
+{
+ lObj = (TList*)FindObject("Objects");
+ if(lObj) {
+ fLhists = (TList*)lObj->FindObject("HistsOfEmcal");
+ }
}
void AliEMCALFolder::DrawQA(const int nsm)
AliEMCALSuperModule* sm = GetSuperModule(0);
if(sm) printf(" TNtuple was created ! sm0 %s \n", sm->GetName());
for(int eta=0; eta<48; eta++) { // eta row
- TDataSet *setEta = sm->FindByName(Form("ETA%2.2i",eta));
+ TFolder *setEta = dynamic_cast<TFolder*>(sm->FindObject(Form("ETA%2.2i",eta)));
if(setEta) {
printf(" ***** eta row %s ******\n", setEta->GetName());
- for(int phi=0; phi<setEta->GetListSize(); phi++) { // cycle on cells (phi directions)
- cell = (AliEMCALCell*)setEta->At(phi);
+ TList* l = (TList*)setEta->GetListOfFolders();
+ for(int phi=0; phi<l->GetSize(); phi++) { // cycle on cells (phi directions)
+ cell = (AliEMCALCell*)l->At(phi);
if(cell) {
cell->FillCellNtuple(fCellNtuple);
//printf(" fill cell %s : %s \n", cell->GetName(), cell->GetTitle());
}
fCellNtuple->Print();
lObj->Add(fCellNtuple);
- // --
- CreateAndFillAdditionalHists();
}
void AliEMCALFolder::CreateAndFillAdditionalHists()
{
+ CreateCellNtuple();
TH1::AddDirectory(0);
fLhists->Add(new TH1F("02_CCoutOnEdge", "cc out on edge of calorimeter (in MeV)", 70, 12., 19.));
fLhists->Add(new TH1F("03_CCoutInside", "cc out inside of calorimeter (in MeV)", 70, 12., 19.));
+ fLhists->Add(new TH1F("04_CCoutOnEdge2", "cc out on edge of calorimeter(2) (in MeV)", 70, 12., 19.));
// Fill
Float_t* args;
for(Int_t i=0; i<(Int_t)fCellNtuple->GetEntries(); i++){
Int_t phi = (Int_t)args[5];
Int_t eta = (Int_t)args[6];
Double_t cc = (Double_t)args[8]*1000.;
- if((phi==0||phi==23) || (eta==0||eta==47)) u::FillH1(fLhists, 2, cc);
- else u::FillH1(fLhists, 3, cc);
+ if ((phi==0||phi==23) || (eta==0||eta==47)) u::FillH1(fLhists, 2, cc);
+ else if((phi==1||phi==22) || (eta==1||eta==46)) u::FillH1(fLhists, 4, cc); // next to edge
+ else u::FillH1(fLhists, 3, cc);
}
// Drawing
Int_t wh=530, ww=750;
h1->SetTitle("CC distribution after #pi^{0} calibration");
h1->SetXTitle(" MeV ");
h1->SetYTitle(" N ");
- TLatex *lat1 = u::lat(Form("rel.width = %4.2f%%",
+ TLatex *lat1 = u::Lat(Form("rel.width = %4.2f%%",
100.*h1->GetRMS()/ h1->GetMean()), 16.5, 100., 12, 0.045);
- TLatex *lat2 = u::lat(Form("rel.width = %4.2f%% (from fit)",
+ TLatex *lat2 = u::Lat(Form("rel.width = %4.2f%% (from fit)",
100.*g->GetParameter(2)/ g->GetParameter(1)), 16.5, 70., 12, 0.045);
if(0) {
u::DrawHist(hccFirst,2,1,"same",3);
+ // Ideal calibration - Jul 18, 2007
+ Double_t ADCCHANNELEC = 0.0153, ccIdeal = ADCCHANNELEC*1.e+3;
+ Double_t ym = h1->GetMaximum();
+ TLine *l = new TLine(ccIdeal,-ym*0.05, ccIdeal,ym*1.05);
+ l->SetLineColor(kBlue);
+ l->Draw();
+
TLegend *leg = new TLegend(0.1,0.6, 0.45,0.85);
leg->AddEntry(hccFirst, "Initial cc ", "L");
leg->AddEntry(h1, "Final cc", "L");
+ TLegendEntry *le=0;
+ if(l) {
+ le = leg->AddEntry(l, "Ideal calibration", "L");
+ le->SetTextColor(l->GetLineColor());
+ }
+
+ TH1 *hCcEdge = (TH1*)fLhists->At(2);
+ TH1 *hCcEdge2 = (TH1*)fLhists->At(4);
+ if(1) {
+ u::DrawHist(hCcEdge,2,kGreen,"same",1);
+ le = leg->AddEntry(hCcEdge , "Edge cell", "L");
+ le->SetTextColor(hCcEdge->GetLineColor());
+
+ u::DrawHist(hCcEdge2,2, 28,"same",1); // 28 - like brown
+ le = leg->AddEntry(hCcEdge2 , "Edge cell (2)", "L");
+ le->SetTextColor(hCcEdge2->GetLineColor());
+ u::DrawHist(h1,2,1,"same");
+ }
+
leg->Draw();
c->Update();
for(int m=0; m<12; m++) {
AliEMCALSuperModule *sm = new AliEMCALSuperModule(m);
Add(sm);
+ sm->SetParent(this);
}
}
/* Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
* See cxx source for full Copyright notice */
-/* $log:$ */
+/* $Log$ */
//_________________________________________________________________________
-// Top EMCAL folder - keep everyrhing for calibration task
+// Top EMCAL folder - keep everyrhing for calibration task
+// Initial version was created with TDataSet staf
+// TObjectSet -> TFolder; Sep 5, 2007
//
//*-- Author: Aleksei Pavlinov (WSU, Detroit, USA)
// --- ROOT system ---
-#include <TObjectSet.h>
+#include <TFolder.h>
#include <TString.h>
class AliEMCALGeometry;
class TList;
class TNtuple;
-class AliEMCALFolder : public TObjectSet {
+class AliEMCALFolder : public TFolder {
public:
void FillPi0Candidate(const Double_t mgg, Int_t absIdMax, Int_t nm);
// Define CC
void FitAllSMs(); // SM0 now
- // Service routine
+ // Service routine
AliEMCALCalibCoefs* GetCCTable(const char* name);
AliEMCALCalibCoefs* GetCCFirst() {return GetCCTable(fgkCCFirstName.Data());}
AliEMCALCalibCoefs* GetCCIn() {return GetCCTable(fgkCCinName.Data());}
AliEMCALCalibCoefs* GetCCOut(){return GetCCTable(fgkCCoutName.Data());}
Int_t GetSMNumber(AliESDCaloCluster* cl);
// Recalibration staf - Jun 18,2007
- static AliEMCALRecPoint *GetRecPoint(AliESDCaloCluster *cl,AliEMCALCalibCoefs *tOld,AliEMCALCalibCoefs *tNew, TList *l=0);
+ static AliEMCALRecPoint *GetRecPoint(AliESDCaloCluster *cl,AliEMCALCalibCoefs *tOld,AliEMCALCalibCoefs *tNew,
+ TList *l=0, Double_t deff=-1., Double_t w0=-1., Double_t phiSlope=0.0);
// MENU
void Save(const char *fn = "EMCALFOLDER.root", const char *opt="RECREATE"); // *MENU*
static AliEMCALFolder* Read(const char *fn = "EMCALFOLDER.root", const char *opt="READ");
protected:
TList* BookHists();
- Int_t fCounter; // Coonter of iteration
+ Int_t fCounter; // Counter of iteration
//
AliEMCALGeometry *fGeometry; //
//
static const TString fgkCCFirstName; // name of first calib.table
static const TString fgkCCinName; // name of initial calib.coefs. table
static const TString fgkCCoutName; // name of out calib.coefs. table
+ static const TString fgkDirOfRootFiles; // name of directory for saving EMCAL folder
- void TestSMStruct();
+ void TestSMStruct(); // *MENU*
ClassDef(AliEMCALFolder,2) // EMCAL folder
#include <TH1.h>
#include <TH2.h>
#include <TF1.h>
+#include <TGraph.h>
+#include <TGraphErrors.h>
#include <TLatex.h>
#include <TChain.h>
#include <TList.h>
#include "AliESDCaloCluster.h"
#include "AliEMCALRecPoint.h"
+#include "AliRunLoader.h"
using namespace std;
//fill 1d histogram
static TH1* hid=0;
if(l == 0) return;
- if(ind < l->GetSize()){
+ if(ind>=0 && ind < l->GetSize()){
hid = (TH1*)l->At(ind);
hid->Fill(x,w);
}
//fill 2d histogram
static TH2* hid=0;
if(l == 0) return;
- if(ind < l->GetSize()){
+ if(ind>=0 && ind < l->GetSize()){
hid = (TH2*)l->At(ind);
hid->Fill(x,y,w);
}
}
}
-TLatex *AliEMCALHistoUtilities::lat(const char *text, Float_t x,Float_t y, Int_t align, Float_t tsize, short tcolor)
+void AliEMCALHistoUtilities::ResetListOfHists(TList *l)
+{
+ if(l==0) return;
+
+ for(int i=0; i<l->GetSize(); i++) {
+ TH1F* h = (TH1F*)l->At(i);
+ h->Reset();
+ }
+}
+
+TLatex *AliEMCALHistoUtilities::Lat(const char *text, Float_t x,Float_t y, Int_t align, Float_t tsize, short tcolor)
{
double y1=y;
TLatex *latex = new TLatex;
return latex;
}
+TGraph *AliEMCALHistoUtilities::DrawGraph(Int_t n, Double_t *x, Double_t *y, Int_t markerColor,
+Int_t markerStyle, const char* opt, const char* tit, const char* xTit,const char* yTit, Int_t ifun,
+const char *optFit, const char *fun)
+{
+ /* Drawing options
+ chopt='L' : A simple polyline between every points is drawn
+ chopt='F' : A fill area is drawn ('CF' draw a smooth fill area)
+ chopt='A' : Axis are drawn around the graph
+ chopt='C' : A smooth Curve is drawn
+ chopt='*' : A Star is plotted at each point
+ chopt='P' : Idem with the current marker
+ chopt='B' : A Bar chart is drawn at each point
+ chopt='1' : ylow=rwymin
+ chopt='X+' : The X-axis is drawn on the top side of the plot.
+ chopt='Y+' : The Y-axis is drawn on the right side of the plot.
+
+ Fitting options
+ The list of fit options is given in parameter option.
+ option = "W" Set all errors to 1
+ = "U" Use a User specified fitting algorithm (via SetFCN)
+ = "Q" Quiet mode (minimum printing)
+ = "V" Verbose mode (default is between Q and V)
+ = "B" Use this option when you want to fix one or more parameters
+ and the fitting function is like "gaus","expo","poln","landau".
+ = "R" Use the Range specified in the function range
+ = "N" Do not store the graphics function, do not draw
+ = "0" Do not plot the result of the fit. By default the fitted function
+ is drawn unless the option"N" above is specified.
+ = "+" Add this new fitted function to the list of fitted functions
+ (by default, any previous function is deleted)
+ = "C" In case of linear fitting, not calculate the chisquare
+ (saves time)
+ = "F" If fitting a polN, switch to minuit fitter
+ = "ROB" In case of linear fitting, compute the LTS regression
+ coefficients (robust(resistant) regression), using
+ the default fraction of good points
+ "ROB=0.x" - compute the LTS regression coefficients, using
+ 0.x as a fraction of good points
+
+ */
+ printf("AliEMCALHistoUtilities::drawGraph started \n");
+ printf("Drawing opt |%s| : Fitting opt |%s|\n", opt, optFit);
+
+ TGraph *gr = new TGraph(n, x, y);
+ gr->SetMarkerColor(markerColor);
+ gr->SetLineColor(markerColor);
+ gr->SetMarkerStyle(markerStyle);
+ gr->SetTitle(tit);
+ gr->Draw(opt);
+
+ TString ctmp(opt);
+ if(ctmp.Contains("A")) {
+ gr->GetHistogram()->SetXTitle(xTit);
+ gr->GetHistogram()->SetYTitle(yTit);
+ }
+ ctmp = optFit;
+ if(ifun>=0) {
+ TString sf("pol"); sf += ifun;
+ gr->Fit(sf.Data(),optFit);
+ printf("\n ** Fit by Polynomial of degree %i : %s **\n", ifun, sf.Data());
+ } else if(!ctmp.Contains("no",TString::kIgnoreCase)){
+ printf("\n ** ifun %i : %s **\n", ifun, fun);
+ gr->Fit(fun, optFit);
+ }
+
+ return gr;
+}
+
+TGraphErrors *AliEMCALHistoUtilities::DrawGraphErrors(const Int_t n,Double_t *x,Double_t *y,Double_t *ex,
+Double_t *ey, Int_t markerColor, Int_t markerStyle, const char* opt, const char* tit,
+const char* xTit,char* yTit, Int_t ifun, const char *optFit, const char *fun)
+{
+ printf("AliEMCALHistoUtilities::drawGraphErrors started \n");
+ printf("Drawing opt |%s| : ifun %i: Fitting opt |%s|, fun |%s|\n",
+ opt, ifun, optFit, fun);
+
+ TGraphErrors *gr = new TGraphErrors(n, x,y,ex,ey);
+ gr->SetMarkerColor(markerColor);
+ gr->SetLineColor(markerColor);
+ gr->SetMarkerStyle(markerStyle);
+ if(tit&&strlen(tit)>0) gr->SetTitle(tit);
+
+ TString ctmp(opt);
+ if(ctmp.Contains("A")) {
+ gr->GetHistogram()->SetXTitle(xTit);
+ gr->GetHistogram()->SetYTitle(yTit);
+ }
+ if(ifun>0) {
+ if(ifun != 999) {
+ TString sf("pol"); sf += ifun;
+ gr->Fit(sf.Data(),optFit);
+ printf("\n ** Fit by Polynomial of degree %i : %s **\n", ifun, sf.Data());
+ } else {
+ gr->Fit(fun, optFit);
+ printf("\n ** Fit by %s **\n", fun);
+ }
+ } else {
+ if(strlen(optFit)) {
+ printf("\n ** ifun %i : %s **\n", ifun, fun);
+ gr->Fit(fun, optFit);
+ }
+ }
+
+ gr->Draw(opt);
+
+ return gr;
+}
+
+TF1* AliEMCALHistoUtilities::GetResolutionFunction(const char *opt, TString &latexName)
+{
+ TString OPT(opt);
+ OPT.ToUpper();
+ TF1 *fres=0;
+ if (OPT.Contains("FRES1")) {
+ fres = new TF1("fres","[0]+[1]/sqrt(x)", 0.0, 101.);
+ latexName = "#frac{#sigma_{E}}{E} = A+#frac{B}{#sqrt{E}}";
+ } else if(OPT.Contains("FRES2")) {
+ fres = new TF1("fres","sqrt([0]*[0]+[1]*[1]/x)", 0.0, 101.);
+ latexName = "#sqrt{A^{2}+#frac{B^{2}}{E}}";
+ }
+ if(fres) {
+ fres->SetParName(0,"A");
+ fres->SetParName(1,"B");
+
+ fres->SetParameter(0, 2.0);
+ fres->SetParameter(1, 6.6);
+ fres->SetLineWidth(2);
+ fres->SetLineColor(kRed);
+ }
+ return fres;
+}
+
void AliEMCALHistoUtilities::InitChain(TChain *chain, const char* nameListOfFiles, Int_t nFileMax)
{
// Read name of files from text file with nameListOfFiles and added to the chain
// chain->Lookup();
}
+AliRunLoader* AliEMCALHistoUtilities::InitKinematics(const Int_t nev, const char* galiceName)
+{
+ static AliRunLoader *RL = 0;
+ if(RL == 0 || nev%1000==0) {
+ if(RL) {
+ RL->UnloadgAlice();
+ delete RL;
+ }
+ RL = AliRunLoader::Open(galiceName,AliConfig::GetDefaultEventFolderName(),"read");
+ RL->LoadgAlice(); // obligatory
+ }
+ if(RL) {
+ RL->GetEvent(nev%1000);
+ RL->LoadKinematics();
+ /* Get what you need after that
+ RL->LoadHits();
+ AliStack *stack=RL->Stack();
+ AliESDCaloCluster * clus = esd->GetCaloCluster(n);
+ Int_t label = clus->GetLabel(); // what is this
+ TParticle *primary=stack->Particle(label);
+ */
+ }
+ return RL;
+}
+
+Double_t AliEMCALHistoUtilities::GetMomentum(const char* nameListOfFiles)
+{
+ // Get momentum value from string - like /....100GEV/..
+ Double_t p=-1; // negative if undefined
+ TString st(nameListOfFiles);
+ if(st.Length()>=5) {
+ st.ToUpper();
+ Ssiz_t ind1 = st.Index("GEV"), ind2=0;
+ if(ind1>0) {
+ for(Int_t i=2; i<=10; i++) {
+ ind2 = st.Index("/",ind1-i);
+ if(ind2>0 && ind2<ind1) break;
+ }
+ TString mom = st(ind2+1, ind1-ind2-1);
+ if(mom.IsFloat()) p = mom.Atof();
+ printf(" dir |%s| : mom |%s| : p %f : ind2,1 %i->%i\n", st.Data(), mom.Data(), p, ind2, ind1);
+ }
+ }
+ return p;
+}
+
int AliEMCALHistoUtilities::ParseString(const TString &topt, TObjArray &Opt)
{
// Moved from AliEMCALGeometry
}
return Opt.GetEntries();
}
+
// Analysis utilites
Bool_t AliEMCALHistoUtilities::GetLorentzVectorFromESDCluster(TLorentzVector &v, const AliESDCaloCluster* cl)
{
return f;
}
+// Calibration stuff
+Double_t AliEMCALHistoUtilities::GetCorrectionCoefficientForGamma_1(const Double_t eRec)
+{
+ // Correction to rec.energy - Jul 15, 2007
+ // E(gamma) = corr * E(eRec);
+ /* corr = p0*(eRec-7.5)+p1*(eRec-7.5)*(eRec-7.5)+p2; 0.0<eRec<10.0
+ 1 p0 6.07157e-05 1.15179e-04 -0.00000e+00 1.20997e-03
+ 2 p1 1.50019e-04 3.13566e-05 -0.00000e+00 1.22531e-02
+ 3 p2 9.99019e-01 4.08251e-04 -0.00000e+00 1.40606e-03
+
+ corr = p3 + p4*eRec + p5*eRec*eRec
+ 1 p3 9.97135e-01 5.31970e-04 1.37962e-09 1.68120e-08
+ 2 p4 3.15740e-04 2.53371e-05 1.11475e-11 1.74192e-04
+ 3 p5 -1.35383e-06 2.19495e-07 -5.82864e-13 4.52277e-02
+ */
+ static Double_t p0=6.07157e-05, p1=1.50019e-04, p2=9.99019e-01;
+ static Double_t p3=9.97135e-01, p4=3.15740e-04, p5=-1.35383e-06;
+ static Double_t corr=1.0;
+ if(eRec>=0.0 && eRec <=10.0) {
+ corr = p0*(eRec-7.5) + p1*(eRec-7.5)*(eRec-7.5) + p2;
+ } else if(eRec>10.0 && eRec <=105.0) {
+ corr = p3 + p4*eRec + p5*eRec*eRec;
+ }
+ //printf(" eRec %f | corr %f \n", eRec, corr);
+ return corr;
+}
+Double_t AliEMCALHistoUtilities::GetCorrectedEnergyForGamma_1(const Double_t eRec)
+{
+ return GetCorrectionCoefficientForGamma_1(eRec) * eRec;
+}
-#ifndef AliEMCALHistoUtilities_H
-#define AliEMCALHistoUtilities_H
-/* Copyright(c) 1998-2002, ALICE Experiment at CERN, All rights reserved. *
+#ifndef ALIEMCALHISTOUTILITIES_H
+#define ALIEMCALHISTOUTILITIES_H
+/* Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
* See cxx source for full Copyright notice */
/* $Id$ */
#include <TNamed.h>
class TList;
+class TString;
class TH1;
+class TGraph;
+class TGraphErrors;
class TF1;
class TLatex;
class TChain;
class AliESDCaloCluster;
class AliEMCALRecPoint;
+class AliRunLoader;
class AliEMCALHistoUtilities: public TNamed {
- public:
- AliEMCALHistoUtilities(const char *name="emcalUtilitiesRoutines",
+ public: AliEMCALHistoUtilities(const char *name="emcalUtilitiesRoutines",
const char *tit="EMCAL utility routines");
AliEMCALHistoUtilities(const AliEMCALHistoUtilities &) : TNamed("", ""){
Fatal("cpy ctor", "not implemented") ; }
const char* opt="RECREATE");
static void AddToNameAndTitle(TH1 *h=0, const char *name=0, const char *title=0);
static void AddToNameAndTitleToList(TList *l=0, const char *name=0, const char *title=0);
- static TLatex *lat(const char *text="", Float_t x=0.0,Float_t y=0.0, Int_t align=12, Float_t tsize=0.05, short tcolor = 1);
+ static void ResetListOfHists(TList *l);
+
+ static TLatex *Lat(const char *text="", Float_t x=0.0,Float_t y=0.0, Int_t align=12, Float_t tsize=0.05, short tcolor = 1);
+ static TGraph *DrawGraph(Int_t n=4, Double_t *x=0, Double_t *y=0, Int_t markerColor=4,
+ Int_t markerStyle=4, const char* opt="", const char* tit="", const char* xTit=" jet E_{t} ",
+ const char* yTit="", Int_t ifun=0, const char *optFit="W+", const char *fun="");
+ static TGraphErrors *DrawGraphErrors(const Int_t n=4,Double_t *x=0,Double_t *y=0,Double_t *ex=0,
+ Double_t *ey=0, Int_t markerColor=4,Int_t markerStyle=4, const char* opt="",
+ const char* tit="", const char* xTit=" jet E_{t} ",
+ char* yTit="", Int_t ifun=0, const char *optFit="W+", const char *fun="");
// TChain
static void InitChain(TChain *chain=0, const char* nameListOfFiles=0, Int_t nFileMax=0);
- //
+ static AliRunLoader* InitKinematics(const Int_t nev=0, const char* galiceName="galice.root");
+ //
+ static Double_t GetMomentum(const char* nameListOfFiles);
static int ParseString(const TString &topt, TObjArray &Opt);
// Analysis utilites
static Bool_t GetLorentzVectorFromESDCluster(TLorentzVector &v, const AliESDCaloCluster *cl);
//
static Double_t Gi(Double_t *x, Double_t *par);
static Double_t GiPol2(Double_t *x, Double_t *par);
+ // Calibration stuff
+ static Double_t GetCorrectionCoefficientForGamma_1(const Double_t eRec);
+ static Double_t GetCorrectedEnergyForGamma_1(const Double_t eRec);
+ static TF1* GetResolutionFunction(const char *opt, TString &latexName);
AliEMCALHistoUtilities & operator = (const AliEMCALHistoUtilities &) {
Fatal("operator =", "not implemented") ; return *this ; }
ClassDef(AliEMCALHistoUtilities,1) // EMCAL utility routines
};
-#endif // AliEMCALHistoUtilities_H
+#endif // ALIEMCALHISTOUTILITIES_H
* provided "as is" without express or implied warranty. *
**************************************************************************/
-/* $Log$ */
+/* $Log$co: warning: `/* $Log' is obsolescent; use ` * $Log'.
+
+ * Revision 1.1 2007/08/08 15:58:01 hristov
+ * New calibration classes. They depend on TTable, so libTable.so is added to the list of Root libraries. (Aleksei)
+ * */
//_________________________________________________________________________
// Set of parameters for pi0 selection
#include "AliEMCALPi0SelectionParam.h"
-TableClassImpl(AliEMCALPi0SelectionParam,pi0SelectionParam)
+ClassImp(pi0SelectionParam)
+pi0SelectionParam::pi0SelectionParam() :
+eOfRpMin(0.3), eOfRpMax(30.), massGGMin(0.03), massGGMax(0.28), momPi0Min(1.8), momPi0Max(12.)
+{ }
+//_________________________________________________________________________
+
+ClassImp(AliEMCALPi0SelectionParam)
+AliEMCALPi0SelectionParam::AliEMCALPi0SelectionParam() : TNamed("",""), fTable(0), fCurrentInd(0)
+{
+}
+
+AliEMCALPi0SelectionParam::AliEMCALPi0SelectionParam(const char* name, const Int_t nrow) : TNamed(name,"table of cell information") , fTable(0), fCurrentInd(0)
+{
+ fTable = new TObjArray(nrow);
+}
+
+void AliEMCALPi0SelectionParam::AddAt(pi0SelectionParam* r)
+{
+ (*fTable)[fCurrentInd] = new pi0SelectionParam(*r);
+ fCurrentInd++;
+}
+
+AliEMCALPi0SelectionParam::~AliEMCALPi0SelectionParam()
+{
+ if(fTable) {
+ fTable->Delete();
+ delete fTable;
+ }
+}
+
+pi0SelectionParam* AliEMCALPi0SelectionParam::GetTable(Int_t i) const
+{
+ return (pi0SelectionParam*)fTable->At(i);
+}
void AliEMCALPi0SelectionParam::PrintTable()
{
/* $Id$ */
//_________________________________________________________________________
-// Set of parameters for pi0 selection
+// Set of parameters for pi0 selection
//
//*-- Author: Aleksei Pavlinov (WSU, Detroit, USA)
// --- ROOT system ---
-#include <TTable.h>
+#include <TNamed.h>
+#include <TObjArray.h>
// unit is GeV
-struct pi0SelectionParam {
+class pi0SelectionParam : public TObject{
+ public:
+ pi0SelectionParam();
+ virtual ~pi0SelectionParam() {};
+ virtual const char* GetName() const {return "Pi0Par";}
+
double eOfRpMin; // minimal energy of em.cluster (rec point)
double eOfRpMax; // maximal energy of em.cluster (rec point)
double massGGMin; // minimal mass of gamma,gamma
double massGGMax; // maximal mass of gamma,gamma
double momPi0Min; // minimal pi0 momentum
double momPi0Max; // maximal pi0 momentum
+
+ ClassDef(pi0SelectionParam,1);
};
-class AliEMCALPi0SelectionParam : public TTable {
+class AliEMCALPi0SelectionParam : public TNamed {
public:
- ClassDefTable(AliEMCALPi0SelectionParam , pi0SelectionParam)
+ AliEMCALPi0SelectionParam(); // default constractor
+ AliEMCALPi0SelectionParam(const char* name, const Int_t nrow);
+ virtual ~AliEMCALPi0SelectionParam();
+
+ AliEMCALPi0SelectionParam & operator = (const AliEMCALPi0SelectionParam & /*rvalue*/) {
+ // assignement operator requested by coding convention but not needed
+ Fatal("operator =", "not implemented");
+ return *this;
+ };
+ //
+ void AddAt(pi0SelectionParam* r);
+ pi0SelectionParam* GetTable(Int_t i) const;
+ Int_t GetSize() const {return fTable->GetSize();}
+ Int_t GetNRows() const {return fCurrentInd;}
// Menu
void PrintTable(); // *MENU*
// Set of parameter(s)
static AliEMCALPi0SelectionParam* Set1();
+ //
+ protected:
+ TObjArray *fTable;
+ Int_t fCurrentInd;
- ClassDef(AliEMCALPi0SelectionParam,1) // Set of Parameters For Pi0 Selection
+ ClassDef(AliEMCALPi0SelectionParam, 2) // Set of Parameters For Pi0 Selection
};
#endif // ALIEMCALPI0SELECTIONPARAM_H
/**************************************************************************
- * Copyright(c) 1998-2002, ALICE Experiment at CERN, All rights reserved. *
+ * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
* *
* Author: The ALICE Off-line Project. *
* Contributors are mentioned in the code where appropriate. *
* provided "as is" without express or implied warranty. *
**************************************************************************/
-/* $Log$ */
-
+/* $Log$co: warning: `/* $Log' is obsolescent; use ` * $Log'.
+ * Revision 1.1 2007/08/08 15:58:01 hristov
+ * New calibration classes. (Aleksei)
+ * */
//_________________________________________________________________________
// This is selector for making a set of histogramms from ESD for rec.points
//*-- Authors: Aleksei Pavlinov (WSU)
-//*-- Feb 17, 2007 - May 23, 2007
+//*-- Feb 17, 2007 - Sep 11, 2007
//*-- Pi0 calibration
+//*-- Recalibration, linearity, shower profile
//_________________________________________________________________________
#include "AliEMCALRecPointsQaESDSelector.h"
#include "AliEMCALFolder.h"
#include "AliEMCALSuperModule.h"
#include "AliEMCALCell.h"
+#include "AliEMCALCellInfo.h"
#include "AliEMCALPi0SelectionParam.h"
#include "AliEMCALHistoUtilities.h"
#include "AliEMCALCalibCoefs.h"
#include "AliLog.h"
#include "AliESD.h"
#include "AliEMCALRecPoint.h"
-
+#include "AliRunLoader.h"
+#include "AliStack.h"
+#include "jetfinder/AliEMCALJetMicroDst.h" // Sgpdge
+#include "AliEMCALDigit.h"
+
+#include <iostream>
+#include <iomanip>
+#include <fstream>
#include <assert.h>
#include <map>
#include <string>
#include <TCanvas.h>
#include <TVector3.h>
#include <TLorentzVector.h>
+#include <TParticle.h>
#include <TChain.h>
#include <TFile.h>
#include <TH1F.h>
+#include <TH1D.h>
+#include <TH2F.h>
#include <TF1.h>
#include <TMath.h>
-#include <TArrayS.h>
#include <TBrowser.h>
-#include <TDataSet.h>
+#include <TList.h>
+#include <TCanvas.h>
+#include <TLine.h>
+#include <TObjString.h>
+#include <TArrayI.h>
+#include <TGraphErrors.h>
+#include <TLegend.h>
+#include <TLegendEntry.h>
+#include <TNtuple.h>
using namespace std;
typedef AliEMCALHistoUtilities u;
+Double_t AliEMCALRecPointsQaESDSelector::fDistEff = 3.;
+Double_t AliEMCALRecPointsQaESDSelector::fW0 = 4.5;
+Double_t AliEMCALRecPointsQaESDSelector::fSlopePhiShift = 0.01; // ?? just guess
+
+AliEMCALFolder* AliEMCALRecPointsQaESDSelector::fEMCAL = 0;
+AliEMCALFolder* AliEMCALRecPointsQaESDSelector::fEMCALOld = 0;
+
+AliEMCALGeometry *geo=0;
+
+char *anaOpt[]={
+ "CORR1", // GetCorrectedEnergyForGamma_1(Double_t eRec);
+ "RECALIB",
+ "IDEAL",
+ "PI0",
+ "GAMMA",
+ "KINE", // reading kine file
+ "PROF", // Shower profile: phi direction now
+ "FIT" // define parameters : deff, w0 and phislope
+};
+enum keyOpt{
+ kCORR1,
+ kRECALIB,
+ kIDEAL,
+ kPI0,
+ kGAMMA,
+ kKINE,
+ kPROF,
+ kFIT,
+ //--
+ kEND
+};
+
+int nAnaOpt = sizeof(anaOpt) / sizeof(char*);
+// Enumeration variables
+// enum {kUndefined=-1, kLed=0, kBeam=1};
+// enum {kPatchResponse=0, kClusterResponse=1};
+// enum {kSkip=0, kSaveToMemory=1};
+
ClassImp(AliEMCALRecPointsQaESDSelector)
Int_t nmaxCell = 4*12*24*11;
AliEMCALRecPointsQaESDSelector::AliEMCALRecPointsQaESDSelector() :
AliSelector(),
+ fPmom(0),
fChain(0),
fLofHistsPC(0),
fLofHistsRP(0),
+ fLKineVsRP(0),
+ fLShowerProfile(0),
+ fCellsInfo(0),
fEmcalPool(0),
- fEMCAL(0),
- fEMCALOld(0)
+ fRunOpts(0),
+ fArrOpts(0),
+ fKeyOpts(0)
{
//
// Constructor. Initialization of pointers
AliSelector::Begin(tree);
}
-void AliEMCALRecPointsQaESDSelector::Init(TTree* tree)
-{
- // Init function
-
- AliSelector::Init(tree);
-}
-
void AliEMCALRecPointsQaESDSelector::SlaveBegin(TTree* tree)
{
// The SlaveBegin() function is called after the Begin() function.
AliSelector::SlaveBegin(tree);
}
+void AliEMCALRecPointsQaESDSelector::Init(TTree* tree)
+{
+ // Init function
+
+ AliSelector::Init(tree);
+}
+Bool_t AliEMCALRecPointsQaESDSelector::Notify()
+{
+ // The Notify() function is called when a new file is opened. This
+ // can be either for a new TTree in a TChain or when when a new TTree
+ // is started when using PROOF. Typically here the branch pointers
+ // will be retrieved. It is normaly not necessary to make changes
+ // to the generated code, but the routine can be extended by the
+ // user if needed.
+
+ return AliSelector::Notify();
+}
+
void AliEMCALRecPointsQaESDSelector::InitStructure(Int_t it)
{
- AliEMCALGeometry::GetInstance("SHISH_TRD1_CURRENT_2X2"); // initialize geometry just once
- //AliEMCALGeometry::GetInstance(""); // initialize geometry just once
+ geo = AliEMCALGeometry::GetInstance("SHISH_TRD1_CURRENT_2X2"); // initialize geometry just once
+ fCellsInfo = AliEMCALCellInfo::GetTableForGeometry(geo);
+
+ if(fRunOpts.Length()>0) CheckRunOpts();
- fLofHistsPC = DefineHistsOfRP("PseudoCl");
- fLofHistsRP = DefineHistsOfRP("RP");
+ Int_t key = 0;
+ if(GetKeyOptsValue(kGAMMA)) key = 1;
+ fLofHistsRP = DefineHistsOfRP("RP", fPmom, key);
+ fLofHistsPC = DefineHistsOfRP("PseudoCl", fPmom);
- fEmcalPool = new TObjectSet("PoolOfEMCAL");
- if(it == 1) {
+ fEmcalPool = new TFolder("PoolOfEMCAL","");
+ if(it <= 1) {
fEMCAL = new AliEMCALFolder(it); // folder for first itteration
fEmcalPool->Add(fEMCAL);
}
+ //if(it<=0) SetName("GammaSel"); // For convinience
+ //else SetName("Pi0Sel");
+ if(GetKeyOptsValue(kKINE)) fLKineVsRP = DefineHistsOfKineVsRP("KineVsRP",fPmom, key);
+ if(GetKeyOptsValue(kPROF)) fLShowerProfile = DefineHistsForShowerProfile("ProfY", fPmom);
+}
+
+void AliEMCALRecPointsQaESDSelector::CheckRunOpts()
+{
+ fRunOpts.ToUpper();
+ int nopt = u::ParseString(fRunOpts, fArrOpts);
+ printf("<I> AliEMCALRecPointsQaESDSelector::CheckRunOpts() analyze %i(%i) options : nAnaOpt %i\n",
+ nopt, fArrOpts.GetEntries(), nAnaOpt);
+ if(nopt <=0) return;
+
+ fKeyOpts = new TArrayI(nAnaOpt);
+ for(Int_t i=0; i<nAnaOpt; i++) (*fKeyOpts)[i] = 0;
+
+ for(Int_t i=0; i<fArrOpts.GetEntries(); i++ ) {
+ TObjString *o = (TObjString*)fArrOpts.At(i);
+
+ TString runOpt = o->String();
+ Int_t indj=-1;
+
+ for(Int_t j=0; j<nAnaOpt; j++) {
+ TString opt = anaOpt[j];
+ if(runOpt.Contains(opt,TString::kIgnoreCase)) {
+ indj = j;
+ break;
+ }
+ }
+ if(indj<0) {
+ printf("<E> option |%s| unavailable **\n",
+ runOpt.Data());
+ assert(0);
+ } else {
+ (*fKeyOpts)[indj] = 1;
+ printf("<I> option |%s| is valid : number %i : |%s| \n",
+ runOpt.Data(), indj, anaOpt[indj]);
+ }
+ }
+}
+
+Int_t AliEMCALRecPointsQaESDSelector::GetKeyOptsValue(Int_t key)
+{
+ static Int_t val=0;
+ val = 0;
+ if(fKeyOpts && key>=0 && key<fKeyOpts->GetSize()) {
+ val = fKeyOpts->At(key);
+ }
+ // printf(" key %i : val %i : opt %s \n", key, val, anaOpt[key]);
+ return val;
}
Bool_t AliEMCALRecPointsQaESDSelector::Process(Long64_t entry)
// Assuming that fTree is the pointer to the TChain being processed,
// use fTree->GetTree()->GetEntry(entry).
- if (AliSelector::Process(entry) == kFALSE)
+ if (AliSelector::Process(entry) == kFALSE) {
+ AliDebug(AliLog::kError, Form("Entry %i is not available ",entry));
return kFALSE;
+ }
+
+ // esd->ReadFromTree(tree);
+ // AliESDEvent * esd = new AliESDEvent();
// Check prerequisites
if (!fESD)
AliDebug(AliLog::kError, "ESD branch not available");
return kFALSE;
}
+
static pi0SelectionParam* rPar = GetEmcalFolder()->GetPi0SelectionParRow(0);
static Int_t nEmcalClusters, indOfFirstEmcalRP, nEmcalRP,nEmcalPseudoClusters;
indOfFirstEmcalRP = fESD->GetFirstEMCALCluster();
u::FillH1(fLofHistsRP, 1, double(indOfFirstEmcalRP));
+ static AliRunLoader* RL = 0;
+ static Int_t nev = 0; // Temporary - 0nly for reading one file now !!
+
static AliESDCaloCluster *cl = 0;
nEmcalRP = nEmcalPseudoClusters = 0;
TList *l=0;
static TLorentzVector v, vcl;
int nrp = 0; // # of RP for gg analysis
+ static Double_t erec=0., ecorr=0.0;
for(int i=indOfFirstEmcalRP; i<indOfFirstEmcalRP+nEmcalClusters; i++) {
cl = fESD->GetCaloCluster(i);
if(cl->GetClusterType() == AliESDCaloCluster::kPseudoCluster) {
l = fLofHistsPC;
} else if(cl->GetClusterType() == AliESDCaloCluster::kClusterv1){
nEmcalRP++;
- if(fEMCAL->GetIterationNumber()>1) {
- AliEMCALRecPoint *rp = AliEMCALFolder::GetRecPoint(cl, fEMCAL->GetCCFirst(), fEMCAL->GetCCIn(), fLofHistsRP);
- // if(rp->GetPointEnergy()>=rPar->eOfRpMin && u::GetLorentzVectorFromRecPoint(v, rp)) {
+ if(fEMCAL->GetIterationNumber()>1||GetKeyOptsValue(kIDEAL)||GetKeyOptsValue(kRECALIB)||GetKeyOptsValue(kFIT)) {
+ AliEMCALRecPoint *rp=0;
+ if(GetKeyOptsValue(kFIT) == kFALSE) fDistEff = -1.; // No fitting ; Sep 4, 2007
+ if(GetKeyOptsValue(kIDEAL)) {
+ rp = AliEMCALFolder::GetRecPoint(cl, fEMCAL->GetCCFirst(), 0, fLofHistsRP, fDistEff, fW0, fSlopePhiShift);
+ } else {
+ rp = AliEMCALFolder::GetRecPoint(cl, fEMCAL->GetCCFirst(), fEMCAL->GetCCIn(), fLofHistsRP, fDistEff, fW0, fSlopePhiShift);
+ }
+ if(GetKeyOptsValue(kPROF)) {
+ FillHistsForShowerProfile( GetListShowerProfile(), rp, GetCellsInfo());
+ }
+ //if(rp->GetPointEnergy()>=rPar->eOfRpMin && u::GetLorentzVectorFromRecPoint(v, rp)) {
if(u::GetLorentzVectorFromRecPoint(v, rp)) { // comparing with RP
+ if(GetKeyOptsValue(kCORR1)) {
+ erec = v.Rho();
+ ecorr = u::GetCorrectedEnergyForGamma_1(erec);
+ v.SetRho(ecorr);
+ v.SetE(ecorr); // This is gamma
+ //printf("<1> erec %f | ecorr %f \n", erec, ecorr);
+ }
new(lvM1[nrp]) TLorentzVector(v);
indLv[nrp] = i;
nrp++;
u::FillH1(fLofHistsRP, 11, vcl.P()-v.P());
u::FillH1(fLofHistsRP, 12, TMath::RadToDeg()*(vcl.Theta()-v.Theta()));
u::FillH1(fLofHistsRP, 13, TMath::RadToDeg()*(vcl.Phi()-v.Phi()));
+
+ u::FillH1(fLofHistsRP, 16, v.P());
+ u::FillH2(fLofHistsRP, 17, vcl.P(), vcl.P()-v.P());
l = 0; // no filling
+ if(GetKeyOptsValue(kIDEAL) || GetKeyOptsValue(kRECALIB)) l = fLofHistsRP;
}
- delete rp;
+ if(rp) delete rp;
} else { // first iteration
// if(cl->E()>=rPar->eOfRpMin && u::GetLorentzVectorFromESDCluster(v, cl)) {
if(u::GetLorentzVectorFromESDCluster(v, cl)) { // comparing with RP
// cut 0.4 GeV may be high !
+ if(GetKeyOptsValue(kCORR1)) {
+ erec = v.Rho();
+ ecorr = u::GetCorrectedEnergyForGamma_1(erec);
+ v.SetRho(ecorr);
+ v.SetE(ecorr); // This is gamma now
+ // printf("<2> erec %f | ecorr %f \n", erec, ecorr);
+ // printf(" v.Rho() %f \n", v.Rho());
+ }
new(lvM1[nrp]) TLorentzVector(v);
indLv[nrp] = i;
nrp++;
u::FillH1(fLofHistsPC, 0, double(nEmcalPseudoClusters));
static TLorentzVector *lv1=0, *lv2=0, lgg;
+ for(int i1=0; i1<nrp; i1++){
+ lv1 = (TLorentzVector*)lvM1.At(i1);
+ u::FillH1(fLofHistsRP, 18, lv1->P());
+ }
static Double_t mgg, pgg;
mgg = pgg = 0.;
nrp = lvM1.GetEntriesFast();
if((mgg>=rPar->massGGMin && mgg<=rPar->massGGMax)) {// pi0 candidates
if((pgg>=rPar->momPi0Min && pgg>=rPar->momPi0Min)) {
- fEMCAL->FillPi0Candidate(mgg,fESD->GetCaloCluster(indLv[i1]),fESD->GetCaloCluster(indLv[i2]));
- u::FillH1(fLofHistsRP, 9, pgg);
- u::FillH1(fLofHistsRP,10, lv1->P());
- u::FillH1(fLofHistsRP,10, lv2->P());
+ if(fEMCAL && fEMCAL->GetIterationNumber()>=1) {
+ fEMCAL->FillPi0Candidate(mgg,fESD->GetCaloCluster(indLv[i1]),fESD->GetCaloCluster(indLv[i2]));
+ u::FillH1(fLofHistsRP, 9, pgg);
+ u::FillH1(fLofHistsRP,10, lv1->P());
+ u::FillH1(fLofHistsRP,10, lv2->P());
+ }
}
}
}
}
}
+ // static Int_t fileNumber = 0;
+ static TString curFileName;
+ // if(GetKeyOptsValue(kKINE) && nev<fChain->GetEntries()) {
+ if(GetKeyOptsValue(kKINE)) {
+ // Get galice.root file in current directory
+ if(nev%1000==0) {
+ printf(" current file |%s|\n", fChain->GetCurrentFile()->GetName());
+ curFileName = fChain->GetCurrentFile()->GetName();
+ curFileName.ReplaceAll("AliESDs.","galice.");
+ }
+ RL = u::InitKinematics(nev, curFileName.Data());
+ // Compare kineamtics vs EMCal clusters
+ FillHistsOfKineVsRP(fLKineVsRP, RL, lvM1);
+ }
+
lvM1.Delete();
if(nEmcalClusters != (nEmcalRP+nEmcalPseudoClusters))
Info("Process","nEmcalClusters %i : RP %i + PC %i ",nEmcalClusters, nEmcalRP, nEmcalPseudoClusters);
+ nev++;
+
return kTRUE;
}
// the results graphically or save the results to file.
AliSelector::Terminate();
- /*
- fMultiplicity = dynamic_cast<TH1F*> (fOutput->FindObject("fMultiplicity"));
- if (!fMultiplicity)
- {
- AliDebug(AliLog::kError, Form("ERROR: Histogram not available %p", (void*) fMultiplicity));
- return;
- }
-
- new TCanvas;
- fMultiplicity->DrawCopy();
- */
PrintInfo();
}
+
//
-TList *AliEMCALRecPointsQaESDSelector::DefineHistsOfRP(const char *name)
+TList *AliEMCALRecPointsQaESDSelector::DefineHistsOfRP(const char *name,Double_t p,Int_t keyOpt)
{
- double ADCchannelEC = 0.00305; // ~3mev per adc count
+ printf("<I> DefineHistsOfRP :%s : p %f : keyOpt %i \n", name, p, keyOpt);
+ Double_t ADCchannelEC = 0.0153; // ~15mev per adc count
+ Double_t xma = p*1.4, xmi=0.0, step=0.0, xmic=xmi, xmac = xma;
+ if(xma<0) xma = 20.;
+ Int_t nmax=1000, scale=4, nmaxc = nmax;
gROOT->cd();
TH1::AddDirectory(1);
new TH1F("02_NumberOf", "number of ", 6, -0.5, 5.5);
new TH1F("03_NumberOfDigitsIn", "number of digits(towers) in rec.points ", 101, -0.5, 100.5);
- new TH1F("04_EnergyOf", "energy of ", 1000, 0.0, 100.);
- int nmax=10000;
- double xmi = ADCchannelEC/2., xma = xmi + ADCchannelEC*nmax;
+ if(keyOpt==1 && p>0.1) {
+ if (p>=100.) scale=12;
+ else if(p>=50.) scale=10;
+ else if(p>=30.) scale=8;
+ xma = p + scale*(0.15*TMath::Sqrt(p));
+ xmi = p - scale*(0.15*TMath::Sqrt(p));
+ step = (xma-xmi)/nmax;
+ }
+
+ if(step < 0.0153) {
+ nmax = int((xma-xmi) / ADCchannelEC)+1;
+ xma = xmi + ADCchannelEC*nmax;
+ }
+ new TH1F("04_EnergyOf", "energy of ", nmax, xmi, xma);
+ nmaxc = nmax; xmic=xmi; xmac = xma;
+
+ nmax = 10000;
+ xmi = ADCchannelEC/2.; xma = xmi + ADCchannelEC*nmax;
// All energy(momentum) unit is GeV if don't notice
- new TH1F("05_DigitEnergyIn", "digit energy in ", nmax, xmi, xma);
+ new TH1F("05_DigitEnergyIn", "digit energy in ", nmaxc, xmic, xmac);
new TH1F("06_DigitTimeIn", "digit time in 10ps(0.01ns) ", 1000, 0.0, 3.e+3); // ns/100 = 10 ps
new TH1F("07_DigitAbsIdIn", "digit abs id in ", nmaxCell, -0.5, double(nmaxCell)-0.5);
new TH1F("08_EffMass", "effective mass of #gamma,#gamma(m_{#pi^{0}}=134.9766 MeV)", 100, 0.0, 0.5);
new TH1F("09_MomOfPi0Candidate", "momentum of #pi^{0} candidates (0.085 <mgg<0.185)", 600, 0.0, 30.0);
new TH1F("10_MomOfRpPi0Candidate", "momentum of RP for #pi^{0} candidates (0.085 <mgg<0.185)", 600, 0.0, 30.0);
// Recalibration staf
- new TH1F("11_MomClESD-RpRecalib", "difference of momentum cl(ESD) - rp(Recalib)", 200, -1., +1.);
+ Double_t pmax=11., dpmax=0.20;
+ if(p>0.1) {
+ pmax=p*1.2;
+ dpmax *= p;
+ }
+ new TH1F("11_MomClESD-RpRecalib", "difference of momentum cl(ESD) - rp(Recalib)", 100, -.01, dpmax);
new TH1F("12_ThetaClESD-RpRecalib", "difference of #theta cl(ESD) - rp(Recalib) in degree", 100, -0.05, +0.05);
new TH1F("13_PhiClESD-RpRecalib", "difference of #phi cl(ESD) - rp(Recalib) in degree ", 100, -0.05, +0.05);
// Digi
new TH1F("14_EDigiRecalib", "energy of digits after recalibration", 2000, 0.0, 20.);
- AliEMCALGeometry* g = AliEMCALGeometry::GetInstance();
- new TH1F("15_AbsIdRecalib", "abs Id of digits after recalibration", g->GetNCells(),-0.5,Double_t(g->GetNCells())-0.5);
+ // AliEMCALGeometry* g = AliEMCALGeometry::GetInstance();
+ new TH1F("15_AbsIdRecalib", "abs Id of digits after recalibration", geo->GetNCells(),-0.5,Double_t(geo->GetNCells())-0.5);
+ new TH1F("16_EnergyOfRecalibRp_", "energy of recalibrated rec.points", nmaxc, xmic, xmac); // Jul 12, 2007
+ new TH2F("17_ShiftRecalib_", "E(clESD) - E(recalib)", 110,0.0, pmax, 50,0.0,dpmax); // Jul 13, 2007
- TList *l = u::MoveHistsToList(Form("ListOfHists%s",name), kFALSE);
- u::AddToNameAndTitleToList(l, name, name);
+ // Corrected staff
+ new TH1F("18_EnergyOfCorrectedLV", "energy of corrected LV", nmaxc, xmic, xmac);
+
+ TString st = Form("ListOfHists%s_P=%5.1f",name,p);
+ st.ReplaceAll(" ","");
+ TList *l = u::MoveHistsToList(st.Data(), kFALSE);
+ st = Form("%s_P=%5.1f",name,p);
+ st.ReplaceAll(" ","");
+ u::AddToNameAndTitleToList(l, st.Data(), st.Data());
+
+ return l;
+}
+
+TList* AliEMCALRecPointsQaESDSelector::DefineHistsOfKineVsRP(const char *name, Double_t p, Int_t keyOpt)
+{
+ printf("<I> DefineHistsOfKineVsRP :%s : p %f : keyOpt %i \n", name, p, keyOpt);
+
+ gROOT->cd();
+ TH1::AddDirectory(1);
+ new TH1F("00_hVx",Form("Vx of primary vertex"), 100, -5., +5.); // 00
+ new TH1F("01_hVy",Form("Vy of primary vertex"), 100, -5., +5.); // 01
+ new TH1F("02_hVz",Form("Vz of primary vertex"), 100, -50., +50.); // 02
+
+ // Double_t ADCchannelEC = 0.0153; // ~15mev per adc count
+ Double_t xma = p*1.4, xmi=0.0, sig=0.15*TMath::Sqrt(p);
+ // Double_t step=0.0, xmic=xmi, xmac = xma;
+ if(xma<0) xma = 20.;
+ //Int_t nmax=1000;
+ // scale=4, nmaxc = nmax;
+
+ new TH1F("03_hGidPrimar", "Geant Id of primary particle ", 101, 0.5, 100.5); // 03
+ new TH1F("04_hPmomPrim","momentum of primary particle", 100, xmi, xma); // 04
+ new TH1F("05_hEtaPrim","#eta of primary particle ", 200, -1., 1.); // 05
+ new TH1F("06_hPhiPrim","#phi of primary particle ", 63*2, 0.0,TMath::Pi()*2.); // 06
+ new TH1F("07_hThetaPrim","#theta of primary particle", 63, 0.0,TMath::Pi()); // 07
+ new TH1F("08_hPhiPrimInDegree","#phi of primary particle in degree", 120, 75.0, 195.0); // 08
+ new TH1F("09_hThetaPrimInDegree","#theta of primary particle in degree", 90, 45., 135.); // 09
+
+ // Particle vs cluster
+ new TH1F("10_hE(P)-E(RP)"," E(p) - E(RP) ", 100, -10.*sig, 10.*sig); // 10
+
+ new TH1F("11_hPvsRpAngleInDegree","angle between P and RP (in degree)", 110, 0.0, 1.0); // 11
+ double dphi=0.5; // for fitting
+ new TH1F("12_hPvsRpDPhiInDegree","dif(#phi) between P and RP (in degree)", 100, -dphi, dphi); // 12
+ new TH1F("13_hPvsRpDThetaInDegree","dif(#theta) between P and RP (in degree)", 100, -0.5, 0.5); // 13
+
+ new TH1F("14_hPvsRpDPhiInDegree","dif(#phi) between P and RP (in degree) right part of SM", 100, -dphi, +dphi); // 14
+ new TH1F("15_hPvsRpDPhiInDegree","dif(#phi) between P and RP (in degree) left part of SM", 100, -dphi, dphi); // 15
+
+ new TH1F("16_hPvsRpDThetaInDegree","dif(#theta) between P and RP (even index)", 100, -0.5, 0.5); // 16
+ new TH1F("17_hPvsRpDThetaInDegree","dif(#theta) between P and RP (odd index)", 100, -0.5, 0.5); // 17
+
+ TString st = Form("ListOfHists%s_P=%5.1f",name,p);
+ st.ReplaceAll(" ","");
+ TList *l = u::MoveHistsToList(st.Data(), kFALSE);
+ st = Form("%s_P=%5.1f",name,p);
+ st.ReplaceAll(" ","");
+ u::AddToNameAndTitleToList(l, st.Data(), st.Data());
+ return l;
+}
+
+TList *AliEMCALRecPointsQaESDSelector::DefineHistsForShowerProfile(const char *name, Double_t p)
+{
+ // Aug 1, 2007 - shower profile business as in testbeam analysis
+ // Phi direction here is Y direction in testbeam
+ printf("<I> DefineHistsForShowerProfile: %s : p %f \n",
+ name, p);
+
+ gROOT->cd();
+ TH1::AddDirectory(1);
+
+ // Cell unit
+ new TH1F("00_hPHiMean"," mean in #phi direction ", 24*10, 0.0, 24.);
+ new TH1F("01_hPHiLog"," mean with log in #phi direction ", 24*10, 0.0, 24.);
+ new TH2F("02_XmeanVsXlog", "xmean vs xlog", 240,0.0, 24.,240,0.0, 24.);
+ new TH2F("03_XmeanVsXlog", "xmean vs xlog (system of cell with max energy, first half on #eta)",
+ 50,-0.5,+0.5, 100,-1.,+1.);
+ new TH2F("04_XmeanVsXlog", "xmean vs xlog (system of cell with max energy, second half on #eta)",
+ 50,-0.5,+0.5, 100,-1.,+1.);
+
+ new TH2F("05_PhiShowerProfile", " #phi shower profile - cumulative distribution",
+ 6*25,-3.,3., 201, -0.0025, 1.0025); // 00; x direction in cell unit
+
+ TString st = Form("L%s_P=%5.1f",name,p);
+ st.ReplaceAll(" ","");
+ TList *l = u::MoveHistsToList(st.Data(), kFALSE);
+ st = Form("%s_P=%5.1f",name,p);
+ st.ReplaceAll(" ","");
+ u::AddToNameAndTitleToList(l, st.Data(), st.Data());
return l;
}
+void AliEMCALRecPointsQaESDSelector::FillHistsOfKineVsRP(TList *l, AliRunLoader* RL, TClonesArray &lvM)
+{
+ //
+ // lvM - array of TLorentzVector's which was cretaef from AliESDCaloCluster's
+ //
+
+ if(l==0 || RL==0) return;
+
+ // TNtuple for qucik analysis
+ static TNtuple *nt=0;
+ if(nt==0) {
+ gROOT->cd();
+ Int_t bsize = int(1.e+7);
+ nt = new TNtuple("angle","angle ntuple for quick analysis",
+ "dPhi:dTheta:phiP:thetaP", bsize);
+ gROOT->GetListOfBrowsables()->Add(nt);
+ }
+ static AliStack* st=0;
+ static TParticle *p=0;
+ static Int_t gid=0, ic=0, pdg=0, i=0;
+ gid = ic = pdg = 0;
+
+ st = RL->Stack();
+ if(st == 0) return;
+ // first primary particle
+ p = st->Particle(0);
+ if(p == 0) return;
+
+ u::FillH1(l, ic++, p->Vx());
+ u::FillH1(l, ic++, p->Vy());
+ u::FillH1(l, ic++, p->Vz());
+
+ pdg = p->GetPdgCode();
+ //gid = gMC->IdFromPDG(pdg); // gMc should be defined
+ AliEMCALJetMicroDst::Sgpdge(pdg,gid); // temporary
+
+ u::FillH1(l, ic++, Double_t(gid));
+ u::FillH1(l, ic++, p->P());
+ u::FillH1(l, ic++, p->Eta());
+ u::FillH1(l, ic++, TVector2::Phi_0_2pi(p->Phi()) );
+ u::FillH1(l, ic++, p->Theta());
+ // In degree
+ u::FillH1(l, ic++, TVector2::Phi_0_2pi(p->Phi())*TMath::RadToDeg());
+ u::FillH1(l, ic++, p->Theta()*TMath::RadToDeg()); // 09
+
+ if(lvM.GetEntriesFast() == 0) return;
+ //
+ // Initial kinematics vs Calo clusters
+ //
+ static TLorentzVector *lv = 0, lvp;
+ lvp.SetPxPyPzE(p->Px(),p->Py(),p->Pz(),p->Energy()); // Particle
+
+ Double_t angle = 0.0, angleMin = 180.0, eDiff=0.0, dPhi=0., dTheta=0., phiP=0.0, thetaP=0.0;
+ Int_t indMin = 0;
+ phiP = lvp.Phi() * TMath::RadToDeg();
+ if (phiP<81. || phiP>99.) return; // cut phi boundaries
+ thetaP = lvp.Theta() * TMath::RadToDeg();
+ if (thetaP<56. || thetaP>89.) return; // cut theta boundaries
+
+ for(i=0; i<lvM.GetEntriesFast(); i++) {
+ lv = (TLorentzVector*)lvM.At(i);
+ angle = lv->Angle(lvp.Vect())*TMath::RadToDeg();
+ if(angleMin > angle) {
+ angleMin = angle;
+ indMin = i;
+ }
+ }
+ lv = (TLorentzVector*)lvM.At(indMin);
+ eDiff = lvp.E() - lv->E();
+ u::FillH1(l, ic++, eDiff);
+
+ dPhi = TVector2::Phi_mpi_pi(lvp.Phi()-lv->Phi())*TMath::RadToDeg();
+ dTheta = TVector2::Phi_mpi_pi(lvp.Theta()-lv->Theta())*TMath::RadToDeg();
+ u::FillH1(l, ic++, angleMin);
+ u::FillH1(l, ic++, dPhi);
+ u::FillH1(l, ic++, dTheta);
+
+ if (phiP>=81. && phiP<=90.) {
+ u::FillH1(l, 14, dPhi);
+ } else if(phiP> 90. && phiP<=99.) {
+ u::FillH1(l, 15, dPhi);
+ }
+ // Unfinished - have to get absid of digit with max energy
+ u::FillH1(l, 16, dTheta);
+ u::FillH1(l, 17, dTheta);
+ if(nt) {
+ nt->Fill(dPhi,dTheta,phiP,thetaP);
+ /*
+ tv__tree = (TTree *) gROOT->FindObject("angle");
+ tv__tree->Draw("dPhi:phiP","abs(dPhi)<0.5","", 9182, 0);
+ tv__tree->Draw("dTheta:thetaP","abs(dTheta)<0.5","", 9182, 0);
+ */
+ }
+}
+
+void AliEMCALRecPointsQaESDSelector::FillHistsForShowerProfile
+(TList *l, AliEMCALRecPoint *rp, AliEMCALCellInfo * t)
+{
+ // Aug 1, 2007
+ if(l==0 || rp==0 || t==0) return;
+ // --
+ static Double_t xmean=0., xlog=0.;
+ static cellInfo rMax;
+ static Int_t phiSize;
+
+ if(rp->GetPointEnergy() < 1.0) return;
+ //if(rp->GetPointEnergy() > 8.0) return; // discard merged clusters for pi0 case
+
+ EvalLocalPhiPosition(1., rp, t, xmean, phiSize, rMax);
+ if(phiSize == 0) return; // just one row in cell directions
+
+ EvalLocalPhiPosition(5.5, rp, t, xlog, phiSize, rMax);
+ if(rMax.iPhi>1.5&&rMax.iPhi<21.5 && rMax.iEta>1.5&&rMax.iEta<46.5) {
+ u::FillH1(l, 0, xmean);
+ u::FillH1(l, 1, xlog);
+ u::FillH2(l, 2, xlog, xmean);
+ // Select two central modules
+ if((rMax.iPhim==5 || rMax.iPhim==6)){
+ // Transition to system of cell with max energy
+ xmean -= (double(rMax.iPhi)+0.5);
+ xlog -= (double(rMax.iPhi)+0.5);
+ if(rMax.iEtam>=2 && rMax.iEtam<=12){ // approximatively first half on eta
+ u::FillH2(l, 3, xlog, xmean);
+ } else {// approximatively second half on eta
+ u::FillH2(l, 4, xlog, xmean);
+ }
+ }
+ }
+}
+
+void AliEMCALRecPointsQaESDSelector::EvalLocalPhiPosition(const Double_t wlog, const AliEMCALRecPoint *rp, const AliEMCALCellInfo* t, Double_t &xcog, Int_t &phiSize, cellInfo &rMax)
+{
+ // wlog = 1 - usual center of gravity; >1 - with logarithmic weight.
+ // digits - array of digits
+ // t -
+ //==
+ // xcog - position of center of gravity in cell unit
+ if(wlog<1.0 || rp==0 || t==0) return;
+ xcog = 0;
+
+ static Double_t wtot=0., w=0., edigi=0., e=0., edigiMax=0.;
+ static Int_t absid = 0, phiMin=0, phiMax=0;
+ static cellInfo* r=0;
+
+ e = rp->GetPointEnergy();
+ wtot = 0.0;
+ edigiMax=0.;
+
+ phiMin = 23; phiMax = 0;
+ for(Int_t iDigit=0; iDigit<rp->GetMultiplicity(); iDigit++) {
+ absid = rp->GetAbsId()[iDigit];
+ edigi = rp->GetEnergiesList()[iDigit];
+ if(wlog > 1.0) w = TMath::Max( 0., wlog + TMath::Log(edigi/e));
+ else w = edigi; // just energy
+
+ r = t->GetTable(absid);
+ xcog += w*(Double_t(r->iPhi) + 0.5);
+ wtot += w;
+ if(edigi > edigiMax) {
+ edigiMax = edigi;
+ rMax = (*r);
+ }
+ if(phiMin > r->iPhi) phiMin = r->iPhi;
+ if(phiMax < r->iPhi) phiMax = r->iPhi;
+ }
+ xcog /= wtot;
+ phiSize = phiMax - phiMin;
+ // printf("phiSize %i \n", phiSize);
+}
/* unused now
TList *AliEMCALRecPointsQaESDSelector::DefineHistsOfTowers(const char *name)
{
void AliEMCALRecPointsQaESDSelector::PrintInfo()
{
// Service routine
+ printf("\n %i Entrie(s) | Option(s) |%s| \n", GetOptsArray().GetEntries(), fRunOpts.Data());
+ for(int i=0; i<nAnaOpt; i++) {
+ if(GetKeyOptsValue(i)) printf(" %i |%s| \n", i, anaOpt[i]);
+ }
+
TList *l[2] = {fLofHistsPC, fLofHistsRP};
printf("\n");
for(int i=0; i<2; i++){
TH1F *h = (TH1F*)l[i]->At(2);
printf(" %s \t: %i \n", h->GetTitle(), int(h->GetEntries()));
}
+ printf(" fDistEff %f fW0 %f fSlopePhiShift %f \n", fDistEff, fW0, fSlopePhiShift);
+}
+
+void AliEMCALRecPointsQaESDSelector::SetMomentum(Double_t p)
+{ // Jul 9, 2007
+ fPmom = p;
+ // SetName(Form("%s_p_%f5.1", GetName(), p));
}
+
AliEMCALFolder* AliEMCALRecPointsQaESDSelector::CreateEmcalFolder(const Int_t it)
{
AliEMCALFolder* newFolder = new AliEMCALFolder(it); // folder for iteration #it
AliEMCALFolder* AliEMCALRecPointsQaESDSelector::GetEmcalOldFolder(const Int_t nsm)
{
AliEMCALFolder* folder=0;
- if(fEmcalPool) folder = (AliEMCALFolder*)fEmcalPool->FindByName(Form("EMCAL_%2.2i",nsm));
+ if(fEmcalPool) folder = (AliEMCALFolder*)fEmcalPool->FindObject(Form("EMCAL_%2.2i",nsm));
return folder;
}
void AliEMCALRecPointsQaESDSelector::Browse(TBrowser* b)
{
if(fESD) b->Add(fESD);
- if(fLofHistsPC) b->Add(fLofHistsPC);
- if(fLofHistsRP) b->Add(fLofHistsRP);
if(fChain) b->Add(fChain);
if(fEmcalPool) b->Add(fEmcalPool);
+ if(geo) b->Add(geo);
+ if(fCellsInfo) b->Add(fCellsInfo);
+ //
+ if(fLofHistsPC) b->Add(fLofHistsPC);
+ if(fLofHistsRP) b->Add(fLofHistsRP);
+ if(fLKineVsRP) b->Add(fLKineVsRP);
+ if(fLShowerProfile) b->Add(fLShowerProfile);
// if(u) b->Add(u);
}
return kFALSE;
}
+void AliEMCALRecPointsQaESDSelector::Save(Int_t ver, const char *optIO)
+{ // Aug 3, 2007
+ TString dir("/home/pavlinov/ALICE/SHISHKEBAB/RF/CALIB/"); // Root directory for saving
+ TString nf=dir;
+ if(GetKeyOptsValue(kPROF)) {
+ nf += "PROF/PROFILE_";
+ nf += ver;
+ nf += ".root";
+ TFile f(nf.Data(), optIO);
+ if(f.IsOpen()) {
+ this->Write();
+ f.Close();
+ printf("<I> Save selectort to file |%s| : optIO %s \n",nf.Data(), optIO);
+ } else {
+ printf("<W> File %s exits ! Increase version : ver %i !\n", nf.Data(), ver);
+ }
+ } else {
+ printf("No PROF option -> no saving now !\n");
+ }
+}
+
+AliEMCALRecPointsQaESDSelector* AliEMCALRecPointsQaESDSelector::ReadSelector(const char* nf)
+{
+ AliEMCALRecPointsQaESDSelector* selector=0;
+
+ TH1::AddDirectory(0);
+ TFile f(nf,"READ");
+ if(f.IsOpen()) {
+ TObject *o = f.Get("AliEMCALRecPointsQaESDSelector");
+ if(o) selector = dynamic_cast<AliEMCALRecPointsQaESDSelector *>(o);
+ }
+ printf("<I> read selector %p : file |%s| \n", selector, nf);
+ return selector;
+}
+
void AliEMCALRecPointsQaESDSelector::ReadAllEmcalFolders()
{
if(fEmcalPool==0) {
- fEmcalPool = new TObjectSet("PoolOfEMCAL");
+ fEmcalPool = new TFolder("PoolOfEMCAL","");
for(Int_t it=1; it<=10; it++){
AliEMCALFolder* fold = AliEMCALFolder::Read(Form("EMCALFOLDER_It%i_fit.root",it), "READ");
if(fold) fEmcalPool->Add(fold);
TH1::AddDirectory(0);
Double_t content, error;
- for(Int_t i=0; i<fEmcalPool->GetListSize(); i++) { // cycle on EMCAL folders
- AliEMCALFolder* folder = static_cast<AliEMCALFolder*>(fEmcalPool->At(i));
+ TList* col = (TList*)fEmcalPool->GetListOfFolders();
+ for(Int_t i=0; i<col->GetSize(); i++) { // cycle on EMCAL folders
+ AliEMCALFolder* folder = static_cast<AliEMCALFolder*>(col->At(i));
if(folder==0) continue;
it = folder->GetIterationNumber();
hout->SetMinimum(0.0);
}
+
+TH1F* AliEMCALRecPointsQaESDSelector::FitHistOfRecPointEnergy(const char *opt)
+{
+ TH1::AddDirectory(0);
+
+ TString OPT(opt);
+ OPT.ToUpper();
+
+ Int_t ind = 4, ind2 = 16;
+ if(GetKeyOptsValue(kIDEAL)) {
+ ind = 16; // Jul 12, 2007
+ ind2 = 4;
+ } else if(GetKeyOptsValue(kCORR1)) {
+ ind = 18;
+ ind2 = 4;
+ }
+ TH1F *hold = (TH1F*)fLofHistsRP->At(ind), *h=0;
+ if(hold == 0) return 0;
+ if(hold->GetEntries() <10.) return hold;
+
+ if(OPT.Contains("CLONE")) {
+ TString newName(Form("C_%s",hold->GetName()));
+ h = (TH1F*)hold->Clone(newName.Data());
+ printf(" Clone hist %s -> |%s|%s| \n",hold->GetName(),h->GetName(),h->GetTitle());
+ } else {
+ h = hold;
+ }
+ TH1F* h2 = (TH1F*)fLofHistsRP->At(ind2);
+
+ Double_t xmax = h->GetXaxis()->GetXmax(), xmin = 0.4*xmax;;
+ TF1 *g = u::Gausi("fRecPointE", xmin, xmax, h);
+ g->SetLineColor(kRed);
+ gStyle->SetOptFit(111);
+
+ h->Fit(g,"N","", xmin, xmax);
+ printf(" (1) xmin %f : xmax %f \n", xmin, xmax);
+
+ xmin = g->GetParameter(1) - 4.*g->GetParameter(2);
+ xmax = g->GetParameter(1) + 4.*g->GetParameter(2);
+ h->Fit(g,"Q+","", xmin, xmax);
+ u::DrawHist(h2,1, 1, "same",2);
+ printf(" (2) xmin %f : xmax %f \n", xmin, xmax);
+
+ return h;
+}
+
+TCanvas *AliEMCALRecPointsQaESDSelector::Linearity(TList *l, int ifun)
+{ //Jul 10, 2007
+ if(l==0) {
+ printf("<E> AliEMCALRecPointsQaESDSelector::Linearity :TList is zero ! Bye ! \n");
+ return 0;
+ }
+ Double_t p[9]={0.5, 1., 3., 5., 10., 20., 30., 50., 100.}, ep[9];
+ // see macro rHits.C ->defineSampleFraction
+ TH1F* hErecOverEin = new TH1F("hErecOverEin","Ratio E_{rec}/E_{#gamma} vs E_{#gamma}", 101,0.5,101.5);
+ // Fill hist
+ TArrayD invRat(9), eInvRat(9), erec(9),derec(9), residual(9);
+ TArrayD res(9), eres(9); // resolution
+ Int_t np=9;
+ for(Int_t var=0; var<9; var++){
+ TH1F *h = (TH1F*)l->At(var);
+ TF1 *f = (TF1*)h->GetListOfFunctions()->At(0);
+ Double_t mean = f->GetParameter(1);
+ Double_t emean = f->GetParError(1);
+
+ Int_t bin = hErecOverEin->FindBin(p[var]);
+ hErecOverEin->SetBinContent(bin, mean/p[var]);
+ hErecOverEin->SetBinError(bin, emean/p[var]);
+ //
+ invRat[var] = p[var]/mean;
+ eInvRat[var] = emean/p[var]*invRat[var];
+ erec[var] = mean;
+ derec[var] = emean;
+ // Resolution in %
+ res[var] = 100.*f->GetParameter(2)/p[var];
+ eres[var] = 100.*f->GetParError(2)/p[var];
+ ep[var] = 0.0;
+ }
+
+ TCanvas *c = new TCanvas("Linearity","Linearity", 20,20, 700, 500);
+ gStyle->SetOptStat(0);
+ if(0) { // E(rec) / E(gamma)
+ u::DrawHist(hErecOverEin,2);
+ hErecOverEin->SetMinimum(0.984);
+ hErecOverEin->SetMaximum(1.001);
+ }
+ Int_t markerColor=1;
+ char *fun="", *optFit="";
+ TF1 *f = 0;
+ if(0) {
+ if(ifun==-5) {
+ fun = "fun5"; optFit="R+";
+ f = new TF1(fun,"([0]*(x-7.5)*(x-7.5))*exp([1]+[2]*x+[3]*x*x)+1.", 0.0, 100.);
+ f->FixParameter(0, 1.95380e-05);
+ // f->SetParameter(0, 1.95380e-05);
+ // f->SetParameter(1, 1.0);
+ } else if(ifun==-4) {
+ fun = "fun4"; optFit="R+";
+ f = new TF1(fun,"([0]*(x-7.5)*(x-7.5))*([1]+[2]*abs(x)+[3]*(x-[4])*(x-[4]))+1.", 0.0, 100.);
+ f->FixParameter(0, 1.95380e-05);
+ f->SetParameter(4, 50.);
+ // f = new TF1(fun,"exp(([0]+[1]*x)*(x-[2])*(x-[3]))", 0.0, 100.);
+ //f = new TF1(fun,"([0]+[1]*x)*(x-[2])*(x-[3])", 0.0, 100.);
+ //f->SetParameter(0, 1.);
+ //f->SetParameter(1, 1.);
+
+ // f->SetParameter(2, 5.);
+ //f->SetParLimits(2, 3.,8.);
+ //f->SetParameter(3, 10.);
+ //f->SetParLimits(3, 9.,15.);
+ //f->SetParameter(2, 2.e-4);
+ } else if(ifun==-3) {
+ fun = "fun3"; optFit="R+";
+ f = new TF1(fun,"[0]+[1]*x+[2]*x*x+[3]*x*x*x", 0.0, 11.);
+ } else if(ifun==-2) {
+ fun = "fun2"; optFit="R+";
+
+ /*
+ f = new TF1(fun,"[0]*(x-7.5)+[1]*(x-7.5)*(x-7.5)+[2]", 0.0, 10.1);
+ f->SetParameter(0, 5.98727e-04);
+ f->SetParameter(1, 2.12045e-04);
+ f->SetParameter(2, 1.);
+ */
+ f = new TF1(fun,"[0]+[1]*x+[2]*x*x", 9.0, 100.1);
+ } else if(ifun==-1) {
+ fun = "fun0"; optFit="R+";
+ f = new TF1(fun,"[0]", 0.0, 100.1);
+ } else if(ifun>=2) {
+ fun = Form("pol%i",ifun);
+ optFit="+";
+ }
+ TGraphErrors *gr = u::DrawGraphErrors(np, erec.GetArray(),invRat.GetArray(),derec.GetArray(),eInvRat.GetArray(),
+ // markerColor,21+markerColor,"AP", " Ratio E_{#gamma}/E_{Rec} ", "E_{Rec} "," E_{#gamma}/E_{Rec}",
+ markerColor,21+markerColor,"AP", " Ratio E_{#gamma}/E_{Rec}^{corr} ", "E_{Rec}^{corr} "," E_{#gamma}/E_{Rec}^{corr}",
+ ifun, optFit, fun);
+ gr->GetHistogram()->SetAxisRange(0.0,100.);
+ // double xmi=0.999, xma=1.017;
+ double xmi=0.995, xma=1.005;
+ gr->GetHistogram()->SetMaximum(xma);
+ gr->GetHistogram()->SetMinimum(xmi);
+ gr->GetHistogram()->SetTitleOffset(1.4,"y");
+ if(ifun==0) {
+ f = new TF1("fres", "AliEMCALHistoUtilities::EnergyCorrectionForGamma_1(x)", 0., 101.);
+ f->Draw("same");
+ }
+ }
+ TLine *line = new TLine(0.0,1.0, 100.,1.0);
+ line->Draw();
+
+ if(0) {
+ c->Clear();
+ for(int i=0; i<9; i++) {
+ residual[i] = 100.*(invRat[i] - u::GetCorrectionCoefficientForGamma_1(erec[i])); // in percent
+ printf(" erec %f : residual %5.3f \n", erec[i], residual[i]);
+ }
+ markerColor = 2;
+ TGraphErrors *gr2=u::DrawGraphErrors(np, erec.GetArray(),residual.GetArray(),derec.GetArray(),eInvRat.GetArray(),
+ markerColor,21+markerColor,"AP"," residual in %, rec.point level", "E_{Rec} "," residual in %",
+ -1, "", 0);
+ gr2->GetHistogram()->SetAxisRange(0.0,100.);
+ gr2->GetHistogram()->SetMaximum(0.2);
+ gr2->GetHistogram()->SetMinimum(-0.1);
+ line = new TLine(0.0,0.0, 101.,0.0);
+ line->Draw();
+ TLatex *lat = u::Lat("linearity better 0.2% after correction",20., 0.15, 12, 0.06, 1);
+ if(lat);
+ }
+ if(1) {
+ TString latexName;
+ c->Clear();
+ gStyle->SetOptFit(111);
+ markerColor = 1;
+ ifun = -11;
+ f = u::GetResolutionFunction("FRES2", latexName);
+ f->SetNpx(10000);
+ f->SetLineColor(kBlack);
+
+ TGraphErrors *gr3=u::DrawGraphErrors(np, p,res.GetArray(), ep, eres.GetArray(),
+ markerColor,21+markerColor,"AP"," resolution in %, rec.point level","E_{#gamma} "," resolution in %",
+ ifun, "+", f->GetName());
+ gr3->GetHistogram()->SetAxisRange(0.0,101.);
+ gr3->GetHistogram()->SetMaximum(14.);
+ gr3->GetHistogram()->SetMinimum(0.0);
+ gr3->SetMarkerSize(1.5);
+
+ TLatex *lat = u::Lat(latexName.Data(),82., 11., 12, 0.06, 1);
+ if(lat);
+ // Exp. data
+ TF1 *fexp = new TF1(*f);
+ fexp->SetName("fexp");
+ fexp->SetParameter(0, 2.168);
+ fexp->SetParameter(1, 8.818);
+ fexp->SetLineWidth(1);
+ fexp->SetLineColor(kRed);
+ fexp->SetLineStyle(1);
+ fexp->Draw("same");
+
+ TLegend *leg = new TLegend(0.21,0.36, 0.68,0.82);
+ leg->AddEntry(gr3, "MC data", "P");
+ leg->AddEntry(f, "fit of MC data", "L");
+ TLegendEntry *ent3 = leg->AddEntry(fexp, "#splitline{fit of exp.data}{FNAL, Nov 2005}", "L");
+ ent3->SetTextColor(fexp->GetLineColor());
+ leg->Draw();
+ }
+ c->Update();
+
+ return c;
+}
+
+TCanvas *AliEMCALRecPointsQaESDSelector::DrawKineVsRP(TList *l)
+{ //Jul 25, 2007
+ if(l==0) {
+ printf("<W> AliEMCALRecPointsQaESDSelector::DrawKineVsRP : TList is zero ! \n");
+ return 0;
+ }
+ TCanvas *c = new TCanvas("KineVsRP","KineVsRP", 20,20, 700, 500);
+ gStyle->SetOptStat(1110);
+ c->Divide(2,2);
+
+ c->cd(1);
+ TH1F* h1 = (TH1F*)l->At(10);
+ u::DrawHist(h1,2);
+
+ c->cd(2);
+ TH1F* h2 = (TH1F*)l->At(11);
+ u::DrawHist(h2,2);
+
+ c->cd(3);
+ TH1F* h3 = (TH1F*)l->At(12);
+ u::DrawHist(h3,2);
+ u::DrawHist((TH1F*)l->At(14), 1,kRed,"same");
+ u::DrawHist((TH1F*)l->At(15), 1,kGreen,"same");
+
+ c->cd(4);
+ TH1F* h4 = (TH1F*)l->At(13);
+ u::DrawHist(h4, 2);
+ u::DrawHist((TH1F*)l->At(16), 1,kRed,"same");
+ u::DrawHist((TH1F*)l->At(17), 1,kGreen,"same");
+
+ /*
+ TH1F* h2 = (TH1F*)l->At(11);
+ u::DrawHist(h2,2);
+ */
+
+ c->Update();
+ return c;
+}
+
+TCanvas* AliEMCALRecPointsQaESDSelector::DrawMeanVsLog(TH2F *h2)
+{
+ // h - input histogramm : mean vds log coordinates
+ if(h2==0) return 0;
+
+ TCanvas *c = new TCanvas("KineVsRP","KineVsRP", 20,20, 700, 500);
+
+ TH1F *hid1 = new TH1F("h1","h1 title ", h2->GetNbinsX(),
+ h2->GetXaxis()->GetXmin(), h2->GetXaxis()->GetXmax());
+
+ gROOT->cd();
+ TH1::AddDirectory(1);
+ TString newName;
+ for(int ix=1; ix<=h2->GetNbinsX();ix++) {
+ newName = "hd_"; newName += ix;
+ TH1D *hd = h2->ProjectionY(newName.Data(),ix,ix);
+ if(hd->Integral()>=4.) {
+ // lhd->Add(hd);
+ hid1->SetBinContent(ix, hd->GetMean());
+ hid1->SetBinError(ix, hd->GetRMS());
+ }
+ }
+ TF1 *f1 = new TF1("fcg", "0.5*TMath::SinH(x/[0])/TMath::SinH(0.5/[0])", -0.5, 0.5);
+ f1->SetParameter(0,0.13);
+ f1->SetLineColor(kRed);
+ f1->SetLineWidth(1);
+
+ gStyle->SetOptStat(0);
+ gStyle->SetOptFit(111);
+ hid1->Fit(f1,"R+");
+
+ u::DrawHist(hid1,2);
+
+//u::DrawHist(h2,2);
+
+ c->Update();
+ return c;
+}
+
+TCanvas* AliEMCALRecPointsQaESDSelector::DrawPhiEtaAnglesDistribution(const char* gn)
+{
+ // Aug 6, 2007
+ // Proper geometry should be defined already
+ // dTheta distibution has two peaks which coresponding different cells inside module !
+
+ TCanvas *c = new TCanvas("Geometry","Geometry", 20,20, 700, 500);
+ c->Divide(2,2);
+
+ if(geo==0) geo = AliEMCALGeometry::GetInstance(gn);
+
+ gROOT->cd();
+ TH1::AddDirectory(1);
+ TH1F *hDtheta = new TH1F("hDtheta","#Delta#theta in one SM", 60, -2.0, +1.0); // in degree
+ TH2F *hDtheta2 = new TH2F("hDtheta2","#Delta#theta vs iEta of cell", 48, -0.5, 47.5, 60,-2.0,+1.0);
+
+ TH1F *hDphi = new TH1F("hDphi","#Delta#ph in one SM", 2000, -10.0, +10.0); // in degree
+
+ TVector3 vg3;
+ AliEMCALCellInfo* t = GetCellsInfo();
+ Double_t thetaCell=0., thetaModule=0.;
+ Double_t phiCell=0., dphi=0.;
+
+ for(int absid=0; absid<12*24*4; absid++){
+ cellInfo *r = t->GetTable(absid);
+ geo->GetGlobal(absid, vg3);
+
+ thetaCell = vg3.Theta()*TMath::RadToDeg();
+ thetaModule = 90. - 1.5*r->iEtam;
+ hDtheta->Fill(thetaCell - thetaModule);
+ hDtheta2->Fill(double(r->iEta), thetaCell - thetaModule);
+
+ phiCell = vg3.Phi()*TMath::RadToDeg();
+ dphi = phiCell - 90.;
+ hDphi->Fill(dphi);
+ }
+
+ c->cd(1);
+ u::DrawHist(hDtheta,2);
+ c->cd(2);
+ u::DrawHist(hDtheta2,2);
+ c->cd(3);
+ u::DrawHist(hDphi,2);
+
+ c->Update();
+ return c;
+}
+
+TCanvas* AliEMCALRecPointsQaESDSelector::DrawDeffVsEnergy()
+{
+ // Aug 31, 2007 - obsolete
+ Double_t p[]={0.5, 1.0, 3., 5., 10., 20.,30.,50.,100.};
+ Double_t ep[]={0., 0., 0., 0., 0., 0., 0., 0., 0.};
+ // 2 pars
+ Double_t deff[]={9.07515, 10.0835, 11.2032, 11.4229, 12.3578, 13.0332, 13.3281, 13.7910, 14.3319};
+ Double_t edeff[]={2.69645e-03, 3.52180e-03, 4.19236e-02, 1.21201e-01, 1.58886e-01, 3.96680e-01, 3.29985e-02, 1.17113e-01, 5.22763e-01};
+ //
+ Double_t w0[]={3.44679e+00, 3.82714e+00, 4.15035e+0, 4.36650e+00, 4.51511e+00, 4.65590, 4.63289e+00, 4.66568, 4.68125};
+ Double_t ew0[]={7.58982e-01, 1.26420e-02, 2.36129e-10, 1.21201e-01, 2.12999e-01, 7.95650e-02, 6.15307e-03, 1.88803e-01, 5.18022e-05};
+ int np=sizeof(p)/sizeof(Double_t);
+ printf("<I> AliEMCALRecPointsQaESDSelector::DrawDeffVsEnergy() | np %i \n", np);
+
+ TCanvas *c = new TCanvas("Deff","Deff", 20,20, 700, 500);
+ c->Divide(2,1);
+
+ Int_t markerColor=1;
+ c->cd(1);
+ TF1 *fdeff= new TF1("fdeff","[0]+[1]*log(x)",0.4, 100.4);
+ if(fdeff);
+ TGraphErrors *gr = u::DrawGraphErrors(np, p,deff, ep, edeff,
+ markerColor,21+markerColor,"AP"," D_{eff} vs E_{#gamma} ","E_{#gamma} "," D_{eff} in cm ",
+ -1, "", 0);
+ // -1, "+", fdeff->GetName());
+ gr->GetHistogram()->SetMaximum(15.);
+ gPad->SetLogx(1);
+
+ c->cd(2);
+ TGraphErrors *gr2 = u::DrawGraphErrors(np, p,w0, ep, ew0,
+ markerColor,21+markerColor,"AP"," w_{0} vs E_{#gamma} ","E_{#gamma} "," w_{0} ",
+ -1, "", 0);
+
+ gr2->GetHistogram()->SetMaximum(5.);
+ gPad->SetLogx(1);
+
+ c->Update();
+ return c;
+}
+
+TCanvas* AliEMCALRecPointsQaESDSelector::DrawDeffVsEnergy2(const char *opt)
+{
+ // Aug 28, 2008 - read pars and pars errors from files
+ Double_t p[]={0.5, 1.0, 3., 5., 10., 20.,30.,50.,100.};
+ Double_t ep[]={0., 0., 0., 0., 0., 0., 0., 0., 0.};
+ // 2 pars
+ Double_t deff[9], edeff[9], w0[9], ew0[9]; // max size now
+ TString OPT(opt);
+
+ int np = sizeof(p)/sizeof(Double_t);
+ printf("<I> AliEMCALRecPointsQaESDSelector::DrawDeffVsEnergy2() | np %i \n", np);
+ ReadParsDeffAndW0("/data/r22b/ALICE/CALIB/FIT/", deff, edeff, w0, ew0, 1);
+
+ TCanvas *c = new TCanvas("Deff","Deff", 20,20, 700, 500);
+ c->Divide(2,1);
+
+ TF1 *fdeff = 0, *fw0 = 0;
+ TString optFit(""), funName("");
+ if(OPT.Contains("fit1")) {
+ fdeff= new TF1("fdeff","[0]+[1]*log(x)",0.1, 101.);
+ fdeff->SetLineColor(kRed);
+ fdeff->SetLineWidth(1);
+
+ /* good description - "[0]/(1.+exp([1]*x))"
+ 1 p0 4.82208e+00 9.93617e-04 -0.00000e+00 7.16648e-06
+ 2 p1 -7.79655e-01 4.07500e-03 -0.00000e+00 2.01009e-03
+ better description - "[0]/(1.+exp([1]*(x-[2])))" (like the Woods-Saxon potential)
+ 1 p0 4.83713e+00 1.00437e-03 -6.98984e-10 4.90371e-07
+ 2 p1 -2.77970e-01 3.35587e-03 -1.79239e-09 2.67312e-07
+ 3 p2 4.41116e+00 8.72191e-02 1.82791e-04 1.55643e-05
+ */
+ fw0= new TF1("fw0","[0]/(1.+exp([1]*(x+[2])))",0.1, 101.);
+ fw0->SetLineColor(kRed);
+ fw0->SetLineWidth(1);
+ fw0->SetParameter(0, 4.8);
+ fw0->SetParameter(1, -2.77970e-01);
+ fw0->SetParameter(2, 4.41116);
+ optFit = "+";
+ }
+ if(fdeff) funName = fdeff->GetName();
+
+ Int_t markerColor=1;
+ c->cd(1);
+ gStyle->SetOptFit(111);
+ TGraphErrors *gr = u::DrawGraphErrors(np, p,deff, ep, edeff,
+ markerColor,21+markerColor,"AP"," D_{eff} vs E_{#gamma} ","E_{#gamma} "," D_{eff} in cm ",
+ -1, optFit.Data(), funName.Data());
+ gr->GetHistogram()->SetMaximum(15.);
+ gPad->SetLogx(1);
+ TLegend *leg1 = new TLegend(0.12,0.76, 0.70,0.90);
+ TLegendEntry *le1 = leg1->AddEntry(fdeff, Form("%s",fdeff->GetTitle()), "lp");
+ //TLegendEntry *le1 = leg1->AddEntry(fdeff, Form("%4.2f+%4.2f*log(E_{#gamma})",
+ //fdeff->GetParameter(0),fdeff->GetParameter(1)), "lp");
+ le1->SetTextColor(fdeff->GetLineColor());
+ leg1->Draw();
+
+ c->cd(2);
+ gStyle->SetOptFit(111);
+ if(fw0) funName = fw0->GetName();
+ TGraphErrors *gr2 = u::DrawGraphErrors(np, p,w0, ep, ew0,
+ markerColor,21+markerColor,"AP"," w_{0} vs E_{#gamma} ","E_{#gamma} "," w_{0} ",
+ -1, optFit.Data(), funName.Data());
+
+ gr2->GetHistogram()->SetMaximum(5.);
+
+ TLegend *leg2 = new TLegend(0.17,0.6, 0.99,0.72);
+ TLegendEntry *le2 = leg2->AddEntry(fw0, Form("%s",fw0->GetTitle()), "lp");
+ //TLegendEntry *le2 = leg2->AddEntry(fw0, Form("#frac{%4.2f}{1.+exp(%4.2f*(x+%4.2f)}",
+ //fw0->GetParameter(0),fw0->GetParameter(1),fw0->GetParameter(2)), "lp");
+ le2->SetTextColor(fw0->GetLineColor());
+ leg2->Draw();
+ //gPad->SetLogx(1);
+
+ c->Update();
+ return c;
+}
+
+void AliEMCALRecPointsQaESDSelector::ReadParsDeffAndW0
+(const char *dirName, double *deff, double *edeff, double *w0, double *ew0, const Int_t pri)
+{
+ int strategy = 0, itmp=0;
+ char line[100];
+ for(int var=11; var<=19; var++){
+ int ind = var -11;
+ ifstream fin;
+ TString fname = Form("%s/fitVar%iStrategy%i.txt", dirName, var, strategy);
+ //printf(" open file %s \n", fname.Data());
+ fin.open(fname.Data());
+
+ for(int i=1; i<=2; i++) { // skip to lines
+ fin.getline(line,100).eof();
+ if(pri>=2) printf("%s \n", line);
+ }
+
+ fin >> itmp >> deff[ind] >> edeff[ind];
+ fin >> itmp >> w0[ind] >> ew0[ind];
+ if(pri>=1)
+ printf(" %i def %f+/-%f : def %f+/-%f\n", ind, deff[ind],edeff[ind],w0[ind],ew0[ind]);
+ fin.getline(line,100).eof(); // skip last line
+ fin.close();
+ }
+}
+
+TCanvas* AliEMCALRecPointsQaESDSelector::DrawSpaceResolution()
+{
+ // Sep 4, 2007;
+ // Space resolution vs optimised space resolution
+ Double_t p[]={0.5, 1.0, 3., 5., 10., 20.,30.,50.,100.};
+ Double_t ep[]={0., 0., 0., 0., 0., 0., 0., 0., 0.};
+ Double_t spAng[]={0.1877, 0.1351, 0.08144, 0.06331, 0.05384, 0.04876, 0.04706, 0.04656, 0.04726};
+ Double_t rmsSpAng[]={0.1033, 0.07685, 0.05235, 0.03992, 0.04012, 0.04257, 0.03472, 0.02814, 0.02784};
+ Double_t spAngOpt[]={0.1733, 0.1311, 0.07961, 0.06401, 0.05347, 0.04618, 0.04288, 0.04, 0.03802};
+ Double_t rmsSpAngOpt[]={0.09476, 0.07472, 0.05011, 0.04242, 0.04075, 0.04304, 0.03545, 0.02744, 0.02593};
+
+ int np=sizeof(p)/sizeof(Double_t);
+ printf("<I> AliEMCALRecPointsQaESDSelector::DrawSpaceResolution() | np %i \n", np);
+
+ Double_t* eSpAng = new Double_t[np];
+ Double_t* eSpAngOpt = new Double_t[np];
+ Double_t C=TMath::Sqrt(8000.), C2 = 1000.*TMath::DegToRad();
+ for(int i=0; i<np; i++){
+ spAng[i] *= C2;
+ rmsSpAng[i] *= C2;
+ spAngOpt[i] *= C2;
+ rmsSpAngOpt[i] *= C2;
+
+ eSpAng[i] = spAng[i]/C;
+ eSpAngOpt[i] = spAngOpt[i]/C;
+ }
+
+ TCanvas *c = new TCanvas("Deff","Deff", 20,20, 700, 500);
+ //c->Divide(2,1);
+
+ c->cd(1);
+ gStyle->SetOptFit(111);
+
+ TGraphErrors *gr = u::DrawGraphErrors(np, p,spAng, ep, eSpAng,
+ kBlack, 21,"AP","Angle resolution (mrad) vs E_{#gamma} ","E_{#gamma} "," anlgle resolution (mrad) ",
+ -1, "", 0);
+ gr->GetHistogram()->SetMaximum(4.);
+ gr->GetHistogram()->SetMinimum(0.6);
+ gPad->SetLogx(1);
+ gPad->SetLogy(1);
+
+ TF1 *fang = new TF1("fang","[0]+[1]/sqrt(x)",0.1, 101.);
+ fang->SetLineColor(kRed);
+ fang->SetLineWidth(1);
+
+ TGraphErrors *gr2 = u::DrawGraphErrors(np, p,spAngOpt, ep, eSpAngOpt,
+ kRed,22,"P"," space vs E_{#gamma} ","E_{#gamma} "," anlgle resolution (mrad) ",
+ -1, "+", fang->GetName());
+ TLegend *leg1 = new TLegend(0.40,0.57, 0.92,0.87);
+ TLegendEntry *le1 = leg1->AddEntry(gr, "Initial angle resolution", "p");
+ le1->SetTextColor(gr->GetLineColor());
+
+ TLegendEntry *le2 = leg1->AddEntry(gr2, "Optimized angle resolution", "p");
+ le2->SetTextColor(gr2->GetLineColor());
+
+ TLegendEntry *le3 = leg1->AddEntry(fang,
+ Form("%3.2f+#frac{%4.2f}{#sqrt{E}}",fang->GetParameter(0),fang->GetParameter(1)), "lp");
+ le3->SetTextColor(fang->GetLineColor());
+
+ leg1->Draw();
+
+ c->Update();
+
+ return c;
+}
+
+void AliEMCALRecPointsQaESDSelector::ResetAllListOfHists()
+{
+ u::ResetListOfHists(fLofHistsPC);
+ u::ResetListOfHists(fLofHistsRP);
+ u::ResetListOfHists(fLKineVsRP);
+ u::ResetListOfHists(fLShowerProfile);
+}
+
+void AliEMCALRecPointsQaESDSelector::ReloadChain(Long64_t entry)
+{ // unused now
+ if(fChain) {
+ fChain->LoadTree(entry);
+ }
+}
+
+void AliEMCALRecPointsQaESDSelector::GetInitialParsForFit(const Int_t var, Double_t &deff, Double_t &w0, Double_t &phislope, const int phiCase)
+{
+ // int phiCase=0; // 0 or 1
+ phislope = 0.001;
+ if(var==11) { // stay here
+ deff = 9.07515;
+ w0 = 3.44679;
+ } else if(var==12) {
+ deff = 1.00835e+01;
+ w0 = 3.82714e+00;
+ } else if(var==13) {
+ deff = 1.12032e+01;
+ w0 = 4.15035e+00;
+ } else if(var==14) {
+ deff = 1.14229e+01;
+ w0 = 4.36650;
+ } else if(var==15) {
+ deff = 1.21361e+01;
+ w0 = 4.72875;
+ } else if(var==16) {
+ deff = 1.28096e+01;
+ w0 = 4.81753e+00;
+ } else if(var==17) {
+ deff = 1.33281e+01;
+ w0 = 4.63289e+00;
+ } else if(var==18) {
+ deff = 1.37910e+01;
+ w0 = 4.66568;
+ } else if(var==19) {// Aug 15, 2007
+ switch (phiCase) {
+ case 0: // phislope was defined than deff and w0 were fixed
+ deff = 1.43319e+01;
+ w0 = 4.69279;
+ phislope = 2.41419e-04;
+ break;
+ case 1: // all parameters free
+ deff = 1.46327e+01;
+ w0 = 4.68125;
+ phislope = 8.95559e-04;
+ break;
+ }
+ } else {
+ printf("<E> var % i -> no definition ! \n", var);
+ assert(0);
+ }
+ printf("<I> AliEMCALRecPointsQaESDSelector::GetInitialParForFit()\n deff %f w0 %f phislope %f \n",
+ deff,w0, phislope);
+}
#include "AliSelector.h"
+#include <TObjArray.h>
+
class AliEMCALFolder;
+class AliRunLoader;
+class AliEMCALRecPoint;
+class AliEMCALCellInfo;
+class cellInfo;
class TList;
+class TCanvas;
+class TH1F;
+class TH2F;
class TBrowser;
class TChain;
-class TObjectSet;
+class TFolder;
+class TArrayI;
class AliEMCALRecPointsQaESDSelector : public AliSelector {
public:
AliEMCALRecPointsQaESDSelector();
virtual ~AliEMCALRecPointsQaESDSelector();
- virtual void Begin(TTree* tree);
- virtual void Init(TTree* tree);
- virtual void SlaveBegin(TTree *tree);
+ virtual void Begin(TTree*);
+ virtual void SlaveBegin(TTree* tree);
+ virtual void Init(TTree *tree);
+ virtual Bool_t Notify();
virtual Bool_t Process(Long64_t entry);
virtual void SlaveTerminate();
virtual void Terminate();
//
void InitStructure(Int_t it);
- static TList *DefineHistsOfRP(const char *name="RP");
- // static TList *DefineHistsOfTowers(const char *name="towers");
+ static TList *DefineHistsOfRP(const char *name="RP", Double_t p=110.0, Int_t keyOpt=0);
+ static TList *DefineHistsOfKineVsRP(const char *name="KineVsRP", Double_t p=110.0, Int_t keyOpt=0);
+ static TList *DefineHistsForShowerProfile(const char *name="ProfY", Double_t p=1.);
+ static void FillHistsOfKineVsRP(TList *l, AliRunLoader* RL, TClonesArray &lvM);
+ static void FillHistsForShowerProfile(TList *l, AliEMCALRecPoint *rp, AliEMCALCellInfo* t);
+ static void EvalLocalPhiPosition(const Double_t wlog, const AliEMCALRecPoint *rp, const AliEMCALCellInfo* t, Double_t &xcog, Int_t &phiSize, cellInfo &rMax);
+ //
+ TList *GetListKineVsRP() {return fLKineVsRP;}
+ TList *GetListShowerProfile() {return fLShowerProfile;}
+ // static TList *DefineHistsOfTowers(const char *name="towers");
//
- void FitEffMassHist(); //*MENU*
+ void FitEffMassHist(); // *MENU*
void PrintInfo(); // *MENU*
//
void SetChain(TChain *chain) {fChain = chain;}
TChain* GetChain() {return fChain;}
+ void SetMomentum(Double_t p);
+ Double_t GetMomentum() {return fPmom;}
AliEMCALFolder* CreateEmcalFolder(const Int_t it);
void SetEmcalFolder(AliEMCALFolder* folder);
void SetEmcalOldFolder(AliEMCALFolder* folder);
- AliEMCALFolder* GetEmcalFolder() {return fEMCAL;}
- AliEMCALFolder* GetEmcalOldFolder() {return fEMCALOld;}
AliEMCALFolder* GetEmcalOldFolder(const Int_t nsm);
- //
+ AliEMCALCellInfo* GetCellsInfo() {return fCellsInfo;}
+ //
+ void SetStringOfRunOpts(const char *st) {fRunOpts = st;}
+ TObjArray GetOptsArray() const {return fArrOpts;}
+ Int_t GetKeyOptsValue(Int_t key); // *MENU*
+ void CheckRunOpts();
+ //
virtual void Browse(TBrowser* b);
virtual Bool_t IsFolder() const;
+ //
+ void Save(Int_t ver=0, const char *optIO="NEW"); // *MENU*
+ static AliEMCALRecPointsQaESDSelector* ReadSelector(const char* nf = "/home/pavlinov/ALICE/SHISHKEBAB/RF/CALIB/PROF/PROFILE_0.root");
+
//
//// Pictures staf - Jun 26, 2007
//
void ReadAllEmcalFolders();
void PictVsIterNumber(const Int_t ind=0, const Int_t nsm=0);
-
+ // Gamma
+ TH1F* FitHistOfRecPointEnergy(const char *opt="CLONE");
+ static TCanvas *Linearity(TList *l, Int_t ifun=3);
+ static TCanvas *DrawKineVsRP(TList *l);
+ // Profile
+ static TCanvas *DrawMeanVsLog(TH2F *h2);
+ // Geometry staff
+ TCanvas *DrawPhiEtaAnglesDistribution(const char *gn="SHISH_TRD1_CURRENT_2X2"); // *MENU*
+ // Geometry constants
+ TCanvas *DrawDeffVsEnergy(); // *MENU*
+ TCanvas *DrawDeffVsEnergy2(const char *opt="fit1"); // *MENU*
+ void ReadParsDeffAndW0(const char *dirName="/data/r22b/ALICE/CALIB/FIT/",
+ double *deff=0, double *edeff=0, double *w0=0, double *ew0=0, const Int_t pri=0);
+ TCanvas *DrawSpaceResolution();
protected:
+ Double_t fPmom; // positive if defined
//
TChain* fChain; //! chain if ESD files
- TList* fLofHistsPC; //! list of histograms of pseudo clusters
- TList* fLofHistsRP; //! list of histograms of rec.points
+ TList* fLofHistsPC; // list of histograms of pseudo clusters
+ TList* fLofHistsRP; // list of histograms of rec.points
+ TList* fLKineVsRP; // list of histograms kinematics vs rec.points
+ TList* fLShowerProfile; // list of histograms for shower profile business
+ //
+ AliEMCALCellInfo *fCellsInfo; //
+ TFolder* fEmcalPool; //
+ //
+ // Options - Jul 10, 2007
//
- TObjectSet* fEmcalPool;
- AliEMCALFolder* fEMCAL; //! current EMCAL object
- AliEMCALFolder* fEMCALOld; //! previous EMCAL object
+ TString fRunOpts; // String of running options
+ TObjArray fArrOpts; // Array of options
+ // Options keys
+ TArrayI *fKeyOpts;
+ // Static parameters
+ public:
+ static AliEMCALFolder* GetEmcalFolder() {return fEMCAL;}
+ static AliEMCALFolder* GetEmcalOldFolder() {return fEMCALOld;}
+ static AliEMCALFolder* fEMCAL; // current EMCAL object
+ static AliEMCALFolder* fEMCALOld; // previous EMCAL object
+ //
+ static Double_t fDistEff; // effective depth of electromagnetic shower
+ static Double_t fW0; //
+ static Double_t fSlopePhiShift; // phi shift of cluster = fSlopePhiShift * phi
+ static void SetFitParameters(Double_t deff, Double_t w0, Double_t slope)
+ {
+ fDistEff = deff; fW0 = w0; fSlopePhiShift = slope;
+ }
+ static void GetFitParameters(Double_t &deff, Double_t &w0, Double_t &slope)
+ {
+ deff = fDistEff; w0 = fW0; slope = fSlopePhiShift;
+ }
+ void ResetAllListOfHists();
+ void ReloadChain(Long64_t entry=0);
+ void GetInitialParsForFit(const Int_t var, Double_t &deff, Double_t &w0, Double_t &phislope, const int phiCase=0);
private:
AliEMCALRecPointsQaESDSelector(const AliEMCALRecPointsQaESDSelector&);
ClassDef(AliEMCALRecPointsQaESDSelector, 1);
};
-
#endif
* provided "as is" without express or implied warranty. *
**************************************************************************/
-/* $Id$ */
+/*
+$Log$
+*/
//_________________________________________________________________________
// Top EMCAL folder which will keep all information about EMCAL itself,
// super Modules (SM), modules, towers, set of hists and so on.
+// TObjectSet -> TFolder; Sep 6, 2007
//
//*-- Author: Aleksei Pavlinov (WSU, Detroit, USA)
ClassImp(AliEMCALSuperModule)
-AliEMCALSuperModule::AliEMCALSuperModule() : TObjectSet(),
-fSMNumber(0)
+AliEMCALSuperModule::AliEMCALSuperModule() : TFolder()
+, fParent(0),fLh(0),fSMNumber(0)
{
}
AliEMCALSuperModule::AliEMCALSuperModule(const Int_t m, const char* title) :
-TObjectSet(Form("SM%2.2i",m)),
-fSMNumber(m)
-{
- SetTitle(title);
+TFolder(Form("SM%2.2i",m), title)
+, fParent(0),fLh(0),fSMNumber(m)
+{
}
AliEMCALSuperModule::~AliEMCALSuperModule()
void AliEMCALSuperModule::Init()
{
- if(GetHists()==0) this->AddObject(BookHists(), kTRUE);
+ if(GetHists()==0) {
+ fLh = BookHists();
+ Add(fLh);
+ }
}
void AliEMCALSuperModule::AddCellToEtaRow(AliEMCALCell *cell, const Int_t etaRow)
{
- static TDataSet *set;
- set = FindByName(Form("ETA%2.2i",etaRow));
+ static TFolder *set;
+ set = dynamic_cast<TFolder*>(FindObject(Form("ETA%2.2i",etaRow)));
if(set==0) {
- set = new TDataSet(Form("ETA%2.2i",etaRow));
+ set = new TFolder(Form("ETA%2.2i",etaRow),"eta row");
Add(set);
}
set->Add(cell);
{
Int_t ncells=0;
for(int eta=0; eta<48; eta++) { // eta row
- TDataSet *setEta = FindByName(Form("ETA%2.2i",eta));
+ TFolder *setEta = dynamic_cast<TFolder*>(FindObject(Form("ETA%2.2i",eta)));
if(setEta) {
- printf(" eta %i : %s : cells %i \n", eta, setEta->GetName(), setEta->GetListSize());
- for(int phi=0; phi<setEta->GetListSize(); phi++) { // cycle on cells (phi directions)
- AliEMCALCell* cell = (AliEMCALCell*)setEta->At(phi);
+ TList* l = (TList*)setEta->GetListOfFolders();
+ printf(" eta %i : %s : cells %i \n", eta, setEta->GetName(), l->GetSize());
+ for(int phi=0; phi<l->GetSize(); phi++) { // cycle on cells (phi directions)
+ AliEMCALCell* cell = dynamic_cast<AliEMCALCell*>(l->At(phi));
+ if(cell == 0) continue;
+
cell->FitEffMassHist();
u::FillH1(GetHists(), 1, cell->GetCcIn()*1.e+3);
TH1 *hCCIn = (TH1*)GetHists()->At(1);
TH1 *hCCOut = (TH1*)GetHists()->At(2);
+ if(hCCIn == 0) return;
+ if(hCCIn->GetEntries()<10.) return;
+
hCCIn->SetStats(kFALSE);
hCCOut->SetStats(kFALSE);
- hCCOut->SetTitle("CC in and out");
+ hCCOut->SetTitle("CC in and out; Iter 1; Jul 26; All Statistics");
hCCOut->SetXTitle("cc in MeV");
hCCOut->SetYTitle(" N ");
TLegend *leg = new TLegend(0.5,0.36, 0.97,0.80);
TLegendEntry *leIn = leg->AddEntry(hCCIn, Form("input cc : %6.2f #pm %6.2f", hCCIn->GetMean(),hCCIn->GetRMS()), "L");
leIn->SetTextColor(hCCIn->GetLineColor());
+
+ if(hCCOut->GetEntries()>10.)
leg->AddEntry(hCCOut, Form("output cc : %6.2f #pm %6.2f", hCCOut->GetMean(),hCCOut->GetRMS()), "L");
+
leg->Draw();
if(c) c->Update();
Int_t AliEMCALSuperModule::GetNumberOfCells()
{
Int_t ncells=0;
- for(int eta=0; eta<GetListSize(); eta++) { // cycle on eta row
- TDataSet *set = At(eta);
- ncells += set->GetListSize();
+ TList* l = (TList*)GetListOfFolders();
+ for(int eta=0; eta<l->GetSize(); eta++) { // cycle on eta row
+ TFolder *setEta = dynamic_cast<TFolder*>(l->At(eta));
+ if(setEta==0) continue;
+
+ TList* le = (TList*)setEta->GetListOfFolders();
+ ncells += le->GetSize();
}
return ncells;
}
//_________________________________________________________________________
// Emcal Super Module
//
-//*-- Author: Aleksei Pavlinov (WSU, Detroit, USA)
+//*-- Author: Aleksei Pavlinov (WSU, Detroit, USA)
+// Super Module folder
+// Initial version was created with TDataSet staf
+// TObjectSet -> TFolder; Sep 6, 2007
-// --- ROOT system ---
-#include <TObjectSet.h>
+#include <TFolder.h>
class TList;
class AliEMCALCell;
-class AliEMCALSuperModule : public TObjectSet {
+class AliEMCALSuperModule : public TFolder {
public:
void Init();
void AddCellToEtaRow(AliEMCALCell *cell, const Int_t etaRow);
- TList* GetHists() {return (TList*)fObj;}
+ TList* GetHists() {return fLh;}
+ TObject* GetParent() {return fParent;}
+ void SetParent(TObject *parent) {fParent=parent;}
// MENU
void FitForAllCells(); //*MENU*
void FitEffMassHist(); //*MENU*
protected:
TList* BookHists();
//
- Int_t fSMNumber;
+ TObject* fParent; // parent
+ TList* fLh; // List of hists
+ Int_t fSMNumber;
- ClassDef(AliEMCALSuperModule,1) // EMCAL SuperModule
+ ClassDef(AliEMCALSuperModule,2) // EMCAL SuperModule
};
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff