AliPhysics  29d4213 (29d4213)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliAnalysisTaskEventSelectionRef.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 <TArrayD.h>
16 #include <TClonesArray.h>
17 #include <THashList.h>
18 #include <TMath.h>
19 #include <TObjArray.h>
20 #include <TString.h>
21 
22 #include "AliAnalysisUtils.h"
23 #include "AliAODTrack.h"
24 #include "AliESDEvent.h"
25 #include "AliESDtrack.h"
26 #include "AliESDtrackCuts.h"
27 #include "AliEmcalTriggerPatchInfoAP.h"
28 #include "AliEMCalHistoContainer.h"
29 #include "AliEMCALGeometry.h"
30 #include "AliEMCALRecoUtils.h"
31 #include "AliInputEventHandler.h"
32 #include "AliVCluster.h"
33 #include "AliVEvent.h"
34 
36 
37 #if __cplusplus < 201103L
38 #ifndef nullptr
39 #define nullptr NULL
40 #endif
41 #endif
42 
44 
45 namespace EMCalTriggerPtAnalysis {
46 
47 AliAnalysisTaskEventSelectionRef::AliAnalysisTaskEventSelectionRef():
48  AliAnalysisTaskSE(),
49  fClusterContainerName(""),
50  fAnalysisUtils(nullptr),
51  fHistos(nullptr),
52  fTrackCuts(nullptr),
53  fGeometry(nullptr),
54  fTriggerPatchContainer(nullptr),
55  fClusterContainer(nullptr),
56  fTrackContainer(nullptr)
57 {
58  for(int itrg = 0; itrg < kCPRntrig; itrg++){
59  fOfflineEnergyThreshold[itrg] = -1;
60  }
61 }
62 
63 AliAnalysisTaskEventSelectionRef::AliAnalysisTaskEventSelectionRef(const char *name):
64  AliAnalysisTaskSE(name),
65  fClusterContainerName(""),
66  fAnalysisUtils(nullptr),
67  fHistos(nullptr),
68  fTrackCuts(nullptr),
69  fGeometry(nullptr),
70  fTriggerPatchContainer(nullptr),
71  fClusterContainer(nullptr),
72  fTrackContainer(nullptr)
73 {
74  for(int itrg = 0; itrg < kCPRntrig; itrg++){
75  fOfflineEnergyThreshold[itrg] = -1;
76  }
77  DefineOutput(1, TList::Class());
78 }
79 
80 AliAnalysisTaskEventSelectionRef::~AliAnalysisTaskEventSelectionRef() {
81  if(fAnalysisUtils) delete fAnalysisUtils;
82  if(fTrackCuts) delete fTrackCuts;
83  if(fTrackContainer) delete fTrackContainer;
84 }
85 
86 void AliAnalysisTaskEventSelectionRef::UserCreateOutputObjects(){
87  fAnalysisUtils = new AliAnalysisUtils;
88 
89  fTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(true, 1);
90  fTrackCuts->SetName("Standard Track cuts");
91  fTrackCuts->SetMinNCrossedRowsTPC(120);
92  fTrackCuts->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
93 
94  fHistos = new AliEMCalHistoContainer("Ref");
95 
96  TArrayD ptbinning, energybinning;
97  CreatePtBinning(ptbinning);
98  CreateEnergyBinning(energybinning);
99 
100  TString triggers[6] = {"MB", "EMC7", "EJ1", "EJ2", "EG1", "EG2"};
101  TString selectionStatus[3] = {"All","Accepted","Rejected"};
102  for(TString *trgit = triggers; trgit < triggers + sizeof(triggers)/sizeof(TString); ++trgit){
103  fHistos->CreateTH1(Form("hEventCount%sBeforeEventSelection", trgit->Data()), Form("Event count for trigger %s before event selection", trgit->Data()), 1, 0.5, 1.5);
104  fHistos->CreateTH1(Form("hEventCount%sBeforeOfflineTrigger", trgit->Data()), Form("Event count for trigger %s before offline selection", trgit->Data()), 1, 0.5, 1.5);
105  fHistos->CreateTH1(Form("hEventCount%sAfterOfflineTrigger", trgit->Data()), Form("Event count for trigger %s before offline selection", trgit->Data()), 1, 0.5, 1.5);
106  fHistos->CreateTH1(Form("hVertexTrigger%sBeforeEventSelection", trgit->Data()), Form("Vertex Distribution for trigger %s before event selection", trgit->Data()), 400, -40, 40);
107  fHistos->CreateTH1(Form("hVertexTrigger%sBeforeOfflineTrigger", trgit->Data()), Form("Vertex Distribution for trigger %s before offline trigger", trgit->Data()), 400, -40, 40);
108  fHistos->CreateTH1(Form("hVertexTrigger%sAfterOfflineTrigger", trgit->Data()), Form("Vertex Distribution for trigger %s after offline trigger", trgit->Data()), 400, -40, 40);
109  for(TString *statit = selectionStatus; statit < selectionStatus + sizeof(selectionStatus)/sizeof(TString); ++statit){
110  fHistos->CreateTH1(Form("hTrackPt%s%s", statit->Data(), trgit->Data()), Form("Pt spectrum of tracks in %s events of trigger %s", statit->Data(), trgit->Data()), ptbinning);
111  fHistos->CreateTH1(Form("hClusterEnergy%s%s", statit->Data(), trgit->Data()), Form("Cluster energy spectrum in %s events of trigger %s", statit->Data(), trgit->Data()), energybinning);
112  if(trgit->CompareTo("MB"))
113  fHistos->CreateTH1(Form("hPatchEnergy%s%s", statit->Data(), trgit->Data()), Form("Patch energy spectrum in %s events of trigger %s", statit->Data(), trgit->Data()), energybinning);
114  }
115  // Make histos for patch energy spectra for Min. Bias for different triggers
116  fHistos->CreateTH1(Form("hPatchEnergy%sMinBias", trgit->Data()), Form("Patch energy spectrum for %s patches found in MB events", trgit->Data()), energybinning);
117  }
118 
119 
120  fTrackContainer = new TObjArray;
121  fTrackContainer->SetOwner(kFALSE);
122 
123  PostData(1, fHistos->GetListOfHistograms());
124 }
125 
126 void AliAnalysisTaskEventSelectionRef::UserExec(Option_t *){
127  // Select event
128  if(!fGeometry){
129  fGeometry = AliEMCALGeometry::GetInstanceFromRunNumber(fInputEvent->GetRunNumber());
130  }
131  fTriggerPatchContainer = static_cast<TClonesArray *>(fInputEvent->FindListObject("EmcalTriggers"));
132  fClusterContainer = static_cast<TClonesArray *>(fInputEvent->FindListObject(fClusterContainerName.Data()));
133  TString triggerstring = fInputEvent->GetFiredTriggerClasses();
134  UInt_t selectionstatus = fInputHandler->IsEventSelected();
135  Bool_t isMinBias = selectionstatus & AliVEvent::kINT7,
136  isEJ1 = (selectionstatus & AliVEvent::kEMCEJE) && triggerstring.Contains("EJ1"),
137  isEJ2 = (selectionstatus & AliVEvent::kEMCEJE) && triggerstring.Contains("EJ2"),
138  isEG1 = (selectionstatus & AliVEvent::kEMCEGA) && triggerstring.Contains("EG1"),
139  isEG2 = (selectionstatus & AliVEvent::kEMCEGA) && triggerstring.Contains("EG2"),
140  isEMC7 = (selectionstatus & AliVEvent::kEMC7) && triggerstring.Contains("CEMC7");
141  if(!(isMinBias || isEMC7 || isEG1 || isEG2 || isEJ1 || isEJ2)) return;
142  const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
143  //if(!fInputEvent->IsPileupFromSPD(3, 0.8, 3., 2., 5.)) return; // reject pileup event
144  if(vtx->GetNContributors() < 1) return;
145  if(fInputEvent->IsA() == AliESDEvent::Class() && fAnalysisUtils->IsFirstEventInChunk(fInputEvent)) return;
146  bool isSelected = kTRUE;
147  if(!fAnalysisUtils->IsVertexSelected2013pA(fInputEvent)) isSelected = kFALSE; // Apply new vertex cut
148  if(fAnalysisUtils->IsPileUpEvent(fInputEvent)) isSelected = kFALSE; // Apply new vertex cut
149  // Apply vertex z cut
150  if(vtx->GetZ() < -10. || vtx->GetZ() > 10.) isSelected = kFALSE;
151 
152  // prepare tracks: Apply track selection and extrapolate to EMCAL surface
153  fTrackContainer->Clear();
154  AliVTrack *trk(nullptr); Double_t etaEMCAL(0.), phiEMCAL(0.); Int_t supermoduleID(-1);
155  for(int itrk = 0; itrk < fInputEvent->GetNumberOfTracks(); itrk++){
156  trk = static_cast<AliESDtrack *>(fInputEvent->GetTrack(itrk));
157  if(TMath::Abs(trk->Eta()) > 0.6) continue;
158  if(trk->IsA() == AliESDtrack::Class()){
159  AliESDtrack *origtrack = static_cast<AliESDtrack *>(trk);
160  if(!TrackSelectionESD(origtrack)) continue;
161  AliESDtrack copytrack(*origtrack);
162  AliEMCALRecoUtils::ExtrapolateTrackToEMCalSurface(&copytrack);
163  etaEMCAL = copytrack.GetTrackEtaOnEMCal();
164  phiEMCAL = copytrack.GetTrackPhiOnEMCal();
165  } else {
166  AliAODTrack *origtrack = static_cast<AliAODTrack *>(trk);
167  if(!TrackSelectionAOD(origtrack)) continue;
168  AliAODTrack copytrack(*origtrack);
169  AliEMCALRecoUtils::ExtrapolateTrackToEMCalSurface(&copytrack);
170  etaEMCAL = copytrack.GetTrackEtaOnEMCal();
171  phiEMCAL = copytrack.GetTrackPhiOnEMCal();
172  }
173  if(!fGeometry->SuperModuleNumberFromEtaPhi(etaEMCAL, phiEMCAL, supermoduleID)) continue;
174  if(supermoduleID < 0 || supermoduleID >= 10) continue;
175  fTrackContainer->Add(trk);
176  }
177 
178  // Fill Event counter and reference vertex distributions for the different trigger classes
179  if(isMinBias){
180  FillEventCounterHists("MB", vtx->GetZ(), isSelected, true);
181  }
182  if(isEMC7){
183  FillEventCounterHists("EMC7", vtx->GetZ(), isSelected, IsOfflineSelected(kCPREL0, fTriggerPatchContainer));
184  }
185  if(isEJ2){
186  FillEventCounterHists("EJ2", vtx->GetZ(), isSelected, IsOfflineSelected(kCPREJ2, fTriggerPatchContainer));
187  }
188  if(isEJ1){
189  FillEventCounterHists("EJ1", vtx->GetZ(), isSelected, IsOfflineSelected(kCPREJ1, fTriggerPatchContainer));
190  }
191  if(isEG2){
192  FillEventCounterHists("EG2", vtx->GetZ(), isSelected, IsOfflineSelected(kCPREG2, fTriggerPatchContainer));
193  }
194  if(isEG1){
195  FillEventCounterHists("EG1", vtx->GetZ(), isSelected, IsOfflineSelected(kCPREG1, fTriggerPatchContainer));
196  }
197 
198  PostData(1, fHistos->GetListOfHistograms());
199 }
200 
201 void AliAnalysisTaskEventSelectionRef::FillEventCounterHists(
202  const char *triggerclass,
203  double vtxz,
204  bool isSelected,
205  bool isOfflineSelected
206 )
207 {
208  // Fill reference distribution for the primary vertex before any z-cut
209  fHistos->FillTH1(Form("hVertexTrigger%sBeforeEventSelection", triggerclass), vtxz);
210  fHistos->FillTH1(Form("hEventCount%sBeforeEventSelection", triggerclass), 1.);
211  if(isSelected){
212  // Fill Event counter and reference vertex distributions after event selection
213  fHistos->FillTH1(Form("hEventCount%sBeforeOfflineTrigger", triggerclass), 1);
214  fHistos->FillTH1(Form("hVertexTrigger%sBeforeOfflineTrigger", triggerclass), vtxz);
215  if(isOfflineSelected){
216  fHistos->FillTH1(Form("hEventCount%sAfterOfflineTrigger", triggerclass), 1);
217  fHistos->FillTH1(Form("hVertexTrigger%sAfterOfflineTrigger", triggerclass), vtxz);
218  }
219 
220  // Now make some distributions of quantities in selected and rejected event
221  //... tracks
222  for(TIter trackIter = TIter(fTrackContainer).Begin(); trackIter != TIter::End(); ++trackIter){
223  ProcessTrack(triggerclass, static_cast<AliVTrack *>(*trackIter), isOfflineSelected);
224  }
225 
226  // ... clusters
227  if(fClusterContainer){
228  for(TIter clsit = TIter(fClusterContainer).Begin(); clsit != TIter::End(); ++clsit){
229  ProcessCluster(triggerclass, static_cast<AliVCluster *>(*clsit), isOfflineSelected);
230  }
231  }
232 
233  // ... patches
234  for(TIter patchiter = TIter(fTriggerPatchContainer).Begin(); patchiter != TIter::End(); ++patchiter){
235  ProcessOfflinePatch(triggerclass, static_cast<AliEmcalTriggerPatchInfo *>(*patchiter), isOfflineSelected);
236  }
237  }
238 }
239 
240 void AliAnalysisTaskEventSelectionRef::ProcessTrack(
241  const char *triggerclass,
242  const AliVTrack * track,
243  bool isOfflineSelected
244 )
245 {
246  fHistos->FillTH1(Form("hTrackPtAll%s", triggerclass), TMath::Abs(track->Pt()));
247  if(isOfflineSelected){
248  fHistos->FillTH1(Form("hTrackPtAccepted%s", triggerclass), TMath::Abs(track->Pt()));
249  } else {
250  fHistos->FillTH1(Form("hTrackPtRejected%s", triggerclass), TMath::Abs(track->Pt()));
251  }
252 }
253 
254 void AliAnalysisTaskEventSelectionRef::ProcessCluster(
255  const char *triggerclass,
256  const AliVCluster *clust,
257  bool isOfflineSelected
258 )
259 {
260  if(!clust->IsEMCAL()) return;
261  fHistos->FillTH1(Form("hClusterEnergyAll%s", triggerclass), clust->E());
262  if(isOfflineSelected){
263  fHistos->FillTH1(Form("hClusterEnergyAccepted%s", triggerclass), clust->E());
264  } else {
265  fHistos->FillTH1(Form("hClusterEnergyRejected%s", triggerclass), clust->E());
266  }
267 }
268 
269 void AliAnalysisTaskEventSelectionRef::ProcessOfflinePatch(
270  const char * triggerclass,
271  const AliEmcalTriggerPatchInfo * patch,
272  bool isOfflineSelected
273 )
274 {
275  bool isSingleShower = patch->IsGammaLowSimple();
276  bool isJetPatch = patch->IsJetLowSimple();
277  if(!(isSingleShower || isJetPatch)) return;
278  TString triggerstring(triggerclass);
279  if(!triggerstring.CompareTo("MB")){
280  if(isSingleShower){
281  fHistos->FillTH1("hPatchEnergyEMC7MinBias", patch->GetPatchE());
282  fHistos->FillTH1("hPatchEnergyEG1MinBias", patch->GetPatchE());
283  fHistos->FillTH1("hPatchEnergyEG2MinBias", patch->GetPatchE());
284  } else {
285  fHistos->FillTH1("hPatchEnergyEJ1MinBias", patch->GetPatchE());
286  fHistos->FillTH1("hPatchEnergyEJ2MinBias", patch->GetPatchE());
287  }
288  } else {
289  bool singleShowerTrigger = !triggerstring.CompareTo("EMC7") || !triggerstring.CompareTo("EG1") || !triggerstring.CompareTo("EG2");
290  if((isSingleShower && singleShowerTrigger) || (isJetPatch && !singleShowerTrigger)){
291  fHistos->FillTH1(Form("hPatchEnergyAll%s", triggerclass), patch->GetPatchE());
292  if(isOfflineSelected)
293  fHistos->FillTH1(Form("hPatchEnergyAccepted%s", triggerclass), patch->GetPatchE());
294  else
295  fHistos->FillTH1(Form("hPatchEnergyRejected%s", triggerclass), patch->GetPatchE());
296  }
297  }
298 }
299 
300 Bool_t AliAnalysisTaskEventSelectionRef::IsOfflineSelected(EmcalTriggerClass trgcls, const TClonesArray * const triggerpatches) const {
301  if(fOfflineEnergyThreshold[trgcls] < 0) return true;
302  bool isSingleShower = ((trgcls == kCPREL0) || (trgcls == kCPREG1) || (trgcls == kCPREG2));
303  int nfound = 0;
304  AliEmcalTriggerPatchInfo *patch = NULL;
305  for(TIter patchIter = TIter(triggerpatches).Begin(); patchIter != TIter::End(); ++patchIter){
306  patch = static_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
307  if(!patch->IsOfflineSimple()) continue;
308  if(isSingleShower){
309  if(!patch->IsGammaLowSimple()) continue;
310  } else {
311  if(!patch->IsJetLowSimple()) continue;
312  }
313  if(patch->GetPatchE() > fOfflineEnergyThreshold[trgcls]) nfound++;
314  }
315  return nfound > 0;
316 }
317 
322 void AliAnalysisTaskEventSelectionRef::CreateEnergyBinning(TArrayD& binning) const {
323  std::vector<double> mybinning;
324  std::map<double,double> definitions;
325  definitions.insert(std::pair<double, double>(1, 0.05));
326  definitions.insert(std::pair<double, double>(2, 0.1));
327  definitions.insert(std::pair<double, double>(4, 0.2));
328  definitions.insert(std::pair<double, double>(7, 0.5));
329  definitions.insert(std::pair<double, double>(16, 1));
330  definitions.insert(std::pair<double, double>(32, 2));
331  definitions.insert(std::pair<double, double>(40, 4));
332  definitions.insert(std::pair<double, double>(50, 5));
333  definitions.insert(std::pair<double, double>(100, 10));
334  definitions.insert(std::pair<double, double>(200, 20));
335  double currentval = 0.;
336  mybinning.push_back(currentval);
337  for(std::map<double,double>::iterator id = definitions.begin(); id != definitions.end(); ++id){
338  double limit = id->first, binwidth = id->second;
339  while(currentval < limit){
340  currentval += binwidth;
341  mybinning.push_back(currentval);
342  }
343  }
344  binning.Set(mybinning.size());
345  int ib = 0;
346  for(std::vector<double>::iterator it = mybinning.begin(); it != mybinning.end(); ++it)
347  binning[ib++] = *it;
348 }
349 
354 void AliAnalysisTaskEventSelectionRef::CreatePtBinning(TArrayD& binning) const {
355  std::vector<double> mybinning;
356  std::map<double,double> definitions;
357  definitions.insert(std::pair<double, double>(1, 0.05));
358  definitions.insert(std::pair<double, double>(2, 0.1));
359  definitions.insert(std::pair<double, double>(4, 0.2));
360  definitions.insert(std::pair<double, double>(7, 0.5));
361  definitions.insert(std::pair<double, double>(16, 1));
362  definitions.insert(std::pair<double, double>(36, 2));
363  definitions.insert(std::pair<double, double>(40, 4));
364  definitions.insert(std::pair<double, double>(50, 5));
365  definitions.insert(std::pair<double, double>(100, 10));
366  definitions.insert(std::pair<double, double>(200, 20));
367  double currentval = 0.;
368  mybinning.push_back(currentval);
369  for(std::map<double,double>::iterator id = definitions.begin(); id != definitions.end(); ++id){
370  double limit = id->first, binwidth = id->second;
371  while(currentval < limit){
372  currentval += binwidth;
373  mybinning.push_back(currentval);
374  }
375  }
376  binning.Set(mybinning.size());
377  int ib = 0;
378  for(std::vector<double>::iterator it = mybinning.begin(); it != mybinning.end(); ++it)
379  binning[ib++] = *it;
380 }
381 
387 Bool_t AliAnalysisTaskEventSelectionRef::TrackSelectionESD(AliESDtrack* track) {
388  return fTrackCuts->AcceptTrack(track);
389 }
390 
396 Bool_t AliAnalysisTaskEventSelectionRef::TrackSelectionAOD(AliAODTrack* track) {
397  if(!track->TestFilterBit(AliAODTrack::kTrkGlobal)) return false;
398  if(track->GetTPCNCrossedRows() < 120) return false;
399  return true;
400 }
401 
402 } /* namespace EMCalTriggerPtAnalysis */
Declarartion of class AliEMCalHistoContainer.
ClassImp(EMCalTriggerPtAnalysis::AliAnalysisTaskEventSelectionRef) namespace EMCalTriggerPtAnalysis