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