AliPhysics  1adf5bd (1adf5bd)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskEmcalClustersRef.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2015, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 #include <array>
16 #include <bitset>
17 #include <iostream>
18 #include <map>
19 #include <sstream>
20 #include <vector>
21 
22 #include <TArrayD.h>
23 #include <TClonesArray.h>
24 #include <TGrid.h>
25 #include <THashList.h>
26 #include <THistManager.h>
27 #include <TLinearBinning.h>
28 #include <TLorentzVector.h>
29 #include <TObjArray.h>
30 #include <TParameter.h>
31 #include <TMath.h>
32 
33 #include "AliAnalysisUtils.h"
34 #include "AliCentrality.h"
35 #include "AliClusterContainer.h"
36 #include "AliEMCALGeometry.h"
37 #include "AliEMCALTriggerPatchInfo.h"
39 #include "AliESDEvent.h"
40 #include "AliInputEventHandler.h"
41 #include "AliLog.h"
42 #include "AliOADBContainer.h"
43 #include "AliVCluster.h"
44 #include "AliVVertex.h"
45 #include "AliMultSelection.h"
46 #include "AliMultEstimator.h"
47 
49 
53 
54 namespace EMCalTriggerPtAnalysis {
55 
59 AliAnalysisTaskEmcalClustersRef::AliAnalysisTaskEmcalClustersRef() :
61  fNameClusterContainer(""),
62  fCentralityRange(-999., 999.),
63  fRequestCentrality(false),
64  fEventCentrality(-1)
65 {
66 }
67 
74  fNameClusterContainer(""),
75  fCentralityRange(-999., 999.),
76  fRequestCentrality(false),
77  fEventCentrality(-1)
78 {
79 }
80 
85 }
86 
91 
92  EnergyBinning energybinning;
93  TLinearBinning smbinning(21, -0.5, 20.5), etabinning(100, -0.7, 0.7);
94 
95  /*
96  * Exclusive classes are defined as triggered events
97  * in a class without lower threshold classes firing.
98  * This is needed to make the triggers statistically
99  * independent.
100  */
101  std::array<Double_t, 5> encuts = {1., 2., 5., 10., 20.};
102  Int_t sectorsWithEMCAL[10] = {4, 5, 6, 7, 8, 9, 13, 14, 15, 16};
103  for(auto trg : GetSupportedTriggers()){
104  fHistos->CreateTH1(Form("hEventCount%s", trg.Data()), Form("Event count for trigger class %s", trg.Data()), 1, 0.5, 1.5);
105  fHistos->CreateTH1(Form("hEventCentrality%s", trg.Data()), Form("Event centrality for trigger class %s", trg.Data()), 103, -2., 101.);
106  fHistos->CreateTH1(Form("hVertexZ%s", trg.Data()), Form("z-position of the primary vertex for trigger class %s", trg.Data()), 200, -40., 40.);
107  fHistos->CreateTH1(Form("hClusterEnergy%s", trg.Data()), Form("Cluster energy for trigger class %s", trg.Data()), energybinning);
108  fHistos->CreateTH1(Form("hClusterET%s", trg.Data()), Form("Cluster transverse energy for trigger class %s", trg.Data()), energybinning);
109  fHistos->CreateTH1(Form("hClusterEnergyFired%s", trg.Data()), Form("Cluster energy for trigger class %s, firing the trigger", trg.Data()), energybinning);
110  fHistos->CreateTH1(Form("hClusterETFired%s", trg.Data()), Form("Cluster transverse energy for trigger class %s, firing the trigger", trg.Data()), energybinning);
111  fHistos->CreateTH2(Form("hClusterEnergySM%s", trg.Data()), Form("Cluster energy versus supermodule for trigger class %s", trg.Data()), smbinning, energybinning);
112  fHistos->CreateTH2(Form("hClusterETSM%s", trg.Data()), Form("Cluster transverse energy versus supermodule for trigger class %s", trg.Data()), smbinning, energybinning);
113  fHistos->CreateTH2(Form("hClusterEnergyFiredSM%s", trg.Data()), Form("Cluster energy versus supermodule for trigger class %s, firing the trigger", trg.Data()), smbinning, energybinning);
114  fHistos->CreateTH2(Form("hClusterETFiredSM%s", trg.Data()), Form("Cluster transverse energy versus supermodule for trigger class %s, firing the trigger", trg.Data()), smbinning, energybinning);
115  fHistos->CreateTH2(Form("hEtaEnergy%s", trg.Data()), Form("Cluster energy vs. eta for trigger class %s", trg.Data()), etabinning, energybinning);
116  fHistos->CreateTH2(Form("hEtaET%s", trg.Data()), Form("Cluster transverse energy vs. eta for trigger class %s", trg.Data()), etabinning, energybinning);
117  fHistos->CreateTH2(Form("hEtaEnergyFired%s", trg.Data()), Form("Cluster energy vs. eta for trigger class %s, firing the trigger", trg.Data()), etabinning, energybinning);
118  fHistos->CreateTH2(Form("hEtaETFired%s", trg.Data()), Form("Cluster transverse energy vs. eta for trigger class %s, firing the trigger", trg.Data()), etabinning, energybinning);
119  for(int ism = 0; ism < 20; ism++){
120  fHistos->CreateTH2(Form("hEtaEnergySM%d%s", ism, trg.Data()), Form("Cluster energy vs. eta in Supermodule %d for trigger %s", ism, trg.Data()), etabinning, energybinning);
121  fHistos->CreateTH2(Form("hEtaETSM%d%s", ism, trg.Data()), Form("Cluster transverse energy vs. eta in Supermodule %d for trigger %s", ism, trg.Data()), etabinning, energybinning);
122  fHistos->CreateTH2(Form("hEtaEnergyFiredSM%d%s", ism, trg.Data()), Form("Cluster energy vs. eta in Supermodule %d for trigger %s, firing the trigger", ism, trg.Data()), etabinning, energybinning);
123  fHistos->CreateTH2(Form("hEtaETFiredSM%d%s", ism, trg.Data()), Form("Cluster transverse energy vs. eta in Supermodule %d for trigger %s, firing the trigger", ism, trg.Data()), etabinning, energybinning);
124  }
125  for(int isec = 0; isec < 10; isec++){
126  fHistos->CreateTH2(Form("hEtaEnergySec%d%s", sectorsWithEMCAL[isec], trg.Data()), Form("Cluster energy vs.eta in tracking sector %d for trigger %s", sectorsWithEMCAL[isec], trg.Data()), etabinning, energybinning);
127  fHistos->CreateTH2(Form("hEtaETSec%d%s", sectorsWithEMCAL[isec], trg.Data()), Form("Cluster transverse energy vs.eta in tracking sector %d for trigger %s", sectorsWithEMCAL[isec], trg.Data()), etabinning, energybinning);
128  fHistos->CreateTH2(Form("hEtaEnergyFiredSec%d%s", sectorsWithEMCAL[isec], trg.Data()), Form("Cluster energy vs.eta in tracking sector %d for trigger %s, firing the trigger", sectorsWithEMCAL[isec], trg.Data()), etabinning, energybinning);
129  fHistos->CreateTH2(Form("hEtaETFiredSec%d%s", sectorsWithEMCAL[isec], trg.Data()), Form("Cluster transverse energy vs.eta in tracking sector %d for trigger %s, firing the trigger", sectorsWithEMCAL[isec], trg.Data()), etabinning, energybinning);
130  }
131  for(auto ien : encuts){
132  fHistos->CreateTH2(Form("hEtaPhi%dG%s", static_cast<int>(ien), trg.Data()), Form("cluster #eta-#phi map for clusters with energy larger than %f GeV/c for trigger class %s", ien, trg.Data()), 100, -0.7, 0.7, 200, 0, 2*TMath::Pi());
133  fHistos->CreateTH2(Form("hEtaPhiFired%dG%s", static_cast<int>(ien), trg.Data()), Form("cluster #eta-#phi map for clusters fired the trigger with energy larger than %f GeV/c for trigger class %s", ien, trg.Data()), 200, -0.7, 0.7, 200, 0, 2*TMath::Pi());
134  }
135  }
136 }
137 
139  fEventCentrality = -1;
140  if(fRequestCentrality){
141  AliMultSelection *mult = dynamic_cast<AliMultSelection *>(InputEvent()->FindListObject("MultSelection"));
142  if(!mult){
143  AliErrorStream() << GetName() << ": Centrality selection enabled but no centrality estimator found" << std::endl;
144  return false;
145  }
146  if(mult->IsEventSelected()) return false;
147  fEventCentrality = mult->GetEstimator("V0M")->GetPercentile();
148  AliDebugStream(1) << GetName() << ": Centrality " << fEventCentrality << std::endl;
150  AliDebugStream(1) << GetName() << ": reject centrality: " << fEventCentrality << std::endl;
151  return false;
152  } else {
153  AliDebugStream(1) << GetName() << ": select centrality " << fEventCentrality << std::endl;
154  }
155  } else {
156  AliDebugStream(1) << GetName() << ": No centrality selection applied" << std::endl;
157  }
158  return true;
159 }
160 
161 
167  AliDebugStream(1) << GetName() << ": UserExec start" << std::endl;
168 
169  TList ej1patches, dj1patches, ej2patches, dj2patches, eg1patches, dg1patches, eg2patches, dg2patches;
170  FindPatchesForTrigger("EJ2", fTriggerPatchInfo, ej2patches);
171  FindPatchesForTrigger("DJ2", fTriggerPatchInfo, dj2patches);
172  FindPatchesForTrigger("EJ1", fTriggerPatchInfo, ej1patches);
173  FindPatchesForTrigger("DJ1", fTriggerPatchInfo, dj1patches);
174  FindPatchesForTrigger("EG2", fTriggerPatchInfo, eg2patches);
175  FindPatchesForTrigger("DG2", fTriggerPatchInfo, dg2patches);
176  FindPatchesForTrigger("EG1", fTriggerPatchInfo, eg1patches);
177  FindPatchesForTrigger("DG1", fTriggerPatchInfo, dg1patches);
178 
179  Double_t energy, et, eta, phi;
180  TList *selpatches(nullptr);
181  for(auto clust : GetClusterContainer(fNameClusterContainer.Data())->all()){
182  //AliVCluster *clust = static_cast<AliVCluster *>(*clustIter);
183  if(!clust->IsEMCAL()) continue;
184  if(clust->GetIsExotic()) continue;
185 
186  TLorentzVector posvec;
187  energy = clust->GetNonLinCorrEnergy();
188  et = posvec.Et();
189  clust->GetMomentum(posvec, fVertex);
190  eta = posvec.Eta();
191  phi = posvec.Phi();
192 
193  // fill histograms allEta
194  for(const auto & trg : fSelectedTriggers){
195  selpatches = nullptr;
196  if(trg.Contains("EJ2") != std::string::npos) selpatches = &ej2patches;
197  if(trg.Contains("DJ2") != std::string::npos) selpatches = &dj2patches;
198  if(trg.Contains("EJ1") != std::string::npos) selpatches = &ej1patches;
199  if(trg.Contains("DJ1") != std::string::npos) selpatches = &dj1patches;
200  if(trg.Contains("EG2") != std::string::npos) selpatches = &eg2patches;
201  if(trg.Contains("DG2") != std::string::npos) selpatches = &dg2patches;
202  if(trg.Contains("EG1") != std::string::npos) selpatches = &eg1patches;
203  if(trg.Contains("DG1") != std::string::npos) selpatches = &dg1patches;
204  FillClusterHistograms(trg.Data(), energy, et, eta, phi, nullptr);
205  }
206  }
207  return true;
208 }
209 
210 void AliAnalysisTaskEmcalClustersRef::FillClusterHistograms(const TString &triggerclass, double energy, double transverseenergy, double eta, double phi, TList *fTriggerPatches){
211  Bool_t hasTriggerPatch = fTriggerPatches ? CorrelateToTrigger(eta, phi, fTriggerPatches) : kFALSE;
212  Int_t supermoduleID = -1, sector = -1;
213  Double_t weight = GetTriggerWeight(triggerclass);
214  AliDebugStream(1) << GetName() << ": Using weight " << weight << " for trigger " << triggerclass << std::endl;
215 
216  fGeom->SuperModuleNumberFromEtaPhi(eta, phi, supermoduleID);
217  fHistos->FillTH1(Form("hClusterEnergy%s", triggerclass.Data()), energy, weight);
218  fHistos->FillTH1(Form("hClusterET%s", triggerclass.Data()), transverseenergy, weight);
219  fHistos->FillTH2(Form("hEtaEnergy%s", triggerclass.Data()), eta, energy, weight);
220  fHistos->FillTH2(Form("hEtaET%s", triggerclass.Data()), eta, transverseenergy, weight);
221  if(supermoduleID >= 0){
222  fHistos->FillTH2(Form("hClusterEnergySM%s", triggerclass.Data()), supermoduleID, energy, weight);
223  fHistos->FillTH2(Form("hClusterETSM%s", triggerclass.Data()), supermoduleID, transverseenergy, weight);
224  fHistos->FillTH2(Form("hEtaEnergySM%d%s", supermoduleID, triggerclass.Data()), eta, energy, weight);
225  fHistos->FillTH2(Form("hEtaETSM%d%s", supermoduleID, triggerclass.Data()), eta, transverseenergy, weight);
226  if(supermoduleID < 12)
227  sector = 4 + int(supermoduleID/2); // EMCAL
228  else
229  sector = 13 + int((supermoduleID-12)/2); // DCAL
230  fHistos->FillTH2(Form("hEtaEnergySec%d%s", sector, triggerclass.Data()), eta, energy, weight);
231  fHistos->FillTH2(Form("hEtaETSec%d%s", sector, triggerclass.Data()), eta, transverseenergy, weight);
232  }
233  if(hasTriggerPatch){
234  fHistos->FillTH1(Form("hClusterEnergyFired%s", triggerclass.Data()), energy, weight);
235  fHistos->FillTH1(Form("hClusterETFired%s", triggerclass.Data()), energy, weight);
236  fHistos->FillTH2(Form("hEtaEnergyFired%s", triggerclass.Data()), eta, energy, weight);
237  fHistos->FillTH2(Form("hEtaETFired%s", triggerclass.Data()), eta, energy, weight);
238  if(supermoduleID >= 0){
239  fHistos->FillTH2(Form("hClusterEnergyFiredSM%s", triggerclass.Data()), supermoduleID, energy, weight);
240  fHistos->FillTH2(Form("hClusterETFiredSM%s", triggerclass.Data()), supermoduleID, transverseenergy, weight);
241  fHistos->FillTH2(Form("hEtaEnergyFiredSM%d%s", supermoduleID, triggerclass.Data()), eta, energy,weight);
242  fHistos->FillTH2(Form("hEtaETFiredSM%d%s", supermoduleID, triggerclass.Data()), eta, transverseenergy, weight);
243  fHistos->FillTH2(Form("hEtaEnergyFiredSec%d%s", sector, triggerclass.Data()), eta, energy, weight);
244  fHistos->FillTH2(Form("hEtaETFiredSec%d%s", sector, triggerclass.Data()), eta, transverseenergy, weight);
245  }
246  }
247  Double_t encuts[5] = {1., 2., 5., 10., 20.};
248  for(int ien = 0; ien < 5; ien++){
249  if(energy > encuts[ien]){
250  fHistos->FillTH2(Form("hEtaPhi%dG%s", static_cast<int>(encuts[ien]), triggerclass.Data()), eta, phi, weight);
251  if(hasTriggerPatch){
252  fHistos->FillTH2(Form("hEtaPhiFired%dG%s", static_cast<int>(encuts[ien]), triggerclass.Data()), eta, phi, weight);
253  }
254  }
255  }
256 }
257 
267  for(const auto &t : fSelectedTriggers){
268  Double_t weight = GetTriggerWeight(t);
269  fHistos->FillTH1(Form("hEventCount%s", t.Data()), 1, weight);
270  fHistos->FillTH1(Form("hEventCentrality%s", t.Data()), fEventCentrality, weight);
271  fHistos->FillTH1(Form("hVertexZ%s", t.Data()), fVertex[2], weight);
272  }
273 }
274 
283  Bool_t hasfound = kFALSE;
284  for(TIter patchIter = TIter(fTriggerPatches).Begin(); patchIter != TIter::End(); ++patchIter){
285  Double_t boundaries[4];
286  GetPatchBoundaries(*patchIter, boundaries);
287  Double_t etamin = TMath::Min(boundaries[0], boundaries[1]),
288  etamax = TMath::Max(boundaries[0], boundaries[1]),
289  phimin = TMath::Min(boundaries[2], boundaries[3]),
290  phimax = TMath::Max(boundaries[2], boundaries[3]);
291  if(etaclust > etamin && etaclust < etamax && phiclust > phimin && phiclust < phimax){
292  hasfound = kTRUE;
293  break;
294  }
295  }
296  return hasfound;
297 }
298 
307 void AliAnalysisTaskEmcalClustersRef::FindPatchesForTrigger(TString triggerclass, const TClonesArray * fTriggerPatches, TList &foundtriggers) const {
308  foundtriggers.Clear();
309  if(!fTriggerPatches) return;
311  if(triggerclass == "EG1") myclass = AliEmcalTriggerOfflineSelection::kTrgEG1;
312  if(triggerclass == "EG2") myclass = AliEmcalTriggerOfflineSelection::kTrgEG2;
313  if(triggerclass == "EJ1") myclass = AliEmcalTriggerOfflineSelection::kTrgEJ1;
314  if(triggerclass == "EJ2") myclass = AliEmcalTriggerOfflineSelection::kTrgEJ2;
315  if(triggerclass == "DMC7") myclass = AliEmcalTriggerOfflineSelection::kTrgDL0;
316  if(triggerclass == "DG1") myclass = AliEmcalTriggerOfflineSelection::kTrgDG1;
317  if(triggerclass == "DG2") myclass = AliEmcalTriggerOfflineSelection::kTrgDG2;
318  if(triggerclass == "DJ1") myclass = AliEmcalTriggerOfflineSelection::kTrgDJ1;
319  if(triggerclass == "DJ2") myclass = AliEmcalTriggerOfflineSelection::kTrgDJ2;
320  for(TIter patchiter = TIter(fTriggerPatches).Begin(); patchiter != TIter::End(); ++patchiter){
321  if(!IsOfflineSimplePatch(*patchiter)) continue;
323  if(!SelectDCALPatch(*patchiter)) continue;
324  } else {
325  if(SelectDCALPatch(*patchiter)) continue;
326  }
328  if(!SelectSingleShowerPatch(*patchiter)) continue;
329  } else {
330  if(!SelectJetPatch(*patchiter)) continue;
331  }
332  double threshold = fTriggerSelection ? fTriggerSelection->GetThresholdForTrigger(myclass) : -1;
333  if(GetPatchEnergy(*patchiter) > threshold) foundtriggers.Add(*patchiter);
334  }
335 }
336 
338  AliEMCALTriggerPatchInfo *patch= dynamic_cast<AliEMCALTriggerPatchInfo *>(o);
339  boundaries[0] = patch->GetEtaMin();
340  boundaries[1] = patch->GetEtaMax();
341  boundaries[2] = patch->GetPhiMin();
342  boundaries[3] = patch->GetPhiMax();
343 }
344 
346  AliEMCALTriggerPatchInfo *patch = dynamic_cast<AliEMCALTriggerPatchInfo *>(o);
347  return patch->IsOfflineSimple();
348 }
349 
351  AliEMCALTriggerPatchInfo *patch = dynamic_cast<AliEMCALTriggerPatchInfo *>(o);
352  return patch->GetRowStart() >= 64;
353 }
354 
356  AliEMCALTriggerPatchInfo *patch = dynamic_cast<AliEMCALTriggerPatchInfo *>(o);
357  return patch->IsGammaLowSimple();
358 }
359 
361  AliEMCALTriggerPatchInfo *patch = dynamic_cast<AliEMCALTriggerPatchInfo *>(o);
362  if(!patch->IsOfflineSimple()) return false;
363  return patch->IsJetLowSimple();
364 }
365 
367  double energy = 0.;
368  AliEMCALTriggerPatchInfo *patch = dynamic_cast<AliEMCALTriggerPatchInfo *>(o);
369  energy = patch->GetPatchE();
370  return energy;
371 }
372 
378 {
379  this->SetMinimum(0.);
380  this->AddStep(1, 0.05);
381  this->AddStep(2, 0.1);
382  this->AddStep(4, 0.2);
383  this->AddStep(7, 0.5);
384  this->AddStep(16, 1);
385  this->AddStep(32, 2);
386  this->AddStep(40, 4);
387  this->AddStep(50, 5);
388  this->AddStep(100, 10);
389  this->AddStep(200, 20);
390 }
391 
392 
393 } /* namespace EMCalTriggerPtAnalysis */
void GetPatchBoundaries(TObject *o, Double_t *boundaries) const
std::vector< TString > fSelectedTriggers
! Triggers selected for given event
double Double_t
Definition: External.C:58
Class creating a linear binning, used in the histogram manager.
Bool_t CorrelateToTrigger(Double_t etaclust, Double_t phiclust, TList *triggerpatches) const
void FillTH2(const char *hname, double x, double y, double weight=1., Option_t *opt="")
AliEmcalTriggerOfflineSelection * fTriggerSelection
Offline trigger selection.
Bool_t fRequestCentrality
Switch on request for centrality range.
static Bool_t IsDCAL(EmcalTriggerClass cls)
void AddStep(Double_t max, Double_t binwidth)
void FillClusterHistograms(const TString &triggerclass, double energy, double transversenergy, double eta, double phi, TList *triggerpatches)
const Double_t etamin
TH2 * CreateTH2(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, Option_t *opt="")
int Int_t
Definition: External.C:63
TString fNameClusterContainer
Name of the cluster container in the event.
AliEMCALGeometry * fGeom
!emcal geometry
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
AliClusterContainer * GetClusterContainer(Int_t i=0) const
Helper class creating user defined custom binning.
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
energy
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
Double_t fVertex[3]
!event vertex
AliCutValueRange< double > fCentralityRange
Selected centrality range.
const Double_t etamax
TClonesArray * fTriggerPatchInfo
!trigger patch info array
bool Bool_t
Definition: External.C:53
static Bool_t IsSingleShower(EmcalTriggerClass cls)
const Double_t phimin
void FindPatchesForTrigger(TString triggerclass, const TClonesArray *triggerpatches, TList &foundpatches) const
void SetMinimum(Double_t min)