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