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