AliPhysics  vAN-20150723 (baea2bf)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliAnalysisTaskChargedParticlesRefMC.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 <vector>
16 #include <map>
17 
18 #include <TArrayD.h>
19 #include <TChain.h>
20 #include <TClonesArray.h>
21 #include <TFile.h>
22 #include <THashList.h>
23 #include <TH1.h>
24 #include <TKey.h>
25 #include <TList.h>
26 #include <TProfile.h>
27 #include <TMath.h>
28 #include <TString.h>
29 #include <TSystem.h>
30 #include <TTree.h>
31 
32 #include "AliAnalysisManager.h"
33 #include "AliAnalysisUtils.h"
34 #include "AliAODMCHeader.h"
35 #include "AliAODMCParticle.h"
36 #include "AliAODTrack.h"
38 #include "AliESDtrackCuts.h"
39 #include "AliESDEvent.h"
40 #include "AliESDtrack.h"
41 #include "AliGenPythiaEventHeader.h"
42 #include "AliInputEventHandler.h"
43 #include "AliMCEvent.h"
44 #include "AliVVertex.h"
45 
46 #include "AliEMCalHistoContainer.h"
48 
52 
53 namespace EMCalTriggerPtAnalysis {
54 
58 AliAnalysisTaskChargedParticlesRefMC::AliAnalysisTaskChargedParticlesRefMC():
59  AliAnalysisTaskSE(),
60  fTrackCuts(NULL),
61  fAnalysisUtil(NULL),
62  fHistos(NULL),
63  fPtHard(0),
64  fPtHardBin(0),
65  fNTrials(0),
66  fXsection(0)
67 {
68 }
69 
75  AliAnalysisTaskSE(name),
76  fTrackCuts(NULL),
77  fAnalysisUtil(NULL),
78  fHistos(NULL),
79  fPtHard(0),
80  fPtHardBin(0),
81  fNTrials(0),
82  fXsection(0)
83 {
84  DefineOutput(1, TList::Class());
85 }
86 
91  if(fTrackCuts) delete fTrackCuts;
92  if(fAnalysisUtil) delete fAnalysisUtil;
93  if(fHistos) delete fHistos;
94 }
99  fAnalysisUtil = new AliAnalysisUtils;
100 
101  fTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(true, 1);
102  fTrackCuts->SetName("Standard Track cuts");
103  fTrackCuts->SetMinNCrossedRowsTPC(120);
104  fTrackCuts->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
105 
106  TArrayD oldbinning, newbinning;
107  CreateOldPtBinning(oldbinning);
108  CreateNewPtBinning(newbinning);
109 
110  fHistos = new AliEMCalHistoContainer("Ref");
111  fHistos->CreateTH1("hNtrials", "Number of trials", 1, 0.5, 1.5);
112  fHistos->CreateTProfile("hCrossSection", "PYTHIA cross section", 1, 0.5, 1.5);
113  fHistos->CreateTH1("hNtrialsEvent", "Number of trials (from header, after event selection)", 1, 0.5, 1.5);
114  fHistos->CreateTProfile("hCrossSectionEvent", "PYTHIA cross section (from header, after event selection)", 1, 0.5, 1.5);
115  fHistos->CreateTH1("hPtHard", "Pt of the hard interaction", 1000, 0., 500);
116  TString triggers[6] = {"True", "MB", "EJ1", "EJ2", "EG1", "EG2"};
117  for(TString *trg = triggers; trg < triggers + sizeof(triggers)/sizeof(TString); trg++){
118  fHistos->CreateTH1(Form("hEventCount%s", trg->Data()), Form("Event Counter for trigger class %s", trg->Data()), 1, 0.5, 1.5);
119  fHistos->CreateTH1(Form("hVertexBefore%s", trg->Data()), Form("Vertex distribution before z-cut for trigger class %s", trg->Data()), 500, -50, 50);
120  fHistos->CreateTH1(Form("hVertexAfter%s", trg->Data()), Form("Vertex distribution after z-cut for trigger class %s", trg->Data()), 100, -10, 10);
121  fHistos->CreateTH1(Form("hPtEtaAllOldBinning%s", trg->Data()), Form("Charged particle pt distribution all eta old binning trigger %s", trg->Data()), oldbinning);
122  fHistos->CreateTH1(Form("hPtEtaCentOldBinning%s", trg->Data()), Form("Charged particle pt distribution central eta old binning trigger %s", trg->Data()), oldbinning);
123  fHistos->CreateTH1(Form("hPtEtaAllNewBinning%s", trg->Data()), Form("Charged particle pt distribution all eta new binning trigger %s", trg->Data()), newbinning);
124  fHistos->CreateTH1(Form("hPtEtaCentNewBinning%s", trg->Data()), Form("Charged particle pt distribution central eta new binning trigger %s", trg->Data()), newbinning);
125  fHistos->CreateTH1(Form("hEtaDistAllPt1%s", trg->Data()), Form("Eta distribution without etacut for tracks with Pt above 1 GeV/c trigger %s", trg->Data()), 100, -1., 1.);
126  fHistos->CreateTH1(Form("hEtaDistAllPt2%s", trg->Data()), Form("Eta distribution without etacut for tracks with Pt above 2 GeV/c trigger %s", trg->Data()), 100, -1., 1.);
127  fHistos->CreateTH1(Form("hEtaDistAllPt5%s", trg->Data()), Form("Eta distribution without etacut for tracks with Pt above 5 GeV/c trigger %s", trg->Data()), 100, -1., 1.);
128  fHistos->CreateTH1(Form("hEtaDistAllPt10%s", trg->Data()), Form("Eta distribution without etacut for tracks with Pt above 10 GeV/c trigger %s", trg->Data()), 100, -1., 1.);
129  fHistos->CreateTH1(Form("hEtaDistAllPt20%s", trg->Data()), Form("Eta distribution without etacut for tracks with Pt above 20 GeV/c trigger %s", trg->Data()), 100, -1., 1.);
130  fHistos->CreateTH1(Form("hEtaDistCutPt1%s", trg->Data()), Form("Eta distribution with etacut for tracks with Pt above 1 GeV/c trigger %s", trg->Data()), 100, -1., 1.);
131  fHistos->CreateTH1(Form("hEtaDistCutPt2%s", trg->Data()), Form("Eta distribution with etacut for tracks with Pt above 2 GeV/c trigger %s", trg->Data()), 100, -1., 1.);
132  fHistos->CreateTH1(Form("hEtaDistCutPt5%s", trg->Data()), Form("Eta distribution with etacut for tracks with Pt above 5 GeV/c trigger %s", trg->Data()), 100, -1., 1.);
133  fHistos->CreateTH1(Form("hEtaDistCutPt10%s", trg->Data()), Form("Eta distribution with etacut for tracks with Pt above 10 GeV/c trigger %s", trg->Data()), 100, -1., 1.);
134  fHistos->CreateTH1(Form("hEtaDistCutPt20%s", trg->Data()), Form("Eta distribution with etacut for tracks with Pt above 20 GeV/c trigger %s", trg->Data()), 100, -1., 1.);
135  }
136  PostData(1, fHistos->GetListOfHistograms());
137 }
138 
139 void AliAnalysisTaskChargedParticlesRefMC::UserExec(Option_t*) { // Select event
140  if(!fMCEvent) return;
141  TClonesArray *fTriggerPatches = dynamic_cast<TClonesArray *>(fInputEvent->FindListObject("EmcalTriggers"));
142  if(!fTriggerPatches) return;
143 
144  TString triggerstring = GetFiredTriggerClasses(fTriggerPatches);
145  Bool_t isMinBias = fInputHandler->IsEventSelected() & AliVEvent::kINT7,
146  isEJ1 = triggerstring.Contains("EJ1"),
147  isEJ2 = triggerstring.Contains("EJ2"),
148  isEG1 = triggerstring.Contains("EG1"),
149  isEG2 = triggerstring.Contains("EG2");
150  if(!(isMinBias || isEG1 || isEG2 || isEJ1 || isEJ2)) return;
151  const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
152  //if(!fInputEvent->IsPileupFromSPD(3, 0.8, 3., 2., 5.)) return; // reject pileup event
153  if(vtx->GetNContributors() < 1) return;
154  // Fill reference distribution for the primary vertex before any z-cut
155  fHistos->FillTH1("hVertexBeforeTrue", vtx->GetZ());
156  if(isMinBias) fHistos->FillTH1("hVertexBeforeMB", vtx->GetZ());
157  if(isEJ1) fHistos->FillTH1("hVertexBeforeEJ1", vtx->GetZ());
158  if(isEJ2) fHistos->FillTH1("hVertexBeforeEJ2", vtx->GetZ());
159  if(isEG1) fHistos->FillTH1("hVertexBeforeEG1", vtx->GetZ());
160  if(isEG2) fHistos->FillTH1("hVertexBeforeEG2", vtx->GetZ());
161  if(!fAnalysisUtil->IsVertexSelected2013pA(fInputEvent)) return; // Apply new vertex cut
162  if(fAnalysisUtil->IsPileUpEvent(fInputEvent)) return; // Apply new vertex cut
163  // Apply vertex z cut
164  if(vtx->GetZ() < -10. || vtx->GetZ() > 10.) return;
165 
166  // Fill Event counter and reference vertex distributions for the different trigger classes
167  fHistos->FillTH1("hEventCountTrue", 1);
168  fHistos->FillTH1("hVertexAfterTrue", vtx->GetZ());
169  if(isMinBias){
170  fHistos->FillTH1("hEventCountMB", 1);
171  fHistos->FillTH1("hVertexAfterMB", vtx->GetZ());
172  }
173  if(isEJ1){
174  fHistos->FillTH1("hEventCountEJ1", 1);
175  fHistos->FillTH1("hVertexAfterEJ1", vtx->GetZ());
176  }
177  if(isEJ2){
178  fHistos->FillTH1("hEventCountEJ2", 1);
179  fHistos->FillTH1("hVertexAfterEJ2", vtx->GetZ());
180  }
181  if(isEG1){
182  fHistos->FillTH1("hEventCountEG1", 1);
183  fHistos->FillTH1("hVertexAfterEG1", vtx->GetZ());
184  }
185  if(isEG2){
186  fHistos->FillTH1("hEventCountEG2", 1);
187  fHistos->FillTH1("hVertexAfterEG2", vtx->GetZ());
188  }
189 
190  // Fill PYTHIA histograms from event header
191  AliGenPythiaEventHeader *pyheader = GetPythiaHeader();
192  if(pyheader){
193  fHistos->FillTH1("hNtrialsEvent", 1., pyheader->Trials());
194  fHistos->FillProfile("hCrossSectionEvent", 1., pyheader->GetXsection());
195  fHistos->FillTH1("hPtHard", pyheader->GetPtHard());
196  }
197 
198  // MonteCarlo Loop
199  // Histograms
200  // - Full eta (-0.8, 0.8), new binning
201  // - Full eta (-0.8, 0.8), old binning
202  // - Eta distribution for tracks above 1, 2, 5, 10 GeV/c without eta cut
203  // - Central eta (-0.8, -0.2), new binning,
204  // - Central eta (-0.8, -0.2), old binning,
205  // - Eta distribution for tracks above 1, 2, 5, 10 GeV/c
206  // - Eta distribution for tracks above 1, 2, 5, 10 GeV/c with eta cut
207  int ptmin[5] = {1,2,5,10,20}; // for eta distributions
208  AliVParticle *truepart = NULL;
209  for(int ipart = 0; ipart < fMCEvent->GetNumberOfTracks(); ipart++){
210  truepart = fMCEvent->GetTrack(ipart);
211 
212  // Select only particles within ALICE acceptance
213  if(TMath::Abs(truepart->Eta()) > 0.8) continue;
214  if(TMath::Abs(truepart->Pt()) < 0.1) continue;
215  if(!truepart->Charge()) continue;
216 
217  if(!IsPhysicalPrimary(truepart, fMCEvent)) continue;
218 
219  // Particle selected
220  fHistos->FillTH1("hPtEtaAllNewBinningTrue", TMath::Abs(truepart->Pt()));
221  fHistos->FillTH1("hPtEtaAllOldBinningTrue", TMath::Abs(truepart->Pt()));
222 
223  for(int icut = 0; icut < 5; icut++){
224  if(TMath::Abs(truepart->Pt()) > static_cast<double>(ptmin[icut])){
225  fHistos->FillTH1(Form("hEtaDistAllPt%dTrue", ptmin[icut]), truepart->Eta());
226  }
227  }
228 
229  if(truepart->Eta() > -0.8 && truepart->Eta() < -0.2){
230  fHistos->FillTH1("hPtEtaCentNewBinningTrue", TMath::Abs(truepart->Pt()));
231  fHistos->FillTH1("hPtEtaCentOldBinningTrue", TMath::Abs(truepart->Pt()));
232  for(int icut = 0; icut < 5; icut++){
233  if(TMath::Abs(truepart->Pt()) > static_cast<double>(ptmin[icut])){
234  fHistos->FillTH1(Form("hEtaDistCutPt%dTrue", ptmin[icut]), truepart->Eta());
235  }
236  }
237  }
238  }
239 
240  // Loop over tracks, fill select particles
241  // Histograms
242  // - Full eta (-0.8, 0.8), new binning
243  // - Full eta (-0.8, 0.8), old binning
244  // - Eta distribution for tracks above 1, 2, 5, 10 GeV/c without eta cut
245  // - Central eta (-0.8, -0.2), new binning,
246  // - Central eta (-0.8, -0.2), old binning,
247  // - Eta distribution for tracks above 1, 2, 5, 10 GeV/c
248  // - Eta distribution for tracks above 1, 2, 5, 10 GeV/c with eta cut
249  AliVTrack *checktrack(NULL);
250  AliVParticle *assocMC(NULL);
251  double ptparticle(-1.), etaparticle(-100.);
252  for(int itrk = 0; itrk < fInputEvent->GetNumberOfTracks(); ++itrk){
253  checktrack = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(itrk));
254  if(!checktrack) continue;
255  // Find associated particle
256  assocMC = fMCEvent->GetTrack(TMath::Abs(checktrack->GetLabel()));
257  if(!assocMC) continue; // Fake track
258  if(!IsPhysicalPrimary(assocMC, fMCEvent)) continue;
259 
260  // Select only particles within ALICE acceptance
261  if(TMath::Abs(checktrack->Eta()) > 0.8) continue;
262  if(TMath::Abs(checktrack->Pt()) < 0.1) continue;
263 
264  // Distinguish track selection for ESD and AOD tracks
265  AliESDtrack *esdtrack(NULL);
266  AliAODTrack *aodtrack(NULL);
267  if((esdtrack = dynamic_cast<AliESDtrack *>(checktrack))){
268  if(!TrackSelectionESD(esdtrack)) continue;
269  } else if((aodtrack = dynamic_cast<AliAODTrack *>(checktrack))){
270  if(!TrackSelectionAOD(aodtrack)) continue;
271  } else {
272  continue;
273  }
274 
275  ptparticle = TMath::Abs(assocMC->Pt());
276  etaparticle = assocMC->Eta();
277 
278  // fill histograms allEta
279  if(isMinBias){
280  fHistos->FillTH1("hPtEtaAllNewBinningMB", ptparticle);
281  fHistos->FillTH1("hPtEtaAllOldBinningMB", ptparticle);
282  for(int icut = 0; icut < 5; icut++){
283  if(TMath::Abs(checktrack->Pt()) > static_cast<double>(ptmin[icut])){
284  fHistos->FillTH1(Form("hEtaDistAllPt%dMB", ptmin[icut]), etaparticle);
285  }
286  }
287  }
288  if(isEJ1){
289  fHistos->FillTH1("hPtEtaAllNewBinningEJ1", ptparticle);
290  fHistos->FillTH1("hPtEtaAllOldBinningEJ1", ptparticle);
291  for(int icut = 0; icut < 5; icut++){
292  if(TMath::Abs(checktrack->Pt()) > static_cast<double>(ptmin[icut])){
293  fHistos->FillTH1(Form("hEtaDistAllPt%dEJ1", ptmin[icut]), etaparticle);
294  }
295  }
296  }
297  if(isEJ2){
298  fHistos->FillTH1("hPtEtaAllNewBinningEJ2", ptparticle);
299  fHistos->FillTH1("hPtEtaAllOldBinningEJ2", ptparticle);
300  for(int icut = 0; icut < 5; icut++){
301  if(TMath::Abs(checktrack->Pt()) > static_cast<double>(ptmin[icut])){
302  fHistos->FillTH1(Form("hEtaDistAllPt%dEJ2", ptmin[icut]), etaparticle);
303  }
304  }
305  }
306  if(isEG1){
307  fHistos->FillTH1("hPtEtaAllNewBinningEG1", ptparticle);
308  fHistos->FillTH1("hPtEtaAllOldBinningEG1", ptparticle);
309  for(int icut = 0; icut < 5; icut++){
310  if(TMath::Abs(checktrack->Pt()) > static_cast<double>(ptmin[icut])){
311  fHistos->FillTH1(Form("hEtaDistAllPt%dEG1", ptmin[icut]), etaparticle);
312  }
313  }
314  }
315  if(isEG2){
316  fHistos->FillTH1("hPtEtaAllNewBinningEG2", ptparticle);
317  fHistos->FillTH1("hPtEtaAllOldBinningEG2", ptparticle);
318  for(int icut = 0; icut < 5; icut++){
319  if(TMath::Abs(checktrack->Pt()) > static_cast<double>(ptmin[icut])){
320  fHistos->FillTH1(Form("hEtaDistAllPt%dEG2", ptmin[icut]), etaparticle);
321  }
322  }
323  }
324 
325  if(checktrack->Eta() > -0.8 && checktrack->Eta() < -0.2){
326  // Fill Histograms in central eta
327  if(isMinBias){
328  fHistos->FillTH1("hPtEtaCentNewBinningMB", ptparticle);
329  fHistos->FillTH1("hPtEtaCentOldBinningMB", ptparticle);
330  for(int icut = 0; icut < 5; icut++){
331  if(TMath::Abs(checktrack->Pt()) > static_cast<double>(ptmin[icut])){
332  fHistos->FillTH1(Form("hEtaDistCutPt%dMB", ptmin[icut]), etaparticle);
333  }
334  }
335  }
336  if(isEJ1){
337  fHistos->FillTH1("hPtEtaCentNewBinningEJ1", ptparticle);
338  fHistos->FillTH1("hPtEtaCentOldBinningEJ1", ptparticle);
339  for(int icut = 0; icut < 5; icut++){
340  if(TMath::Abs(checktrack->Pt()) > static_cast<double>(ptmin[icut])){
341  fHistos->FillTH1(Form("hEtaDistCutPt%dEJ1", ptmin[icut]), etaparticle);
342  }
343  }
344  }
345  if(isEJ2){
346  fHistos->FillTH1("hPtEtaCentNewBinningEJ2", ptparticle);
347  fHistos->FillTH1("hPtEtaCentOldBinningEJ2", ptparticle);
348  for(int icut = 0; icut < 5; icut++){
349  if(TMath::Abs(checktrack->Pt()) > static_cast<double>(ptmin[icut])){
350  fHistos->FillTH1(Form("hEtaDistCutPt%dEJ2", ptmin[icut]), etaparticle);
351  }
352  }
353  }
354  if(isEG1){
355  fHistos->FillTH1("hPtEtaCentNewBinningEG1", ptparticle);
356  fHistos->FillTH1("hPtEtaCentOldBinningEG1", ptparticle);
357  for(int icut = 0; icut < 5; icut++){
358  if(TMath::Abs(checktrack->Pt()) > static_cast<double>(ptmin[icut])){
359  fHistos->FillTH1(Form("hEtaDistCutPt%dEG1", ptmin[icut]), etaparticle);
360  }
361  }
362  }
363  if(isEG2){
364  fHistos->FillTH1("hPtEtaCentNewBinningEG2", ptparticle);
365  fHistos->FillTH1("hPtEtaCentOldBinningEG2", ptparticle);
366  for(int icut = 0; icut < 5; icut++){
367  if(TMath::Abs(checktrack->Pt()) > static_cast<double>(ptmin[icut])){
368  fHistos->FillTH1(Form("hEtaDistCutPt%dEG2", ptmin[icut]), etaparticle);
369  }
370  }
371  }
372  }
373  }
374  PostData(1, fHistos->GetListOfHistograms());
375 
376 }
377 
383 {
384  // Called when file changes.
385 
386  TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
387  if (!tree) {
388  AliError(Form("%s - UserNotify: No current tree!",GetName()));
389  return kFALSE;
390  }
391 
392  Float_t xsection = 0;
393  Float_t trials = 0;
394  Int_t pthardbin = 0;
395 
396  TFile *curfile = tree->GetCurrentFile();
397  if (!curfile) {
398  AliError(Form("%s - UserNotify: No current file!",GetName()));
399  return kFALSE;
400  }
401 
402  TChain *chain = dynamic_cast<TChain*>(tree);
403  if (chain) tree = chain->GetTree();
404 
405  Int_t nevents = tree->GetEntriesFast();
406 
407  PythiaInfoFromFile(curfile->GetName(), xsection, trials, pthardbin);
408 
409  fHistos->FillTH1("hNtrials", 1., trials);
410  fHistos->FillProfile("hCrossSection", 1., xsection);
411  //fHistos->FillTH1(pthardbin, nevents);
412 
413  return kTRUE;
414 }
415 
428 Bool_t AliAnalysisTaskChargedParticlesRefMC::PythiaInfoFromFile(const char* currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard) const {
429 
430  TString file(currFile);
431  fXsec = 0;
432  fTrials = 1;
433 
434  if (file.Contains(".zip#")) {
435  Ssiz_t pos1 = file.Index("root_archive",12,0,TString::kExact);
436  Ssiz_t pos = file.Index("#",1,pos1,TString::kExact);
437  Ssiz_t pos2 = file.Index(".root",5,TString::kExact);
438  file.Replace(pos+1,pos2-pos1,"");
439  } else {
440  // not an archive take the basename....
441  file.ReplaceAll(gSystem->BaseName(file.Data()),"");
442  }
443  AliDebug(1,Form("File name: %s",file.Data()));
444 
445  // Get the pt hard bin
446  TString strPthard(file);
447 
448  strPthard.Remove(strPthard.Last('/'));
449  strPthard.Remove(strPthard.Last('/'));
450  if (strPthard.Contains("AOD")) strPthard.Remove(strPthard.Last('/'));
451  strPthard.Remove(0,strPthard.Last('/')+1);
452  if (strPthard.IsDec())
453  pthard = strPthard.Atoi();
454  else
455  AliWarning(Form("Could not extract file number from path %s", strPthard.Data()));
456 
457  // problem that we cannot really test the existance of a file in a archive so we have to live with open error message from root
458  TFile *fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec.root"));
459 
460  if (!fxsec) {
461  // next trial fetch the histgram file
462  fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
463  if (!fxsec) {
464  // not a severe condition but inciate that we have no information
465  return kFALSE;
466  } else {
467  // find the tlist we want to be independtent of the name so use the Tkey
468  TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
469  if (!key) {
470  fxsec->Close();
471  return kFALSE;
472  }
473  TList *list = dynamic_cast<TList*>(key->ReadObj());
474  if (!list) {
475  fxsec->Close();
476  return kFALSE;
477  }
478  fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
479  fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
480  fxsec->Close();
481  }
482  } else { // no tree pyxsec.root
483  TTree *xtree = (TTree*)fxsec->Get("Xsection");
484  if (!xtree) {
485  fxsec->Close();
486  return kFALSE;
487  }
488  UInt_t ntrials = 0;
489  Double_t xsection = 0;
490  xtree->SetBranchAddress("xsection",&xsection);
491  xtree->SetBranchAddress("ntrials",&ntrials);
492  xtree->GetEntry(0);
493  fTrials = ntrials;
494  fXsec = xsection;
495  fxsec->Close();
496  }
497  return kTRUE;
498 }
499 
505  AliGenPythiaEventHeader *pythiaHeader = dynamic_cast<AliGenPythiaEventHeader*>(MCEvent()->GenEventHeader());
506  if (!pythiaHeader) {
507  // Check if AOD
508  AliAODMCHeader* aodMCH = dynamic_cast<AliAODMCHeader*>(InputEvent()->FindListObject(AliAODMCHeader::StdBranchName()));
509  if (aodMCH) {
510  for (UInt_t i = 0;i<aodMCH->GetNCocktailHeaders();i++) {
511  pythiaHeader = dynamic_cast<AliGenPythiaEventHeader*>(aodMCH->GetCocktailHeader(i));
512  if (pythiaHeader) break;
513  }
514  }
515  }
516  return pythiaHeader;
517 }
518 
524  std::vector<double> mybinning;
525  std::map<double,double> definitions;
526  definitions.insert(std::pair<double,double>(2.5, 0.1));
527  definitions.insert(std::pair<double,double>(7., 0.25));
528  definitions.insert(std::pair<double,double>(15., 0.5));
529  definitions.insert(std::pair<double,double>(25., 1.));
530  definitions.insert(std::pair<double,double>(40., 2.5));
531  definitions.insert(std::pair<double,double>(50., 5.));
532  definitions.insert(std::pair<double,double>(100., 10.));
533  double currentval = 0;
534  for(std::map<double,double>::iterator id = definitions.begin(); id != definitions.end(); ++id){
535  double limit = id->first, binwidth = id->second;
536  while(currentval < limit){
537  currentval += binwidth;
538  mybinning.push_back(currentval);
539  }
540  }
541  binning.Set(mybinning.size());
542  int ib = 0;
543  for(std::vector<double>::iterator it = mybinning.begin(); it != mybinning.end(); ++it)
544  binning[ib++] = *it;
545 }
546 
552  std::vector<double> mybinning;
553  std::map<double,double> definitions;
554  definitions.insert(std::pair<double, double>(1, 0.05));
555  definitions.insert(std::pair<double, double>(2, 0.1));
556  definitions.insert(std::pair<double, double>(4, 0.2));
557  definitions.insert(std::pair<double, double>(7, 0.5));
558  definitions.insert(std::pair<double, double>(16, 1));
559  definitions.insert(std::pair<double, double>(36, 2));
560  definitions.insert(std::pair<double, double>(40, 4));
561  definitions.insert(std::pair<double, double>(50, 5));
562  definitions.insert(std::pair<double, double>(100, 10));
563  definitions.insert(std::pair<double, double>(200, 20));
564  double currentval = 0.;
565  mybinning.push_back(currentval);
566  for(std::map<double,double>::iterator id = definitions.begin(); id != definitions.end(); ++id){
567  double limit = id->first, binwidth = id->second;
568  while(currentval < limit){
569  currentval += binwidth;
570  mybinning.push_back(currentval);
571  }
572  }
573  binning.Set(mybinning.size());
574  int ib = 0;
575  for(std::vector<double>::iterator it = mybinning.begin(); it != mybinning.end(); ++it)
576  binning[ib++] = *it;
577 }
578 
585  return fTrackCuts->AcceptTrack(track);
586 }
587 
594  if(!track->TestFilterBit(AliAODTrack::kTrkGlobal)) return false;
595  if(track->GetTPCNCrossedRows() < 120) return false;
596  return true;
597 }
598 
604 TString AliAnalysisTaskChargedParticlesRefMC::GetFiredTriggerClasses(const TClonesArray* triggerpatches) {
605  TString triggerstring = "";
606  Int_t nEJ1 = 0, nEJ2 = 0, nEG1 = 0, nEG2 = 0;
607  double minADC_EJ1 = 260.,
608  minADC_EJ2 = 127.,
609  minADC_EG1 = 140.,
610  minADC_EG2 = 89.;
611  for(TIter patchIter = TIter(triggerpatches).Begin(); patchIter != TIter::End(); ++patchIter){
612  AliEmcalTriggerPatchInfo *patch = dynamic_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
613  if(!patch->IsOfflineSimple()) continue;
614  if(patch->IsJetHighSimple() && patch->GetADCOfflineAmp() > minADC_EJ1) nEJ1++;
615  if(patch->IsJetLowSimple() && patch->GetADCOfflineAmp() > minADC_EJ2) nEJ2++;
616  if(patch->IsGammaHighSimple() && patch->GetADCOfflineAmp() > minADC_EG1) nEG1++;
617  if(patch->IsGammaLowSimple() && patch->GetADCOfflineAmp() > minADC_EG2) nEG2++;
618  }
619  if(nEJ1) triggerstring += "EJ1";
620  if(nEJ2){
621  if(triggerstring.Length()) triggerstring += ",";
622  triggerstring += "EJ2";
623  }
624  if(nEG1){
625  if(triggerstring.Length()) triggerstring += ",";
626  triggerstring += "EG1";
627  }
628  if(nEG2){
629  if(triggerstring.Length()) triggerstring += ",";
630  triggerstring += "EG2";
631  }
632  return triggerstring;
633 }
634 
643 Bool_t AliAnalysisTaskChargedParticlesRefMC::IsPhysicalPrimary(const AliVParticle* const part, AliMCEvent* const mcevent) {
644  Bool_t physprim = false;
645  const AliAODMCParticle *aodmc = dynamic_cast<const AliAODMCParticle *>(part);
646  if(aodmc){
647  physprim = aodmc->IsPhysicalPrimary();
648  } else {
649  physprim = mcevent->IsPhysicalPrimary(part->GetLabel());
650  }
651  return physprim;
652 }
653 
654 } /* namespace EMCalTriggerPtAnalysis */
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
TList * list
void CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Container class for histograms for the high- charged particle analysis.
Main data structure storing all relevant information of EMCAL/DCAL trigger patches.
Declarartion of class AliEMCalHistoContainer.
Bool_t IsPhysicalPrimary(const AliVParticle *const part, AliMCEvent *const mcevent)
Bool_t PythiaInfoFromFile(const char *currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard) const
void CreateTProfile(const char *name, const char *title, int nbinsX, double xmin, double xmax, Option_t *opt="")
Unit test class for charged particle distributions (MC case)
TFile * file
Class to make array of trigger patch objects in AOD/ESD events.
void FillTH1(const char *hname, double x, double weight=1.)
void FillProfile(const char *name, double x, double y, double weight=1.)