]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/JetTasks/AliAnalysisTaskJetServices.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Some clean up in addTask macros, added centrality selection and event classes to...
[u/mrichter/AliRoot.git] / PWG4 / JetTasks / AliAnalysisTaskJetServices.cxx
CommitLineData
bf7b8731 1// **************************************
2// Task used for the correction of determiantion of reconstructed jet spectra
3// Compares input (gen) and output (rec) jets
4// *******************************************
5
6
7/**************************************************************************
8 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
9 * *
10 * Author: The ALICE Off-line Project. *
11 * Contributors are mentioned in the code where appropriate. *
12 * *
13 * Permission to use, copy, modify and distribute this software and its *
14 * documentation strictly for non-commercial purposes is hereby granted *
15 * without fee, provided that the above copyright notice appears in all *
16 * copies and that both the copyright notice and this permission notice *
17 * appear in the supporting documentation. The authors make no claims *
18 * about the suitability of this software for any purpose. It is *
19 * provided "as is" without express or implied warranty. *
20 **************************************************************************/
21
22
23#include <TROOT.h>
24#include <TRandom.h>
b1cd90b5 25#include <TString.h>
bf7b8731 26#include <TSystem.h>
27#include <TInterpreter.h>
28#include <TChain.h>
29#include <TFile.h>
30#include <TKey.h>
31#include <TH1F.h>
32#include <TH2F.h>
33#include <TH3F.h>
34#include <TProfile.h>
35#include <TList.h>
36#include <TLorentzVector.h>
37#include <TClonesArray.h>
c0997643 38#include "TDatabasePDG.h"
bf7b8731 39
40#include "AliAnalysisTaskJetServices.h"
f4132e7d 41#include "AliESDCentrality.h"
fe669ac6 42#include "AliAnalysisDataContainer.h"
43#include "AliAnalysisDataSlot.h"
bf7b8731 44#include "AliAnalysisManager.h"
45#include "AliJetFinder.h"
46#include "AliJetHeader.h"
47#include "AliJetReader.h"
48#include "AliJetReaderHeader.h"
49#include "AliUA1JetHeaderV1.h"
50#include "AliESDEvent.h"
51#include "AliAODEvent.h"
52#include "AliAODHandler.h"
53#include "AliAODTrack.h"
54#include "AliAODJet.h"
55#include "AliAODMCParticle.h"
56#include "AliMCEventHandler.h"
57#include "AliMCEvent.h"
58#include "AliStack.h"
59#include "AliGenPythiaEventHeader.h"
60#include "AliJetKineReaderHeader.h"
61#include "AliGenCocktailEventHeader.h"
62#include "AliInputEventHandler.h"
fe669ac6 63#include "AliPhysicsSelection.h"
3493c3a9 64#include "AliTriggerAnalysis.h"
bf7b8731 65
66#include "AliAnalysisHelperJetTasks.h"
67
68ClassImp(AliAnalysisTaskJetServices)
69
f4132e7d 70AliAODHeader* AliAnalysisTaskJetServices::fgAODHeader = NULL;
71
bf7b8731 72AliAnalysisTaskJetServices::AliAnalysisTaskJetServices(): AliAnalysisTaskSE(),
73 fUseAODInput(kFALSE),
cc0649e4 74 fUsePhysicsSelection(kFALSE),
3493c3a9 75 fMC(kFALSE),
f4132e7d 76 fFilterAODCollisions(kFALSE),
39e7e8ab 77 fPhysicsSelectionFlag(AliVEvent::kMB),
f2dd0695 78 fSelectionInfoESD(0),
6d597ca2 79 fEventCutInfoESD(0),
bf7b8731 80 fAvgTrials(1),
c0997643 81 fVtxXMean(0),
82 fVtxYMean(0),
83 fVtxZMean(0),
84 fVtxRCut(1.),
85 fVtxZCut(8.),
df7848fc 86 fPtMinCosmic(5.),
87 fRIsolMinCosmic(3.),
88 fMaxCosmicAngle(0.01),
f4132e7d 89 fNonStdFile(""),
bf7b8731 90 fh1Xsec(0x0),
91 fh1Trials(0x0),
92 fh1PtHard(0x0),
93 fh1PtHardTrials(0x0),
f2dd0695 94 fh1SelectionInfoESD(0x0),
6d597ca2 95 fh1EventCutInfoESD(0),
bf7b8731 96 fh2TriggerCount(0x0),
97 fh2ESDTriggerCount(0x0),
98 fh2TriggerVtx(0x0),
99 fh2ESDTriggerVtx(0x0),
b5a3f310 100 fh2ESDTriggerRun(0x0),
101 fh2VtxXY(0x0),
df7848fc 102 fh1NCosmicsPerEvent(0x0),
3493c3a9 103 fTriggerAnalysis(0x0),
bf7b8731 104 fHistList(0x0)
105{
b5a3f310 106 fRunRange[0] = fRunRange[1] = 0;
bf7b8731 107}
108
109AliAnalysisTaskJetServices::AliAnalysisTaskJetServices(const char* name):
110 AliAnalysisTaskSE(name),
111 fUseAODInput(kFALSE),
cc0649e4 112 fUsePhysicsSelection(kFALSE),
3493c3a9 113 fMC(kFALSE),
f4132e7d 114 fFilterAODCollisions(kFALSE),
39e7e8ab 115 fPhysicsSelectionFlag(AliVEvent::kMB),
f2dd0695 116 fSelectionInfoESD(0),
6d597ca2 117 fEventCutInfoESD(0),
bf7b8731 118 fAvgTrials(1),
c0997643 119 fVtxXMean(0),
120 fVtxYMean(0),
121 fVtxZMean(0),
122 fVtxRCut(1.),
123 fVtxZCut(8.),
df7848fc 124 fPtMinCosmic(5.),
125 fRIsolMinCosmic(3.),
126 fMaxCosmicAngle(0.01),
f4132e7d 127 fNonStdFile(""),
bf7b8731 128 fh1Xsec(0x0),
129 fh1Trials(0x0),
130 fh1PtHard(0x0),
131 fh1PtHardTrials(0x0),
f2dd0695 132 fh1SelectionInfoESD(0x0),
6d597ca2 133 fh1EventCutInfoESD(0),
bf7b8731 134 fh2TriggerCount(0x0),
135 fh2ESDTriggerCount(0x0),
136 fh2TriggerVtx(0x0),
137 fh2ESDTriggerVtx(0x0),
b5a3f310 138 fh2ESDTriggerRun(0x0),
139 fh2VtxXY(0x0),
df7848fc 140 fh1NCosmicsPerEvent(0x0),
3493c3a9 141 fTriggerAnalysis(0x0),
bf7b8731 142 fHistList(0x0)
143{
b5a3f310 144 fRunRange[0] = fRunRange[1] = 0;
bf7b8731 145 DefineOutput(1,TList::Class());
146}
147
148
149
150Bool_t AliAnalysisTaskJetServices::Notify()
151{
152 //
153 // Implemented Notify() to read the cross sections
154 // and number of trials from pyxsec.root
155 //
156
157 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
158 Float_t xsection = 0;
159 Float_t ftrials = 1;
160
161 fAvgTrials = 1;
162 if(tree){
163 TFile *curfile = tree->GetCurrentFile();
164 if (!curfile) {
165 Error("Notify","No current file");
166 return kFALSE;
167 }
168 if(!fh1Xsec||!fh1Trials){
169 Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
170 return kFALSE;
171 }
172 AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
173 fh1Xsec->Fill("<#sigma>",xsection);
174 // construct a poor man average trials
175 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
176 if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
177 }
178 return kTRUE;
179}
180
181void AliAnalysisTaskJetServices::UserCreateOutputObjects()
182{
183
184 //
185 // Create the output container
186 //
187
188
bf7b8731 189 if (fDebug > 1) printf("AnalysisTaskJetServices::UserCreateOutputObjects() \n");
190
191 OpenFile(1);
192 if(!fHistList)fHistList = new TList();
193
194 Bool_t oldStatus = TH1::AddDirectoryStatus();
195 TH1::AddDirectory(kFALSE);
bf7b8731 196 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
197 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
198 fHistList->Add(fh1Xsec);
199
200 fh1Trials = new TH1F("fh1Trials","trials root file",1,0,1);
201 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
202 fHistList->Add(fh1Trials);
203
204 const Int_t nBinPt = 100;
205 Double_t binLimitsPt[nBinPt+1];
206 for(Int_t iPt = 0;iPt <= nBinPt;iPt++){
207 if(iPt == 0){
208 binLimitsPt[iPt] = 0.0;
209 }
210 else {// 1.0
211 binLimitsPt[iPt] = binLimitsPt[iPt-1] + 2.5;
212 }
213 }
214
215 fh2TriggerCount = new TH2F("fh2TriggerCount",";Trigger No.;constrained;Count",AliAnalysisHelperJetTasks::kTrigger,-0.5,AliAnalysisHelperJetTasks::kTrigger-0.5,kConstraints,-0.5,kConstraints-0.5);
216 fHistList->Add(fh2TriggerCount);
217
218 fh2ESDTriggerCount = new TH2F("fh2ESDTriggerCount",";Trigger No.;constrained;Count",AliAnalysisHelperJetTasks::kTrigger,-0.5,AliAnalysisHelperJetTasks::kTrigger-0.5,kConstraints,-0.5,kConstraints-0.5);
219 fHistList->Add(fh2ESDTriggerCount);
57ca1193 220 const Int_t nBins = AliAnalysisHelperJetTasks::kTrigger*kConstraints;
221 fh2TriggerVtx = new TH2F("fh2TriggerVtx",";Constraint No. * (trig no+1);Vtx (cm);Count",nBins,-0.5,nBins-0.5,400,-20.,20.);
bf7b8731 222 fHistList->Add(fh2TriggerVtx);
223
57ca1193 224 fh2ESDTriggerVtx = new TH2F("fh2ESDTriggerVtx",";Constraint No.* (trg no+1);Vtx (cm);Count",nBins,-0.5,nBins-0.5,400,-20.,20.);
bf7b8731 225 fHistList->Add(fh2ESDTriggerVtx);
226
227
228 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",nBinPt,binLimitsPt);
229 fHistList->Add(fh1PtHard);
230 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",nBinPt,binLimitsPt);
231 fHistList->Add(fh1PtHardTrials);
f2dd0695 232 fh1SelectionInfoESD = new TH1F("fh1SelectionInfoESD","Bit Masks that satisfy the selection info",
233 AliAnalysisHelperJetTasks::kTotalSelections,0.5,AliAnalysisHelperJetTasks::kTotalSelections+0.5);
f2dd0695 234 fHistList->Add(fh1SelectionInfoESD);
6d597ca2 235
236 fh1EventCutInfoESD = new TH1F("fh1EventCutInfoESD","Bit Masks that for the events after each step of cuts",
237 kTotalEventCuts,0.5,kTotalEventCuts+0.5);
238 fHistList->Add(fh1EventCutInfoESD);
239
b5a3f310 240 // 3 decisions, 0 trigger X, X + SPD vertex, X + SPD vertex in range
241 // 3 triggers BB BE/EB EE
242
f4132e7d 243 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 244 fHistList->Add(fh2ESDTriggerRun);
245
246 fh2VtxXY = new TH2F("fh2VtxXY","Beam Spot all INT triggered events;x (cm);y (cm)",160,-10,10,160,-10,10);
247 fHistList->Add(fh2VtxXY);
bf7b8731 248 // =========== Switch on Sumw2 for all histos ===========
249 for (Int_t i=0; i<fHistList->GetEntries(); ++i) {
250 TH1 *h1 = dynamic_cast<TH1*>(fHistList->At(i));
251 if (h1){
252 h1->Sumw2();
253 continue;
254 }
255 THnSparse *hn = dynamic_cast<THnSparse*>(fHistList->At(i));
256 if(hn)hn->Sumw2();
257 }
258
df7848fc 259 fh1NCosmicsPerEvent = new TH1F("fh1NCosmicsPerEvent","Number of cosmic candidates per event",10,0.,10.);
260 fHistList->Add(fh1NCosmicsPerEvent),
261
bf7b8731 262
263 TH1::AddDirectory(oldStatus);
f4132e7d 264
265 // Add an AOD branch for replication
266 if(fNonStdFile.Length()){
267 if (fDebug > 1) AliInfo("Replicating header");
268 fgAODHeader = new AliAODHeader;
269 AddAODBranch("AliAODHeader",&fgAODHeader,fNonStdFile.Data());
270 }
bf7b8731 271}
272
273void AliAnalysisTaskJetServices::Init()
274{
275 //
276 // Initialization
277 //
bf7b8731 278 if (fDebug > 1) printf("AnalysisTaskJetServices::Init() \n");
279
280}
281
282void AliAnalysisTaskJetServices::UserExec(Option_t */*option*/)
283{
284
285 //
286 // Execute analysis for current event
287 //
288
289 AliAODEvent *aod = 0;
290 AliESDEvent *esd = 0;
f4132e7d 291
292
3493c3a9 293
bf7b8731 294 AliAnalysisHelperJetTasks::Selected(kTRUE,kFALSE); // set slection to false
f4132e7d 295 AliAnalysisHelperJetTasks::EventClass(kTRUE,0);
f2dd0695 296 fSelectionInfoESD = 0; // reset
6d597ca2 297 fEventCutInfoESD = 0; // reset
f2dd0695 298 AliAnalysisHelperJetTasks::SelectInfo(kTRUE,fSelectionInfoESD); // set slection to false
299
bf7b8731 300
f4132e7d 301 static AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
302
303
304
bf7b8731 305 if(fUseAODInput){
306 aod = dynamic_cast<AliAODEvent*>(InputEvent());
307 if(!aod){
308 Printf("%s:%d AODEvent not found in Input Manager %d",(char*)__FILE__,__LINE__,fUseAODInput);
309 return;
310 }
311 // fethc the header
312 }
313 else{
314 // assume that the AOD is in the general output...
315 aod = AODEvent();
316 if(!aod){
317 Printf("%s:%d AODEvent not found in the Output",(char*)__FILE__,__LINE__);
318 return;
319 }
320 esd = dynamic_cast<AliESDEvent*>(InputEvent());
321 }
322
3493c3a9 323 fSelectionInfoESD |= kNoEventCut;
324 fEventCutInfoESD |= kNoEventCut;
6d597ca2 325
326 Bool_t esdVtxValid = false;
327 Bool_t esdVtxIn = false;
328 Bool_t aodVtxValid = false;
329 Bool_t aodVtxIn = false;
f2dd0695 330
b5a3f310 331 if(esd){
3493c3a9 332 // trigger analyisis
333 if(!fTriggerAnalysis){
334 fTriggerAnalysis = new AliTriggerAnalysis;
335 fTriggerAnalysis->SetAnalyzeMC(fMC);
336 fTriggerAnalysis->SetSPDGFOThreshhold(1);
337 }
338 // fTriggerAnalysis->FillTriggerClasses(esd);
339 Bool_t v0A = fTriggerAnalysis->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0A);
340 Bool_t v0C = fTriggerAnalysis->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0C);
341 Bool_t v0ABG = fTriggerAnalysis->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0ABG);
342 Bool_t v0CBG = fTriggerAnalysis->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0CBG);
343 Bool_t spdFO = fTriggerAnalysis->SPDFiredChips(esd, 0);;
344 if(v0A)fSelectionInfoESD |= AliAnalysisHelperJetTasks::kV0A;
345 if(v0C)fSelectionInfoESD |= AliAnalysisHelperJetTasks::kV0C;
346 if(!(v0ABG||v0CBG))fSelectionInfoESD |= AliAnalysisHelperJetTasks::kNoV0BG;
347 if(spdFO)fSelectionInfoESD |= AliAnalysisHelperJetTasks::kSPDFO;
b5a3f310 348 }
f4132e7d 349
c0997643 350 // Apply additional constraints
6d597ca2 351 Bool_t esdEventSelected = IsEventSelected(esd);
352 Bool_t esdEventPileUp = IsEventPileUp(esd);
353 Bool_t esdEventCosmic = IsEventCosmic(esd);
f2dd0695 354
6d597ca2 355 Bool_t aodEventSelected = IsEventSelected(aod);
f2dd0695 356
39e7e8ab 357 Bool_t physicsSelection = ((fInputHandler->IsEventSelected())&fPhysicsSelectionFlag);
f4132e7d 358 if(aodH&&physicsSelection&&fFilterAODCollisions){
359 aodH->SetFillAOD(kTRUE);
360 }
361
362
6d597ca2 363 fEventCutInfoESD |= kPhysicsSelectionCut; // other alreay set via IsEventSelected
364 fh1EventCutInfoESD->Fill(fEventCutInfoESD);
365
f2dd0695 366 if(esdEventSelected) fSelectionInfoESD |= AliAnalysisHelperJetTasks::kVertexIn;
367 if(esdEventPileUp) fSelectionInfoESD |= AliAnalysisHelperJetTasks::kIsPileUp;
368 if(esdEventCosmic) fSelectionInfoESD |= AliAnalysisHelperJetTasks::kIsCosmic;
369 if(physicsSelection) fSelectionInfoESD |= AliAnalysisHelperJetTasks::kPhysicsSelection;
370
371
372 // here we have all selection information, fill histogram
373 for(unsigned int i = 1;i<(UInt_t)fh1SelectionInfoESD->GetNbinsX();i++){
374 if((i&fSelectionInfoESD)==i)fh1SelectionInfoESD->Fill(i);
375 }
6d597ca2 376 AliAnalysisHelperJetTasks::SelectInfo(kTRUE,fSelectionInfoESD);
c0997643 377
378 if(esd&&aod&&fDebug){
379 if(esdEventSelected&&!aodEventSelected){
380 Printf("%s:%d Different Selection for ESD and AOD",(char*)__FILE__,__LINE__);
381 const AliESDVertex *vtxESD = esd->GetPrimaryVertex();
bf7b8731 382 const AliAODVertex *vtxAOD = aod->GetPrimaryVertex();
c0997643 383 Printf("ESD Vtx %s %s %d",vtxESD->GetName(),vtxESD->GetTitle(),vtxESD->GetNContributors());
384 vtxESD->Print();
385 Printf("AOD Vtx %s %s %d",vtxAOD->GetName(),vtxAOD->GetTitle(),vtxAOD->GetNContributors());
386 vtxAOD->Print();
bf7b8731 387 }
c0997643 388 }
389
390 // loop over all possible triggers for esd
391
392 if(esd){
393 const AliESDVertex *vtxESD = esd->GetPrimaryVertex();
f4132e7d 394 esdVtxValid = IsVertexValid(vtxESD);
395 esdVtxIn = IsVertexIn(vtxESD);
57ca1193 396 Float_t zvtx = vtxESD->GetZ();
f4132e7d 397 Int_t iCl = GetEventClass(esd);
398 AliAnalysisHelperJetTasks::EventClass(kTRUE,iCl);
399 Bool_t cand = physicsSelection;
400
401 Printf("%s:%d %d %d %d Icl %d",(char*)__FILE__,__LINE__,esdVtxValid,esdVtxIn,cand,iCl);
402
403 if(cand){
404 fh2ESDTriggerCount->Fill(0.,kSelectedALICE);
405 fh2ESDTriggerCount->Fill(iCl,kSelectedALICE);
406 fh2ESDTriggerVtx->Fill(kSelectedALICE*(iCl+1),zvtx);
407 }
408 // if(!fUsePhysicsSelection)cand = AliAnalysisHelperJetTasks::IsTriggerFired(esd,AliAnalysisHelperJetTasks::kMB1);
409 if(esdVtxValid){
410 fh2ESDTriggerCount->Fill(0.,kTriggeredVertex);
411 fh2ESDTriggerCount->Fill(iCl,kTriggeredVertex);
412 fh2ESDTriggerVtx->Fill(iCl,zvtx);
413 if(esdVtxIn){
414 fh2ESDTriggerCount->Fill(0.,kTriggeredVertexIn);
415 fh2ESDTriggerCount->Fill(iCl,kTriggeredVertexIn);
416 fh2ESDTriggerVtx->Fill(kTriggeredVertexIn*(iCl+1),zvtx);
bf7b8731 417 }
f4132e7d 418 if(cand){
419 fh2ESDTriggerCount->Fill(0.,kSelectedALICEVertexValid);
420 fh2ESDTriggerCount->Fill(iCl,kSelectedALICEVertexValid);
421 fh2ESDTriggerVtx->Fill(kSelectedALICEVertexValid*(iCl+1),zvtx);
c0997643 422 }
bf7b8731 423 }
f4132e7d 424 if(cand&&esdVtxIn){
425 fh2ESDTriggerCount->Fill(0.,kSelectedALICEVertexIn);
426 fh2ESDTriggerCount->Fill(iCl,kSelectedALICEVertexIn);
427 fh2ESDTriggerVtx->Fill(kSelectedALICEVertexIn*(iCl+1),zvtx);
428 fh2ESDTriggerVtx->Fill(kSelected*(iCl+1),zvtx);
429 fh2ESDTriggerCount->Fill(iCl,kSelected);
430 fh2ESDTriggerCount->Fill(0.,kSelected);
431 AliAnalysisHelperJetTasks::Selected(kTRUE,kTRUE);// select this event
432 }
bf7b8731 433 }
434
c0997643 435 if(aod){
436 const AliAODVertex *vtxAOD = aod->GetPrimaryVertex();
6d597ca2 437 aodVtxValid = IsVertexValid(vtxAOD);
438 aodVtxIn = IsVertexIn(vtxAOD);
439
c0997643 440 for(int it = AliAnalysisHelperJetTasks::kAcceptAll;it < AliAnalysisHelperJetTasks::kTrigger;it++){
441 Bool_t aodTrig = kFALSE;
442 aodTrig = AliAnalysisHelperJetTasks::IsTriggerFired(aod,(AliAnalysisHelperJetTasks::Trigger)it);
6d597ca2 443 Bool_t cand = aod->GetHeader()->GetOfflineTrigger()&AliVEvent::kMB;
c0997643 444 if(aodTrig)fh2TriggerCount->Fill(it,kAllTriggered);
6d597ca2 445 if(aodVtxValid){
c0997643 446 Float_t zvtx = vtxAOD->GetZ();
57ca1193 447 if(aodTrig){
448 fh2TriggerCount->Fill(it,kTriggeredVertex);
449 fh2TriggerVtx->Fill(kTriggeredVertex*(it+1),zvtx);
450 }
07c3e78c 451 if(cand&&aodEventSelected){
452 fh2TriggerCount->Fill(it,kSelected);
453 }
c0997643 454 if(aodTrig&&aodEventSelected){
57ca1193 455 fh2TriggerVtx->Fill(kTriggeredVertexIn*(it+1),zvtx);
c0997643 456 fh2TriggerCount->Fill(it,kTriggeredVertexIn);
457 }
e0c120d9 458 if(cand){
459 fh2TriggerCount->Fill(it,kSelectedALICEVertexValid);
57ca1193 460 fh2TriggerVtx->Fill(kSelectedALICEVertexValid*(it+1),zvtx);
461 if(aodEventSelected){
462 fh2TriggerCount->Fill(it,kSelectedALICEVertexIn);
463 fh2TriggerVtx->Fill(kSelectedALICEVertexIn*(it+1),zvtx);
464 }
e0c120d9 465 }
c0997643 466 }
467 }
468 }
bf7b8731 469
470 if (fDebug > 1)printf(" AliAnalysisTaskJetServices: Analysing event # %5d\n", (Int_t) fEntry);
471
472
473 Double_t ptHard = 0;
474 Double_t nTrials = 1; // Trials for MC trigger
475
476 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
477 AliMCEvent* mcEvent = MCEvent();
478 // AliStack *pStack = 0;
479 if(mcEvent){
480 AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
481 if(pythiaGenHeader){
482 nTrials = pythiaGenHeader->Trials();
483 ptHard = pythiaGenHeader->GetPtHard();
484 int iProcessType = pythiaGenHeader->ProcessType();
485 // 11 f+f -> f+f
486 // 12 f+barf -> f+barf
487 // 13 f+barf -> g+g
488 // 28 f+g -> f+g
489 // 53 g+g -> f+barf
490 // 68 g+g -> g+g
491 if (fDebug > 10)Printf("%d iProcessType %d",__LINE__, iProcessType);
492 if(fDebug>20)AliAnalysisHelperJetTasks::PrintStack(mcEvent);
493 fh1PtHard->Fill(ptHard);
494 fh1PtHardTrials->Fill(ptHard,nTrials);
495
496 }// if pythia gen header
497 }
498
499 // trigger selection
500
f4132e7d 501 // replication of
502 if(fNonStdFile.Length()&&aod){
503 if (fgAODHeader){
504 *fgAODHeader = *(dynamic_cast<AliAODHeader*>(aod->GetHeader()));
505 }
506 }
507
bf7b8731 508 PostData(1, fHistList);
f4132e7d 509 }
bf7b8731 510
6d597ca2 511Bool_t AliAnalysisTaskJetServices::IsEventSelected(const AliESDEvent* esd){
c0997643 512 if(!esd)return kFALSE;
513 const AliESDVertex *vtx = esd->GetPrimaryVertex();
6d597ca2 514 return IsVertexIn(vtx); // vertex in calls vertex valid
515}
516
517
518Bool_t AliAnalysisTaskJetServices::IsEventSelected(const AliAODEvent* aod) const {
519 if(!aod)return kFALSE;
520 const AliAODVertex *vtx = aod->GetPrimaryVertex();
521 return IsVertexIn(vtx); // VertexIn calls VertexValid
522}
523
524Bool_t AliAnalysisTaskJetServices::IsVertexValid ( const AliESDVertex* vtx) {
525
c0997643 526 // We can check the type of the vertex by its name and title
527 // if vertexer failed title is empty (default c'tor) and ncontributors is 0
528 // vtx name title
529 // Tracks PrimaryVertex VertexerTracksNoConstraint
6d597ca2 530 // Tracks PrimaryVertex VertexerTracksWithConstraint
c0997643 531 // TPC TPCVertex VertexerTracksNoConstraint
6d597ca2 532 // TPC TPCVertex VertexerTracksWithConstraint
c0997643 533 // SPD SPDVertex vertexer: 3D
534 // SPD SPDVertex vertexer: Z
535
6d597ca2 536 Int_t nCont = vtx->GetNContributors();
537
538 if(nCont>=1){
539 fEventCutInfoESD |= kContributorsCut1;
540 if(nCont>=2){
541 fEventCutInfoESD |= kContributorsCut2;
542 if(nCont>=3){
543 fEventCutInfoESD |= kContributorsCut3;
544 }
545 }
546 }
547
548 if(nCont<3)return kFALSE;
c0997643 549
6d597ca2 550 // do not want tpc only primary vertex
551 TString vtxName(vtx->GetName());
552 if(vtxName.Contains("TPCVertex")){
553 fEventCutInfoESD |= kVertexTPC;
554 return kFALSE;
555 }
556 if(vtxName.Contains("SPDVertex"))fEventCutInfoESD |= kVertexSPD;
557 if(vtxName.Contains("PrimaryVertex"))fEventCutInfoESD |= kVertexGlobal;
c0997643 558
c0997643 559
6d597ca2 560 TString vtxTitle(vtx->GetTitle());
561 if(vtxTitle.Contains("vertexer: Z")){
562 if(vtx->GetDispersion()>0.02)return kFALSE;
563 }
564 fEventCutInfoESD |= kSPDDispersionCut;
565 return kTRUE;
566}
567
568
569Bool_t AliAnalysisTaskJetServices::IsVertexValid ( const AliAODVertex* vtx) const {
570
571 // We can check the type of the vertex by its name and title
572 // if vertexer failed title is empty (default c'tor) and ncontributors is 0
573 // vtx name title
574 // Tracks PrimaryVertex VertexerTracksNoConstraint
575 // TPC TPCVertex VertexerTracksNoConstraint
576 // SPD SPDVertex vertexer: 3D
577 // SPD SPDVertex vertexer: Z
578
579 if(vtx->GetNContributors()<3)return kFALSE;
c0997643 580 // do not want tpc only primary vertex
581 TString vtxName(vtx->GetName());
582 if(vtxName.Contains("TPCVertex"))return kFALSE;
6d597ca2 583
584 // no dispersion yet...
585 /*
586 TString vtxTitle(vtx->GetTitle());
587 if(vtxTitle.Contains("vertexer: Z")){
588 if(vtx->GetDispersion()>0.02)return kFALSE;
589 }
590 */
591 return kTRUE;
592}
593
594
595Bool_t AliAnalysisTaskJetServices::IsVertexIn (const AliESDVertex* vtx) {
596
597 if(!IsVertexValid(vtx))return kFALSE;
598
c0997643 599 Float_t zvtx = vtx->GetZ();
600 Float_t yvtx = vtx->GetY();
601 Float_t xvtx = vtx->GetX();
602
6d597ca2 603 xvtx -= fVtxXMean;
604 yvtx -= fVtxYMean;
605 zvtx -= fVtxZMean;
606
607
608
609 if(TMath::Abs(zvtx)>fVtxZCut){
610 return kFALSE;
611 }
612 fEventCutInfoESD |= kVertexZCut;
613 Float_t r2 = yvtx*yvtx+xvtx*xvtx;
614 if(r2>(fVtxRCut*fVtxRCut)){
615 return kFALSE;
616 }
617 fEventCutInfoESD |= kVertexRCut;
618 return kTRUE;
619}
620
621
622Bool_t AliAnalysisTaskJetServices::IsVertexIn (const AliAODVertex* vtx) const {
623
624 if(!IsVertexValid(vtx))return kFALSE;
625
626 Float_t zvtx = vtx->GetZ();
627 Float_t yvtx = vtx->GetY();
628 Float_t xvtx = vtx->GetX();
c0997643 629
630 xvtx -= fVtxXMean;
631 yvtx -= fVtxYMean;
632 zvtx -= fVtxZMean;
633
634 Float_t r2 = yvtx*yvtx+xvtx*xvtx;
635
6d597ca2 636 Bool_t vertexIn = TMath::Abs(zvtx)<fVtxZCut&&r2<(fVtxRCut*fVtxRCut);
637 return vertexIn;
c0997643 638}
639
6d597ca2 640Bool_t AliAnalysisTaskJetServices::IsEventPileUp(const AliESDEvent* esd) const{
f2dd0695 641 if(!esd)return kFALSE;
642 return esd->IsPileupFromSPD();
643}
644
6d597ca2 645Bool_t AliAnalysisTaskJetServices::IsEventCosmic(const AliESDEvent* esd) const {
f2dd0695 646 if(!esd)return kFALSE;
647 // add track cuts for which we look for cosmics...
df7848fc 648
649 Bool_t isCosmic = kFALSE;
650 Int_t nTracks = esd->GetNumberOfTracks();
651 Int_t nCosmicCandidates = 0;
652
653 for (Int_t iTrack1 = 0; iTrack1 < nTracks; iTrack1++) {
654 AliESDtrack* track1 = (AliESDtrack*)esd->GetTrack(iTrack1);
655 if (!track1) continue;
75946d5d 656 UInt_t status1 = track1->GetStatus();
657 //If track is ITS stand alone track, skip the track
658 if (((status1 & AliESDtrack::kITSin) == 0 || (status1 & AliESDtrack::kTPCin))) continue;
df7848fc 659 if(track1->Pt()<fPtMinCosmic) continue;
660 //Start 2nd track loop to look for correlations
661 for (Int_t iTrack2 = iTrack1+1; iTrack2 < nTracks; iTrack2++) {
662 AliESDtrack* track2 = (AliESDtrack*)esd->GetTrack(iTrack2);
663 if(!track2) continue;
75946d5d 664 UInt_t status2 = track2->GetStatus();
665 //If track is ITS stand alone track, skip the track
666 if (((status2 & AliESDtrack::kITSin) == 0 || (status2 & AliESDtrack::kTPCin))) continue;
df7848fc 667 if(track2->Pt()<fPtMinCosmic) continue;
668 //Check if back-to-back
669 Double_t mom1[3],mom2[3];
670 track1->GetPxPyPz(mom1);
671 track2->GetPxPyPz(mom2);
672 TVector3 momv1(mom1[0],mom1[1],mom1[2]);
673 TVector3 momv2(mom2[0],mom2[1],mom2[2]);
674 Float_t theta = (float)(momv1.Phi()-momv2.Phi());
675 if(theta<-0.5*TMath::Pi()) theta+=2.*TMath::Pi();
676
677 Float_t deltaPhi = track1->Phi()-track2->Phi();
678 if(deltaPhi<-0.5*TMath::Pi()) deltaPhi+=2.*TMath::Pi();
679
680 Float_t rIsol = (float)(TMath::Sqrt( deltaPhi*deltaPhi+(track1->Eta()-track2->Eta())*(track1->Eta()-track2->Eta()) ));
681 if(rIsol<fRIsolMinCosmic) continue;
682
683 if(TMath::Abs(TMath::Pi()-theta)<fMaxCosmicAngle) {
75946d5d 684 nCosmicCandidates+=1;
df7848fc 685 isCosmic = kTRUE;
686 }
687
688 }
689 }
690
691 fh1NCosmicsPerEvent->Fill((float)nCosmicCandidates);
692
693 return isCosmic;
f2dd0695 694}
695
696
f4132e7d 697Int_t AliAnalysisTaskJetServices::GetEventClass(AliESDEvent *esd){
698
699 Float_t cent = 999;
700 if(esd->GetCentrality()){
701 cent = esd->GetCentrality()->GetCentralityPercentile("V0M");
702 }
703 if(cent>50)return 4;
704 if(cent>30)return 3;
705 if(cent>10)return 2;
706 return 1;
bf7b8731 707
f4132e7d 708}
f2dd0695 709
bf7b8731 710void AliAnalysisTaskJetServices::Terminate(Option_t */*option*/)
711{
712 // Terminate analysis
713 //
bf7b8731 714}
f4132e7d 715