AliPhysics  cdeda5a (cdeda5a)
 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 
5 
6 #include <TList.h>
7 
8 #include "AliClusterContainer.h"
9 
13 
14 // Actually registers the class with the base class
16 
17 //________________________________________________________________________
19  AliEmcalCorrectionComponent("AliEmcalCorrectionClusterNonLinearity"),
20  fEnergyDistBefore(0),
21  fEnergyTimeHistBefore(0),
22  fEnergyDistAfter(0),
23  fEnergyTimeHistAfter(0)
24 
25 {
26  // Default constructor
27  AliDebug(3, Form("%s", __PRETTY_FUNCTION__));
28 
29 }
30 
31 //________________________________________________________________________
33 {
34  // Destructor
35 }
36 
37 //________________________________________________________________________
39 {
40  // Initialization
41  AliDebug(3, Form("%s", __PRETTY_FUNCTION__));
43  // Do base class initializations and if it fails -> bail out
44  //AliAnalysisTaskEmcal::ExecOnce();
45  //if (!fInitialized) return;
46 
47  GetProperty("createHistos", fCreateHisto);
48 
49  std::string nonLinFunctStr = "";
50  GetProperty("nonLinFunct", nonLinFunctStr);
51  UInt_t nonLinFunct = nonlinearityFunctionMap[nonLinFunctStr];
52 
53  // init reco utils
54  if (!fRecoUtils)
55  fRecoUtils = new AliEMCALRecoUtils;
56  fRecoUtils->SetNonLinearityFunction(nonLinFunct);
57 
58  if (fRecoUtils) {
59  fRecoUtils->InitNonLinearityParam();
60  fRecoUtils->Print("");
61  }
62 
63  return kTRUE;
64 }
65 
66 //________________________________________________________________________
68 {
69  AliDebug(3, Form("%s", __PRETTY_FUNCTION__));
71 
72  // Create my user objects.
73  if (fCreateHisto){
74  fEnergyDistBefore = new TH1F("hEnergyDistBefore","hEnergyDistBefore;E_{clus} (GeV)",1500,0,150);
76  fEnergyTimeHistBefore = new TH2F("hEnergyTimeDistBefore","hEnergyTimeDistBefore;E_{clus} (GeV);time",1500,0,150,500,0,1e-6);
78  fEnergyDistAfter = new TH1F("hEnergyDistAfter","hEnergyDistAfter;E_{clus} (GeV)",1500,0,150);
80  fEnergyTimeHistAfter = new TH2F("hEnergyTimeDistAfter","hEnergyTimeDistAfter;E_{clus} (GeV);time",1500,0,150,500,0,1e-6);
82 
83  // Take ownership of output list
84  fOutput->SetOwner(kTRUE);
85  }
86 }
87 
88 //________________________________________________________________________
90 {
91  // Run
92  AliDebug(3, Form("%s", __PRETTY_FUNCTION__));
94 
95  // loop over clusters
96  AliVCluster *clus = 0;
97  AliClusterContainer * clusCont = 0;
98  TIter nextClusCont(&fClusterCollArray);
99  while ((clusCont = static_cast<AliClusterContainer*>(nextClusCont()))) {
100 
101  if (!clusCont) return kFALSE;
102  auto clusItCont = clusCont->all_momentum();
103 
104  for (AliClusterIterableMomentumContainer::iterator clusIterator = clusItCont.begin(); clusIterator != clusItCont.end(); ++clusIterator) {
105  clus = static_cast<AliVCluster *>(clusIterator->second);
106 
107  if (!clus->IsEMCAL()) continue;
108 
109  if (fCreateHisto) {
110  fEnergyDistBefore->Fill(clus->E());
111  fEnergyTimeHistBefore->Fill(clus->E(), clus->GetTOF());
112  }
113 
114  if (fRecoUtils) {
115  if (fRecoUtils->GetNonLinearityFunction() != AliEMCALRecoUtils::kNoCorrection) {
116  Double_t energy = fRecoUtils->CorrectClusterEnergyLinearity(clus);
117  clus->SetNonLinCorrEnergy(energy);
118  }
119  }
120 
121  // Fill histograms only if cluster is not exotic, as in ClusterMaker (the clusters are flagged, not removed)
122  if (fCreateHisto && !clus->GetIsExotic()) {
123  fEnergyDistAfter->Fill(clus->GetNonLinCorrEnergy());
124  fEnergyTimeHistAfter->Fill(clus->GetNonLinCorrEnergy(), clus->GetTOF());
125  }
126  }
127  }
128 
129  return kTRUE;
130 }
double Double_t
Definition: External.C:58
Definition: External.C:236
bidirectional stl iterator over the EMCAL iterable container
TH2F * fEnergyTimeHistBefore
!energy/time distribution before
AliEMCALRecoUtils * fRecoUtils
Pointer to RecoUtils.
TObjArray fClusterCollArray
Cluster collection array.
unsigned int UInt_t
Definition: External.C:33
Base class for correction components in the EMCal correction framework.
TList * fOutput
! List of output histograms
energy
Bool_t fCreateHisto
Flag to make some basic histograms.
const AliClusterIterableMomentumContainer all_momentum() const
Cluster energy non-linearity correction component in the EMCal correction framework.
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
std::map< std::string, AliEMCALRecoUtils::NonlinearityFunctions > nonlinearityFunctionMap
bool Bool_t
Definition: External.C:53
Container structure for EMCAL clusters.
static RegisterCorrectionComponent< AliEmcalCorrectionClusterNonLinearity > reg
TH2F * fEnergyTimeHistAfter
!energy/time distribution after
bool GetProperty(std::string propertyName, T &property, bool requiredProperty=true, std::string correctionName="")
Retrieve property.