AliAODEvent::GetHeader() returns AliVHeader
[u/mrichter/AliRoot.git] / PWGJE / AliAnalysisTaskJetServices.cxx
CommitLineData
d3a3f33d 1
bf7b8731 2// **************************************
3// Task used for the correction of determiantion of reconstructed jet spectra
4// Compares input (gen) and output (rec) jets
5// *******************************************
6
7
8/**************************************************************************
9 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
10 * *
11 * Author: The ALICE Off-line Project. *
12 * Contributors are mentioned in the code where appropriate. *
13 * *
14 * Permission to use, copy, modify and distribute this software and its *
15 * documentation strictly for non-commercial purposes is hereby granted *
16 * without fee, provided that the above copyright notice appears in all *
17 * copies and that both the copyright notice and this permission notice *
18 * appear in the supporting documentation. The authors make no claims *
19 * about the suitability of this software for any purpose. It is *
20 * provided "as is" without express or implied warranty. *
21 **************************************************************************/
22
23
24#include <TROOT.h>
25#include <TRandom.h>
b1cd90b5 26#include <TString.h>
bf7b8731 27#include <TSystem.h>
28#include <TInterpreter.h>
29#include <TChain.h>
30#include <TFile.h>
31#include <TKey.h>
32#include <TH1F.h>
33#include <TH2F.h>
34#include <TH3F.h>
35#include <TProfile.h>
36#include <TList.h>
37#include <TLorentzVector.h>
40445651 38#include <TRandom3.h>
bf7b8731 39#include <TClonesArray.h>
c0997643 40#include "TDatabasePDG.h"
bf7b8731 41
42#include "AliAnalysisTaskJetServices.h"
b6dd6ad2 43#include "AliCentrality.h"
fe669ac6 44#include "AliAnalysisDataContainer.h"
45#include "AliAnalysisDataSlot.h"
bf7b8731 46#include "AliAnalysisManager.h"
47#include "AliJetFinder.h"
48#include "AliJetHeader.h"
49#include "AliJetReader.h"
50#include "AliJetReaderHeader.h"
51#include "AliUA1JetHeaderV1.h"
52#include "AliESDEvent.h"
53#include "AliAODEvent.h"
54#include "AliAODHandler.h"
55#include "AliAODTrack.h"
2eded560 56#include "AliAODVZERO.h"
bf7b8731 57#include "AliAODJet.h"
58#include "AliAODMCParticle.h"
59#include "AliMCEventHandler.h"
60#include "AliMCEvent.h"
61#include "AliStack.h"
62#include "AliGenPythiaEventHeader.h"
63#include "AliJetKineReaderHeader.h"
64#include "AliGenCocktailEventHeader.h"
65#include "AliInputEventHandler.h"
fe669ac6 66#include "AliPhysicsSelection.h"
3493c3a9 67#include "AliTriggerAnalysis.h"
bf7b8731 68
69#include "AliAnalysisHelperJetTasks.h"
70
71ClassImp(AliAnalysisTaskJetServices)
72
2eded560 73AliAODVZERO* AliAnalysisTaskJetServices::fgAODVZERO = NULL;
f4132e7d 74AliAODHeader* AliAnalysisTaskJetServices::fgAODHeader = NULL;
80a790fd 75TClonesArray* AliAnalysisTaskJetServices::fgAODVertices = NULL;
f4132e7d 76
40445651 77AliAnalysisTaskJetServices::AliAnalysisTaskJetServices():
78 AliAnalysisTaskSE(),
bf7b8731 79 fUseAODInput(kFALSE),
cc0649e4 80 fUsePhysicsSelection(kFALSE),
3493c3a9 81 fMC(kFALSE),
f4132e7d 82 fFilterAODCollisions(kFALSE),
39e7e8ab 83 fPhysicsSelectionFlag(AliVEvent::kMB),
f2dd0695 84 fSelectionInfoESD(0),
6d597ca2 85 fEventCutInfoESD(0),
40445651 86 fFilterMask(0),
2eded560 87 fRPMethod(0),
21df3cc7 88 fCollisionType(kPbPb),
3d1dba29 89 fNTrigger(0),
bf7b8731 90 fAvgTrials(1),
c0997643 91 fVtxXMean(0),
92 fVtxYMean(0),
93 fVtxZMean(0),
94 fVtxRCut(1.),
95 fVtxZCut(8.),
e5ee1c19 96 fVtxDeltaZCut(0.1),
df7848fc 97 fPtMinCosmic(5.),
98 fRIsolMinCosmic(3.),
99 fMaxCosmicAngle(0.01),
40445651 100 fCentrality(101),
101 fTrackRecEtaWindow(0.9),
102 fMinTrackPt(0.15),
103 fRPAngle(0),
2eded560 104 fPsiVZEROA(0),
105 fPsiVZEROC(0),
3d1dba29 106 fTriggerBit(0x0),
40445651 107 fRandomizer(0),
f4132e7d 108 fNonStdFile(""),
3d1dba29 109 fTriggerName(0x0),
bf7b8731 110 fh1Xsec(0x0),
111 fh1Trials(0x0),
432abb2b 112 fh1AvgTrials(0x0),
bf7b8731 113 fh1PtHard(0x0),
114 fh1PtHardTrials(0x0),
f2dd0695 115 fh1SelectionInfoESD(0x0),
6d597ca2 116 fh1EventCutInfoESD(0),
742ee86c 117 fh1CentralityESD(0),
118 fh1Centrality(0),
40445651 119 fh1RP(0),
e5697f4d 120 fh2ReducedTrigger(0),
3d1dba29 121 fh2CentralityTriggerESD(0),
122 fh2CentralityTrigger(0),
bf7b8731 123 fh2TriggerCount(0x0),
124 fh2ESDTriggerCount(0x0),
125 fh2TriggerVtx(0x0),
126 fh2ESDTriggerVtx(0x0),
b5a3f310 127 fh2ESDTriggerRun(0x0),
128 fh2VtxXY(0x0),
df7848fc 129 fh1NCosmicsPerEvent(0x0),
2eded560 130 fp1RPXA(0x0),
131 fp1RPYA(0x0),
132 fp1RPXC(0x0),
133 fp1RPYC(0x0),
b360e4be 134 fp1CalibRPXA(0x0),
135 fp1CalibRPYA(0x0),
136 fp1CalibRPXC(0x0),
137 fp1CalibRPYC(0x0),
2eded560 138 fh2RPAC(0x0),
139 fh2RPAT(0x0),
140 fh2RPCT(0x0),
141 fh2XYA(0x0),
142 fh2XYC(0x0),
40445651 143 fh2RPCentrality(0x0),
2eded560 144 fh2RPACentrality(0x0),
145 fh2RPCCentrality(0x0),
3493c3a9 146 fTriggerAnalysis(0x0),
bf7b8731 147 fHistList(0x0)
148{
b5a3f310 149 fRunRange[0] = fRunRange[1] = 0;
bf7b8731 150}
151
152AliAnalysisTaskJetServices::AliAnalysisTaskJetServices(const char* name):
153 AliAnalysisTaskSE(name),
154 fUseAODInput(kFALSE),
cc0649e4 155 fUsePhysicsSelection(kFALSE),
3493c3a9 156 fMC(kFALSE),
f4132e7d 157 fFilterAODCollisions(kFALSE),
39e7e8ab 158 fPhysicsSelectionFlag(AliVEvent::kMB),
f2dd0695 159 fSelectionInfoESD(0),
6d597ca2 160 fEventCutInfoESD(0),
40445651 161 fFilterMask(0),
2eded560 162 fRPMethod(0),
21df3cc7 163 fCollisionType(kPbPb),
3d1dba29 164 fNTrigger(0),
bf7b8731 165 fAvgTrials(1),
c0997643 166 fVtxXMean(0),
167 fVtxYMean(0),
168 fVtxZMean(0),
169 fVtxRCut(1.),
170 fVtxZCut(8.),
e5ee1c19 171 fVtxDeltaZCut(0.1),
df7848fc 172 fPtMinCosmic(5.),
173 fRIsolMinCosmic(3.),
174 fMaxCosmicAngle(0.01),
40445651 175 fCentrality(101),
176 fTrackRecEtaWindow(0.9),
177 fMinTrackPt(0.15),
178 fRPAngle(0),
2eded560 179 fPsiVZEROA(0),
180 fPsiVZEROC(0),
3d1dba29 181 fTriggerBit(0x0),
40445651 182 fRandomizer(0),
f4132e7d 183 fNonStdFile(""),
3d1dba29 184 fTriggerName(0x0),
bf7b8731 185 fh1Xsec(0x0),
186 fh1Trials(0x0),
432abb2b 187 fh1AvgTrials(0x0),
bf7b8731 188 fh1PtHard(0x0),
189 fh1PtHardTrials(0x0),
f2dd0695 190 fh1SelectionInfoESD(0x0),
6d597ca2 191 fh1EventCutInfoESD(0),
742ee86c 192 fh1CentralityESD(0),
193 fh1Centrality(0),
40445651 194 fh1RP(0),
e5697f4d 195 fh2ReducedTrigger(0),
3d1dba29 196 fh2CentralityTriggerESD(0),
197 fh2CentralityTrigger(0),
bf7b8731 198 fh2TriggerCount(0x0),
199 fh2ESDTriggerCount(0x0),
200 fh2TriggerVtx(0x0),
201 fh2ESDTriggerVtx(0x0),
b5a3f310 202 fh2ESDTriggerRun(0x0),
203 fh2VtxXY(0x0),
df7848fc 204 fh1NCosmicsPerEvent(0x0),
2eded560 205 fp1RPXA(0x0),
206 fp1RPYA(0x0),
207 fp1RPXC(0x0),
208 fp1RPYC(0x0),
b360e4be 209 fp1CalibRPXA(0x0),
210 fp1CalibRPYA(0x0),
211 fp1CalibRPXC(0x0),
212 fp1CalibRPYC(0x0),
2eded560 213 fh2RPAC(0x0),
214 fh2RPAT(0x0),
215 fh2RPCT(0x0),
216 fh2XYA(0x0),
217 fh2XYC(0x0),
40445651 218 fh2RPCentrality(0x0),
2eded560 219 fh2RPACentrality(0x0),
220 fh2RPCCentrality(0x0),
3493c3a9 221 fTriggerAnalysis(0x0),
bf7b8731 222 fHistList(0x0)
223{
b5a3f310 224 fRunRange[0] = fRunRange[1] = 0;
bf7b8731 225 DefineOutput(1,TList::Class());
226}
227
228
229
230Bool_t AliAnalysisTaskJetServices::Notify()
231{
232 //
233 // Implemented Notify() to read the cross sections
234 // and number of trials from pyxsec.root
235 //
236
237 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
238 Float_t xsection = 0;
239 Float_t ftrials = 1;
240
241 fAvgTrials = 1;
242 if(tree){
243 TFile *curfile = tree->GetCurrentFile();
244 if (!curfile) {
245 Error("Notify","No current file");
246 return kFALSE;
247 }
248 if(!fh1Xsec||!fh1Trials){
249 Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
250 return kFALSE;
251 }
252 AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
253 fh1Xsec->Fill("<#sigma>",xsection);
254 // construct a poor man average trials
255 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
256 if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
432abb2b 257 fh1Trials->Fill("#sum{avg ntrials}",ftrials);
bf7b8731 258 }
259 return kTRUE;
260}
261
262void AliAnalysisTaskJetServices::UserCreateOutputObjects()
263{
264
265 //
266 // Create the output container
267 //
268
40445651 269 fRandomizer = new TRandom3(0);
bf7b8731 270 if (fDebug > 1) printf("AnalysisTaskJetServices::UserCreateOutputObjects() \n");
271
272 OpenFile(1);
273 if(!fHistList)fHistList = new TList();
40440b28 274 fHistList->SetOwner();
275 PostData(1, fHistList); // post data in any case once
bf7b8731 276
277 Bool_t oldStatus = TH1::AddDirectoryStatus();
278 TH1::AddDirectory(kFALSE);
bf7b8731 279 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
280 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
281 fHistList->Add(fh1Xsec);
282
283 fh1Trials = new TH1F("fh1Trials","trials root file",1,0,1);
284 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
285 fHistList->Add(fh1Trials);
286
432abb2b 287 fh1AvgTrials = new TH1F("fh1AvgTrials","avg trials root file",1,0,1);
288 fh1AvgTrials->GetXaxis()->SetBinLabel(1,"#sum{avg ntrials}");
289 fHistList->Add(fh1AvgTrials);
290
66fd4b2b 291 const Int_t nBinPt = 125;
bf7b8731 292 Double_t binLimitsPt[nBinPt+1];
293 for(Int_t iPt = 0;iPt <= nBinPt;iPt++){
294 if(iPt == 0){
295 binLimitsPt[iPt] = 0.0;
296 }
297 else {// 1.0
66fd4b2b 298 binLimitsPt[iPt] = binLimitsPt[iPt-1] + 2.;
bf7b8731 299 }
300 }
301
742ee86c 302
eb1739b7 303 fh1CentralityESD = new TH1F("fh1CentralityESD","cent",103,-1,102);
742ee86c 304 fHistList->Add(fh1CentralityESD);
305
eb1739b7 306 fh1Centrality = new TH1F("fh1Centrality","cent",103,-1,102);
742ee86c 307 fHistList->Add(fh1Centrality);
308
40445651 309 fh1RP = new TH1F("fh1RP","RP;#Psi",440, -1.*TMath::Pi(), 2.*TMath::Pi());
310 fHistList->Add(fh1RP);
311
e5697f4d 312 fh2ReducedTrigger = new TH2F("fh2ReducedTrigger","red trigger;red trigger;centrality",1<<fNTrigger,-0.5,(1<<fNTrigger)-0.5,100,-0.5,99.5);
313 fHistList->Add(fh2ReducedTrigger);
3d1dba29 314
315 fh2CentralityTriggerESD = new TH2F("fh2CentralityTriggerESD",";cent;trigger no",103,-1,102,fNTrigger,-0.5,fNTrigger-0.5);
316 fHistList->Add(fh2CentralityTriggerESD);
317
318 fh2CentralityTrigger = new TH2F("fh2CentralityTrigger",";cent;trigger no",103,-1,102,fNTrigger,-0.5,fNTrigger-0.5);
319 fHistList->Add(fh2CentralityTrigger);
320
321 for(int i = 0;i<fNTrigger;++i){
5f5ae725 322 fh2CentralityTriggerESD->GetYaxis()->SetBinLabel(i+1,fTriggerName[i].Data());
323 fh2CentralityTrigger->GetYaxis()->SetBinLabel(i+1,fTriggerName[i].Data());
3d1dba29 324 }
325
326
d95fc15a 327 fh2TriggerCount = new TH2F("fh2TriggerCount",";Trigger No.;constrained;Count",6,-0.5,5.5,kConstraints,-0.5,kConstraints-0.5);
bf7b8731 328 fHistList->Add(fh2TriggerCount);
329
d95fc15a 330 fh2ESDTriggerCount = new TH2F("fh2ESDTriggerCount",";Trigger No.;constrained;Count",6,-0.5,5.5,kConstraints,-0.5,kConstraints-0.5);
bf7b8731 331 fHistList->Add(fh2ESDTriggerCount);
d95fc15a 332 const Int_t nBins = 6*kConstraints;
57ca1193 333 fh2TriggerVtx = new TH2F("fh2TriggerVtx",";Constraint No. * (trig no+1);Vtx (cm);Count",nBins,-0.5,nBins-0.5,400,-20.,20.);
bf7b8731 334 fHistList->Add(fh2TriggerVtx);
335
57ca1193 336 fh2ESDTriggerVtx = new TH2F("fh2ESDTriggerVtx",";Constraint No.* (trg no+1);Vtx (cm);Count",nBins,-0.5,nBins-0.5,400,-20.,20.);
bf7b8731 337 fHistList->Add(fh2ESDTriggerVtx);
338
339
340 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",nBinPt,binLimitsPt);
341 fHistList->Add(fh1PtHard);
342 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",nBinPt,binLimitsPt);
343 fHistList->Add(fh1PtHardTrials);
f2dd0695 344 fh1SelectionInfoESD = new TH1F("fh1SelectionInfoESD","Bit Masks that satisfy the selection info",
345 AliAnalysisHelperJetTasks::kTotalSelections,0.5,AliAnalysisHelperJetTasks::kTotalSelections+0.5);
f2dd0695 346 fHistList->Add(fh1SelectionInfoESD);
6d597ca2 347
348 fh1EventCutInfoESD = new TH1F("fh1EventCutInfoESD","Bit Masks that for the events after each step of cuts",
349 kTotalEventCuts,0.5,kTotalEventCuts+0.5);
350 fHistList->Add(fh1EventCutInfoESD);
351
b5a3f310 352 // 3 decisions, 0 trigger X, X + SPD vertex, X + SPD vertex in range
353 // 3 triggers BB BE/EB EE
354
f4132e7d 355 fh2ESDTriggerRun = new TH2F("fh2ESDTriggerRun","Eventclass vs run number:run;trigger",(Int_t)(1+fRunRange[1]-fRunRange[0]),fRunRange[0]-0.5,fRunRange[1]+0.5,10,-0.5,9.5);
b5a3f310 356 fHistList->Add(fh2ESDTriggerRun);
357
358 fh2VtxXY = new TH2F("fh2VtxXY","Beam Spot all INT triggered events;x (cm);y (cm)",160,-10,10,160,-10,10);
359 fHistList->Add(fh2VtxXY);
bf7b8731 360 // =========== Switch on Sumw2 for all histos ===========
361 for (Int_t i=0; i<fHistList->GetEntries(); ++i) {
362 TH1 *h1 = dynamic_cast<TH1*>(fHistList->At(i));
363 if (h1){
364 h1->Sumw2();
365 continue;
366 }
367 THnSparse *hn = dynamic_cast<THnSparse*>(fHistList->At(i));
368 if(hn)hn->Sumw2();
369 }
370
df7848fc 371 fh1NCosmicsPerEvent = new TH1F("fh1NCosmicsPerEvent","Number of cosmic candidates per event",10,0.,10.);
2eded560 372 fHistList->Add(fh1NCosmicsPerEvent),
40445651 373
40445651 374
2eded560 375 fh2RPCentrality = new TH2F("fh2RPCentrality" ,"Reaction Plane Angle from tracks" , 20, 0.,100., 180, 0, TMath::Pi());
40445651 376 fHistList->Add(fh2RPCentrality);
377
40445651 378
2eded560 379 fh2RPACentrality = new TH2F("fh2RPACentrality" ,"Reaction Plane Angle from vzero A" , 20, 0.,100., 180, 0, TMath::Pi());
380 fHistList->Add(fh2RPACentrality);
40445651 381
2eded560 382 fh2RPCCentrality = new TH2F("fh2RPCCentrality" ,"Reaction Plane Angle from vzero C" , 20, 0.,100., 180, 0, TMath::Pi());
383 fHistList->Add(fh2RPCCentrality);
40445651 384
2eded560 385 fh2RPAC = new TH2F("fh2RPAC" ,"Reaction Plane Angle vzero a vs c" , 180, 0, TMath::Pi(), 180, 0, TMath::Pi());
386 fHistList->Add( fh2RPAC);
40445651 387
2eded560 388 fh2RPAT = new TH2F("fh2RPAT" ,"Reaction Plane Angle vzero a vs tracks" , 180, 0, TMath::Pi(), 180, 0, TMath::Pi());
389 fHistList->Add( fh2RPAT);
390
391 fh2RPCT = new TH2F("fh2RPCT" ,"Reaction Plane Angle vzero c vs tracks" , 180, 0, TMath::Pi(), 180, 0, TMath::Pi());
392 fHistList->Add( fh2RPCT);
393
394 fh2XYA = new TH2F("fh2XYA" ,"XY vzeroa ;X;Y;" ,100,-0.3,0.3,100,-0.3,0.3);
395 fHistList->Add(fh2XYA);
396 fh2XYC = new TH2F("fh2XYC" ,"XY vzeroc ;X;Y;" ,100,-0.3,0.3,100,-0.3,0.3);
397 fHistList->Add(fh2XYC);
398
399 // profiles for mean
400 fp1RPXA = new TProfile("fp1RPXA","mean vzeroa x vs run number;run;x",(Int_t)(1+fRunRange[1]-fRunRange[0]),fRunRange[0]-0.5,fRunRange[1]+0.5);
401 fHistList->Add(fp1RPXA);
402
403 fp1RPYA = new TProfile("fp1RPYA","mean vzeroa y vs run number;run;y",(Int_t)(1+fRunRange[1]-fRunRange[0]),fRunRange[0]-0.5,fRunRange[1]+0.5);
404 fHistList->Add(fp1RPYA);
405
406
407 fp1RPXC = new TProfile("fp1RPXC","mean vzeroc x vs run number;run;x",(Int_t)(1+fRunRange[1]-fRunRange[0]),fRunRange[0]-0.5,fRunRange[1]+0.5);
408 fHistList->Add(fp1RPXC);
409
410 fp1RPYC = new TProfile("fp1RPYC","mean vzeroa y vs run number;run;y",(Int_t)(1+fRunRange[1]-fRunRange[0]),fRunRange[0]-0.5,fRunRange[1]+0.5);
411 fHistList->Add(fp1RPYC);
df7848fc 412
bf7b8731 413
414 TH1::AddDirectory(oldStatus);
f4132e7d 415
416 // Add an AOD branch for replication
417 if(fNonStdFile.Length()){
418 if (fDebug > 1) AliInfo("Replicating header");
419 fgAODHeader = new AliAODHeader;
420 AddAODBranch("AliAODHeader",&fgAODHeader,fNonStdFile.Data());
2eded560 421 if (fDebug > 1) AliInfo("Replicating vzeros");
422 fgAODVZERO = new AliAODVZERO;
423 AddAODBranch("AliAODVZERO",&fgAODVZERO,fNonStdFile.Data());
070f06d2 424 if (fDebug > 1) AliInfo("Replicating primary vertices");
070f06d2 425 fgAODVertices = new TClonesArray("AliAODVertex",3);
80a790fd 426 fgAODVertices->SetName("vertices");
070f06d2 427 AddAODBranch("TClonesArray",&fgAODVertices,fNonStdFile.Data());
f4132e7d 428 }
bf7b8731 429}
430
431void AliAnalysisTaskJetServices::Init()
432{
433 //
434 // Initialization
435 //
bf7b8731 436 if (fDebug > 1) printf("AnalysisTaskJetServices::Init() \n");
437
438}
439
440void AliAnalysisTaskJetServices::UserExec(Option_t */*option*/)
441{
442
443 //
444 // Execute analysis for current event
445 //
446
447 AliAODEvent *aod = 0;
448 AliESDEvent *esd = 0;
f4132e7d 449
450
3493c3a9 451
bf7b8731 452 AliAnalysisHelperJetTasks::Selected(kTRUE,kFALSE); // set slection to false
f4132e7d 453 AliAnalysisHelperJetTasks::EventClass(kTRUE,0);
40445651 454 AliAnalysisHelperJetTasks::ReactionPlane(kTRUE,0); // set slection to false
f2dd0695 455 fSelectionInfoESD = 0; // reset
6d597ca2 456 fEventCutInfoESD = 0; // reset
f2dd0695 457 AliAnalysisHelperJetTasks::SelectInfo(kTRUE,fSelectionInfoESD); // set slection to false
458
bf7b8731 459
f4132e7d 460 static AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
461
462
463
bf7b8731 464 if(fUseAODInput){
465 aod = dynamic_cast<AliAODEvent*>(InputEvent());
466 if(!aod){
467 Printf("%s:%d AODEvent not found in Input Manager %d",(char*)__FILE__,__LINE__,fUseAODInput);
468 return;
469 }
470 // fethc the header
471 }
472 else{
473 // assume that the AOD is in the general output...
474 aod = AODEvent();
475 if(!aod){
476 Printf("%s:%d AODEvent not found in the Output",(char*)__FILE__,__LINE__);
477 return;
478 }
479 esd = dynamic_cast<AliESDEvent*>(InputEvent());
480 }
ffab794c 481 if(aod&&fDebug>2){
482 aod->Print();
483 Printf("Vertices %d",aod->GetNumberOfVertices());
484 Printf("tracks %d",aod->GetNumberOfTracks());
485 Printf("jets %d",aod->GetNJets());
486 }
3493c3a9 487 fSelectionInfoESD |= kNoEventCut;
488 fEventCutInfoESD |= kNoEventCut;
6d597ca2 489
490 Bool_t esdVtxValid = false;
491 Bool_t esdVtxIn = false;
492 Bool_t aodVtxValid = false;
493 Bool_t aodVtxIn = false;
f2dd0695 494
c0997643 495 // Apply additional constraints
6d597ca2 496 Bool_t esdEventSelected = IsEventSelected(esd);
497 Bool_t esdEventPileUp = IsEventPileUp(esd);
498 Bool_t esdEventCosmic = IsEventCosmic(esd);
f2dd0695 499
6d597ca2 500 Bool_t aodEventSelected = IsEventSelected(aod);
f2dd0695 501
39e7e8ab 502 Bool_t physicsSelection = ((fInputHandler->IsEventSelected())&fPhysicsSelectionFlag);
20b170eb 503
f4132e7d 504
6d597ca2 505 fEventCutInfoESD |= kPhysicsSelectionCut; // other alreay set via IsEventSelected
e5ee1c19 506 fh1EventCutInfoESD->Fill (fEventCutInfoESD);
6d597ca2 507
f2dd0695 508 if(esdEventSelected) fSelectionInfoESD |= AliAnalysisHelperJetTasks::kVertexIn;
509 if(esdEventPileUp) fSelectionInfoESD |= AliAnalysisHelperJetTasks::kIsPileUp;
510 if(esdEventCosmic) fSelectionInfoESD |= AliAnalysisHelperJetTasks::kIsCosmic;
511 if(physicsSelection) fSelectionInfoESD |= AliAnalysisHelperJetTasks::kPhysicsSelection;
512
513
3d1dba29 514
515
f2dd0695 516 // here we have all selection information, fill histogram
517 for(unsigned int i = 1;i<(UInt_t)fh1SelectionInfoESD->GetNbinsX();i++){
518 if((i&fSelectionInfoESD)==i)fh1SelectionInfoESD->Fill(i);
519 }
6d597ca2 520 AliAnalysisHelperJetTasks::SelectInfo(kTRUE,fSelectionInfoESD);
c0997643 521
522 if(esd&&aod&&fDebug){
523 if(esdEventSelected&&!aodEventSelected){
524 Printf("%s:%d Different Selection for ESD and AOD",(char*)__FILE__,__LINE__);
525 const AliESDVertex *vtxESD = esd->GetPrimaryVertex();
bf7b8731 526 const AliAODVertex *vtxAOD = aod->GetPrimaryVertex();
c0997643 527 Printf("ESD Vtx %s %s %d",vtxESD->GetName(),vtxESD->GetTitle(),vtxESD->GetNContributors());
528 vtxESD->Print();
529 Printf("AOD Vtx %s %s %d",vtxAOD->GetName(),vtxAOD->GetTitle(),vtxAOD->GetNContributors());
530 vtxAOD->Print();
bf7b8731 531 }
c0997643 532 }
533
534 // loop over all possible triggers for esd
535
aa3ba8d2 536 Float_t cent = 0;
21df3cc7 537 if(fCollisionType==kPbPb){
0a918d8d 538 if(aod)cent = ((AliVAODHeader*)aod->GetHeader())->GetCentrality();
21df3cc7 539 if(fDebug)Printf("%s:%d %3.3f",(char*)__FILE__,__LINE__,cent);
540 if(cent<0)cent = 101;
541 }
40445651 542 fCentrality = cent;
543 fRPAngle = 0;
2eded560 544 fPsiVZEROA = 0;
545 fPsiVZEROC = 0;
546
40445651 547
c0997643 548 if(esd){
549 const AliESDVertex *vtxESD = esd->GetPrimaryVertex();
e5ee1c19 550 const AliESDVertex *vtxESDSPD = esd->GetPrimaryVertexSPD();
f4132e7d 551 esdVtxValid = IsVertexValid(vtxESD);
e5ee1c19 552 esdVtxIn = IsVertexIn(vtxESD,vtxESDSPD);
d8f21f85 553 if(aodH&&physicsSelection&&fFilterAODCollisions&&aod){
40440b28 554 if(fDebug)Printf("%s:%d Centrality %3.3f vtxin %d",(char*)__FILE__,__LINE__,cent,esdVtxIn);
d8f21f85 555 if(cent<=80&&esdVtxIn){
556 aodH->SetFillAOD(kTRUE);
c9931fda 557 aodH->SetFillExtension(kTRUE);
d8f21f85 558 }
559 }
560
561
57ca1193 562 Float_t zvtx = vtxESD->GetZ();
f4132e7d 563 Int_t iCl = GetEventClass(esd);
564 AliAnalysisHelperJetTasks::EventClass(kTRUE,iCl);
565 Bool_t cand = physicsSelection;
566
ffab794c 567 if(fDebug)Printf("%s:%d %d %d %d Icl %d",(char*)__FILE__,__LINE__,esdVtxValid,esdVtxIn,cand,iCl);
d8f21f85 568 fh2ESDTriggerCount->Fill(0.,kAllTriggered);
569 fh2ESDTriggerCount->Fill(iCl,kAllTriggered);
f4132e7d 570 if(cand){
571 fh2ESDTriggerCount->Fill(0.,kSelectedALICE);
572 fh2ESDTriggerCount->Fill(iCl,kSelectedALICE);
573 fh2ESDTriggerVtx->Fill(kSelectedALICE*(iCl+1),zvtx);
574 }
575 // if(!fUsePhysicsSelection)cand = AliAnalysisHelperJetTasks::IsTriggerFired(esd,AliAnalysisHelperJetTasks::kMB1);
576 if(esdVtxValid){
577 fh2ESDTriggerCount->Fill(0.,kTriggeredVertex);
578 fh2ESDTriggerCount->Fill(iCl,kTriggeredVertex);
579 fh2ESDTriggerVtx->Fill(iCl,zvtx);
580 if(esdVtxIn){
581 fh2ESDTriggerCount->Fill(0.,kTriggeredVertexIn);
582 fh2ESDTriggerCount->Fill(iCl,kTriggeredVertexIn);
583 fh2ESDTriggerVtx->Fill(kTriggeredVertexIn*(iCl+1),zvtx);
bf7b8731 584 }
f4132e7d 585 if(cand){
586 fh2ESDTriggerCount->Fill(0.,kSelectedALICEVertexValid);
587 fh2ESDTriggerCount->Fill(iCl,kSelectedALICEVertexValid);
588 fh2ESDTriggerVtx->Fill(kSelectedALICEVertexValid*(iCl+1),zvtx);
c0997643 589 }
bf7b8731 590 }
742ee86c 591
c08c3ad2 592 if(cand&&esdVtxIn&&iCl<5){
f4132e7d 593 fh2ESDTriggerCount->Fill(0.,kSelectedALICEVertexIn);
594 fh2ESDTriggerCount->Fill(iCl,kSelectedALICEVertexIn);
595 fh2ESDTriggerVtx->Fill(kSelectedALICEVertexIn*(iCl+1),zvtx);
596 fh2ESDTriggerVtx->Fill(kSelected*(iCl+1),zvtx);
597 fh2ESDTriggerCount->Fill(iCl,kSelected);
598 fh2ESDTriggerCount->Fill(0.,kSelected);
599 AliAnalysisHelperJetTasks::Selected(kTRUE,kTRUE);// select this event
742ee86c 600 if(esd->GetCentrality()){
601 Float_t tmpCent = 100;
602 tmpCent = esd->GetCentrality()->GetCentralityPercentile("V0M");
f5c02503 603 if(tmpCent<0)tmpCent = 101;
742ee86c 604 fh1CentralityESD->Fill(tmpCent);
3d1dba29 605 UInt_t ir = 0;
606 for(int it = 0;it<fNTrigger;++it){
607 if(fInputHandler->IsEventSelected()&fTriggerBit[it]){
608 fh2CentralityTriggerESD->Fill(tmpCent,it);
609 ir |= (1<<it);
610 }
611 }
e5697f4d 612 fh2ReducedTrigger->Fill(ir,cent);
742ee86c 613 }
f4132e7d 614 }
bf7b8731 615 }
616
ffab794c 617
618
c0997643 619 if(aod){
620 const AliAODVertex *vtxAOD = aod->GetPrimaryVertex();
e5ee1c19 621 const AliAODVertex *vtxAODSPD = aod->GetPrimaryVertexSPD();
6d597ca2 622 aodVtxValid = IsVertexValid(vtxAOD);
e5ee1c19 623 aodVtxIn = IsVertexIn(vtxAOD,vtxAODSPD);
ffab794c 624 Float_t zvtx = vtxAOD->GetZ();
625 Int_t iCl = GetEventClass(aod);
626 AliAnalysisHelperJetTasks::EventClass(kTRUE,iCl);
0a918d8d 627 Bool_t cand = ((AliVAODHeader*)aod->GetHeader())->GetOfflineTrigger()&fPhysicsSelectionFlag;
628 if(fDebug)Printf("%s:%d AOD selection %d %d",(char*)__FILE__,__LINE__,cand,((AliVAODHeader*)aod->GetHeader())->GetOfflineTrigger());
d8f21f85 629 fh2TriggerCount->Fill(0.,kAllTriggered);
630 fh2TriggerCount->Fill(iCl,kAllTriggered);
ffab794c 631 if(cand){
632 fh2TriggerCount->Fill(0.,kSelectedALICE);
633 fh2TriggerCount->Fill(iCl,kSelectedALICE);
634 fh2TriggerVtx->Fill(kSelectedALICE*(iCl+1),zvtx);
635 }
636 if(aodVtxValid){
637 fh2TriggerCount->Fill(0.,kTriggeredVertex);
638 fh2TriggerCount->Fill(iCl,kTriggeredVertex);
639 fh2TriggerVtx->Fill(iCl,zvtx);
640 if(aodVtxIn){
641 fh2TriggerCount->Fill(0.,kTriggeredVertexIn);
642 fh2TriggerCount->Fill(iCl,kTriggeredVertexIn);
643 fh2TriggerVtx->Fill(kTriggeredVertexIn*(iCl+1),zvtx);
c0997643 644 }
ffab794c 645 if(cand){
646 fh2TriggerCount->Fill(0.,kSelectedALICEVertexValid);
647 fh2TriggerCount->Fill(iCl,kSelectedALICEVertexValid);
648 fh2TriggerVtx->Fill(kSelectedALICEVertexValid*(iCl+1),zvtx);
649 }
650 }
90a006c2 651
c08c3ad2 652 if(cand&&aodVtxIn&&iCl<5){
ffab794c 653 fh2TriggerCount->Fill(0.,kSelectedALICEVertexIn);
654 fh2TriggerCount->Fill(iCl,kSelectedALICEVertexIn);
655 fh2TriggerVtx->Fill(kSelectedALICEVertexIn*(iCl+1),zvtx);
656 fh2TriggerVtx->Fill(kSelected*(iCl+1),zvtx);
657 fh2TriggerCount->Fill(iCl,kSelected);
658 fh2TriggerCount->Fill(0.,kSelected);
742ee86c 659 fh1Centrality->Fill(cent);
3d1dba29 660 for(int it = 0;it<fNTrigger;++it){
661 if(fInputHandler->IsEventSelected()&fTriggerBit[it])fh2CentralityTrigger->Fill(cent,it);
662 }
742ee86c 663 AliAnalysisHelperJetTasks::Selected(kTRUE,kTRUE);// select this event
90a006c2 664 if(aodH&&cand&&fFilterAODCollisions&&!esd){
21df3cc7 665 if(fCentrality<=80&&aodVtxIn){
90a006c2 666 aodH->SetFillAOD(kTRUE);
667 aodH->SetFillExtension(kTRUE);
668 }
669 }
670
40445651 671 TList recTracks;
672 GetListOfTracks(&recTracks);
2eded560 673 CalculateReactionPlaneAngleVZERO(aod);
0a918d8d 674 fRPAngle = ((AliVAODHeader*)aod->GetHeader())->GetEventplane();
d3a3f33d 675 fh1RP->Fill(fRPAngle);
676 fh2RPCentrality->Fill(fCentrality,fRPAngle);
2eded560 677 fh2RPACentrality->Fill(fCentrality,fPsiVZEROA);
678 fh2RPCCentrality->Fill(fCentrality,fPsiVZEROC);
679 fh2RPAC->Fill(fPsiVZEROA,fPsiVZEROC);
680 fh2RPAT->Fill(fPsiVZEROA,fRPAngle);
681 fh2RPCT->Fill(fPsiVZEROC,fRPAngle);
682 if(fRPMethod==kRPTracks)AliAnalysisHelperJetTasks::ReactionPlane(kTRUE,fRPAngle); // set slection to false
683 else if(fRPMethod==kRPVZEROA)AliAnalysisHelperJetTasks::ReactionPlane(kTRUE,fPsiVZEROA); // set slection to false
684 else if(fRPMethod==kRPVZEROC)AliAnalysisHelperJetTasks::ReactionPlane(kTRUE,fPsiVZEROA); // set slection to false
d3a3f33d 685
21df3cc7 686 if(fUseAODInput&&fCentrality<=80){
221d8aa9 687 if(fFilterAODCollisions&&aod){
742ee86c 688 aodH->SetFillAOD(kTRUE);
eedaa08e 689 aodH->SetFillExtension(kTRUE);
221d8aa9 690 }
691 }
221d8aa9 692 }
c0997643 693 }
bf7b8731 694
695 if (fDebug > 1)printf(" AliAnalysisTaskJetServices: Analysing event # %5d\n", (Int_t) fEntry);
696
697
698 Double_t ptHard = 0;
699 Double_t nTrials = 1; // Trials for MC trigger
700
432abb2b 701 fh1AvgTrials->Fill("#sum{avg ntrials}",fAvgTrials);
bf7b8731 702 AliMCEvent* mcEvent = MCEvent();
703 // AliStack *pStack = 0;
704 if(mcEvent){
705 AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
706 if(pythiaGenHeader){
707 nTrials = pythiaGenHeader->Trials();
708 ptHard = pythiaGenHeader->GetPtHard();
709 int iProcessType = pythiaGenHeader->ProcessType();
710 // 11 f+f -> f+f
711 // 12 f+barf -> f+barf
712 // 13 f+barf -> g+g
713 // 28 f+g -> f+g
714 // 53 g+g -> f+barf
715 // 68 g+g -> g+g
716 if (fDebug > 10)Printf("%d iProcessType %d",__LINE__, iProcessType);
717 if(fDebug>20)AliAnalysisHelperJetTasks::PrintStack(mcEvent);
718 fh1PtHard->Fill(ptHard);
719 fh1PtHardTrials->Fill(ptHard,nTrials);
720
721 }// if pythia gen header
722 }
723
724 // trigger selection
725
f4132e7d 726 // replication of
727 if(fNonStdFile.Length()&&aod){
728 if (fgAODHeader){
729 *fgAODHeader = *(dynamic_cast<AliAODHeader*>(aod->GetHeader()));
2eded560 730 Double_t q[kRPMethods] = {fRPAngle,fPsiVZEROA,fPsiVZEROC};
731 fgAODHeader->SetQTheta(q,kRPMethods);
732 }
733 if (fgAODVZERO){
734 *fgAODVZERO = *(dynamic_cast<AliAODVZERO*>(aod->GetVZEROData()));
f4132e7d 735 }
80a790fd 736 if(fgAODVertices){
80a790fd 737 fgAODVertices->Delete();
738 TClonesArray &vertices = *fgAODVertices;
80a790fd 739 const AliAODVertex *vtxAOD = aod->GetPrimaryVertex();
070f06d2 740 // we only use some basic information,
741
742
743 Double_t pos[3];
744 Double_t covVtx[6];
745
746 vtxAOD->GetXYZ(pos); // position
747 vtxAOD->GetCovMatrix(covVtx); //covariance matrix
748 Int_t jVertices = 0;
749 AliAODVertex * vtx = new(vertices[jVertices++])
750 AliAODVertex(pos, covVtx, vtxAOD->GetChi2perNDF(), NULL, -1, AliAODVertex::kPrimary);
751 vtx->SetName(vtxAOD->GetName());
752 vtx->SetTitle(vtxAOD->GetTitle());
070f06d2 753 TString vtitle = vtxAOD->GetTitle();
d8f21f85 754 vtx->SetNContributors(vtxAOD->GetNContributors());
070f06d2 755
756 // Add SPD "main" vertex
757 const AliAODVertex *vtxS = aod->GetPrimaryVertexSPD();
758 vtxS->GetXYZ(pos); // position
759 vtxS->GetCovMatrix(covVtx); //covariance matrix
760 AliAODVertex * mVSPD = new(vertices[jVertices++])
761 AliAODVertex(pos, covVtx, vtxS->GetChi2perNDF(), NULL, -1, AliAODVertex::kMainSPD);
762 mVSPD->SetName(vtxS->GetName());
763 mVSPD->SetTitle(vtxS->GetTitle());
764 mVSPD->SetNContributors(vtxS->GetNContributors());
765
766 // Add tpc only vertex
767 if(esd){
768 const AliESDVertex *vtxT = esd->GetPrimaryVertexTPC();
769 vtxT->GetXYZ(pos); // position
770 vtxT->GetCovMatrix(covVtx); //covariance matrix
771 AliAODVertex * mVTPC = new(vertices[jVertices++])
772 AliAODVertex(pos, covVtx, vtxT->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainTPC);
773 mVTPC->SetName(vtxT->GetName());
774 mVTPC->SetTitle(vtxT->GetTitle());
775 mVTPC->SetNContributors(vtxT->GetNContributors());
776 }
80a790fd 777 }
f4132e7d 778 }
779
bf7b8731 780 PostData(1, fHistList);
80a790fd 781}
bf7b8731 782
6d597ca2 783Bool_t AliAnalysisTaskJetServices::IsEventSelected(const AliESDEvent* esd){
c0997643 784 if(!esd)return kFALSE;
785 const AliESDVertex *vtx = esd->GetPrimaryVertex();
e5ee1c19 786 const AliESDVertex *vtxSPD = esd->GetPrimaryVertex();
787 return IsVertexIn(vtx,vtxSPD); // vertex in calls vertex valid
6d597ca2 788}
789
80a790fd 790AliAnalysisTaskJetServices::~AliAnalysisTaskJetServices(){
791 if(fgAODVertices){
792 fgAODVertices->Delete();
793 delete fgAODVertices;
794 }
40445651 795 delete fRandomizer;
070f06d2 796 if(fgAODHeader)delete fgAODHeader;
b360e4be 797 if(fgAODVZERO)delete fgAODVZERO;
798 delete fp1CalibRPXA;
799 delete fp1CalibRPYA;
800 delete fp1CalibRPXC;
801 delete fp1CalibRPYC;
3d1dba29 802 delete [] fTriggerBit;
803 delete [] fTriggerName;
b360e4be 804
80a790fd 805}
806
6d597ca2 807
b360e4be 808void AliAnalysisTaskJetServices::SetV0Centroids(TProfile *xa,TProfile *ya,TProfile *xc, TProfile *yc){
809
810 if(xa){
811 if(fp1CalibRPXA)delete fp1CalibRPXA;
812 fp1CalibRPXA = (TProfile*)xa->Clone(Form("%sCalib",xa->GetName()));
813 }
814 else{
815 Printf("%s:%d centroid histogram is 0x0",(char*)__FILE__,__LINE__);
816 }
817 if(ya){
818 if(fp1CalibRPYA)delete fp1CalibRPYA;
819 fp1CalibRPYA = (TProfile*)ya->Clone(Form("%sCalib",ya->GetName()));
820 }
821 else{
822 Printf("%s:%d centroid histogram is 0x0",(char*)__FILE__,__LINE__);
823 }
824 if(xc){
825 if(fp1CalibRPXC)delete fp1CalibRPXC;
826 fp1CalibRPXC = (TProfile*)xc->Clone(Form("%sCalib",xc->GetName()));
827 }
828 else{
829 Printf("%s:%d centroid histogram is 0x0",(char*)__FILE__,__LINE__);
830 }
831 if(ya){
832 if(fp1CalibRPYC)delete fp1CalibRPYC;
833 fp1CalibRPYC = (TProfile*)yc->Clone(Form("%sCalib",yc->GetName()));
834 }
835 else{
836 Printf("%s:%d centroid histogram is 0x0",(char*)__FILE__,__LINE__);
837 }
838}
839
6d597ca2 840Bool_t AliAnalysisTaskJetServices::IsEventSelected(const AliAODEvent* aod) const {
841 if(!aod)return kFALSE;
842 const AliAODVertex *vtx = aod->GetPrimaryVertex();
e5ee1c19 843 const AliAODVertex *vtxSPD = aod->GetPrimaryVertexSPD();
844 return IsVertexIn(vtx,vtxSPD); // VertexIn calls VertexValid
6d597ca2 845}
846
847Bool_t AliAnalysisTaskJetServices::IsVertexValid ( const AliESDVertex* vtx) {
848
c0997643 849 // We can check the type of the vertex by its name and title
850 // if vertexer failed title is empty (default c'tor) and ncontributors is 0
851 // vtx name title
852 // Tracks PrimaryVertex VertexerTracksNoConstraint
6d597ca2 853 // Tracks PrimaryVertex VertexerTracksWithConstraint
c0997643 854 // TPC TPCVertex VertexerTracksNoConstraint
6d597ca2 855 // TPC TPCVertex VertexerTracksWithConstraint
c0997643 856 // SPD SPDVertex vertexer: 3D
857 // SPD SPDVertex vertexer: Z
858
d3a3f33d 859 if(!vtx){
860 Printf("%s:%d No ESD vertex found",(char*)__FILE__,__LINE__);
861 return kFALSE;
862 }
6d597ca2 863 Int_t nCont = vtx->GetNContributors();
6d597ca2 864 if(nCont>=1){
865 fEventCutInfoESD |= kContributorsCut1;
866 if(nCont>=2){
867 fEventCutInfoESD |= kContributorsCut2;
868 if(nCont>=3){
869 fEventCutInfoESD |= kContributorsCut3;
870 }
871 }
872 }
873
874 if(nCont<3)return kFALSE;
c0997643 875
6d597ca2 876 // do not want tpc only primary vertex
877 TString vtxName(vtx->GetName());
878 if(vtxName.Contains("TPCVertex")){
879 fEventCutInfoESD |= kVertexTPC;
880 return kFALSE;
881 }
882 if(vtxName.Contains("SPDVertex"))fEventCutInfoESD |= kVertexSPD;
883 if(vtxName.Contains("PrimaryVertex"))fEventCutInfoESD |= kVertexGlobal;
c0997643 884
c0997643 885
6d597ca2 886 TString vtxTitle(vtx->GetTitle());
887 if(vtxTitle.Contains("vertexer: Z")){
888 if(vtx->GetDispersion()>0.02)return kFALSE;
889 }
890 fEventCutInfoESD |= kSPDDispersionCut;
891 return kTRUE;
892}
893
894
895Bool_t AliAnalysisTaskJetServices::IsVertexValid ( const AliAODVertex* vtx) const {
896
897 // We can check the type of the vertex by its name and title
898 // if vertexer failed title is empty (default c'tor) and ncontributors is 0
899 // vtx name title
900 // Tracks PrimaryVertex VertexerTracksNoConstraint
901 // TPC TPCVertex VertexerTracksNoConstraint
902 // SPD SPDVertex vertexer: 3D
903 // SPD SPDVertex vertexer: Z
d8f21f85 904
d3a3f33d 905 if(!vtx){
906 Printf("%s:%d No AOD vertex found",(char*)__FILE__,__LINE__);
907 return kFALSE;
908 }
909
9410e6a2 910 TString vtxName(vtx->GetName());
d8f21f85 911 if(fDebug){
912 Printf(" n contrib %d",vtx->GetNContributors());
9410e6a2 913 Printf("vtxName: %s",vtxName.Data());
d8f21f85 914 vtx->Print();
915 }
6d597ca2 916
d8f21f85 917 // if(vtx->GetNContributors()<3)return kFALSE;
c0997643 918 // do not want tpc only primary vertex
9410e6a2 919
c0997643 920 if(vtxName.Contains("TPCVertex"))return kFALSE;
6d597ca2 921
e5ee1c19 922 // no dispersion yet...
6d597ca2 923 /*
924 TString vtxTitle(vtx->GetTitle());
925 if(vtxTitle.Contains("vertexer: Z")){
926 if(vtx->GetDispersion()>0.02)return kFALSE;
927 }
928 */
929 return kTRUE;
930}
931
932
e5ee1c19 933Bool_t AliAnalysisTaskJetServices::IsVertexIn (const AliESDVertex* vtx,const AliESDVertex* vtxSPD) {
6d597ca2 934
935 if(!IsVertexValid(vtx))return kFALSE;
936
c0997643 937 Float_t zvtx = vtx->GetZ();
938 Float_t yvtx = vtx->GetY();
939 Float_t xvtx = vtx->GetX();
e5ee1c19 940 Float_t deltaZ = zvtx - vtxSPD->GetZ();
941 if(fDebug)Printf("%s:%d deltaz %3.3f ",__FILE__,__LINE__,deltaZ);
c0997643 942
6d597ca2 943 xvtx -= fVtxXMean;
944 yvtx -= fVtxYMean;
945 zvtx -= fVtxZMean;
946
947
948
949 if(TMath::Abs(zvtx)>fVtxZCut){
950 return kFALSE;
951 }
952 fEventCutInfoESD |= kVertexZCut;
953 Float_t r2 = yvtx*yvtx+xvtx*xvtx;
954 if(r2>(fVtxRCut*fVtxRCut)){
955 return kFALSE;
956 }
957 fEventCutInfoESD |= kVertexRCut;
e5ee1c19 958
959
960
961
962 if(TMath::Abs(deltaZ)>fVtxDeltaZCut)return kFALSE;
963 fEventCutInfoESD |= kVertexDeltaZCut;
6d597ca2 964 return kTRUE;
965}
966
967
e5ee1c19 968Bool_t AliAnalysisTaskJetServices::IsVertexIn (const AliAODVertex* vtx,const AliAODVertex* vtxSPD) const {
6d597ca2 969
970 if(!IsVertexValid(vtx))return kFALSE;
971
972 Float_t zvtx = vtx->GetZ();
973 Float_t yvtx = vtx->GetY();
974 Float_t xvtx = vtx->GetX();
e5ee1c19 975 Float_t deltaZ = zvtx - vtxSPD->GetZ();
976 if(fDebug)Printf("%s:%d deltaz %3.3f ",__FILE__,__LINE__,deltaZ);
c0997643 977
978 xvtx -= fVtxXMean;
979 yvtx -= fVtxYMean;
980 zvtx -= fVtxZMean;
981
982 Float_t r2 = yvtx*yvtx+xvtx*xvtx;
983
6d597ca2 984 Bool_t vertexIn = TMath::Abs(zvtx)<fVtxZCut&&r2<(fVtxRCut*fVtxRCut);
e5ee1c19 985
986
987 if(TMath::Abs(deltaZ)>fVtxDeltaZCut)vertexIn = kFALSE;
988 if(fDebug)Printf("%s:%d vertexIn %d ",__FILE__,__LINE__,vertexIn);
989
990
6d597ca2 991 return vertexIn;
c0997643 992}
993
6d597ca2 994Bool_t AliAnalysisTaskJetServices::IsEventPileUp(const AliESDEvent* esd) const{
f2dd0695 995 if(!esd)return kFALSE;
996 return esd->IsPileupFromSPD();
997}
998
6d597ca2 999Bool_t AliAnalysisTaskJetServices::IsEventCosmic(const AliESDEvent* esd) const {
f2dd0695 1000 if(!esd)return kFALSE;
1001 // add track cuts for which we look for cosmics...
df7848fc 1002
1003 Bool_t isCosmic = kFALSE;
1004 Int_t nTracks = esd->GetNumberOfTracks();
1005 Int_t nCosmicCandidates = 0;
1006
1007 for (Int_t iTrack1 = 0; iTrack1 < nTracks; iTrack1++) {
1008 AliESDtrack* track1 = (AliESDtrack*)esd->GetTrack(iTrack1);
1009 if (!track1) continue;
75946d5d 1010 UInt_t status1 = track1->GetStatus();
1011 //If track is ITS stand alone track, skip the track
1012 if (((status1 & AliESDtrack::kITSin) == 0 || (status1 & AliESDtrack::kTPCin))) continue;
df7848fc 1013 if(track1->Pt()<fPtMinCosmic) continue;
1014 //Start 2nd track loop to look for correlations
1015 for (Int_t iTrack2 = iTrack1+1; iTrack2 < nTracks; iTrack2++) {
1016 AliESDtrack* track2 = (AliESDtrack*)esd->GetTrack(iTrack2);
1017 if(!track2) continue;
75946d5d 1018 UInt_t status2 = track2->GetStatus();
1019 //If track is ITS stand alone track, skip the track
1020 if (((status2 & AliESDtrack::kITSin) == 0 || (status2 & AliESDtrack::kTPCin))) continue;
df7848fc 1021 if(track2->Pt()<fPtMinCosmic) continue;
1022 //Check if back-to-back
1023 Double_t mom1[3],mom2[3];
1024 track1->GetPxPyPz(mom1);
1025 track2->GetPxPyPz(mom2);
1026 TVector3 momv1(mom1[0],mom1[1],mom1[2]);
1027 TVector3 momv2(mom2[0],mom2[1],mom2[2]);
1028 Float_t theta = (float)(momv1.Phi()-momv2.Phi());
1029 if(theta<-0.5*TMath::Pi()) theta+=2.*TMath::Pi();
1030
1031 Float_t deltaPhi = track1->Phi()-track2->Phi();
1032 if(deltaPhi<-0.5*TMath::Pi()) deltaPhi+=2.*TMath::Pi();
1033
1034 Float_t rIsol = (float)(TMath::Sqrt( deltaPhi*deltaPhi+(track1->Eta()-track2->Eta())*(track1->Eta()-track2->Eta()) ));
1035 if(rIsol<fRIsolMinCosmic) continue;
1036
1037 if(TMath::Abs(TMath::Pi()-theta)<fMaxCosmicAngle) {
75946d5d 1038 nCosmicCandidates+=1;
df7848fc 1039 isCosmic = kTRUE;
1040 }
1041
1042 }
1043 }
1044
1045 fh1NCosmicsPerEvent->Fill((float)nCosmicCandidates);
1046
1047 return isCosmic;
f2dd0695 1048}
1049
1050
f4132e7d 1051Int_t AliAnalysisTaskJetServices::GetEventClass(AliESDEvent *esd){
1052
aa3ba8d2 1053 Float_t cent = 0;
21df3cc7 1054 if(fCollisionType==kPbPb){
1055 if(esd->GetCentrality()){
1056 cent = esd->GetCentrality()->GetCentralityPercentile("V0M");
1057 }
1058 if(cent<0)cent = 101;
1059 if(cent>80||cent<0)return 5;
1060 if(cent>50)return 4;
1061 if(cent>30)return 3;
1062 if(cent>10)return 2;
1063 return 1;
f4132e7d 1064 }
f4132e7d 1065 return 1;
f4132e7d 1066}
f2dd0695 1067
ffab794c 1068
1069Int_t AliAnalysisTaskJetServices::GetEventClass(AliAODEvent *aod){
1070
21df3cc7 1071 if(fCollisionType==kPbPb){
0a918d8d 1072 Float_t cent = ((AliVAODHeader*)aod->GetHeader())->GetCentrality();
21df3cc7 1073 if(cent>80||cent<0)return 5;
1074 if(cent>50)return 4;
1075 if(cent>30)return 3;
1076 if(cent>10)return 2;
1077 return 1;
1078 }
ffab794c 1079 return 1;
1080
1081}
1082
1083
bf7b8731 1084void AliAnalysisTaskJetServices::Terminate(Option_t */*option*/)
1085{
1086 // Terminate analysis
1087 //
bf7b8731 1088}
f4132e7d 1089
2eded560 1090Bool_t AliAnalysisTaskJetServices::CalculateReactionPlaneAngleVZERO(AliAODEvent *aod){
40445651 1091
2eded560 1092 // const Double_t arr_eta[11]={-3.7, -3.2, -2.7, -2.2, -1.7,0, 2.8, 3.4, 3.9, 4.5,5.1};
1093 if(!aod)return kFALSE;
1094 static bool bFirst = true;
1095 static Double_t v0phi[64] = {0,};
1096
1097 if(bFirst){
1098 int is=0;
1099 for(int iArm = 0; iArm<2; iArm++){
1100 for(int iRing = 0; iRing<4 ; iRing++){
1101 for(int iSec = 0; iSec<8 ; iSec++){
1102 v0phi[is] = 22.5 + 45. * iSec;
1103 v0phi[is] *= TMath::Pi()/180;
1104 // cout<< is <<" "<< v0phi[is]<<endl;
1105 is++;
1106 }
1107 }
40445651 1108 }
2eded560 1109 bFirst = false;
40445651 1110 }
1111
2eded560 1112 //
1113 const AliAODVZERO *aodVZERO = aod->GetVZEROData();
1114 Double_t numYZNA = 0,numXZNA = 0,sumZNA = 0;
1115 Double_t meanXA = 0,meanYA = 0;
1116
1117 Double_t numYZNC = 0,numXZNC = 0,sumZNC = 0;
1118 Double_t meanXC = 0,meanYC = 0;
1119
b360e4be 1120
1121
1122 static Int_t iOldRun = -1;
1123 static Int_t iFoundBin = -1;
2e90089b 1124
1125 if(aod->GetRunNumber()!=iOldRun&&(fp1CalibRPYA)){
b360e4be 1126 // search only or the bin in case of new runs
1127 iFoundBin = -1;
1128 Int_t ib = fp1CalibRPYA->FindBin(aod->GetRunNumber());
1129 Float_t err = fp1CalibRPYA->GetBinError(ib);
1130 if(err>0){// value can be zero...
1131 iFoundBin = ib;
1132 }
1133 else{
1134 Int_t ibLo = ib-1;
1135 Int_t ibUp = ib+1;
1136 while(iFoundBin<0&&(ibLo>0||ibUp<=fp1CalibRPYA->GetNbinsX())){
1137 err = fp1CalibRPYA->GetBinError(ibLo);
1138 if(err>0){
1139 iFoundBin = ibLo;
1140 }
1141 else{
1142 err = fp1CalibRPYA->GetBinError(ibUp);
1143 if(err>0)iFoundBin = ibUp;
1144 }
1145 ibUp++;
1146 ibLo--;
1147 }
1148 }
1149 iOldRun = aod->GetRunNumber();
1150 }
1151
2e90089b 1152 if(fDebug)Printf("%s:%d iFoundBin %d",(char*)__FILE__,__LINE__,iFoundBin);
b360e4be 1153
d889ce29 1154 if(iFoundBin>0&&(fp1CalibRPYA)){
b360e4be 1155 meanXA = fp1CalibRPXA->GetBinContent(iFoundBin);
1156 meanYA = fp1CalibRPYA->GetBinContent(iFoundBin);
1157 meanXC = fp1CalibRPXC->GetBinContent(iFoundBin);
1158 meanYC = fp1CalibRPYC->GetBinContent(iFoundBin);
1159 }
1160
2e90089b 1161 if(fDebug)Printf("%s:%d iFoundBin %1.3E %1.3E %1.3E %1.3E",(char*)__FILE__,__LINE__,meanXA,meanYA,meanXC,meanYC);
b360e4be 1162
2eded560 1163 for (int i=0; i<64; i++) {
1164 Double_t mult = aodVZERO->GetMultiplicity(i);
1165 Double_t phi = v0phi[i];
1166 if (mult>0) {
1167 if (i<32) { //C-side
1168 Double_t wZNC= mult;
1169 numYZNC += sin(2.*phi)*wZNC;
1170 numXZNC += cos(2.*phi)*wZNC;
1171 sumZNC+=wZNC;
1172 }
1173 else if(i>31){ //A-side
1174 Double_t wZNA=mult;
1175 numYZNA += sin(2.*phi)*wZNA;
1176 numXZNA += cos(2.*phi)*wZNA;
1177 sumZNA+=wZNA;
1178 }
1179 }// mult>0
1180 }// 64 sectors
1181
40445651 1182
a374d56f 1183 Double_t XC = 0;
1184 Double_t YC = 0;
1185 Double_t XA = 0;
1186 Double_t YA = 0;
1187
1188 if(sumZNC!=0){
2eded560 1189
a374d56f 1190 XC = numXZNC/sumZNC;
1191 YC = numYZNC/sumZNC;
6622701d 1192 fPsiVZEROC = 0.5*TMath::ATan2(YC-meanYC, XC-meanXC);
a374d56f 1193 if(fPsiVZEROC>TMath::Pi()){fPsiVZEROC-=TMath::Pi();}
1194 if(fPsiVZEROC<0){fPsiVZEROC+=TMath::Pi();}
1195
1196
1197 fh2XYC->Fill(XC-meanXC,YC-meanYC); // control
1198 fp1RPXC->Fill(aod->GetRunNumber(),XC);
1199 fp1RPYC->Fill(aod->GetRunNumber(),YC);
1200
2eded560 1201
a374d56f 1202 }
40445651 1203
a374d56f 1204 if(sumZNA!=0){
1205 XA = numXZNA/sumZNA;
1206 YA = numYZNA/sumZNA;
1207 fPsiVZEROA = 0.5*TMath::ATan2(YA-meanYA, XA-meanXA);
1208 if(fPsiVZEROA>TMath::Pi()){fPsiVZEROA-=TMath::Pi();}
1209 if(fPsiVZEROA<0){fPsiVZEROA+=TMath::Pi();}
1210
1211 fh2XYA->Fill(XA-meanXA,YA-meanYA); // control
1212 fp1RPXA->Fill(aod->GetRunNumber(),XA);
1213 fp1RPYA->Fill(aod->GetRunNumber(),YA);
1214
1215 }
40445651 1216 return kTRUE;
40445651 1217
40445651 1218}
1219
40445651 1220Int_t AliAnalysisTaskJetServices::GetListOfTracks(TList *list){
1221 Int_t iCount = 0;
1222 AliAODEvent *aod = 0;
1223 if(fUseAODInput)aod = dynamic_cast<AliAODEvent*>(InputEvent());
1224 else aod = AODEvent();
1225 if(!aod){
1226 return iCount;
1227 }
1228 for(int it = 0;it < aod->GetNumberOfTracks();++it){
1229 AliAODTrack *tr = aod->GetTrack(it);
1230 if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
1231 if(TMath::Abs(tr->Eta())>fTrackRecEtaWindow)continue;
1232 if(tr->Pt()<fMinTrackPt)continue;
1233 list->Add(tr);
1234 iCount++;
1235 }
1236 list->Sort();
1237 return iCount;
1238
1239}
1240
3d1dba29 1241void AliAnalysisTaskJetServices::SetNTrigger(Int_t n){
1242 if(n>0){
1243 fNTrigger = n;
1244 delete [] fTriggerName;
1245 fTriggerName = new TString [fNTrigger];
1246 delete [] fTriggerBit;fTriggerBit = 0;
1247 fTriggerBit = new UInt_t [fNTrigger];
1248 }
1249 else{
1250 fNTrigger = 0;
1251 }
1252}
1253
1254void AliAnalysisTaskJetServices::SetTrigger(Int_t i,UInt_t it,const char* c){
1255 if(i<fNTrigger){
1256 Printf("%d",it);
1257 Printf("%p",c);
1258 Printf("%s",c);
1259 Printf("%p",&fTriggerName[i]);
1260
1261 fTriggerBit[i] = it;
1262 fTriggerName[i] = c; // placement new
1263 // new(&fTriggerName[i]) TString(c); // placement new
1264 Printf("%s",fTriggerName[i].Data());
1265
1266 }
1267}