AliPhysics  vAN-20151012 (2287573)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
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 <map>
16 #include <vector>
17 
18 #include <TArrayD.h>
19 #include <TClonesArray.h>
20 #include <THashList.h>
21 #include <TString.h>
22 
23 #include "AliAnalysisUtils.h"
24 #include "AliEMCalHistoContainer.h"
26 #include "AliInputEventHandler.h"
27 #include "AliLog.h"
28 
30 
32 
33 namespace EMCalTriggerPtAnalysis {
34 
38 AliAnalysisTaskEmcalPatchesRef::AliAnalysisTaskEmcalPatchesRef() :
39  AliAnalysisTaskSE(),
40  fAnalysisUtil(NULL),
41  fHistos(NULL),
42  fTriggerStringFromPatches(kFALSE)
43 {
44  for(int itrg = 0; itrg < kEPRntrig; itrg++){
45  fOfflineEnergyThreshold[itrg] = -1;
46  }
47 }
48 
53 AliAnalysisTaskEmcalPatchesRef::AliAnalysisTaskEmcalPatchesRef(const char *name):
54  AliAnalysisTaskSE(name),
55  fAnalysisUtil(NULL),
56  fHistos(NULL),
57  fTriggerStringFromPatches(kFALSE)
58 {
59  for(int itrg = 0; itrg < kEPRntrig; itrg++){
60  fOfflineEnergyThreshold[itrg] = -1;
61  }
62  DefineOutput(1, TList::Class());
63 }
64 
68 AliAnalysisTaskEmcalPatchesRef::~AliAnalysisTaskEmcalPatchesRef() {
69 }
70 
76 void AliAnalysisTaskEmcalPatchesRef::UserCreateOutputObjects(){
77  fAnalysisUtil = new AliAnalysisUtils;
78 
79  TArrayD energybinning;
80  CreateEnergyBinning(energybinning);
81  fHistos = new AliEMCalHistoContainer("Ref");
82  TString triggers[15] = {"MB", "EMC7", "EJ1", "EJ2", "EG1", "EG2", "EG2excl", "EJ2excl", "MBexcl", "E1combined", "E1Jonly", "E1Gonly", "E2combined", "E2Jonly", "E2Gonly"};
83  TString patchtype[5] = {"EG1", "EG2", "EJ1", "EJ2", "EMC7"};
84  Double_t encuts[5] = {1., 2., 5., 10., 20.};
85  for(TString *trg = triggers; trg < triggers+15; trg++){
86  fHistos->CreateTH1(Form("hEventCount%s", trg->Data()), Form("Event count for trigger class %s", trg->Data()), 1, 0.5, 1.5);
87  for(int ipatch = 0; ipatch < 5; ipatch++){
88  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);
89  for(int ien = 0; ien < 5; ien++){
90  fHistos->CreateTH2(Form("h%sEtaPhi%dG%s", patchtype[ipatch].Data(), static_cast<int>(encuts[ien]), trg->Data()), Form("%s-patch #eta-#phi map for clusters with energy larger than %f GeV/c for trigger class %s", patchtype[ipatch].Data(), encuts[ien], trg->Data()), 100, -0.7, 0.7, 100, 1.4, 3.2);
91  }
92  }
93  }
94  PostData(1, fHistos->GetListOfHistograms());
95 
96 }
97 
102 void AliAnalysisTaskEmcalPatchesRef::UserExec(Option_t *){
103  TClonesArray *patches = dynamic_cast<TClonesArray *>(fInputEvent->FindListObject("EmcalTriggers"));
104  TString triggerstring = "";
105  if(fTriggerStringFromPatches){
106  triggerstring = GetFiredTriggerClassesFromPatches(patches);
107  } else {
108  triggerstring = fInputEvent->GetFiredTriggerClasses();
109  }
110  UInt_t selectionstatus = fInputHandler->IsEventSelected();
111  Bool_t isMinBias = selectionstatus & AliVEvent::kINT7,
112  isEMC7 = (selectionstatus & AliVEvent::kEMC7) && triggerstring.Contains("EMC7") && IsOfflineSelected(kEPREL0, patches),
113  isEJ1 = (selectionstatus & AliVEvent::kEMCEJE) && triggerstring.Contains("EJ1") && IsOfflineSelected(kEPREJ1, patches),
114  isEJ2 = (selectionstatus & AliVEvent::kEMCEJE) && triggerstring.Contains("EJ2") && IsOfflineSelected(kEPREJ2, patches),
115  isEG1 = (selectionstatus & AliVEvent::kEMCEGA) && triggerstring.Contains("EG1") && IsOfflineSelected(kEPREG1, patches),
116  isEG2 = (selectionstatus & AliVEvent::kEMCEGA) && triggerstring.Contains("EG2") && IsOfflineSelected(kEPREG2, patches);
117  if(!(isMinBias || isEMC7 || isEG1 || isEG2 || isEJ1 || isEJ2)) return;
118  const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
119  //if(!fInputEvent->IsPileupFromSPD(3, 0.8, 3., 2., 5.)) return; // reject pileup event
120  if(vtx->GetNContributors() < 1) return;
121  // Fill reference distribution for the primary vertex before any z-cut
122  if(!fAnalysisUtil->IsVertexSelected2013pA(fInputEvent)) return; // Apply new vertex cut
123  if(fAnalysisUtil->IsPileUpEvent(fInputEvent)) return; // Apply new vertex cut
124  // Apply vertex z cut
125  if(vtx->GetZ() < -10. || vtx->GetZ() > 10.) return;
126 
127  // Fill Event counter and reference vertex distributions for the different trigger classes
128  if(isMinBias){
129  fHistos->FillTH1("hEventCountMB", 1);
130  // Check for exclusive classes
131  if(!(isEG1 || isEG2 || isEJ1 || isEJ2)){
132  fHistos->FillTH1("hEventCountMBexcl", 1);
133  }
134  }
135  if(isEMC7){
136  fHistos->FillTH1("hEventCountEMC7", 1);
137  }
138  if(isEJ1){
139  fHistos->FillTH1("hEventCountEJ1", 1);
140  if(isEG1 || isEG2){
141  fHistos->FillTH1("hEventCountE1combined", 1);
142  } else {
143  fHistos->FillTH1("hEventCountE1Jonly", 1);
144  }
145 
146  }
147  if(isEJ2){
148  fHistos->FillTH1("hEventCountEJ2", 1);
149  // Check for exclusive classes
150  if(!isEJ1){
151  fHistos->FillTH1("hEventCountEJ2excl", 1);
152  }
153  if(isEG1 || isEG2){
154  fHistos->FillTH1("hEventCountE2combined", 1);
155  } else {
156  fHistos->FillTH1("hEventCountE2Jonly", 1);
157  }
158 
159  }
160  if(isEG1){
161  fHistos->FillTH1("hEventCountEG1", 1);
162  }
163  if(isEG2){
164  fHistos->FillTH1("hEventCountEG2", 1);
165  // Check for exclusive classes
166  if(!isEG1){
167  fHistos->FillTH1("hEventCountEG2excl", 1);
168  }
169  }
170 
171  if(!patches){
172  AliError("Trigger patch container not available");
173  return;
174  }
175 
176  Double_t vertexpos[3];
177  fInputEvent->GetPrimaryVertex()->GetXYZ(vertexpos);
178 
179  Double_t energy, eta, phi;
180  for(TIter patchIter = TIter(patches).Begin(); patchIter != TIter::End(); ++patchIter){
181  AliEmcalTriggerPatchInfo *patch = static_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
182  if(!patch->IsOfflineSimple()) continue;
183 
184  std::vector<TString> patchnames;
185  if(patch->IsJetHighSimple())
186  patchnames.push_back("EJ1");
187  if(patch->IsJetLowSimple())
188  patchnames.push_back("EJ2");
189  if(patch->IsGammaHighSimple())
190  patchnames.push_back("EG1");
191  if(patch->IsGammaLowSimple()){
192  patchnames.push_back("EG2");
193  patchnames.push_back("EMC7"); // Treat EG2 patch also as EMC7 patch (single shower);
194  }
195  if(!patchnames.size()){
196  // Undefined patch type - ignore
197  continue;
198  }
199 
200  TLorentzVector posvec;
201  energy = patch->GetPatchE();
202  eta = patch->GetEtaGeo();
203  phi = patch->GetPhiGeo();
204 
205  // fill histograms allEta
206  for(std::vector<TString>::iterator nameit = patchnames.begin(); nameit != patchnames.end(); ++nameit){
207  if(isMinBias){
208  FillPatchHistograms("MB", *nameit, energy, eta, phi);
209  // check for exclusive classes
210  if(!(isEG1 || isEG2 || isEJ1 || isEJ2)){
211  FillPatchHistograms("MBexcl", *nameit, energy, eta, phi);
212  }
213  }
214  if(isEMC7){
215  FillPatchHistograms("EMC7", *nameit, energy, eta, phi);
216  }
217  if(isEJ1){
218  FillPatchHistograms("EJ1", *nameit, energy, eta, phi);
219  if(isEG1 || isEG2) {
220  FillPatchHistograms("E1combined", *nameit, energy, eta, phi);
221  } else {
222  FillPatchHistograms("E1Jonly", *nameit, energy, eta, phi);
223  }
224  }
225  if(isEJ2){
226  FillPatchHistograms("EJ2", *nameit, energy, eta, phi);
227  // check for exclusive classes
228  if(!isEJ1){
229  FillPatchHistograms("EJ2excl", *nameit, energy, eta, phi);
230  }
231  if(isEG1 || isEG2){
232  FillPatchHistograms("E2combined", *nameit, energy, eta, phi);
233  } else {
234  FillPatchHistograms("E2Jonly", *nameit, energy, eta, phi);
235  }
236  }
237  if(isEG1){
238  FillPatchHistograms("EG1", *nameit, energy, eta, phi);
239  if(!(isEJ1 || isEJ2)){
240  FillPatchHistograms("E1Gonly", *nameit, energy, eta, phi);
241  }
242  }
243  if(isEG2){
244  FillPatchHistograms("EG2", *nameit, energy, eta, phi);
245  // check for exclusive classes
246  if(!isEG1){
247  FillPatchHistograms("EG2excl", *nameit, energy, eta, phi);
248  }
249  if(!(isEJ2 || isEJ1)){
250  FillPatchHistograms("E2Gonly", *nameit, energy, eta, phi);
251  }
252  }
253  }
254  }
255  PostData(1, fHistos->GetListOfHistograms());
256 }
257 
266 void AliAnalysisTaskEmcalPatchesRef::FillPatchHistograms(TString triggerclass, TString patchname, double energy, double eta, double phi){
267  fHistos->FillTH1(Form("h%sPatchEnergy%s", patchname.Data(), triggerclass.Data()), energy);
268  Double_t encuts[5] = {1., 2., 5., 10., 20.};
269  for(int ien = 0; ien < 5; ien++){
270  if(energy > encuts[ien]){
271  fHistos->FillTH2(Form("h%sEtaPhi%dG%s", patchname.Data(), static_cast<int>(encuts[ien]), triggerclass.Data()), eta, phi);
272  }
273  }
274 }
275 
280 void AliAnalysisTaskEmcalPatchesRef::CreateEnergyBinning(TArrayD& binning) const {
281  std::vector<double> mybinning;
282  std::map<double,double> definitions;
283  definitions.insert(std::pair<double, double>(1, 0.05));
284  definitions.insert(std::pair<double, double>(2, 0.1));
285  definitions.insert(std::pair<double, double>(4, 0.2));
286  definitions.insert(std::pair<double, double>(7, 0.5));
287  definitions.insert(std::pair<double, double>(16, 1));
288  definitions.insert(std::pair<double, double>(32, 2));
289  definitions.insert(std::pair<double, double>(40, 4));
290  definitions.insert(std::pair<double, double>(50, 5));
291  definitions.insert(std::pair<double, double>(100, 10));
292  definitions.insert(std::pair<double, double>(200, 20));
293  double currentval = 0.;
294  mybinning.push_back(currentval);
295  for(std::map<double,double>::iterator id = definitions.begin(); id != definitions.end(); ++id){
296  double limit = id->first, binwidth = id->second;
297  while(currentval < limit){
298  currentval += binwidth;
299  mybinning.push_back(currentval);
300  }
301  }
302  binning.Set(mybinning.size());
303  int ib = 0;
304  for(std::vector<double>::iterator it = mybinning.begin(); it != mybinning.end(); ++it)
305  binning[ib++] = *it;
306 }
307 
316 Bool_t AliAnalysisTaskEmcalPatchesRef::IsOfflineSelected(EmcalTriggerClass trgcls, const TClonesArray * const triggerpatches) const {
317  if(fOfflineEnergyThreshold[trgcls] < 0) return true;
318  bool isSingleShower = ((trgcls == kEPREL0) || (trgcls == kEPREG1) || (trgcls == kEPREG2));
319  int nfound = 0;
320  AliEmcalTriggerPatchInfo *patch = NULL;
321  for(TIter patchIter = TIter(triggerpatches).Begin(); patchIter != TIter::End(); ++patchIter){
322  patch = static_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
323  if(!patch->IsOfflineSimple()) continue;
324  if(isSingleShower){
325  if(!patch->IsGammaLowSimple()) continue;
326  } else {
327  if(!patch->IsJetLowSimple()) continue;
328  }
329  if(patch->GetPatchE() > fOfflineEnergyThreshold[trgcls]) nfound++;
330  }
331  return nfound > 0;
332 }
333 
334 
340 TString AliAnalysisTaskEmcalPatchesRef::GetFiredTriggerClassesFromPatches(const TClonesArray* triggerpatches) const {
341  TString triggerstring = "";
342  if(!triggerpatches) return triggerstring;
343  Int_t nEJ1 = 0, nEJ2 = 0, nEG1 = 0, nEG2 = 0;
344  double minADC_EJ1 = 260.,
345  minADC_EJ2 = 127.,
346  minADC_EG1 = 140.,
347  minADC_EG2 = 89.;
348  for(TIter patchIter = TIter(triggerpatches).Begin(); patchIter != TIter::End(); ++patchIter){
349  AliEmcalTriggerPatchInfo *patch = dynamic_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
350  if(!patch->IsOfflineSimple()) continue;
351  if(patch->IsJetHighSimple() && patch->GetADCOfflineAmp() > minADC_EJ1) nEJ1++;
352  if(patch->IsJetLowSimple() && patch->GetADCOfflineAmp() > minADC_EJ2) nEJ2++;
353  if(patch->IsGammaHighSimple() && patch->GetADCOfflineAmp() > minADC_EG1) nEG1++;
354  if(patch->IsGammaLowSimple() && patch->GetADCOfflineAmp() > minADC_EG2) nEG2++;
355  }
356  if(nEJ1) triggerstring += "EJ1";
357  if(nEJ2){
358  if(triggerstring.Length()) triggerstring += ",";
359  triggerstring += "EJ2";
360  }
361  if(nEG1){
362  if(triggerstring.Length()) triggerstring += ",";
363  triggerstring += "EG1";
364  }
365  if(nEG2){
366  if(triggerstring.Length()) triggerstring += ",";
367  triggerstring += "EG2";
368  }
369  return triggerstring;
370 }
371 
372 
373 } /* namespace EMCalTriggerPtAnalysis */
ClassImp(EMCalTriggerPtAnalysis::AliAnalysisTaskEmcalPatchesRef) namespace EMCalTriggerPtAnalysis
Main data structure storing all relevant information of EMCAL/DCAL trigger patches.
Declarartion of class AliEMCalHistoContainer.
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
Class to make array of trigger patch objects in AOD/ESD events.