AliPhysics  b7e5564 (b7e5564)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskEmcalPatchesRef.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 <iostream>
16 #include <map>
17 #include <vector>
18 
19 #include <TArrayD.h>
20 #include <TClonesArray.h>
21 #include <THistManager.h>
22 #include <TLinearBinning.h>
23 #include <THashList.h>
24 #include <TString.h>
25 
26 #include "AliAnalysisUtils.h"
27 #include "AliESDEvent.h"
28 #include "AliEMCALTriggerPatchInfo.h"
30 #include "AliInputEventHandler.h"
31 #include "AliLog.h"
32 #include "AliMultSelection.h"
33 #include "AliMultEstimator.h"
34 
36 
38 
39 namespace EMCalTriggerPtAnalysis {
40 
44 AliAnalysisTaskEmcalPatchesRef::AliAnalysisTaskEmcalPatchesRef() :
45  AliAnalysisTaskSE(),
46  fAnalysisUtil(NULL),
47  fTriggerSelection(NULL),
48  fHistos(NULL),
49  fRequestAnalysisUtil(kTRUE),
50  fTriggerStringFromPatches(kFALSE),
51  fCentralityRange(-999., 999.),
52  fVertexRange(-999., 999.),
53  fRequestCentrality(false)
54 {
55 }
56 
61 AliAnalysisTaskEmcalPatchesRef::AliAnalysisTaskEmcalPatchesRef(const char *name):
62  AliAnalysisTaskSE(name),
63  fAnalysisUtil(NULL),
64  fTriggerSelection(NULL),
65  fHistos(NULL),
66  fRequestAnalysisUtil(kTRUE),
67  fTriggerStringFromPatches(kFALSE),
68  fCentralityRange(-999., 999.),
69  fVertexRange(-999., 999.),
70  fRequestCentrality(false)
71 {
72  DefineOutput(1, TList::Class());
73 }
74 
78 AliAnalysisTaskEmcalPatchesRef::~AliAnalysisTaskEmcalPatchesRef() {
79 }
80 
86 void AliAnalysisTaskEmcalPatchesRef::UserCreateOutputObjects(){
87  AliInfoStream() << "Creating histograms for task " << GetName() << std::endl;
88  fAnalysisUtil = new AliAnalysisUtils;
89 
90  EnergyBinning energybinning;
91  TLinearBinning etabinning(100, -0.7, 0.7);
92  fHistos = new THistManager("Ref");
93  TString triggers[18] = {"MB", "EMC7", "DMC7",
94  "EJ1", "EJ2", "EG1", "EG2", "DJ1", "DJ2", "DG1", "DG2",
95  "MBexcl", "EMC7excl", "DMC7excl", "EG2excl", "EJ2excl", "DG2excl", "DJ2excl"
96  };
97  TString patchtype[10] = {"EG1", "EG2", "EJ1", "EJ2", "EMC7", "DG1", "DG2", "DJ1", "DJ2", "DMC7"};
98  Double_t encuts[5] = {1., 2., 5., 10., 20.};
99  for(TString *trg = triggers; trg < triggers+18; trg++){
100  fHistos->CreateTH1(Form("hEventCount%s", trg->Data()), Form("Event count for trigger class %s", trg->Data()), 1, 0.5, 1.5);
101  fHistos->CreateTH1(Form("hEventCentrality%s", trg->Data()), Form("Event centrality for trigger class %s", trg->Data()), 103, -2., 101.);
102  for(int ipatch = 0; ipatch < 10; ipatch++){
103  fHistos->CreateTH1(Form("h%sPatchEnergy%s", patchtype[ipatch].Data(), trg->Data()), Form("%s-patch energy for trigger class %s", patchtype[ipatch].Data(), trg->Data()), energybinning);
104  fHistos->CreateTH1(Form("h%sPatchET%s", patchtype[ipatch].Data(), trg->Data()), Form("%s-patch transverse energy for trigger class %s", patchtype[ipatch].Data(), trg->Data()), energybinning);
105  fHistos->CreateTH2(Form("h%sPatchEnergyEta%s", patchtype[ipatch].Data(), trg->Data()), Form("%s-patch energy for trigger class %s", patchtype[ipatch].Data(), trg->Data()), energybinning, etabinning);
106  fHistos->CreateTH2(Form("h%sPatchETEta%s", patchtype[ipatch].Data(), trg->Data()), Form("%s-patch transverse energy for trigger class %s", patchtype[ipatch].Data(), trg->Data()), energybinning, etabinning);
107  for(int ien = 0; ien < 5; ien++){
108  fHistos->CreateTH2(Form("h%sEtaPhi%dG%s", patchtype[ipatch].Data(), static_cast<int>(encuts[ien]), trg->Data()), Form("%s-patch #eta-#phi map for patches with energy larger than %f GeV/c for trigger class %s", patchtype[ipatch].Data(), encuts[ien], trg->Data()), 100, -0.7, 0.7, 200, 0, TMath::TwoPi());
109  fHistos->CreateTH2(Form("h%sColRow%dG%s", patchtype[ipatch].Data(), static_cast<int>(encuts[ien]), trg->Data()), Form("%s-patch col-row map for patches with energy larger than %f GeV/c for trigger class %s", patchtype[ipatch].Data(), encuts[ien], trg->Data()), 48, -0.5, 47.5, 104, -0.5, 103.5);
110  }
111  }
112  }
113  PostData(1, fHistos->GetListOfHistograms());
114  AliDebugStream(1) << "Histograms done" << std::endl;
115 }
116 
121 void AliAnalysisTaskEmcalPatchesRef::UserExec(Option_t *){
122  AliDebugStream(1) << GetName() << ": Start function" << std::endl;
123  TClonesArray *patches = dynamic_cast<TClonesArray *>(fInputEvent->FindListObject("EmcalTriggers"));
124  TString triggerstring = "";
125  if(fTriggerStringFromPatches){
126  triggerstring = GetFiredTriggerClassesFromPatches(patches);
127  } else {
128  triggerstring = fInputEvent->GetFiredTriggerClasses();
129  }
130  UInt_t selectionstatus = fInputHandler->IsEventSelected();
131  Bool_t isMinBias = selectionstatus & AliVEvent::kINT7,
132  isEMC7 = (selectionstatus & AliVEvent::kEMC7) && triggerstring.Contains("EMC7"),
133  isEJ1 = (selectionstatus & AliVEvent::kEMCEJE) && triggerstring.Contains("EJ1"),
134  isEJ2 = (selectionstatus & AliVEvent::kEMCEJE) && triggerstring.Contains("EJ2"),
135  isEG1 = (selectionstatus & AliVEvent::kEMCEGA) && triggerstring.Contains("EG1"),
136  isEG2 = (selectionstatus & AliVEvent::kEMCEGA) && triggerstring.Contains("EG2"),
137  isDMC7 = (selectionstatus & AliVEvent::kEMC7) && triggerstring.Contains("DMC7"),
138  isDJ1 = (selectionstatus & AliVEvent::kEMCEJE) && triggerstring.Contains("DJ1"),
139  isDJ2 = (selectionstatus & AliVEvent::kEMCEJE) && triggerstring.Contains("DJ2"),
140  isDG1 = (selectionstatus & AliVEvent::kEMCEGA) && triggerstring.Contains("DG1"),
141  isDG2 = (selectionstatus & AliVEvent::kEMCEGA) && triggerstring.Contains("DG2");
142  if(patches && fTriggerSelection){
143  isEMC7 &= fTriggerSelection->IsOfflineSelected(AliEmcalTriggerOfflineSelection::kTrgEL0, patches);
144  isEJ1 &= fTriggerSelection->IsOfflineSelected(AliEmcalTriggerOfflineSelection::kTrgEJ1, patches);
145  isEJ2 &= fTriggerSelection->IsOfflineSelected(AliEmcalTriggerOfflineSelection::kTrgEJ2, patches);
146  isEG1 &= fTriggerSelection->IsOfflineSelected(AliEmcalTriggerOfflineSelection::kTrgEG1, patches);
147  isEG2 &= fTriggerSelection->IsOfflineSelected(AliEmcalTriggerOfflineSelection::kTrgEG2, patches);
148  isDMC7 &= fTriggerSelection->IsOfflineSelected(AliEmcalTriggerOfflineSelection::kTrgDL0, patches);
149  isDJ1 &= fTriggerSelection->IsOfflineSelected(AliEmcalTriggerOfflineSelection::kTrgDJ1, patches);
150  isDJ2 &= fTriggerSelection->IsOfflineSelected(AliEmcalTriggerOfflineSelection::kTrgDJ2, patches);
151  isDG1 &= fTriggerSelection->IsOfflineSelected(AliEmcalTriggerOfflineSelection::kTrgDG1, patches);
152  isDG2 &= fTriggerSelection->IsOfflineSelected(AliEmcalTriggerOfflineSelection::kTrgDG2, patches);
153 
154  }
155  if(!(isMinBias || isEMC7 || isEG1 || isEG2 || isEJ1 || isEJ2 || isDMC7 || isDG1 || isDG2 || isDJ1 || isDJ2)){
156  AliDebugStream(1) << GetName() << ": Reject trigger" << std::endl;
157  return;
158  }
159  AliDebugStream(1) << "Event selected" << std::endl;
160  AliMultSelection *mult = dynamic_cast<AliMultSelection *>(InputEvent()->FindListObject("MultSelection"));
161  // In case a centrality estimator is used, event selection,
162  // otherwise ignore event selection from multiplicity task
163  if(fRequestCentrality){
164  if(mult && !mult->IsEventSelected()) return;
165  }
166  double centrality = mult ? mult->GetEstimator("V0M")->GetPercentile() : -1;
167  AliDebugStream(1) << GetName() << ": Centrality " << centrality << std::endl;
168  if(!fCentralityRange.IsInRange(centrality)){
169  AliDebugStream(1) << GetName() << ": reject centrality: " << centrality << std::endl;
170  return;
171  } else {
172  AliDebugStream(1) << GetName() << ": select centrality " << centrality << std::endl;
173  }
174  const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
175  if(!vtx) vtx = fInputEvent->GetPrimaryVertexSPD();
176  //if(!fInputEvent->IsPileupFromSPD(3, 0.8, 3., 2., 5.)) return; // reject pileup event
177  if(vtx->GetNContributors() < 1){
178  AliDebug(1, Form("%s: Reject contributor\n", GetName()));
179  return;
180  }
181  // Fill reference distribution for the primary vertex before any z-cut
182  if(fRequestAnalysisUtil){
183  AliDebugStream(1) << GetName() << ": Reject analysis util" << std::endl;
184  if(fInputEvent->IsA() == AliESDEvent::Class() && fAnalysisUtil->IsFirstEventInChunk(fInputEvent)) return;
185  if(!fAnalysisUtil->IsVertexSelected2013pA(fInputEvent)) return; // Apply new vertex cut
186  if(fAnalysisUtil->IsPileUpEvent(fInputEvent)) return; // Apply new vertex cut
187  }
188  // Apply vertex z cut
189  if(!fVertexRange.IsInRange(vtx->GetZ())){
190  AliDebugStream(1) << GetName() << ": Reject Z(" << vtx->GetZ() << ")" << std::endl;
191  return;
192  }
193  AliDebugStream(1) << GetName() << ": Event Selected" << std::endl;
194 
195  // Fill Event counter and reference vertex distributions for the different trigger classes
196  if(isMinBias){
197  fHistos->FillTH1("hEventCountMB", 1);
198  fHistos->FillTH1("hEventCentralityMB", centrality);
199  // Check for exclusive classes
200  if(!(isEMC7 || isDMC7 || isEJ1 || isEJ2 || isEG1 || isEG2 || isDJ1 || isDJ2 || isDG1 || isDG2)){
201  fHistos->FillTH1("hEventCountMBexcl", 1);
202  fHistos->FillTH1("hEventCentralityMBexcl", centrality);
203  }
204  }
205 
206  // L0 triggers
207  if(isEMC7){
208  fHistos->FillTH1("hEventCountEMC7", 1);
209  fHistos->FillTH1("hEventCentralityEMC7", centrality);
210  if(!(isEJ1 || isEJ2 || isEG1 || isEG2)){
211  fHistos->FillTH1("hEventCountEMC7excl", 1);
212  fHistos->FillTH1("hEventCentralityEMC7excl", centrality);
213  }
214  }
215  if(isDMC7){
216  fHistos->FillTH1("hEventCountDMC7", 1);
217  fHistos->FillTH1("hEventCentralityDMC7", centrality);
218  if(!(isDJ1 || isDJ2 || isDG1 || isDG2)){
219  fHistos->FillTH1("hEventCountDMC7excl", 1);
220  fHistos->FillTH1("hEventCentralityDMC7excl", centrality);
221  }
222  }
223 
224  // L1 jet triggers
225  if(isEJ1){
226  fHistos->FillTH1("hEventCountEJ1", 1);
227  fHistos->FillTH1("hEventCentralityEJ1", centrality);
228  }
229  if(isDJ1){
230  fHistos->FillTH1("hEventCountEJ1", 1);
231  fHistos->FillTH1("hEventCentralityDJ1", centrality);
232  }
233 
234  if(isEJ2){
235  fHistos->FillTH1("hEventCountEJ2", 1);
236  fHistos->FillTH1("hEventCentralityEJ2", centrality);
237  // Check for exclusive classes
238  if(!isEJ1){
239  fHistos->FillTH1("hEventCountEJ2excl", 1);
240  fHistos->FillTH1("hEventCentralityEJ2excl", centrality);
241  }
242  }
243  if(isDJ2){
244  fHistos->FillTH1("hEventCountDJ2", 1);
245  fHistos->FillTH1("hEventCentralityDJ2", centrality);
246  // Check for exclusive classes
247  if(!isDJ1){
248  fHistos->FillTH1("hEventCountDJ2excl", 1);
249  fHistos->FillTH1("hEventCentralityDJ2excl", centrality);
250  }
251  }
252 
253  // L1 gamma triggers
254  if(isEG1){
255  fHistos->FillTH1("hEventCountEG1", 1);
256  fHistos->FillTH1("hEventCentralityEG1", centrality);
257  }
258  if(isDG1){
259  fHistos->FillTH1("hEventCountDG1", 1);
260  fHistos->FillTH1("hEventCentralityDG1", centrality);
261  }
262  if(isEG2){
263  fHistos->FillTH1("hEventCountEG2", 1);
264  fHistos->FillTH1("hEventCentralityEG2", centrality);
265  // Check for exclusive classes
266  if(!isEG1){
267  fHistos->FillTH1("hEventCountEG2excl", 1);
268  fHistos->FillTH1("hEventCentralityEG2excl", centrality);
269  }
270  }
271  if(isDG2){
272  fHistos->FillTH1("hEventCountDG2", 1);
273  fHistos->FillTH1("hEventCentralityDG2", centrality);
274  // Check for exclusive classes
275  if(!isEG1){
276  fHistos->FillTH1("hEventCountDG2excl", 1);
277  fHistos->FillTH1("hEventCentralityDG2excl", centrality);
278  }
279  }
280 
281  if(!patches){
282  AliErrorStream() << GetName() << ": Trigger patch container not available" << std::endl;
283  return;
284  }
285 
286  AliDebugStream(1) << GetName() << ": Number of trigger patches " << patches->GetEntries() << std::endl;
287 
288  Double_t vertexpos[3];
289  fInputEvent->GetPrimaryVertex()->GetXYZ(vertexpos);
290 
291  Double_t energy, eta, phi, et;
292  Int_t col, row;
293  for(TIter patchIter = TIter(patches).Begin(); patchIter != TIter::End(); ++patchIter){
294  if(!IsOfflineSimplePatch(*patchIter)) continue;
295  AliEMCALTriggerPatchInfo *patch = static_cast<AliEMCALTriggerPatchInfo *>(*patchIter);
296 
297  bool isDCAL = SelectDCALPatch(*patchIter),
298  isSingleShower = SelectSingleShowerPatch(*patchIter),
299  isJetPatch = SelectJetPatch(*patchIter);
300 
301  std::vector<TString> patchnames;
302  if(isJetPatch){
303  if(isDCAL){
304  patchnames.push_back("DJ1");
305  patchnames.push_back("DJ2");
306  } else {
307  patchnames.push_back("EJ1");
308  patchnames.push_back("EJ2");
309  }
310  }
311  if(isSingleShower){
312  if(isDCAL){
313  patchnames.push_back("DMC7");
314  patchnames.push_back("DG1");
315  patchnames.push_back("DG2");
316  } else {
317  patchnames.push_back("EMC7");
318  patchnames.push_back("EG1");
319  patchnames.push_back("EG2");
320  }
321  }
322  if(!patchnames.size()){
323  // Undefined patch type - ignore
324  continue;
325  }
326 
327  TLorentzVector posvec;
328  energy = patch->GetPatchE();
329  eta = patch->GetEtaGeo();
330  phi = patch->GetPhiGeo();
331  col = patch->GetColStart();
332  row = patch->GetRowStart();
333  et = patch->GetLorentzVectorCenterGeo().Et();
334 
335  // fill histograms allEta
336  for(std::vector<TString>::iterator nameit = patchnames.begin(); nameit != patchnames.end(); ++nameit){
337  if(isMinBias){
338  FillPatchHistograms("MB", *nameit, energy, et, eta, phi, col, row);
339  // check for exclusive classes
340  if(!(isEMC7 || isDMC7 || isEJ1 || isEJ2 || isEG1 || isEG2 || isDJ1 || isDJ2 || isDG1 || isDG2)){
341  FillPatchHistograms("MBexcl", *nameit, energy, et, eta, phi, col, row);
342  }
343  }
344 
345  if(isEMC7){
346  FillPatchHistograms("EMC7", *nameit, energy, et, eta, phi, col, row);
347  if(!(isEJ1 || isEJ2 || isEG1 || isEG2)){
348  FillPatchHistograms("EMC7excl", *nameit, energy, et, eta, phi, col, row);
349  }
350  }
351  if(isDMC7){
352  FillPatchHistograms("DMC7", *nameit, energy, et, eta, phi, col, row);
353  if(!(isDJ1 || isDJ2 || isDG1 || isDG2)){
354  FillPatchHistograms("DMC7excl", *nameit, energy, et, eta, phi, col, row);
355  }
356  }
357 
358  if(isEJ1){
359  FillPatchHistograms("EJ1", *nameit, energy, et, eta, phi, col, row);
360  }
361  if(isDJ1){
362  FillPatchHistograms("DJ1", *nameit, energy, et, eta, phi, col, row);
363  }
364 
365 
366  if(isEJ2){
367  FillPatchHistograms("EJ2", *nameit, energy, et, eta, phi, col, row);
368  // check for exclusive classes
369  if(!isEJ1){
370  FillPatchHistograms("EJ2excl", *nameit, energy, et, eta, phi, col, row);
371  }
372  }
373  if(isDJ2){
374  FillPatchHistograms("DJ2", *nameit, energy, et, eta, phi, col, row);
375  // check for exclusive classes
376  if(!isDJ1){
377  FillPatchHistograms("DJ2excl", *nameit, energy, et, eta, phi, col, row);
378  }
379  }
380 
381  if(isEG1){
382  FillPatchHistograms("EG1", *nameit, energy, et, eta, phi, col, row);
383  }
384  if(isDG1){
385  FillPatchHistograms("DG1", *nameit, energy, et, eta, phi, col, row);
386  }
387 
388  if(isEG2){
389  FillPatchHistograms("EG2", *nameit, energy, et, eta, phi, col, row);
390  // check for exclusive classes
391  if(!isEG1){
392  FillPatchHistograms("EG2excl", *nameit, energy, et, eta, phi, col, row);
393  }
394  }
395  if(isDG2){
396  FillPatchHistograms("DG2", *nameit, energy, et, eta, phi, col, row);
397  // check for exclusive classes
398  if(!isDG1){
399  FillPatchHistograms("DG2excl", *nameit, energy, et, eta, phi, col, row);
400  }
401  }
402  }
403  }
404  PostData(1, fHistos->GetListOfHistograms());
405 }
406 
415 void AliAnalysisTaskEmcalPatchesRef::FillPatchHistograms(TString triggerclass, TString patchname, double energy, double transverseenergy, double eta, double phi, int col, int row){
416  fHistos->FillTH1(Form("h%sPatchEnergy%s", patchname.Data(), triggerclass.Data()), energy);
417  fHistos->FillTH1(Form("h%sPatchET%s", patchname.Data(), triggerclass.Data()), transverseenergy);
418  fHistos->FillTH2(Form("h%sPatchEnergyEta%s", patchname.Data(), triggerclass.Data()), energy, eta);
419  fHistos->FillTH2(Form("h%sPatchETEta%s", patchname.Data(), triggerclass.Data()), transverseenergy, eta);
420  Double_t encuts[5] = {1., 2., 5., 10., 20.};
421  for(int ien = 0; ien < 5; ien++){
422  if(energy > encuts[ien]){
423  fHistos->FillTH2(Form("h%sEtaPhi%dG%s", patchname.Data(), static_cast<int>(encuts[ien]), triggerclass.Data()), eta, phi);
424  fHistos->FillTH2(Form("h%sColRow%dG%s", patchname.Data(), static_cast<int>(encuts[ien]), triggerclass.Data()), col, row);
425  }
426  }
427 }
428 
434 TString AliAnalysisTaskEmcalPatchesRef::GetFiredTriggerClassesFromPatches(const TClonesArray* triggerpatches) const {
435  TString triggerstring = "";
436  if(!triggerpatches) return triggerstring;
437  Int_t nEJ1 = 0, nEJ2 = 0, nEG1 = 0, nEG2 = 0, nDJ1 = 0, nDJ2 = 0, nDG1 = 0, nDG2 = 0;
438  double minADC_J1 = 260.,
439  minADC_J2 = 127.,
440  minADC_G1 = 140.,
441  minADC_G2 = 89.;
442  for(TIter patchIter = TIter(triggerpatches).Begin(); patchIter != TIter::End(); ++patchIter){
443  AliEMCALTriggerPatchInfo *patch = dynamic_cast<AliEMCALTriggerPatchInfo *>(*patchIter);
444  if(!patch->IsOfflineSimple()) continue;
445  if(patch->IsJetHighSimple() && patch->GetADCOfflineAmp() > minADC_J1){
446  if(patch->IsDCalPHOS()) nDJ1++;
447  else nEJ1++;
448  }
449  if(patch->IsJetLowSimple() && patch->GetADCOfflineAmp() > minADC_J2){
450  if(patch->IsDCalPHOS()) nDJ2++;
451  else nEJ2++;
452  }
453  if(patch->IsGammaHighSimple() && patch->GetADCOfflineAmp() > minADC_G1){
454  if(patch->IsDCalPHOS()) nDG1++;
455  else nEG1++;
456  }
457  if(patch->IsGammaLowSimple() && patch->GetADCOfflineAmp() > minADC_G2){
458  if(patch->IsDCalPHOS()) nDG2++;
459  else nEG2++;
460  }
461  }
462  if(nEJ1) triggerstring += "EJ1";
463  if(nEJ2){
464  if(triggerstring.Length()) triggerstring += ",";
465  triggerstring += "EJ2";
466  }
467  if(nEG1){
468  if(triggerstring.Length()) triggerstring += ",";
469  triggerstring += "EG1";
470  }
471  if(nEG2){
472  if(triggerstring.Length()) triggerstring += ",";
473  triggerstring += "EG2";
474  }
475  if(nDJ1){
476  if(triggerstring.Length()) triggerstring += ",";
477  triggerstring += "DJ1";
478  }
479  if(nEJ2){
480  if(triggerstring.Length()) triggerstring += ",";
481  triggerstring += "DJ2";
482  }
483  if(nDG1){
484  if(triggerstring.Length()) triggerstring += ",";
485  triggerstring += "DG1";
486  }
487  if(nDG2){
488  if(triggerstring.Length()) triggerstring += ",";
489  triggerstring += "DG2";
490  }
491  return triggerstring;
492 }
493 
494 void AliAnalysisTaskEmcalPatchesRef::GetPatchBoundaries(TObject *o, Double_t *boundaries) const {
495  AliEMCALTriggerPatchInfo *patch = dynamic_cast<AliEMCALTriggerPatchInfo *>(o);
496  boundaries[0] = patch->GetEtaMin();
497  boundaries[1] = patch->GetEtaMax();
498  boundaries[2] = patch->GetPhiMin();
499  boundaries[3] = patch->GetPhiMax();
500 }
501 
502 bool AliAnalysisTaskEmcalPatchesRef::IsOfflineSimplePatch(TObject *o) const {
503  AliEMCALTriggerPatchInfo *patch = dynamic_cast<AliEMCALTriggerPatchInfo *>(o);
504  return patch->IsOfflineSimple();
505 }
506 
507 bool AliAnalysisTaskEmcalPatchesRef::SelectDCALPatch(TObject *o) const {
508  AliEMCALTriggerPatchInfo *patch = dynamic_cast<AliEMCALTriggerPatchInfo *>(o);
509  return patch->GetRowStart() >= 64;
510 }
511 
512 bool AliAnalysisTaskEmcalPatchesRef::SelectSingleShowerPatch(TObject *o) const{
513  AliEMCALTriggerPatchInfo *patch = dynamic_cast<AliEMCALTriggerPatchInfo *>(o);
514  if(!patch->IsOfflineSimple()) return false;
515  return patch->IsGammaLowSimple();
516 }
517 
518 bool AliAnalysisTaskEmcalPatchesRef::SelectJetPatch(TObject *o) const{
519  AliEMCALTriggerPatchInfo *patch = dynamic_cast<AliEMCALTriggerPatchInfo *>(o);
520  if(!patch->IsOfflineSimple()) return false;
521  return patch->IsJetLowSimple();
522 }
523 
524 double AliAnalysisTaskEmcalPatchesRef::GetPatchEnergy(TObject *o) const {
525  double energy = 0.;
526  AliEMCALTriggerPatchInfo *patch = dynamic_cast<AliEMCALTriggerPatchInfo *>(o);
527  energy = patch->GetPatchE();
528  return energy;
529 }
530 
531 AliAnalysisTaskEmcalPatchesRef::EnergyBinning::EnergyBinning():
533 {
534  this->SetMinimum(0.);
535  this->AddStep(1., 0.05);
536  this->AddStep(2., 0.1);
537  this->AddStep(4, 0.2);
538  this->AddStep(7, 0.5);
539  this->AddStep(16, 1);
540  this->AddStep(32, 2);
541  this->AddStep(40, 4);
542  this->AddStep(50, 5);
543  this->AddStep(100, 10);
544  this->AddStep(200, 20);
545 }
546 
547 
548 } /* namespace EMCalTriggerPtAnalysis */
Class creating a linear binning, used in the histogram manager.
ClassImp(EMCalTriggerPtAnalysis::AliAnalysisTaskEmcalPatchesRef) namespace EMCalTriggerPtAnalysis
centrality
Helper class creating user defined custom binning.
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
energy
Container class for histograms for the high- charged particle analysis.
Definition: THistManager.h:43