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