AliPhysics  608b256 (608b256)
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 const std::map <std::string, AliEMCALRecoUtils::NonlinearityFunctions> AliEmcalCorrectionClusterNonLinearity::fgkNonlinearityFunctionMap = {
18  { "kPi0MC", AliEMCALRecoUtils::kPi0MC },
19  { "kPi0GammaGamma", AliEMCALRecoUtils::kPi0GammaGamma },
20  { "kPi0GammaConversion", AliEMCALRecoUtils::kPi0GammaConversion },
21  { "kNoCorrection", AliEMCALRecoUtils::kNoCorrection },
22  { "kBeamTest", AliEMCALRecoUtils::kBeamTest },
23  { "kBeamTestCorrected", AliEMCALRecoUtils::kBeamTestCorrected },
24  { "kPi0MCv2", AliEMCALRecoUtils::kPi0MCv2 },
25  { "kPi0MCv3", AliEMCALRecoUtils::kPi0MCv3 },
26  { "kBeamTestCorrectedv2", AliEMCALRecoUtils::kBeamTestCorrectedv2 },
27  { "kSDMv5", AliEMCALRecoUtils::kSDMv5 },
28  { "kPi0MCv5", AliEMCALRecoUtils::kPi0MCv5 },
29  { "kSDMv6", AliEMCALRecoUtils::kSDMv6 },
30  { "kPi0MCv6", AliEMCALRecoUtils::kPi0MCv6 },
31  { "kBeamTestCorrectedv3", AliEMCALRecoUtils::kBeamTestCorrectedv3 },
32  { "kPCMv1", AliEMCALRecoUtils::kPCMv1 },
33  { "kPCMplusBTCv1", AliEMCALRecoUtils::kPCMplusBTCv1 },
34  { "kPCMsysv1", AliEMCALRecoUtils::kPCMsysv1 },
35  { "kBeamTestCorrectedv4", AliEMCALRecoUtils::kBeamTestCorrectedv4 }
36 };
37 
42  AliEmcalCorrectionComponent("AliEmcalCorrectionClusterNonLinearity"),
43  fEnergyDistBefore(0),
44  fEnergyTimeHistBefore(0),
45  fEnergyDistAfter(0),
46  fEnergyTimeHistAfter(0),
47  fSetForceClusterE(kFALSE)
48 {
49 }
50 
55 {
56 }
57 
62 {
63  // Initialization
65 
66  std::string nonLinFunctStr = "";
67  GetProperty("nonLinFunct", nonLinFunctStr);
68  UInt_t nonLinFunct = fgkNonlinearityFunctionMap.at(nonLinFunctStr);
69 
70  GetProperty("setForceClusterE", fSetForceClusterE);
71 
72  // init reco utils
73  if (!fRecoUtils)
75  fRecoUtils->SetNonLinearityFunction(nonLinFunct);
76 
77  if (fRecoUtils) {
79  fRecoUtils->Print("");
80  }
81 
82  return kTRUE;
83 }
84 
89 {
91 
92  // Create my user objects.
93  if (fCreateHisto){
94  fEnergyDistBefore = new TH1F("hEnergyDistBefore","hEnergyDistBefore;E_{clus} (GeV)",1500,0,150);
96  fEnergyTimeHistBefore = new TH2F("hEnergyTimeDistBefore","hEnergyTimeDistBefore;E_{clus} (GeV);time (s)",1500,0,150,500,-1e-6,1e-6);
98  fEnergyDistAfter = new TH1F("hEnergyDistAfter","hEnergyDistAfter;E_{clus} (GeV)",1500,0,150);
100  fEnergyTimeHistAfter = new TH2F("hEnergyTimeDistAfter","hEnergyTimeDistAfter;E_{clus} (GeV);time (s)",1500,0,150,500,-1e-6,1e-6);
102 
103  // Take ownership of output list
104  fOutput->SetOwner(kTRUE);
105  }
106 }
107 
112 {
114 
115  // loop over clusters
116  AliVCluster *clus = 0;
117  AliClusterContainer * clusCont = 0;
118  TIter nextClusCont(&fClusterCollArray);
119  while ((clusCont = static_cast<AliClusterContainer*>(nextClusCont()))) {
120 
121  if (!clusCont) return kFALSE;
122  auto clusItCont = clusCont->all_momentum();
123 
124  for (AliClusterIterableMomentumContainer::iterator clusIterator = clusItCont.begin(); clusIterator != clusItCont.end(); ++clusIterator) {
125  clus = static_cast<AliVCluster *>(clusIterator->second);
126 
127  if (!clus->IsEMCAL()) continue;
128 
129  if (fCreateHisto) {
130  fEnergyDistBefore->Fill(clus->E());
131  fEnergyTimeHistBefore->Fill(clus->E(), clus->GetTOF());
132  }
133 
134  if (fRecoUtils) {
137  clus->SetNonLinCorrEnergy(energy);
138  if ( fSetForceClusterE ) clus->SetE(energy);
139  }
140  }
141 
142  // Fill histograms only if cluster is not exotic, as in ClusterMaker (the clusters are flagged, not removed)
143  if (fCreateHisto && !clus->GetIsExotic()) {
144  Float_t energy = clus->GetNonLinCorrEnergy();
145  if(fSetForceClusterE) energy = clus->E();
146 
147  fEnergyDistAfter->Fill(energy);
148  fEnergyTimeHistAfter->Fill(energy, clus->GetTOF());
149  }
150  }
151  }
152 
153  return kTRUE;
154 }
double Double_t
Definition: External.C:58
Definition: External.C:236
static const std::map< std::string, AliEMCALRecoUtils::NonlinearityFunctions > fgkNonlinearityFunctionMap
Relates string to the non-linearity function enumeration for YAML configuration.
bidirectional stl iterator over the EMCAL iterable container
energy
Definition: HFPtSpectrum.C:45
void SetNonLinearityFunction(Int_t fun)
TH2F * fEnergyTimeHistBefore
!energy/time distribution before
Bool_t fSetForceClusterE
Only for backwards compatibility, force cluster->E() to be set to the cluster non-linearity corrected...
Some utilities for cluster and cell treatment.
AliEMCALRecoUtils * fRecoUtils
Pointer to RecoUtils.
TObjArray fClusterCollArray
Cluster collection array.
void Print(const Option_t *) const
unsigned int UInt_t
Definition: External.C:33
float Float_t
Definition: External.C:68
Base class for correction components in the EMCal correction framework.
Int_t GetNonLinearityFunction() const
TList * fOutput
! List of output histograms
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.
bool Bool_t
Definition: External.C:53
Float_t CorrectClusterEnergyLinearity(AliVCluster *clu)
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.