AliPhysics  026afea (026afea)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliEmcalCorrectionClusterNonLinearity.cxx
Go to the documentation of this file.
1 // AliEmcalCorrectionClusterNonLinearity
2 //
3 // Author: C.Loizides, S.Aiola
4 
6 
10 
11 // Actually registers the class with the base class
13 
14 //________________________________________________________________________
16  AliEmcalCorrectionComponent("AliEmcalCorrectionClusterNonLinearity"),
17  fEnergyDistBefore(0),
18  fEnergyTimeHistBefore(0),
19  fEnergyDistAfter(0),
20  fEnergyTimeHistAfter(0)
21 
22 {
23  // Default constructor
24  AliDebug(3, Form("%s", __PRETTY_FUNCTION__));
25 
26 }
27 
28 //________________________________________________________________________
30 {
31  // Destructor
32 }
33 
34 //________________________________________________________________________
36 {
37  // Initialization
38  AliDebug(3, Form("%s", __PRETTY_FUNCTION__));
40  // Do base class initializations and if it fails -> bail out
41  //AliAnalysisTaskEmcal::ExecOnce();
42  //if (!fInitialized) return;
43 
44  GetProperty("createHistos", fCreateHisto);
45 
46  std::string nonLinFunctStr = "";
47  GetProperty("nonLinFunct", nonLinFunctStr);
48  UInt_t nonLinFunct = nonlinearityFunctionMap[nonLinFunctStr];
49 
51  Double_t minE = 0.;
52  GetProperty("clusterEMin", minE);
53  Double_t minPt = 0.;
54  GetProperty("clusterPtMin", minPt);
55 
56  // Settings from sample run macro
57  fClusCont->SetClusECut(minE);
58  fClusCont->SetClusPtCut(minPt);
59 
60  // init reco utils
61  if (!fRecoUtils)
62  fRecoUtils = new AliEMCALRecoUtils;
63  fRecoUtils->SetNonLinearityFunction(nonLinFunct);
64 
65  if (fRecoUtils) {
66  fRecoUtils->InitNonLinearityParam();
67  fRecoUtils->Print("");
68  }
69 
70  // Create my user objects.
71  if (fCreateHisto){
72  fEnergyDistBefore = new TH1F("hEnergyDistBefore","hEnergyDistBefore;E_{clus} (GeV)",1500,0,150);
74  fEnergyTimeHistBefore = new TH2F("hEnergyTimeDistBefore","hEnergyTimeDistBefore;E_{clus} (GeV);time",1500,0,150,500,0,1e-6);
76  fEnergyDistAfter = new TH1F("hEnergyDistAfter","hEnergyDistAfter;E_{clus} (GeV)",1500,0,150);
78  fEnergyTimeHistAfter = new TH2F("hEnergyTimeDistAfter","hEnergyTimeDistAfter;E_{clus} (GeV);time",1500,0,150,500,0,1e-6);
80 
81  // Take ownership of output list
82  fOutput->SetOwner(kTRUE);
83  }
84 
85  return kTRUE;
86 }
87 
88 //________________________________________________________________________
90 {
91  // Run
92  AliDebug(3, Form("%s", __PRETTY_FUNCTION__));
94 
95  if (!fClusCont) return kFALSE;
96 
97  // loop over clusters
98  fClusCont->ResetCurrentID();
99  AliVCluster *clus = 0;
100  while ((clus = fClusCont->GetNextCluster())) {
101  if (!clus->IsEMCAL()) continue;
102 
103  if (fCreateHisto) {
104  fEnergyDistBefore->Fill(clus->E());
105  fEnergyTimeHistBefore->Fill(clus->E(), clus->GetTOF());
106  }
107 
108  if (fRecoUtils) {
109  if (fRecoUtils->GetNonLinearityFunction() != AliEMCALRecoUtils::kNoCorrection) {
110  Double_t energy = fRecoUtils->CorrectClusterEnergyLinearity(clus);
111  clus->SetNonLinCorrEnergy(energy);
112  }
113  }
114 
115  // Fill histograms only if cluster is not exotic, as in ClusterMaker (the clusters are flagged, not removed)
116  if (fCreateHisto && !clus->GetIsExotic()) {
117  fEnergyDistAfter->Fill(clus->GetNonLinCorrEnergy());
118  fEnergyTimeHistAfter->Fill(clus->GetNonLinCorrEnergy(), clus->GetTOF());
119  }
120  }
121 
122  return kTRUE;
123 }
void AddContainer(inputObjectType type)
double Double_t
Definition: External.C:58
Definition: External.C:236
TH2F * fEnergyTimeHistBefore
!energy/time distribution before
AliEMCALRecoUtils * fRecoUtils
! pointer to reco utils
AliClusterContainer * fClusCont
! pointer to the cluster container
unsigned int UInt_t
Definition: External.C:33
void GetProperty(std::string propertyName, T &property, bool requiredProperty=true, std::string correctionName="")
Retrieve property.
AliVCluster * GetNextCluster()
TList * fOutput
! list of output histograms
energy
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
std::map< std::string, AliEMCALRecoUtils::NonlinearityFunctions > nonlinearityFunctionMap
void SetClusPtCut(Double_t cut)
bool Bool_t
Definition: External.C:53
void SetClusECut(Double_t cut)
static RegisterCorrectionComponent< AliEmcalCorrectionClusterNonLinearity > reg
TH2F * fEnergyTimeHistAfter
!energy/time distribution after