AliPhysics  5dd2c10 (5dd2c10)
 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 }
45 
50 AliAnalysisTaskEmcalPatchesRef::AliAnalysisTaskEmcalPatchesRef(const char *name):
51  AliAnalysisTaskSE(name),
52  fAnalysisUtil(NULL),
53  fHistos(NULL),
54  fTriggerStringFromPatches(kFALSE)
55 {
56  DefineOutput(1, TList::Class());
57 }
58 
62 AliAnalysisTaskEmcalPatchesRef::~AliAnalysisTaskEmcalPatchesRef() {
63 }
64 
70 void AliAnalysisTaskEmcalPatchesRef::UserCreateOutputObjects(){
71  fAnalysisUtil = new AliAnalysisUtils;
72 
73  TArrayD energybinning;
74  CreateEnergyBinning(energybinning);
75  fHistos = new AliEMCalHistoContainer("Ref");
76  TString triggers[14] = {"MB", "EJ1", "EJ2", "EG1", "EG2", "EG2excl", "EJ2excl", "MBexcl", "E1combined", "E1Jonly", "E1Gonly", "E2combined", "E2Jonly", "E2Gonly"};
77  TString patchtype[4] = {"EG1", "EG2", "EJ1", "EJ2"};
78  Double_t encuts[5] = {1., 2., 5., 10., 20.};
79  for(TString *trg = triggers; trg < triggers+14; trg++){
80  fHistos->CreateTH1(Form("hEventCount%s", trg->Data()), Form("Event count for trigger class %s", trg->Data()), 1, 0.5, 1.5);
81  for(int ipatch = 0; ipatch < 4; ipatch++){
82  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);
83  for(int ien = 0; ien < 5; ien++){
84  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);
85  }
86  }
87  }
88  PostData(1, fHistos->GetListOfHistograms());
89 
90 }
91 
96 void AliAnalysisTaskEmcalPatchesRef::UserExec(Option_t *){
97  TClonesArray *patches = dynamic_cast<TClonesArray *>(fInputEvent->FindListObject("EmcalTriggers"));
98  TString triggerstring = "";
99  if(fTriggerStringFromPatches){
100  triggerstring = GetFiredTriggerClassesFromPatches(patches);
101  } else {
102  triggerstring = fInputEvent->GetFiredTriggerClasses();
103  }
104  UInt_t selectionstatus = fInputHandler->IsEventSelected();
105  Bool_t isMinBias = selectionstatus & AliVEvent::kINT7,
106  isEJ1 = (selectionstatus & AliVEvent::kEMCEJE) && triggerstring.Contains("EJ1"),
107  isEJ2 = (selectionstatus & AliVEvent::kEMCEJE) && triggerstring.Contains("EJ2"),
108  isEG1 = (selectionstatus & AliVEvent::kEMCEGA) && triggerstring.Contains("EG1"),
109  isEG2 = (selectionstatus & AliVEvent::kEMCEGA) && triggerstring.Contains("EG2");
110  if(!(isMinBias || isEG1 || isEG2 || isEJ1 || isEJ2)) return;
111  const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
112  //if(!fInputEvent->IsPileupFromSPD(3, 0.8, 3., 2., 5.)) return; // reject pileup event
113  if(vtx->GetNContributors() < 1) return;
114  // Fill reference distribution for the primary vertex before any z-cut
115  if(!fAnalysisUtil->IsVertexSelected2013pA(fInputEvent)) return; // Apply new vertex cut
116  if(fAnalysisUtil->IsPileUpEvent(fInputEvent)) return; // Apply new vertex cut
117  // Apply vertex z cut
118  if(vtx->GetZ() < -10. || vtx->GetZ() > 10.) return;
119 
120  // Fill Event counter and reference vertex distributions for the different trigger classes
121  if(isMinBias){
122  fHistos->FillTH1("hEventCountMB", 1);
123  // Check for exclusive classes
124  if(!(isEG1 || isEG2 || isEJ1 || isEJ2)){
125  fHistos->FillTH1("hEventCountMBexcl", 1);
126  }
127  }
128  if(isEJ1){
129  fHistos->FillTH1("hEventCountEJ1", 1);
130  if(isEG1 || isEG2){
131  fHistos->FillTH1("hEventCountE1combined", 1);
132  } else {
133  fHistos->FillTH1("hEventCountE1Jonly", 1);
134  }
135 
136  }
137  if(isEJ2){
138  fHistos->FillTH1("hEventCountEJ2", 1);
139  // Check for exclusive classes
140  if(!isEJ1){
141  fHistos->FillTH1("hEventCountEJ2excl", 1);
142  }
143  if(isEG1 || isEG2){
144  fHistos->FillTH1("hEventCountE2combined", 1);
145  } else {
146  fHistos->FillTH1("hEventCountE2Jonly", 1);
147  }
148 
149  }
150  if(isEG1){
151  fHistos->FillTH1("hEventCountEG1", 1);
152  }
153  if(isEG2){
154  fHistos->FillTH1("hEventCountEG2", 1);
155  // Check for exclusive classes
156  if(!isEG1){
157  fHistos->FillTH1("hEventCountEG2excl", 1);
158  }
159  }
160 
161  if(!patches){
162  AliError("Trigger patch container not available");
163  return;
164  }
165 
166  Double_t vertexpos[3];
167  fInputEvent->GetPrimaryVertex()->GetXYZ(vertexpos);
168 
169  Double_t energy, eta, phi;
170  for(TIter patchIter = TIter(patches).Begin(); patchIter != TIter::End(); ++patchIter){
171  AliEmcalTriggerPatchInfo *patch = static_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
172  if(!patch->IsOfflineSimple()) continue;
173 
174  std::vector<TString> patchnames;
175  if(patch->IsJetHighSimple())
176  patchnames.push_back("EJ1");
177  if(patch->IsJetLowSimple())
178  patchnames.push_back("EJ2");
179  if(patch->IsGammaHighSimple())
180  patchnames.push_back("EG1");
181  if(patch->IsGammaLowSimple())
182  patchnames.push_back("EG2");
183  if(!patchnames.size()){
184  // Undefined patch type - ignore
185  continue;
186  }
187 
188  TLorentzVector posvec;
189  energy = patch->GetPatchE();
190  eta = patch->GetEtaGeo();
191  phi = patch->GetPhiGeo();
192 
193  // fill histograms allEta
194  for(std::vector<TString>::iterator nameit = patchnames.begin(); nameit != patchnames.end(); ++nameit){
195  if(isMinBias){
196  FillPatchHistograms("MB", *nameit, energy, eta, phi);
197  // check for exclusive classes
198  if(!(isEG1 || isEG2 || isEJ1 || isEJ2)){
199  FillPatchHistograms("MBexcl", *nameit, energy, eta, phi);
200  }
201  }
202  if(isEJ1){
203  FillPatchHistograms("EJ1", *nameit, energy, eta, phi);
204  if(isEG1 || isEG2) {
205  FillPatchHistograms("E1combined", *nameit, energy, eta, phi);
206  } else {
207  FillPatchHistograms("E1Jonly", *nameit, energy, eta, phi);
208  }
209  }
210  if(isEJ2){
211  FillPatchHistograms("EJ2", *nameit, energy, eta, phi);
212  // check for exclusive classes
213  if(!isEJ1){
214  FillPatchHistograms("EJ2excl", *nameit, energy, eta, phi);
215  }
216  if(isEG1 || isEG2){
217  FillPatchHistograms("E2combined", *nameit, energy, eta, phi);
218  } else {
219  FillPatchHistograms("E2Jonly", *nameit, energy, eta, phi);
220  }
221  }
222  if(isEG1){
223  FillPatchHistograms("EG1", *nameit, energy, eta, phi);
224  if(!(isEJ1 || isEJ2)){
225  FillPatchHistograms("E1Gonly", *nameit, energy, eta, phi);
226  }
227  }
228  if(isEG2){
229  FillPatchHistograms("EG2", *nameit, energy, eta, phi);
230  // check for exclusive classes
231  if(!isEG1){
232  FillPatchHistograms("EG2excl", *nameit, energy, eta, phi);
233  }
234  if(!(isEJ2 || isEJ1)){
235  FillPatchHistograms("E2Gonly", *nameit, energy, eta, phi);
236  }
237  }
238  }
239  }
240  PostData(1, fHistos->GetListOfHistograms());
241 }
242 
251 void AliAnalysisTaskEmcalPatchesRef::FillPatchHistograms(TString triggerclass, TString patchname, double energy, double eta, double phi){
252  fHistos->FillTH1(Form("h%sPatchEnergy%s", patchname.Data(), triggerclass.Data()), energy);
253  Double_t encuts[5] = {1., 2., 5., 10., 20.};
254  for(int ien = 0; ien < 5; ien++){
255  if(energy > encuts[ien]){
256  fHistos->FillTH2(Form("h%sEtaPhi%dG%s", patchname.Data(), static_cast<int>(encuts[ien]), triggerclass.Data()), eta, phi);
257  }
258  }
259 }
260 
265 void AliAnalysisTaskEmcalPatchesRef::CreateEnergyBinning(TArrayD& binning) const {
266  std::vector<double> mybinning;
267  std::map<double,double> definitions;
268  definitions.insert(std::pair<double, double>(1, 0.05));
269  definitions.insert(std::pair<double, double>(2, 0.1));
270  definitions.insert(std::pair<double, double>(4, 0.2));
271  definitions.insert(std::pair<double, double>(7, 0.5));
272  definitions.insert(std::pair<double, double>(16, 1));
273  definitions.insert(std::pair<double, double>(32, 2));
274  definitions.insert(std::pair<double, double>(40, 4));
275  definitions.insert(std::pair<double, double>(50, 5));
276  definitions.insert(std::pair<double, double>(100, 10));
277  definitions.insert(std::pair<double, double>(200, 20));
278  double currentval = 0.;
279  mybinning.push_back(currentval);
280  for(std::map<double,double>::iterator id = definitions.begin(); id != definitions.end(); ++id){
281  double limit = id->first, binwidth = id->second;
282  while(currentval < limit){
283  currentval += binwidth;
284  mybinning.push_back(currentval);
285  }
286  }
287  binning.Set(mybinning.size());
288  int ib = 0;
289  for(std::vector<double>::iterator it = mybinning.begin(); it != mybinning.end(); ++it)
290  binning[ib++] = *it;
291 }
292 
298 TString AliAnalysisTaskEmcalPatchesRef::GetFiredTriggerClassesFromPatches(const TClonesArray* triggerpatches) const {
299  TString triggerstring = "";
300  if(!triggerpatches) return triggerstring;
301  Int_t nEJ1 = 0, nEJ2 = 0, nEG1 = 0, nEG2 = 0;
302  double minADC_EJ1 = 260.,
303  minADC_EJ2 = 127.,
304  minADC_EG1 = 140.,
305  minADC_EG2 = 89.;
306  for(TIter patchIter = TIter(triggerpatches).Begin(); patchIter != TIter::End(); ++patchIter){
307  AliEmcalTriggerPatchInfo *patch = dynamic_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
308  if(!patch->IsOfflineSimple()) continue;
309  if(patch->IsJetHighSimple() && patch->GetADCOfflineAmp() > minADC_EJ1) nEJ1++;
310  if(patch->IsJetLowSimple() && patch->GetADCOfflineAmp() > minADC_EJ2) nEJ2++;
311  if(patch->IsGammaHighSimple() && patch->GetADCOfflineAmp() > minADC_EG1) nEG1++;
312  if(patch->IsGammaLowSimple() && patch->GetADCOfflineAmp() > minADC_EG2) nEG2++;
313  }
314  if(nEJ1) triggerstring += "EJ1";
315  if(nEJ2){
316  if(triggerstring.Length()) triggerstring += ",";
317  triggerstring += "EJ2";
318  }
319  if(nEG1){
320  if(triggerstring.Length()) triggerstring += ",";
321  triggerstring += "EG1";
322  }
323  if(nEG2){
324  if(triggerstring.Length()) triggerstring += ",";
325  triggerstring += "EG2";
326  }
327  return triggerstring;
328 }
329 
330 
331 } /* 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.