AliPhysics  29d4213 (29d4213)
 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"
25 #include "AliEmcalTriggerPatchInfoAP.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, etabinning;
80  CreateEnergyBinning(energybinning);
81  CreateLinearBinning(etabinning, 100, -0.7, 0.7);
82  fHistos = new AliEMCalHistoContainer("Ref");
83  TString triggers[15] = {"MB", "EMC7", "EJ1", "EJ2", "EG1", "EG2", "EG2excl", "EJ2excl", "MBexcl", "E1combined", "E1Jonly", "E1Gonly", "E2combined", "E2Jonly", "E2Gonly"};
84  TString patchtype[5] = {"EG1", "EG2", "EJ1", "EJ2", "EMC7"};
85  Double_t encuts[5] = {1., 2., 5., 10., 20.};
86  for(TString *trg = triggers; trg < triggers+15; trg++){
87  fHistos->CreateTH1(Form("hEventCount%s", trg->Data()), Form("Event count for trigger class %s", trg->Data()), 1, 0.5, 1.5);
88  for(int ipatch = 0; ipatch < 5; ipatch++){
89  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);
90  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);
91  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);
92  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);
93  for(int ien = 0; ien < 5; ien++){
94  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);
95  }
96  }
97  }
98  PostData(1, fHistos->GetListOfHistograms());
99 
100 }
101 
106 void AliAnalysisTaskEmcalPatchesRef::UserExec(Option_t *){
107  TClonesArray *patches = dynamic_cast<TClonesArray *>(fInputEvent->FindListObject("EmcalTriggers"));
108  TString triggerstring = "";
109  if(fTriggerStringFromPatches){
110  triggerstring = GetFiredTriggerClassesFromPatches(patches);
111  } else {
112  triggerstring = fInputEvent->GetFiredTriggerClasses();
113  }
114  UInt_t selectionstatus = fInputHandler->IsEventSelected();
115  Bool_t isMinBias = selectionstatus & AliVEvent::kINT7,
116  isEMC7 = (selectionstatus & AliVEvent::kEMC7) && triggerstring.Contains("EMC7") && IsOfflineSelected(kEPREL0, patches),
117  isEJ1 = (selectionstatus & AliVEvent::kEMCEJE) && triggerstring.Contains("EJ1") && IsOfflineSelected(kEPREJ1, patches),
118  isEJ2 = (selectionstatus & AliVEvent::kEMCEJE) && triggerstring.Contains("EJ2") && IsOfflineSelected(kEPREJ2, patches),
119  isEG1 = (selectionstatus & AliVEvent::kEMCEGA) && triggerstring.Contains("EG1") && IsOfflineSelected(kEPREG1, patches),
120  isEG2 = (selectionstatus & AliVEvent::kEMCEGA) && triggerstring.Contains("EG2") && IsOfflineSelected(kEPREG2, patches);
121  if(!(isMinBias || isEMC7 || isEG1 || isEG2 || isEJ1 || isEJ2)) return;
122  const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
123  //if(!fInputEvent->IsPileupFromSPD(3, 0.8, 3., 2., 5.)) return; // reject pileup event
124  if(vtx->GetNContributors() < 1) return;
125  // Fill reference distribution for the primary vertex before any z-cut
126  if(!fAnalysisUtil->IsVertexSelected2013pA(fInputEvent)) return; // Apply new vertex cut
127  if(fAnalysisUtil->IsPileUpEvent(fInputEvent)) return; // Apply new vertex cut
128  // Apply vertex z cut
129  if(vtx->GetZ() < -10. || vtx->GetZ() > 10.) return;
130 
131  // Fill Event counter and reference vertex distributions for the different trigger classes
132  if(isMinBias){
133  fHistos->FillTH1("hEventCountMB", 1);
134  // Check for exclusive classes
135  if(!(isEG1 || isEG2 || isEJ1 || isEJ2)){
136  fHistos->FillTH1("hEventCountMBexcl", 1);
137  }
138  }
139  if(isEMC7){
140  fHistos->FillTH1("hEventCountEMC7", 1);
141  }
142  if(isEJ1){
143  fHistos->FillTH1("hEventCountEJ1", 1);
144  if(isEG1 || isEG2){
145  fHistos->FillTH1("hEventCountE1combined", 1);
146  } else {
147  fHistos->FillTH1("hEventCountE1Jonly", 1);
148  }
149 
150  }
151  if(isEJ2){
152  fHistos->FillTH1("hEventCountEJ2", 1);
153  // Check for exclusive classes
154  if(!isEJ1){
155  fHistos->FillTH1("hEventCountEJ2excl", 1);
156  }
157  if(isEG1 || isEG2){
158  fHistos->FillTH1("hEventCountE2combined", 1);
159  } else {
160  fHistos->FillTH1("hEventCountE2Jonly", 1);
161  }
162 
163  }
164  if(isEG1){
165  fHistos->FillTH1("hEventCountEG1", 1);
166  }
167  if(isEG2){
168  fHistos->FillTH1("hEventCountEG2", 1);
169  // Check for exclusive classes
170  if(!isEG1){
171  fHistos->FillTH1("hEventCountEG2excl", 1);
172  }
173  }
174 
175  if(!patches){
176  AliError("Trigger patch container not available");
177  return;
178  }
179 
180  Double_t vertexpos[3];
181  fInputEvent->GetPrimaryVertex()->GetXYZ(vertexpos);
182 
183  Double_t energy, eta, phi, et;
184  for(TIter patchIter = TIter(patches).Begin(); patchIter != TIter::End(); ++patchIter){
185  AliEmcalTriggerPatchInfo *patch = static_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
186  if(!patch->IsOfflineSimple()) continue;
187 
188  std::vector<TString> patchnames;
189  if(patch->IsJetHighSimple())
190  patchnames.push_back("EJ1");
191  if(patch->IsJetLowSimple())
192  patchnames.push_back("EJ2");
193  if(patch->IsGammaHighSimple())
194  patchnames.push_back("EG1");
195  if(patch->IsGammaLowSimple()){
196  patchnames.push_back("EG2");
197  patchnames.push_back("EMC7"); // Treat EG2 patch also as EMC7 patch (single shower);
198  }
199  if(!patchnames.size()){
200  // Undefined patch type - ignore
201  continue;
202  }
203 
204  TLorentzVector posvec;
205  energy = patch->GetPatchE();
206  eta = patch->GetEtaGeo();
207  phi = patch->GetPhiGeo();
208  et = patch->GetLorentzVectorCenterGeo().Et();
209 
210  // fill histograms allEta
211  for(std::vector<TString>::iterator nameit = patchnames.begin(); nameit != patchnames.end(); ++nameit){
212  if(isMinBias){
213  FillPatchHistograms("MB", *nameit, energy, et, eta, phi);
214  // check for exclusive classes
215  if(!(isEG1 || isEG2 || isEJ1 || isEJ2)){
216  FillPatchHistograms("MBexcl", *nameit, energy, et, eta, phi);
217  }
218  }
219  if(isEMC7){
220  FillPatchHistograms("EMC7", *nameit, energy, et, eta, phi);
221  }
222  if(isEJ1){
223  FillPatchHistograms("EJ1", *nameit, energy, et, eta, phi);
224  if(isEG1 || isEG2) {
225  FillPatchHistograms("E1combined", *nameit, energy, et, eta, phi);
226  } else {
227  FillPatchHistograms("E1Jonly", *nameit, energy, et, eta, phi);
228  }
229  }
230  if(isEJ2){
231  FillPatchHistograms("EJ2", *nameit, energy, et, eta, phi);
232  // check for exclusive classes
233  if(!isEJ1){
234  FillPatchHistograms("EJ2excl", *nameit, energy, et, eta, phi);
235  }
236  if(isEG1 || isEG2){
237  FillPatchHistograms("E2combined", *nameit, energy, et, eta, phi);
238  } else {
239  FillPatchHistograms("E2Jonly", *nameit, energy, et, eta, phi);
240  }
241  }
242  if(isEG1){
243  FillPatchHistograms("EG1", *nameit, energy, et, eta, phi);
244  if(!(isEJ1 || isEJ2)){
245  FillPatchHistograms("E1Gonly", *nameit, energy, et, eta, phi);
246  }
247  }
248  if(isEG2){
249  FillPatchHistograms("EG2", *nameit, energy, et, eta, phi);
250  // check for exclusive classes
251  if(!isEG1){
252  FillPatchHistograms("EG2excl", *nameit, energy, et, eta, phi);
253  }
254  if(!(isEJ2 || isEJ1)){
255  FillPatchHistograms("E2Gonly", *nameit, energy, et, eta, phi);
256  }
257  }
258  }
259  }
260  PostData(1, fHistos->GetListOfHistograms());
261 }
262 
271 void AliAnalysisTaskEmcalPatchesRef::FillPatchHistograms(TString triggerclass, TString patchname, double energy, double transverseenergy, double eta, double phi){
272  fHistos->FillTH1(Form("h%sPatchEnergy%s", patchname.Data(), triggerclass.Data()), energy);
273  fHistos->FillTH1(Form("h%sPatchET%s", patchname.Data(), triggerclass.Data()), transverseenergy);
274  fHistos->FillTH2(Form("h%sPatchEnergyEta%s", patchname.Data(), triggerclass.Data()), eta, energy);
275  fHistos->FillTH2(Form("h%sPatchETEta%s", patchname.Data(), triggerclass.Data()), eta, transverseenergy);
276  Double_t encuts[5] = {1., 2., 5., 10., 20.};
277  for(int ien = 0; ien < 5; ien++){
278  if(energy > encuts[ien]){
279  fHistos->FillTH2(Form("h%sEtaPhi%dG%s", patchname.Data(), static_cast<int>(encuts[ien]), triggerclass.Data()), eta, phi);
280  }
281  }
282 }
283 
288 void AliAnalysisTaskEmcalPatchesRef::CreateEnergyBinning(TArrayD& binning) const {
289  std::vector<double> mybinning;
290  std::map<double,double> definitions;
291  definitions.insert(std::pair<double, double>(1, 0.05));
292  definitions.insert(std::pair<double, double>(2, 0.1));
293  definitions.insert(std::pair<double, double>(4, 0.2));
294  definitions.insert(std::pair<double, double>(7, 0.5));
295  definitions.insert(std::pair<double, double>(16, 1));
296  definitions.insert(std::pair<double, double>(32, 2));
297  definitions.insert(std::pair<double, double>(40, 4));
298  definitions.insert(std::pair<double, double>(50, 5));
299  definitions.insert(std::pair<double, double>(100, 10));
300  definitions.insert(std::pair<double, double>(200, 20));
301  double currentval = 0.;
302  mybinning.push_back(currentval);
303  for(std::map<double,double>::iterator id = definitions.begin(); id != definitions.end(); ++id){
304  double limit = id->first, binwidth = id->second;
305  while(currentval < limit){
306  currentval += binwidth;
307  mybinning.push_back(currentval);
308  }
309  }
310  binning.Set(mybinning.size());
311  int ib = 0;
312  for(std::vector<double>::iterator it = mybinning.begin(); it != mybinning.end(); ++it)
313  binning[ib++] = *it;
314 }
315 
323 void AliAnalysisTaskEmcalPatchesRef::CreateLinearBinning(TArrayD& binning, int nbins, double min, double max) const {
324  double binwidth = (max-min)/static_cast<double>(nbins);
325  binning.Set(nbins+1);
326  binning[0] = min;
327  double currentlimit = min + binwidth;
328  for(int ibin = 0; ibin < nbins; ibin++){
329  binning[ibin+1] = currentlimit;
330  currentlimit += binwidth;
331  }
332 }
333 
342 Bool_t AliAnalysisTaskEmcalPatchesRef::IsOfflineSelected(EmcalTriggerClass trgcls, const TClonesArray * const triggerpatches) const {
343  if(fOfflineEnergyThreshold[trgcls] < 0) return true;
344  bool isSingleShower = ((trgcls == kEPREL0) || (trgcls == kEPREG1) || (trgcls == kEPREG2));
345  int nfound = 0;
346  AliEmcalTriggerPatchInfo *patch = NULL;
347  for(TIter patchIter = TIter(triggerpatches).Begin(); patchIter != TIter::End(); ++patchIter){
348  patch = static_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
349  if(!patch->IsOfflineSimple()) continue;
350  if(isSingleShower){
351  if(!patch->IsGammaLowSimple()) continue;
352  } else {
353  if(!patch->IsJetLowSimple()) continue;
354  }
355  if(patch->GetPatchE() > fOfflineEnergyThreshold[trgcls]) nfound++;
356  }
357  return nfound > 0;
358 }
359 
360 
366 TString AliAnalysisTaskEmcalPatchesRef::GetFiredTriggerClassesFromPatches(const TClonesArray* triggerpatches) const {
367  TString triggerstring = "";
368  if(!triggerpatches) return triggerstring;
369  Int_t nEJ1 = 0, nEJ2 = 0, nEG1 = 0, nEG2 = 0;
370  double minADC_EJ1 = 260.,
371  minADC_EJ2 = 127.,
372  minADC_EG1 = 140.,
373  minADC_EG2 = 89.;
374  for(TIter patchIter = TIter(triggerpatches).Begin(); patchIter != TIter::End(); ++patchIter){
375  AliEmcalTriggerPatchInfo *patch = dynamic_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
376  if(!patch->IsOfflineSimple()) continue;
377  if(patch->IsJetHighSimple() && patch->GetADCOfflineAmp() > minADC_EJ1) nEJ1++;
378  if(patch->IsJetLowSimple() && patch->GetADCOfflineAmp() > minADC_EJ2) nEJ2++;
379  if(patch->IsGammaHighSimple() && patch->GetADCOfflineAmp() > minADC_EG1) nEG1++;
380  if(patch->IsGammaLowSimple() && patch->GetADCOfflineAmp() > minADC_EG2) nEG2++;
381  }
382  if(nEJ1) triggerstring += "EJ1";
383  if(nEJ2){
384  if(triggerstring.Length()) triggerstring += ",";
385  triggerstring += "EJ2";
386  }
387  if(nEG1){
388  if(triggerstring.Length()) triggerstring += ",";
389  triggerstring += "EG1";
390  }
391  if(nEG2){
392  if(triggerstring.Length()) triggerstring += ",";
393  triggerstring += "EG2";
394  }
395  return triggerstring;
396 }
397 
398 
399 } /* namespace EMCalTriggerPtAnalysis */
ClassImp(EMCalTriggerPtAnalysis::AliAnalysisTaskEmcalPatchesRef) namespace EMCalTriggerPtAnalysis
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
const Int_t nbins