AliPhysics  vAN-20150924 (e816f45)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliAnalysisTaskEmcalClustersRef.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 <bitset>
16 #include <iostream>
17 #include <map>
18 #include <sstream>
19 #include <vector>
20 
21 #include <TArrayD.h>
22 #include <TClonesArray.h>
23 #include <THashList.h>
24 #include <TLorentzVector.h>
25 #include <TString.h>
26 
27 #include "AliAnalysisUtils.h"
28 #include "AliEMCALGeometry.h"
30 #include "AliEMCalHistoContainer.h"
31 #include "AliInputEventHandler.h"
32 #include "AliLog.h"
33 #include "AliVCluster.h"
34 #include "AliVVertex.h"
35 
37 
39 
40 namespace EMCalTriggerPtAnalysis {
41 
45 AliAnalysisTaskEmcalClustersRef::AliAnalysisTaskEmcalClustersRef() :
46  AliAnalysisTaskSE(),
47  fAnalysisUtil(NULL),
48  fHistos(NULL),
49  fGeometry(NULL),
50  fClusterContainer(""),
51  fTriggerStringFromPatches(kFALSE)
52 {
53 
54 }
55 
60 AliAnalysisTaskEmcalClustersRef::AliAnalysisTaskEmcalClustersRef(const char *name) :
61  AliAnalysisTaskSE(name),
62  fAnalysisUtil(),
63  fHistos(NULL),
64  fGeometry(NULL),
65  fClusterContainer(""),
66  fTriggerStringFromPatches(kFALSE)
67 {
68  DefineOutput(1, TList::Class());
69 }
70 
74 AliAnalysisTaskEmcalClustersRef::~AliAnalysisTaskEmcalClustersRef() {
75 }
76 
80 void AliAnalysisTaskEmcalClustersRef::UserCreateOutputObjects(){
81  fAnalysisUtil = new AliAnalysisUtils;
82 
83  TArrayD energybinning;
84  CreateEnergyBinning(energybinning);
85  TArrayD smbinning(14);
86  int ilimit = 0;
87  for(double borders = -0.5; borders < 13.5; borders +=1) smbinning[ilimit++] = borders;
88  fHistos = new AliEMCalHistoContainer("Ref");
89  TString triggers[17] = {
90  "MB", "EMC7",
91  "EJ1", "EJ2", "EG1", "EG2",
92  "EMC7excl","EG1excl", "EG2excl", "EJ1excl", "EJ2excl",
93  "E1combined", "E1Jonly", "E1Gonly", "E2combined", "E2Jonly", "E2Gonly"
94  };
95  Double_t encuts[5] = {1., 2., 5., 10., 20.};
96  for(TString *trg = triggers; trg < triggers + sizeof(triggers)/sizeof(TString); trg++){
97  fHistos->CreateTH1(Form("hEventCount%s", trg->Data()), Form("Event count for trigger class %s", trg->Data()), 1, 0.5, 1.5);
98  fHistos->CreateTH1(Form("hClusterEnergy%s", trg->Data()), Form("Cluster energy for trigger class %s", trg->Data()), energybinning);
99  fHistos->CreateTH1(Form("hClusterEnergyFired%s", trg->Data()), Form("Cluster energy for trigger class %s, firing the trigger", trg->Data()), energybinning);
100  fHistos->CreateTH2(Form("hClusterEnergySM%s", trg->Data()), Form("Cluster energy versus supermodule for trigger class %s", trg->Data()), smbinning, energybinning);
101  fHistos->CreateTH2(Form("hClusterEnergyFiredSM%s", trg->Data()), Form("Cluster energy versus supermodule for trigger class %s, firing the trigger", trg->Data()), smbinning, energybinning);
102  for(int ien = 0; ien < 5; ien++){
103  fHistos->CreateTH2(Form("hEtaPhi%dG%s", static_cast<int>(encuts[ien]), trg->Data()), Form("cluster #eta-#phi map for clusters with energy larger than %f GeV/c for trigger class %s", encuts[ien], trg->Data()), 100, -0.7, 0.7, 100, 1.4, 3.2);
104  fHistos->CreateTH2(Form("hEtaPhiFired%dG%s", static_cast<int>(encuts[ien]), trg->Data()), Form("cluster #eta-#phi map for clusters fired the trigger with energy larger than %f GeV/c for trigger class %s", encuts[ien], trg->Data()), 100, -0.7, 0.7, 100, 1.4, 3.2);
105  }
106  }
107  PostData(1, fHistos->GetListOfHistograms());
108 }
109 
110 
115 void AliAnalysisTaskEmcalClustersRef::UserExec(Option_t *){
116  if(!fGeometry){
117  fGeometry = AliEMCALGeometry::GetInstanceFromRunNumber(InputEvent()->GetRunNumber());
118  }
119  TString triggerstring = "";
120  TClonesArray *triggerpatches = dynamic_cast<TClonesArray *>(fInputEvent->FindListObject("EmcalTriggers"));
121 
122  if(fTriggerStringFromPatches){
123  triggerstring = GetFiredTriggerClassesFromPatches(triggerpatches);
124  } else {
125  triggerstring = fInputEvent->GetFiredTriggerClasses();
126  }
127 
128  UInt_t selectionstatus = fInputHandler->IsEventSelected();
129  std::stringstream triggerdebug;
130  triggerdebug << "Offline bits: " << std::bitset<sizeof(UInt_t) * 8>(selectionstatus);
131  AliDebug(2, triggerdebug.str().c_str());
132  Bool_t isMinBias = selectionstatus & AliVEvent::kINT7,
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  isEMC7 = (selectionstatus & AliVEvent::kEMC7) && triggerstring.Contains("CEMC7");
138  if(!(isMinBias || isEMC7 || isEG1 || isEG2 || isEJ1 || isEJ2)) return;
139  const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
140  //if(!fInputEvent->IsPileupFromSPD(3, 0.8, 3., 2., 5.)) return; // reject pileup event
141  if(vtx->GetNContributors() < 1) return;
142  // Fill reference distribution for the primary vertex before any z-cut
143  if(!fAnalysisUtil->IsVertexSelected2013pA(fInputEvent)) return; // Apply new vertex cut
144  if(fAnalysisUtil->IsPileUpEvent(fInputEvent)) return; // Apply new vertex cut
145  // Apply vertex z cut
146  if(vtx->GetZ() < -10. || vtx->GetZ() > 10.) return;
147 
148  // Fill Event counter and reference vertex distributions for the different trigger classes
149  if(isMinBias){
150  fHistos->FillTH1("hEventCountMB", 1);
151  }
152  if(isEMC7){
153  fHistos->FillTH1("hEventCountEMC7", 1);
154  if(!isMinBias){
155  fHistos->FillTH1("hEventCountEMC7excl", 1);
156  }
157  }
158  if(isEJ2){
159  fHistos->FillTH1("hEventCountEJ2", 1);
160  // Check for exclusive classes
161  if(!(isMinBias)){
162  fHistos->FillTH1("hEventCountEJ2excl", 1);
163  }
164  if(isEG1 || isEG2){
165  fHistos->FillTH1("hEventCountE2combined", 1);
166  } else {
167  fHistos->FillTH1("hEventCountE2Jonly", 1);
168  }
169 
170  }
171  if(isEJ1){
172  fHistos->FillTH1("hEventCountEJ1", 1);
173  // Check for exclusive classes
174  if(!(isMinBias || isEJ2)){
175  fHistos->FillTH1("hEventCountEJ1excl", 1);
176  }
177  if(isEG1 || isEG2){
178  fHistos->FillTH1("hEventCountE1combined", 1);
179  } else {
180  fHistos->FillTH1("hEventCountE1Jonly", 1);
181  }
182  }
183  if(isEG2){
184  fHistos->FillTH1("hEventCountEG2", 1);
185  // Check for exclusive classes
186  if(!(isMinBias)){
187  fHistos->FillTH1("hEventCountEG2excl", 1);
188  }
189  if(!(isEJ1 || isEJ2)){
190  fHistos->FillTH1("hEventCountE2Gonly", 1);
191  }
192  }
193  if(isEG1){
194  fHistos->FillTH1("hEventCountEG1", 1);
195  // Check for exclusive classes
196  if(!(isMinBias || isEG2)){
197  fHistos->FillTH1("hEventCountEG1excl", 1);
198  }
199  if(!(isEJ1 || isEJ2)){
200  fHistos->FillTH1("hEventCountE1Gonly", 1);
201  }
202  }
203 
204  TClonesArray *clusterArray = dynamic_cast<TClonesArray *>(fInputEvent->FindListObject(fClusterContainer.Data()));
205  if(!clusterArray){
206  AliError("Cluster array not found");
207  return;
208  }
209 
210  Double_t vertexpos[3];
211  fInputEvent->GetPrimaryVertex()->GetXYZ(vertexpos);
212 
213  Double_t energy, eta, phi;
214  for(TIter clustIter = TIter(clusterArray).Begin(); clustIter != TIter::End(); ++clustIter){
215  AliVCluster *clust = static_cast<AliVCluster *>(*clustIter);
216  if(!clust->IsEMCAL()) continue;
217 
218  TLorentzVector posvec;
219  energy = clust->E();
220  clust->GetMomentum(posvec, vertexpos);
221  eta = posvec.Eta();
222  phi = posvec.Phi();
223 
224  // fill histograms allEta
225  if(isMinBias){
226  FillClusterHistograms("MB", energy, eta, phi, NULL);
227  }
228  if(isEMC7){
229  FillClusterHistograms("EMC7", energy, eta, phi, NULL);
230  if(!isMinBias){
231  FillClusterHistograms("EMC7excl", energy, eta, phi, NULL);
232  }
233  }
234  if(isEJ2){
235  TList ej2patches;
236  FindPatchesForTrigger("EJ2", triggerpatches, ej2patches);
237  FillClusterHistograms("EJ2", energy, eta, phi, &ej2patches);
238  // check for exclusive classes
239  if(!isMinBias){
240  FillClusterHistograms("EJ2excl", energy, eta, phi, &ej2patches);
241  }
242  if(isEG1 || isEG2){
243  FillClusterHistograms("E2combined", energy, eta, phi, &ej2patches);
244  } else {
245  FillClusterHistograms("E2Jonly", energy, eta, phi, &ej2patches);
246  }
247  }
248  if(isEJ1){
249  TList ej1patches;
250  FindPatchesForTrigger("EJ1", triggerpatches, ej1patches);
251  FillClusterHistograms("EJ1", energy, eta, phi, &ej1patches);
252  // check for exclusive classes
253  if(!(isMinBias || isEJ2)){
254  FillClusterHistograms("EJ1excl", energy, eta, phi, &ej1patches);
255  }
256  if(isEG1 || isEG2) {
257  FillClusterHistograms("E1combined", energy, eta, phi, &ej1patches);
258  } else {
259  FillClusterHistograms("E1Jonly", energy, eta, phi, &ej1patches);
260  }
261  }
262  if(isEG2){
263  TList eg2patches;
264  FindPatchesForTrigger("EG2", triggerpatches, eg2patches);
265  FillClusterHistograms("EG2", energy, eta, phi, &eg2patches);
266  // check for exclusive classes
267  if(!isMinBias){
268  FillClusterHistograms("EG2excl", energy, eta, phi, &eg2patches);
269  }
270  if(!(isEJ2 || isEJ1)){
271  FillClusterHistograms("E2Gonly", energy, eta, phi, &eg2patches);
272  }
273  }
274  if(isEG1){
275  TList eg1patches;
276  FindPatchesForTrigger("EG1", triggerpatches, eg1patches);
277  FillClusterHistograms("EG1", energy, eta, phi, &eg1patches);
278  if(!(isEG2 || isMinBias))
279  FillClusterHistograms("EG1excl", energy, eta, phi, &eg1patches);
280  if(!(isEJ1 || isEJ2)){
281  FillClusterHistograms("E1Gonly", energy, eta, phi, &eg1patches);
282  }
283  }
284  }
285  PostData(1, fHistos->GetListOfHistograms());
286 }
287 
288 void AliAnalysisTaskEmcalClustersRef::FillClusterHistograms(TString triggerclass, double energy, double eta, double phi, TList *triggerpatches){
289  Bool_t hasTriggerPatch = triggerpatches ? CorrelateToTrigger(eta, phi, triggerpatches) : kFALSE;
290  Int_t supermoduleID = -1;
291  fGeometry->SuperModuleNumberFromEtaPhi(eta, phi, supermoduleID);
292  fHistos->FillTH1(Form("hClusterEnergy%s", triggerclass.Data()), energy);
293  if(supermoduleID >= 0)
294  fHistos->FillTH2(Form("hClusterEnergySM%s", triggerclass.Data()), supermoduleID, energy);
295  if(hasTriggerPatch){
296  fHistos->FillTH1(Form("hClusterEnergyFired%s", triggerclass.Data()), energy);
297  if(supermoduleID >= 0)
298  fHistos->FillTH2(Form("hClusterEnergyFiredSM%s", triggerclass.Data()), supermoduleID, energy);
299  }
300  Double_t encuts[5] = {1., 2., 5., 10., 20.};
301  for(int ien = 0; ien < 5; ien++){
302  if(energy > encuts[ien]){
303  fHistos->FillTH2(Form("hEtaPhi%dG%s", static_cast<int>(encuts[ien]), triggerclass.Data()), eta, phi);
304  if(hasTriggerPatch){
305  fHistos->FillTH2(Form("hEtaPhiFired%dG%s", static_cast<int>(encuts[ien]), triggerclass.Data()), eta, phi);
306  }
307  }
308  }
309 }
310 
315 void AliAnalysisTaskEmcalClustersRef::CreateEnergyBinning(TArrayD& binning) const {
316  std::vector<double> mybinning;
317  std::map<double,double> definitions;
318  definitions.insert(std::pair<double, double>(1, 0.05));
319  definitions.insert(std::pair<double, double>(2, 0.1));
320  definitions.insert(std::pair<double, double>(4, 0.2));
321  definitions.insert(std::pair<double, double>(7, 0.5));
322  definitions.insert(std::pair<double, double>(16, 1));
323  definitions.insert(std::pair<double, double>(32, 2));
324  definitions.insert(std::pair<double, double>(40, 4));
325  definitions.insert(std::pair<double, double>(50, 5));
326  definitions.insert(std::pair<double, double>(100, 10));
327  definitions.insert(std::pair<double, double>(200, 20));
328  double currentval = 0.;
329  mybinning.push_back(currentval);
330  for(std::map<double,double>::iterator id = definitions.begin(); id != definitions.end(); ++id){
331  double limit = id->first, binwidth = id->second;
332  while(currentval < limit){
333  currentval += binwidth;
334  mybinning.push_back(currentval);
335  }
336  }
337  binning.Set(mybinning.size());
338  int ib = 0;
339  for(std::vector<double>::iterator it = mybinning.begin(); it != mybinning.end(); ++it)
340  binning[ib++] = *it;
341 }
342 
350 Bool_t AliAnalysisTaskEmcalClustersRef::CorrelateToTrigger(Double_t etaclust, Double_t phiclust, TList *triggerpatches) const {
351  Bool_t hasfound = kFALSE;
352  for(TIter patchIter = TIter(triggerpatches).Begin(); patchIter != TIter::End(); ++patchIter){
353  AliEmcalTriggerPatchInfo *mypatch = static_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
354  Double_t etamin = TMath::Min(mypatch->GetEtaMin(), mypatch->GetEtaMax()),
355  etamax = TMath::Max(mypatch->GetEtaMin(), mypatch->GetEtaMax()),
356  phimin = TMath::Min(mypatch->GetPhiMin(), mypatch->GetPhiMax()),
357  phimax = TMath::Max(mypatch->GetPhiMin(), mypatch->GetPhiMax());
358  if(etaclust > etamin && etaclust < etamax && phiclust > phimin && phiclust < phimax){
359  hasfound = kTRUE;
360  break;
361  }
362  }
363  return hasfound;
364 }
365 
372 void AliAnalysisTaskEmcalClustersRef::FindPatchesForTrigger(TString triggerclass, const TClonesArray * triggerpatches, TList &foundtriggers) const {
373  double minADC_EJ1 = 260.,
374  minADC_EJ2 = 127.,
375  minADC_EG1 = 140.,
376  minADC_EG2 = 89.;
377  if(!triggerpatches) return;
378  Bool_t isEG1 = (triggerclass == "EG1"),
379  isEG2 = (triggerclass == "EG2"),
380  isEJ1 = (triggerclass == "EJ1"),
381  isEJ2 = (triggerclass == "EJ2");
382  AliEmcalTriggerPatchInfo *mypatch = NULL;
383  for(TIter patchiter = TIter(triggerpatches).Begin(); patchiter != TIter::End(); ++patchiter){
384  mypatch = dynamic_cast<AliEmcalTriggerPatchInfo *>(*patchiter);
385  if(!mypatch->IsOfflineSimple()) continue;
386  if(isEG1){
387  if(mypatch->IsGammaHighSimple() && mypatch->GetADCAmp() > minADC_EG1)
388  foundtriggers.Add(mypatch);
389  }
390  if(isEG2){
391  if(mypatch->IsGammaLowSimple() && mypatch->GetADCAmp() > minADC_EG2)
392  foundtriggers.Add(mypatch);
393  }
394  if(isEJ1){
395  if(mypatch->IsJetHighSimple() && mypatch->GetADCAmp() > minADC_EJ1)
396  foundtriggers.Add(mypatch);
397  }
398  if(isEJ2){
399  if(mypatch->IsJetLowSimple() && mypatch->GetADCAmp() > minADC_EJ2)
400  foundtriggers.Add(mypatch);
401  }
402  }
403 }
404 
410 TString AliAnalysisTaskEmcalClustersRef::GetFiredTriggerClassesFromPatches(const TClonesArray* triggerpatches) const {
411  TString triggerstring = "";
412  if(!triggerpatches) return triggerstring;
413  Int_t nEJ1 = 0, nEJ2 = 0, nEG1 = 0, nEG2 = 0;
414  double minADC_EJ1 = 260.,
415  minADC_EJ2 = 127.,
416  minADC_EG1 = 140.,
417  minADC_EG2 = 89.;
418  for(TIter patchIter = TIter(triggerpatches).Begin(); patchIter != TIter::End(); ++patchIter){
419  AliEmcalTriggerPatchInfo *patch = dynamic_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
420  if(!patch->IsOfflineSimple()) continue;
421  if(patch->IsJetHighSimple() && patch->GetADCOfflineAmp() > minADC_EJ1) nEJ1++;
422  if(patch->IsJetLowSimple() && patch->GetADCOfflineAmp() > minADC_EJ2) nEJ2++;
423  if(patch->IsGammaHighSimple() && patch->GetADCOfflineAmp() > minADC_EG1) nEG1++;
424  if(patch->IsGammaLowSimple() && patch->GetADCOfflineAmp() > minADC_EG2) nEG2++;
425  }
426  if(nEJ1) triggerstring += "EJ1";
427  if(nEJ2){
428  if(triggerstring.Length()) triggerstring += ",";
429  triggerstring += "EJ2";
430  }
431  if(nEG1){
432  if(triggerstring.Length()) triggerstring += ",";
433  triggerstring += "EG1";
434  }
435  if(nEG2){
436  if(triggerstring.Length()) triggerstring += ",";
437  triggerstring += "EG2";
438  }
439  return triggerstring;
440 }
441 
442 
443 } /* namespace EMCalTriggerPtAnalysis */
ClassImp(EMCalTriggerPtAnalysis::AliAnalysisTaskEmcalClustersRef) namespace EMCalTriggerPtAnalysis
const Double_t etamin
Main data structure storing all relevant information of EMCAL/DCAL trigger patches.
Declarartion of class AliEMCalHistoContainer.
energy
const Double_t etamax
Class to make array of trigger patch objects in AOD/ESD events.
Int_t GetRunNumber(TString)
Definition: PlotMuonQA.C:2235
const Double_t phimin