]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/JetTasks/AliAnalysisTaskJetChem.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Bastian Bathen:
[u/mrichter/AliRoot.git] / PWG4 / JetTasks / AliAnalysisTaskJetChem.cxx
CommitLineData
9b2de807 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id:$ */
17
18
19#include <TChain.h>
20#include <TFile.h>
21#include <TList.h>
22#include <TH1F.h>
23#include <TH2F.h>
24#include <TProfile.h>
25#include <TVector3.h>
26#include <TLorentzVector.h>
27#include <TMath.h>
28#include <TTree.h>
29#include <TDatabasePDG.h>
30
31#include "AliAnalysisTaskJetChem.h"
32#include "AliAnalysisManager.h"
33#include "AliMCEventHandler.h"
34#include "AliAODInputHandler.h"
35#include "AliAODHandler.h"
36#include "AliAODJet.h"
37#include "AliAODTrack.h"
38#include "AliAODMCParticle.h"
39
40#include "AliGenPythiaEventHeader.h"
41#include "AliAnalysisHelperJetTasks.h"
42#include "AliESDtrack.h"
43#include "AliExternalTrackParam.h"
44
9280dfa6 45//#include <iostream> // TEST!!!
46
47
9b2de807 48//
49// Analysis class for jet chemistry studies
50// based on AliAnalysisTaskUE by Arian Abrahantes Quintana and Enesto Lopez
51// contact: Oliver Busch, o.busch@gsi.de
52
53// This class needs as input AOD with track and Jets
54// the output is a list of histograms
55//
56// AOD can be either connected to the InputEventHandler
57// for a chain of AOD files
58// or
59// to the OutputEventHandler
60// for a chain of ESD files, so this case class should be
61// in the train after the Jet finder
62//
63
64
65
66
67ClassImp(AliAnalysisTaskJetChem)
68
69////////////////////////////////////////////////////////////////////////
70
71
72//____________________________________________________________________
73AliAnalysisTaskJetChem:: AliAnalysisTaskJetChem(const char* name): AliAnalysisTaskSE(name),
74fDebug(kFALSE),
75fDeltaAOD(kFALSE),
76fDeltaAODBranch(""),
77fAODBranch("jets"),
78fDeltaAODBranchMC(""),
79fAODBranchMC("jetsMC"),
80fArrayJetsAOD(0x0),
81fArrayJetsMC(0x0),
82fAOD(0x0),
83fAODjets(0x0),
84fListOfHistos(0x0),
85fJetsOnFly(kFALSE),
86fUseLOConeJets(kFALSE),
87fUseLOConeMCJets(kFALSE),
88fUsePythiaJets(kFALSE),
89fConeRadius(0.),
90fTrackPtCutJF(0),
91fFilterBitJF(0x01),
92fRequireITSRefitJF(kFALSE),
93fRejectK0TracksJF(kFALSE),
94fJetPtCut(5.0),
95fJetEtaCut(0.5),
96fFilterBit(0x10),
97fTrackPtCut(0.),
98fTrackEtaCut(0.9),
99fUseOnFlyV0s(kFALSE),
100fCutnSigdEdx(0.),
101fUseAODMCTracksForUE(kFALSE),
102fAreaReg(1.0),
103fAvgTrials(1),
9280dfa6 104fhPrimVertexNCont(0x0),
9b2de807 105fhPrimVertexRho(0x0),
106fhPrimVertexZ(0x0),
107fhNJets(0x0),
108fhNJetsMC(0x0),
9280dfa6 109fhLeadingEta(0x0),
110fhLeadingNTracksVsEta(0x0),
111fhLeadingPtVsEta(0x0),
9b2de807 112fhLeadingPhi(0x0),
113fhLeadingPt(0x0),
114fhLeadingPtDiffr(0x0),
115fhLeadingEtaMC(0x0),
116fhLeadingPhiMC(0x0),
117fhLeadingPtMC(0x0),
118fhLeadingPtMCDiffr(0x0),
9280dfa6 119fhPhiEtaTracksNoCut(0x0),
120fhPtTracksNoCut(0x0),
121fhPhiEtaTracks(0x0),
122fhPtTracks(0x0),
123fhTrackMult(0x0),
9b2de807 124fhEtaMCTracks(0x0),
125fhPhiMCTracks(0x0),
126fhPtMCTracks(0x0),
127fhnTracksVsPtLeading(0x0),
128fhdNdEtaPhiDist(0x0),
129fhRegionSumPtMaxVsEt(0x0),
130fhRegionMultMaxVsEt(0x0),
131fhRegionSumPtMinVsEt(0x0),
132fhRegionMultMinVsEt(0x0),
133fhNV0s(0x0),
134fhV0onFly(0x0),
135fhV0DCADaughters(0x0),
136fhV0Radius(0x0),
137fhV0DCAToVertex(0x0),
138fhV0DCAToVertexK0(0x0),
139fhV0InvMassK0(0x0),
140fhV0InvMassK0JetEvt(0x0),
141fhV0InvMassLambda(0x0),
142fhV0InvMassAntiLambda(0x0),
143fhV0InvMassLambdaJetEvt(0x0),
144fhV0InvMassAntiLambdaJetEvt(0x0),
145fhdROpanK0VsPt(0x0),
146fhdPhiJetV0(0x0),
147fhdPhiJetK0(0x0),
148fhdRJetK0(0x0),
149fhdNdptV0(0x0),
150fhdNdptK0(0x0),
151fhPtVsEtaK0(0x0),
152fhV0InvMassK0DCA(0x0),
153fhV0InvMassK0DCAdEdx(0x0),
154fhV0InvMassK0DCAPID(0x0),
155fhV0InvMassLambdaDCAdEdx(0x0),
156fhV0InvMassAntiLambdaDCAdEdx(0x0),
157fhdNdptK0DCA(0x0),
158fhdNdptK0DCAdEdx(0x0),
159fhV0InvMassK0Min(0x0),
160fhV0InvMassLambdaMin(0x0),
161fhV0InvMassAntiLambdaMin(0x0),
162fhV0InvMassK0Max(0x0),
163fhV0InvMassLambdaMax(0x0),
164fhV0InvMassAntiLambdaMax(0x0),
165fhV0InvMassK0Jet(0x0),
166fhV0InvMassLambdaJet(0x0),
167fhV0InvMassAntiLambdaJet(0x0),
168fhV0InvMassK0Lambda(0x0),
169fhdNdptK0JetEvt(0x0),
9280dfa6 170fhdNdzK0(0x0),
171fhdNdzK05to10(0x0),
172fhdNdzK010to20(0x0),
173fhdNdzK020to30(0x0),
174fhdNdzK030to40(0x0),
175fhdNdzK040to60(0x0),
9b2de807 176fhdNdxiK0(0x0),
9280dfa6 177fhdNdzLambda(0x0),
178fhdNdzAntiLambda(0x0),
9b2de807 179fhdNdzK0Max(0x0),
180fhdNdxiK0Max(0x0),
181fhdNdzLambdaMax(0x0),
182fhdNdxiLambdaMax(0x0),
183fhdNdptK0Max(0x0),
184fhdNdptLambdaMax(0x0),
185fhdNdzK0Min(0x0),
186fhdNdxiK0Min(0x0),
187fhdNdzLambdaMin(0x0),
188fhdNdxiLambdaMin(0x0),
189fhdNdptK0Min(0x0),
190fhdNdptLambdaMin(0x0),
191fhdNdzK0Jet(0x0),
192fhdNdxiK0Jet(0x0),
193fhdNdzLambdaJet(0x0),
194fhdNdxiLambdaJet(0x0),
195fhdNdptK0Jet(0x0),
196fhdNdptLambdaJet(0x0),
197fhdEdxVsMomV0(0x0),
198fhdEdxVsMomV0pidEdx(0x0),
199fhdEdxVsMomV0piPID(0x0),
200fhdPhiJetK0MC(0x0),
201fhdRJetK0MC(0x0),
9b2de807 202fhdRV0MC(0x0),
9280dfa6 203fhdNdptchPiMCMax(0x0),
204fhdNdptK0MCMax(0x0),
205fhdNdptchKMCMax(0x0),
206fhdNdptpMCMax(0x0),
207fhdNdptpBarMCMax(0x0),
208fhdNdptLambdaMCMax(0x0),
209fhdNdptLambdaBarMCMax(0x0),
210fhdNdptchPiMCMin(0x0),
211fhdNdptK0MCMin(0x0),
212fhdNdptchKMCMin(0x0),
213fhdNdptpMCMin(0x0),
214fhdNdptpBarMCMin(0x0),
215fhdNdptLambdaMCMin(0x0),
216fhdNdptLambdaBarMCMin(0x0),
217fhdNdptOmegaMCMin(0x0),
218fhdNdptOmegaBarMCMin(0x0),
219fhdNdptchPiMCJet(0x0),
220fhdNdptK0MCJet(0x0),
221fhdNdptchKMCJet(0x0),
222fhdNdptpMCJet(0x0),
223fhdNdptpBarMCJet(0x0),
224fhdNdptLambdaMCJet(0x0),
225fhdNdptLambdaBarMCJet(0x0),
9b2de807 226fhPIDMC(0x0),
9280dfa6 227fhPIDMC_quarkEv(0x0),
228fhPIDMC_gluonEv(0x0),
9b2de807 229fhPIDMCAll(0x0),
230fhPIDMCMin(0x0),
231fhPIDMCJet(0x0),
9280dfa6 232fhPIDMCMotherK0(0x0),
233fhPIDMCGrandMotherK0(0x0),
234fhPIDMCMotherChK(0x0),
235fhPIDMCMotherK0Trans(0x0),
236fhPIDMCGrandMotherK0Trans(0x0),
237fhPIDMCMotherChKTrans(0x0),
9b2de807 238fhdNdptgammaMC(0x0),
239fhdNdptchPiMC(0x0),
240fhdNdptpi0MC(0x0),
241fhdNdptK0MC(0x0),
242fhdNdptchKMC(0x0),
243fhdNdptpMC(0x0),
244fhdNdptpBarMC(0x0),
9280dfa6 245fhdNdptLambdaMC(0x0),
246fhdNdptLambdaBarMC(0x0),
9b2de807 247fhdNdptOmegaMC(0x0),
248fhdNdptOmegaBarMC(0x0),
249fhdNdxiMC(0x0),
250fhdNdxiK0MC(0x0),
251fhdNdxiK0MCJet(0x0),
252fhdNdzK0MC(0x0),
253fhdNdzK0MCJet(0x0),
254fhdNdptK0MCJetEvt(0x0),
255fhnJetsAODvsMC(0x0),
256fhLeadingPtAODvsMC(0x0),
257fhLeadingEtaAODvsMC(0x0),
258fhLeadingPhiAODvsMC(0x0),
259fhnTracksLeadingAODvsMC(0x0),
260fhLeadingdRAODMC(0x0),
261fhLeadingPtAODvsMCdRcut(0x0),
262fhdnTracksVsdPtLeadingAODMC(0x0),
263fhnTracksJetVsPtAOD(0x0),
264fhnTracksJetVsPtAODquarkEv(0x0),
265fhRadiusJetVsPtAOD(0x0),
266fhnTracksJetVsPtMC(0x0),
267fhnTracksJetVsPtMCquarkEv(0x0),
268fhRadiusJetVsPtMC(0x0),
269fhnTracksJetVsPtMCK0(0x0),
270fhnTracksJetVsPtMCK0quarkEv(0x0),
271fhRadiusJetVsPtMCK0(0x0),
272fhnTracksJetVsPtAODK0(0x0),
273fhnTracksJetVsPtAODK0quarkEv(0x0),
274fhRadiusJetVsPtAODK0(0x0),
275fhnTracksJetVsPtAODpKch(0x0),
276fhRadiusJetVsPtAODpKch(0x0),
277fhPythiaProcess(0x0),
278fhPythiaProcessK0(0x0),
279fhPythiaProcessKch(0x0),
280fhPythiaProcessp(0x0),
281fhPythiaProcesspbar(0x0),
9280dfa6 282fhdNdzJets5to10(0x0),
283fhdNdzJets10to20(0x0),
284fhdNdzJets20to30(0x0),
285fhdNdzJets30to40(0x0),
286fhdNdzJets40to60(0x0),
287fhdNdxiJets5to10(0x0),
288fhdNdxiJets10to20(0x0),
289fhdNdxiJets20to30(0x0),
290fhdNdxiJets30to40(0x0),
291fhdNdxiJets40to60(0x0),
292fhdNdptTracksJetPt5to10(0x0),
293fhdNdptTracksJetPt10to20(0x0),
294fhdNdptTracksJetPt20to30(0x0),
295fhdNdptTracksJetPt30to40(0x0),
296fhdNdptTracksJetPt40to60(0x0),
9b2de807 297fh1Xsec(0x0),
9f77f0ce 298fh1Trials(0x0),
299fpdgdb(0x0){
9b2de807 300 // Default constructor
301
302 fAreaReg = 2*TMath::Pi()/6.0 * 2*fTrackEtaCut;
303 fpdgdb = TDatabasePDG::Instance();
304
9b2de807 305 // Output slot #1 writes into a TList container, 0 reserved for std AOD output
306 DefineOutput(1, TList::Class());
307}
308
309//______________________________________________________________
310Bool_t AliAnalysisTaskJetChem::UserNotify()
311{
312 //
313 // read the cross sections
314 // and number of trials from pyxsec.root
315 //
316
317 fAvgTrials = 1;
318 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
319 Float_t xsection = 0;
320 Float_t ftrials = 1;
321 if(tree){
322 TFile *curfile = tree->GetCurrentFile();
323 if (!curfile) {
324 Error("Notify","No current file");
325 return kFALSE;
326 }
327 if(!fh1Xsec||!fh1Trials){
328 Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
329 return kFALSE;
330 }
331 AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
332 fh1Xsec->Fill("<#sigma>",xsection);
333 // construct average trials
334 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
335 if(ftrials>=nEntries)fAvgTrials = ftrials/nEntries;
336 }
337 return kTRUE;
338}
339
340//____________________________________________________________________
341void AliAnalysisTaskJetChem::UserCreateOutputObjects()
342{
343 // Create the output container
344 //
345 AliInfo("UserCreateOutPutObjects()");
346 //
347 // Histograms
348
349 OpenFile(1);
350 CreateHistos();
351 PostData(1, fListOfHistos); // PostData at least once
352
353}
354
355//____________________________________________________________________
356void AliAnalysisTaskJetChem::UserExec(Option_t */*option*/)
357{
358 // Execute analysis for current event
359 //
360 if (fDebug > 3) AliInfo( " Processing event..." );
361
362 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
363
364 if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
365 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
366 if (fDebug > 1) AliInfo(" ==== Tracks from AliAODInputHandler");
367 // Case when jets are reconstructed on the fly from AOD tracks
368 // (the Jet Finder is using the AliJetAODReader) of InputEventHandler
369 // and put in the OutputEventHandler AOD. Useful whe you want to reconstruct jets with
370 // different parameters to default ones stored in the AOD or to use a different algorithm
371 if( fJetsOnFly ) {
372 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
373 if( handler && handler->InheritsFrom("AliAODHandler") ) {
374 fAODjets = ((AliAODHandler*)handler)->GetAOD();
375 if (fDebug > 1) AliInfo(" ==== Jets from AliAODHandler");
376 }
377 } else {
378 fAODjets = fAOD;
379 if (fDebug > 1) AliInfo(" ==== Jets from AliAODInputHandler");
380 }
381 } else {
382 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
383 if( handler && handler->InheritsFrom("AliAODHandler") ) {
384 fAOD = ((AliAODHandler*)handler)->GetAOD();
385 fAODjets = fAOD;
386 if (fDebug > 1) AliInfo(" ==== Tracks and Jets from AliAODHandler");
387 } else {
388 AliFatal("I can't get any AOD Event Handler");
389 return;
390 }
391 }
392
393 // -------------
394
395 // fetch the pythia header info and get the trials
396 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
397 Float_t nTrials = 1;
398 if (mcHandler) {
399 AliMCEvent* mcEvent = mcHandler->MCEvent();
400 if (mcEvent) {
401 AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
402 if(pythiaGenHeader){
403 nTrials = pythiaGenHeader->Trials();
404 }
405 }
406 }
407
408 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
409
410 AnalyseEvent();
411
412 // Post the data
413 PostData(1, fListOfHistos);
414}
415
416//____________________________________________________________________
417void AliAnalysisTaskJetChem::AnalyseEvent()
418{
419 // check Trigger
420
9280dfa6 421 // TString firedTriggerClasses = fAOD->GetFiredTriggerClasses();
9b2de807 422 // AliInfo(Form("firedTriggerClasses %s",firedTriggerClasses.Data()));
423 // if(firedTriggerClasses.Length() > 0 && !firedTriggerClasses.Contains("CINT1B")) return;
424
425 AliAODVertex *primVertex = fAOD->GetPrimaryVertex(); // is this from tracks or from SPD ? SPD should only be used if no track vertex
426 if(!primVertex){
9280dfa6 427 AliInfo("no prim Vertex found - skip event");
428 fhPrimVertexNCont->Fill(-1);
9b2de807 429 return;
3568c3d6 430 }
9280dfa6 431
3e3b1bae 432 TString primVtxName(primVertex->GetName());
433
434 if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
435
3568c3d6 436 AliInfo("found TPC prim Vertex - skip event");
437 fhPrimVertexNCont->Fill(-1);
438 return;
439 }
9280dfa6 440
3568c3d6 441 Int_t nVertContributors = primVertex->GetNContributors();
442 fhPrimVertexNCont->Fill(nVertContributors);
443
444 //cout<<" prim vertex name "<<primVertex->GetName()<<" nCont "<<nVertContributors<<endl;
445
446 if(nVertContributors<1){ // eventually check if not SPD vertex ???
447 AliInfo("prim vertex no contributors - skip event");
448 return;
449 }
9280dfa6 450
3568c3d6 451 Double_t vertX = primVertex->GetX();
452 Double_t vertY = primVertex->GetY();
453 Double_t vertZ = primVertex->GetZ();
454
455 Double_t vertRho = TMath::Sqrt(vertX*vertX+vertY*vertY);
456
457 fhPrimVertexRho->Fill(vertRho);
458 fhPrimVertexZ->Fill(vertZ);
459
460 if(TMath::Abs(vertZ)>10){
461 AliInfo(Form("prim vertex z=%f - skip event",vertZ));
462 return;
9b2de807 463 }
464
3568c3d6 465
9b2de807 466 Int_t pythiaPID = GetPythiaProcessID();
9280dfa6 467
9b2de807 468 // ------------------------------------------------
469 // Find Leading Jets 1,2,3
470 // (could be skipped if Jets are sort by Pt...)
471
472 //Int_t index1 = -1;
473 //Int_t nTracksLeading = 0;
474
475 Int_t nJetsAOD = 0;
9280dfa6 476 AliAODJet* leadingJetAOD = NULL; // non-zero if any jet in acc and leading jet survives pt cut
9b2de807 477 Int_t indexLeadingAOD = -1;
478 Double_t ptLeadingAOD = 0.;
479 Int_t indexMaxRegionAOD = 0; // initialize with '0', '+-1' transverse regions
480
481 Int_t nJetsMC = 0;
482 AliAODJet* leadingJetMC = NULL; // non-zero if leading jet survives acc, pt cut
483 Int_t indexLeadingMC = -1;
484 Double_t ptLeadingMC = 0.;
485 Int_t indexMaxRegionMC = 0; // initialize with '0', '+-1' transverse regions
486
9280dfa6 487 Double_t ptLeadingAODAllEta = 0;
488 AliAODJet* leadingJetAODAllEta = NULL;
489
9b2de807 490 if(fUseLOConeJets) fArrayJetsAOD = FindChargedParticleJets();
491 else{ // Use jets on AOD/DeltaAOD
492
493 if(fDeltaAOD){
494 if (fDebug > 1) AliInfo(" ==== Jets From Delta-AODs !");
495 if (fDebug > 1) AliInfo(Form(" ==== Reading Branch: %s ",fDeltaAODBranch.Data()));
496 fArrayJetsAOD = (TClonesArray*)fAODjets->GetList()->FindObject(fDeltaAODBranch.Data());
497 if (!fArrayJetsAOD){
498 AliFatal(" No jet-array! ");
499 return;
500 }
501 }
502 else{
503 if (fDebug > 1) AliInfo(" ==== AOD jets: Read Standard-AODs !");
504 if (fDebug > 1) AliInfo(Form(" ==== Reading Branch: %s ",fAODBranch.Data()));
505
506 fArrayJetsAOD = (TClonesArray*)fAODjets->FindListObject(fAODBranch.Data());
507 }
508 }
509
510
511 if(fUseLOConeMCJets) fArrayJetsMC = FindChargedParticleJetsMC();
512 else if(fUsePythiaJets) fArrayJetsMC = GetPythiaJets();
513 else{
514
515 if(fDeltaAOD){
516 if (fDebug > 1) AliInfo(" ==== MC Jets From Delta-AODs !");
517 if (fDebug > 1) AliInfo(Form(" ==== Reading Branch: %s ",fDeltaAODBranchMC.Data()));
518 fArrayJetsMC = (TClonesArray*)fAODjets->GetList()->FindObject(fDeltaAODBranchMC.Data());
519 if (!fArrayJetsMC){
520 AliFatal(" No jet-array! ");
521 return;
522 }
523 }
524 else{
525 if (fDebug > 1) AliInfo(" ==== MC jets: Read Standard-AODs !");
526 if (fDebug > 1) AliInfo(Form(" ==== Reading Branch: %s ",fAODBranchMC.Data()));
527
528 fArrayJetsMC = (TClonesArray*)fAODjets->FindListObject(fAODBranchMC.Data());
529 }
530 }
531
532 if(fArrayJetsAOD) nJetsAOD = fArrayJetsAOD->GetEntries();
533 if(fArrayJetsMC) nJetsMC = fArrayJetsMC->GetEntries();
534
535 fhNJets->Fill(nJetsAOD);
536 fhNJetsMC->Fill(nJetsMC);
537 fhnJetsAODvsMC->Fill(nJetsMC,nJetsAOD);
538
539 if(fDebug>1) AliInfo(Form("AOD %d jets",nJetsAOD));
540
541
9280dfa6 542 // for xcheck: find leading jet in large eta range
543 for(Int_t i=0; i<nJetsAOD; i++){
544
545 AliAODJet* jet = (AliAODJet*) fArrayJetsAOD->At(i);
546 Double_t jetPt = jet->Pt();
547
548 if(jetPt > ptLeadingAODAllEta){
549 ptLeadingAODAllEta = jetPt;
550 leadingJetAODAllEta = jet;
551 }
552 }
553
554 if(leadingJetAODAllEta){
555 //cout<<" trackRefs entries "<<leadingJetAODAllEta->GetRefTracks()->GetEntriesFast()<<endl;
556 fhLeadingNTracksVsEta->Fill(leadingJetAODAllEta->Eta(),leadingJetAODAllEta->GetRefTracks()->GetEntriesFast());
557 fhLeadingPtVsEta->Fill(leadingJetAODAllEta->Eta(),leadingJetAODAllEta->Pt());
558
559 }
560
9b2de807 561 // find leading jet AOD
562 for(Int_t i=0; i<nJetsAOD; ++i){
563
564 AliAODJet* jet = (AliAODJet*) fArrayJetsAOD->At(i);
565 Double_t jetPt = jet->Pt();
566 Double_t jetEta = jet->Eta();
567
568 if((jetPt > ptLeadingAOD) && (TMath::Abs(jetEta)<fJetEtaCut)){
569 ptLeadingAOD = jetPt;
570 indexLeadingAOD = i;
571 }
572 }
573
574 // find leading jet MC
575 for(Int_t i=0; i<nJetsMC; ++i){
576
577 AliAODJet* jet = (AliAODJet*) fArrayJetsMC->At(i);
578 Double_t jetPt = jet->Pt();
579 Double_t jetEta = jet->Eta();
580
581 if((jetPt > ptLeadingMC) && (TMath::Abs(jetEta)<fJetEtaCut)){
582 ptLeadingMC = jetPt;
583 indexLeadingMC = i;
584 }
585 }
586
587
588 //cout<<" ptLeadingAOD "<<ptLeadingAOD<<" MC "<<ptLeadingMC<<endl;
589
590 if(indexLeadingAOD>=0){ // event with jet
591
592 Double_t etaLeadingAOD = ((AliAODJet*) fArrayJetsAOD->At(indexLeadingAOD))->Eta();
593 Double_t phiLeadingAOD = ((AliAODJet*) fArrayJetsAOD->At(indexLeadingAOD))->Phi();
594 Int_t nTracksLeadingAOD = ((AliAODJet*) fArrayJetsAOD->At(indexLeadingAOD))->GetRefTracks()->GetEntriesFast();
595
596 if(fDebug>1) AliInfo(Form("\n Pt Leading AOD Jet = %6.1f eta=%5.3f, nTracks %d",ptLeadingAOD,etaLeadingAOD,nTracksLeadingAOD));
597
598 fhLeadingEta->Fill(etaLeadingAOD);
9280dfa6 599
9b2de807 600 if(TMath::Abs(etaLeadingAOD)<fJetEtaCut){ // leading jet eta cut
601
602 fhnTracksVsPtLeading->Fill(ptLeadingAOD,nTracksLeadingAOD);
603 fhLeadingPt->Fill(ptLeadingAOD);
604 if(IsDiffractiveEvent(pythiaPID)) fhLeadingPtDiffr->Fill(ptLeadingAOD);
605
606 if(ptLeadingAOD>fJetPtCut){ // leading jet pt cut
607
608 leadingJetAOD = (AliAODJet*) fArrayJetsAOD->At(indexLeadingAOD);
609
9280dfa6 610 fhLeadingPhi->Fill(phiLeadingAOD);
9b2de807 611
612 // ----------------------------------------------
613 // Find max and min regions
614 Double_t sumPtRegionPosit = 0.;
615 Double_t sumPtRegionNegat = 0.;
616 Int_t nTrackRegionPosit = 0;
617 Int_t nTrackRegionNegat = 0;
618
619 Int_t nTracks = fAOD->GetNTracks();
620
621 for (Int_t ipart=0; ipart<nTracks; ++ipart) {
622
623 AliAODTrack* part = fAOD->GetTrack(ipart);
624 if ( !part->TestFilterBit(fFilterBit) ) continue; // track cut selection
625 if (!part->IsPrimaryCandidate()) continue; // reject whatever is not linked to collision point
626 // PID Selection: Reject everything but hadrons
627 Bool_t isHadron = part->GetMostProbablePID()==AliAODTrack::kPion ||
628 part->GetMostProbablePID()==AliAODTrack::kKaon ||
629 part->GetMostProbablePID()==AliAODTrack::kProton;
630 if (!isHadron ) continue;
631 if (!part->Charge() ) continue; //Only charged
632 if (part->Pt() < fTrackPtCut ) continue;
633 if(TMath::Abs(part->Eta()) > fTrackEtaCut ) continue;
634
635 TVector3 partVect(part->Px(), part->Py(), part->Pz());
636
637 Int_t region = IsTrackInsideRegion(leadingJetAOD,&partVect );
638
639 if (region > 0) {
640 sumPtRegionPosit += part->Pt();
641 nTrackRegionPosit++;
642 }
643 if (region < 0) {
644 sumPtRegionNegat += part->Pt();
645 nTrackRegionNegat++;
646 }
647 } // tracks loop
648
649 // fill sumPt and mult density for transverse regions
650
651 if( sumPtRegionPosit > sumPtRegionNegat ) {
652 FillSumPtRegion( ptLeadingAOD, sumPtRegionPosit/fAreaReg, sumPtRegionNegat/fAreaReg );
653 }
654 else {
655 FillSumPtRegion( ptLeadingAOD, sumPtRegionNegat/fAreaReg, sumPtRegionPosit/fAreaReg );
656 }
657 if (nTrackRegionPosit > nTrackRegionNegat ) {
658 FillMultRegion( ptLeadingAOD, nTrackRegionPosit/fAreaReg, nTrackRegionNegat/fAreaReg);
659 }
660 else {
661 FillMultRegion( ptLeadingAOD, nTrackRegionNegat/fAreaReg, nTrackRegionPosit/fAreaReg);
662 }
663
664 indexMaxRegionAOD = (sumPtRegionPosit > sumPtRegionNegat) ? 1 : -1;
665
666 } // leading jet pt cut
667 } // leading jet eta cut
668 } // jet event
669
670
671 if(indexLeadingMC>=0){ // event with MC jet
672
673 Double_t etaLeadingMC = ((AliAODJet*) fArrayJetsMC->At(indexLeadingMC))->Eta();
674 Double_t phiLeadingMC = ((AliAODJet*) fArrayJetsMC->At(indexLeadingMC))->Phi();
675 Int_t nTracksLeadingMC = ((AliAODJet*) fArrayJetsMC->At(indexLeadingMC))->GetRefTracks()->GetEntriesFast();
676
677 if(fDebug>1) AliInfo(Form("\n Pt Leading MC Jet = %6.1f eta=%5.3f, nTracks %d",ptLeadingMC,etaLeadingMC,nTracksLeadingMC));
678
679 fhLeadingEtaMC->Fill(etaLeadingMC);
680
681 if(TMath::Abs(etaLeadingMC)<fJetEtaCut){ // leading jet eta cut
682
683 fhLeadingPtMC->Fill(ptLeadingMC);
684 if(IsDiffractiveEvent(pythiaPID)) fhLeadingPtMCDiffr->Fill(ptLeadingMC);
685
686 if(ptLeadingMC>fJetPtCut){ // leading jet pt cut
687
688 leadingJetMC = (AliAODJet*) fArrayJetsMC->At(indexLeadingMC);
689
9280dfa6 690 fhLeadingPhiMC->Fill(phiLeadingMC); // -pi to pi
9b2de807 691
692 // ----------------------------------------------
693 // Find max and min regions
694 Double_t sumPtRegionPosit = 0;
695 Double_t sumPtRegionNegat = 0;
696 Int_t nTrackRegionPosit = 0;
697 Int_t nTrackRegionNegat = 0;
698
699 TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
700
701 Int_t ntrks = farray->GetEntries();
702 if (fDebug>1) AliInfo(Form("In UE MC analysis tracks %d \n",ntrks));
703
704 for(Int_t i =0 ; i < ntrks; i++){
705
706 AliAODMCParticle* mctrk = (AliAODMCParticle*)farray->At(i);
707 //Cuts
708 if (!(mctrk->IsPhysicalPrimary())) continue;
709 //if (!(mctrk->IsPrimary())) continue;
710
711 if (mctrk->Charge() == 0 || mctrk->Charge()==-99) continue;
712
713 if (mctrk->Pt() < fTrackPtCut ) continue;
714 if( TMath::Abs(mctrk->Eta()) > fTrackEtaCut ) continue;
715
716 Bool_t isHadron = TMath::Abs(mctrk->GetPdgCode())==211 ||
717 TMath::Abs(mctrk->GetPdgCode())==2212 ||
718 TMath::Abs(mctrk->GetPdgCode())==321;
719
720 if (!isHadron) continue;
721
722 TVector3 partVect(mctrk->Px(), mctrk->Py(), mctrk->Pz());
723
724 Int_t region = IsTrackInsideRegion(leadingJetMC,&partVect );
725
726 if (region > 0) {
727 sumPtRegionPosit += mctrk->Pt();
728 nTrackRegionPosit++;
729 }
730 if (region < 0) {
731 sumPtRegionNegat += mctrk->Pt();
732 nTrackRegionNegat++;
733 }
734 } // AliAODMCParticle loop
735
736 indexMaxRegionMC = (sumPtRegionPosit > sumPtRegionNegat) ? 1 : -1;
737
738 } // leading jet pt
739 } // leading jet eta
740 } // found jet
741
742
743 Bool_t foundK0AOD = kFALSE; // event with K0 ?
744 Bool_t foundK0MC = kFALSE; // event with K0 ?
745
746 CheckV0s(leadingJetAOD,indexMaxRegionAOD,foundK0AOD); // here leadingJetAOD/MC nonzero if jet passes eta & pt cut
747 CheckMCParticles(leadingJetMC,indexMaxRegionMC,foundK0MC);
748 CompLeadingJets(leadingJetAOD,leadingJetMC,pythiaPID,foundK0AOD,foundK0MC);
9280dfa6 749 FillReferencePlotsTracks();
750 FillReferenceFF(leadingJetAOD);
751
9b2de807 752
753 if(fUseLOConeJets && fArrayJetsAOD){
754 fArrayJetsAOD->Delete(); // no 'Clear': AliAODjet contains TMomentum and TRefArray
2ba76308 755 delete fArrayJetsAOD;
9b2de807 756 }
757 if(fUseLOConeMCJets && fArrayJetsMC){
758 fArrayJetsMC->Delete(); // no 'Clear': AliAODjet contains TMomentum and TRefArray
2ba76308 759 delete fArrayJetsMC;
9b2de807 760 }
761}
762
763// __________________________________________________________________
764
765
766Double_t AliAnalysisTaskJetChem::GetJetRadius(const AliAODJet* jet, const Double_t energyFrac){
767
768 // calc jet radius containing fraction energyFrac of full jet pt
769
9f77f0ce 770 const Int_t kArraySize = 1000;
771 const Int_t kInitVal = -999;
772
9b2de807 773 Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
774
9f77f0ce 775 if(nTracks>kArraySize){
776 AliError(Form("nTracks in jet %d exceeds max array size",nTracks));
777 return -1;
778 }
779
780
781 Double_t deltaR[kArraySize];
782 Double_t pt[kArraySize];
783 Int_t index[kArraySize];
784 for(int i=0; i<kArraySize; i++) index[i] = kInitVal;
9b2de807 785
9f77f0ce 786
9b2de807 787 Double_t ptTot = 0;
788
789 for(int i=0; i<nTracks; i++){
790
791 AliAODTrack* track = (AliAODTrack*) jet->GetRefTracks()->At(i);
792
793 TLorentzVector *mom4Jet = jet->MomentumVector();
794 TVector3 mom3Jet = mom4Jet->Vect();
795
796 Double_t trackMom[3];
797 track->PxPyPz(trackMom);
798 TVector3 mom3Track(trackMom);
799
800 Double_t dR = mom3Jet.DeltaR(mom3Track);
801
802 deltaR[i] = dR;
803 pt[i] = track->Pt();
804
805 ptTot += pt[i];
806 }
807
808 //cout<<" ptTot "<<ptTot<<" jetPt "<<jet->Pt()<<endl; // Xcheck
809
810 TMath::Sort(nTracks,deltaR,index,kFALSE); // sort in decreasing order
811
812 Double_t ptSum = 0;
813
814 for(int i=0; i<nTracks; i++){
815
816 Int_t ind = index[i];
817
818 Double_t ptTrack = pt[ind];
819 Double_t dR = deltaR[ind];
820
821 ptSum += ptTrack;
822
823 if(ptSum >= ptTot*energyFrac) return dR;
824 }
825
826 return -1;
827
828}
829
830//____________________________________________________________________
831void AliAnalysisTaskJetChem::FillSumPtRegion( Double_t leadingE, Double_t ptMax, Double_t ptMin )
832{
833 // Fill sumPt of control regions
834
835 fhRegionSumPtMaxVsEt->Fill( leadingE, ptMax );
836 fhRegionSumPtMinVsEt->Fill( leadingE, ptMin );
837}
838
839//____________________________________________________________________
840void AliAnalysisTaskJetChem::FillMultRegion(Double_t leadingE, Double_t nTrackPtmax, Double_t nTrackPtmin)
841{
842 // Fill Nch multiplicity of control regions
843
844 fhRegionMultMaxVsEt->Fill( leadingE, nTrackPtmax );
845 fhRegionMultMinVsEt->Fill( leadingE, nTrackPtmin );
846
847}
848
849//____________________________________________________________________
850Int_t AliAnalysisTaskJetChem::IsTrackInsideRegion(const AliAODJet* aodjetVect,const TVector3 *partVect)
851{
852 // return de region in delta phi
853 // -1 negative delta phi
854 // 1 positive delta phi
855 // 0 outside region
856
857 TLorentzVector* jetVectLorentz = aodjetVect->MomentumVector();
858 TVector3 jetVect = jetVectLorentz->Vect();
859
9b2de807 860 static const Double_t k60rad = 60.*TMath::Pi()/180.;
861 static const Double_t k120rad = 120.*TMath::Pi()/180.;
862
863 Int_t region = 0;
864 if( TMath::Abs(partVect->Eta()) > fTrackEtaCut ) return 0;
865 // transverse regions
866 if (jetVect.DeltaPhi(*partVect) < -k60rad && jetVect.DeltaPhi(*partVect) > -k120rad ) region = -1;
867 if (jetVect.DeltaPhi(*partVect) > k60rad && jetVect.DeltaPhi(*partVect) < k120rad ) region = 1;
868
869 return region;
870}
871
872//____________________________________________________________________
873TClonesArray* AliAnalysisTaskJetChem::FindChargedParticleJetsMC()
874{
875 // CDF jet finder:
876 // loop over pt-ordered list of tracks, combine tracks within jet cone, recalc cone axis after each step
877 // based on implementation by Arian Abrahantes Quintana and Enesto Lopez
878
879 TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
880 if(!farray){
881 AliInfo("no mcparticles branch");
882 return 0;
883 }
884
885 Int_t nTracks = farray->GetEntries();
886
887 if( !nTracks ) return 0;
888 TObjArray tracks(nTracks);
889
890 for (Int_t ipart=0; ipart<nTracks; ++ipart) {
891
892 AliAODMCParticle* part = (AliAODMCParticle*)farray->At(ipart);
893
894 if(!part->IsPhysicalPrimary()) continue;
895
896 // exclude neutrinos
897 Int_t pdg = TMath::Abs(part->GetPdgCode());
898 if((pdg == 12 || pdg == 14 || pdg == 16)) return kFALSE;
899
9280dfa6 900 if( !part->Charge() ) continue; // comment / uncomment here
9b2de807 901 fhEtaMCTracks->Fill(part->Eta());
902
903 if(TMath::Abs(part->Eta()) > fTrackEtaCut) continue;
904 fhPtMCTracks->Fill(part->Pt());
905
906 if( part->Pt() < fTrackPtCutJF ) continue;
9280dfa6 907 fhPhiMCTracks->Fill(part->Phi());
9b2de807 908
909 tracks.AddLast(part);
910
911 }
912
913 QSortTracks( tracks, 0, tracks.GetEntriesFast() );
914
915 nTracks = tracks.GetEntriesFast();
916
917 if( !nTracks ) return 0;
918 TObjArray *jets = new TObjArray(nTracks);
919 TIter itrack(&tracks);
920 while( nTracks ) {
921 // 2- Start with the highest pT particle ...
922 Float_t px,py,pz,pt;
923
924 AliAODMCParticle* track = (AliAODMCParticle*)itrack.Next();
925
926 if( !track ) continue;
927 px = track->Px();
928 py = track->Py();
929 pz = track->Pz();
930 pt = track->Pt();
931 jets->AddLast( new AliAODJet(px,py,pz,pt) ); // Use the energy member to store Pt
932 AliAODJet* jet = (AliAODJet*)jets->Last();
933 jet->AddTrack(track);
934 tracks.Remove(track);
935
936
937 TVector2 jetVect2(jet->Px(),jet->Py());
938
939 // 3- Go to the next highest pT particle not already included...
940 AliAODMCParticle* track1;
941 while ( (track1 = (AliAODMCParticle*)(itrack.Next())) ) {
942 TVector2 trackVect2(track1->Px(),track1->Py());
943 Double_t dphi = trackVect2.DeltaPhi(jetVect2);
944
945 Double_t r = TMath::Sqrt( (jet->Eta()-track1->Eta())*(jet->Eta()-track1->Eta()) +
946 dphi*dphi );
947
948 if( r < fConeRadius ) {
949 Double_t fPt = jet->E()+track1->Pt(); // Scalar sum of Pt
950 // recalculating the centroid
951 Double_t eta = jet->Eta()*jet->E()/fPt + track1->Eta()*track1->Pt()/fPt;
952
953 jetVect2.SetMagPhi(jet->E()/fPt,jetVect2.Phi());
954 trackVect2.SetMagPhi(track1->Pt()/fPt,trackVect2.Phi());
955
956 TVector2 sumVect2 = jetVect2+trackVect2;
957 Double_t phi = sumVect2.Phi();
958
959 //jet->SetPtEtaPhiE( 1., eta, phi, fPt );
960 ((TLorentzVector*) jet->MomentumVector())->SetPtEtaPhiE(1,eta,phi,fPt);
961
962 jet->AddTrack(track1);
963 tracks.Remove(track1);
964 }
965 }
966
967 tracks.Compress();
968
969 nTracks = tracks.GetEntries();
970 // 4- Continue until all particles are in a jet.
971 itrack.Reset();
972 } // end while nTracks
973
974 // Convert to AODjets....
975 Int_t njets = jets->GetEntriesFast();
976
977 TClonesArray* aodjets = new TClonesArray("AliAODJet",njets);
978 aodjets->SetOwner(kTRUE);
979
980 Int_t count = 0;
981 for(Int_t ijet=0; ijet<njets; ++ijet) {
982 AliAODJet* jet = (AliAODJet*)jets->At(ijet);
983
984 if (jet->E() < fJetPtCut) continue;
985 Float_t px, py,pz,en; // convert to 4-vector
986 px = jet->E() * TMath::Cos(jet->Phi()); // Pt * cos(phi)
987 py = jet->E() * TMath::Sin(jet->Phi()); // Pt * sin(phi)
988 pz = jet->E() / TMath::Tan(2.0 * TMath::ATan(TMath::Exp(-jet->Eta())));
989 en = TMath::Sqrt(px * px + py * py + pz * pz);
990 jet->SetPxPyPzE(px,py,pz,en);
991
992
993 TClonesArray &tmpaodjets = *aodjets;
994 new(tmpaodjets[count++]) AliAODJet(*jet);
995 //aodjets->AddLast( new AliAODJet(*jet));
996 //aodjets->AddLast( new AliAODJet(px, py, pz, en) );
997 }
998 // Order jets according to their pT .
999 QSortTracks( *aodjets, 0, aodjets->GetEntriesFast() );
1000
1001 // debug
1002 //if (fDebug>3) AliInfo(Form(" %d Charged jets found\n",njets));
1003
1004 jets->Delete(); // OB - should I cleanup or leave it to garbage collection ?
1005 delete jets;
1006
1007 return aodjets;
1008}
1009
1010// ___________________________________________________________________
1011
1012Bool_t AliAnalysisTaskJetChem::IsTrackFromK0(const Int_t indexTrack){
1013
1014 // check wether track with index is from K0 (= V0 with proper inv mass)
1015
1016 TClonesArray* tracks = fAOD->GetTracks();
1017
1018 for(int i=0; i<fAOD->GetNumberOfV0s(); i++){ // loop over V0s
1019
1020 AliAODv0* v0 = fAOD->GetV0(i);
1021
1022 Bool_t isOnFly = v0->GetOnFlyStatus();
1023
1024 if( (fUseOnFlyV0s && !isOnFly) || (!fUseOnFlyV0s && isOnFly) ) continue;
1025
1026 Double_t massK0 = v0->MassK0Short();
1027
1028 // FIXME: here should go cuts (at least V0 quality)
1029
1030 if(IsK0InvMass(massK0)){
1031
1032 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0));
1033 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1));
1034
1035 Int_t indexPos = tracks->IndexOf(trackPos);
1036 Int_t indexNeg = tracks->IndexOf(trackNeg);
1037
1038
1039 if(indexPos == indexTrack){ return kTRUE;}
1040 if(indexNeg == indexTrack){ return kTRUE;}
1041
1042 }
1043 }
1044
1045 return kFALSE;
1046
1047}
1048
1049//____________________________________________________________________
1050TClonesArray* AliAnalysisTaskJetChem::FindChargedParticleJets()
1051{
1052
1053 // CDF jet finder:
1054 // loop over pt-ordered list of tracks, combine tracks within jet cone, recalc cone axis after each step
1055 // based on implementation by Arian Abrahantes Quintana and Enesto Lopez
1056 // ref: PHYSICAL REVIEW D 65 092002, CDF Collaboration
1057
1058
1059 // 1 - Order all charged particles according to their pT .
1060 Int_t nTracks = fAOD->GetNTracks();
1061
1062 if( !nTracks ) return 0;
1063 TObjArray tracks(nTracks);
1064
1065 // standard cuts + ITS refit = stdrd PWG4 cut (compose them for productions with old ESDFilter task)
1066
1067 for (Int_t ipart=0; ipart<nTracks; ++ipart) {
1068 AliAODTrack* part = fAOD->GetTrack( ipart );
1069
1070 UInt_t status = part->GetStatus();
1071
1072 if( !part->TestFilterBit(fFilterBitJF) ) continue; // track cut selection
1073 if(fRequireITSRefitJF && ((status&AliESDtrack::kITSrefit)==0)) continue;
1074
1075 if(!part->Charge() ) continue;
9b2de807 1076
1077 if(TMath::Abs(part->Eta()) > fTrackEtaCut) continue;
1078
1079 if(fRejectK0TracksJF && IsTrackFromK0(ipart)) continue;
1080
9b2de807 1081 if( part->Pt() < fTrackPtCutJF ) continue;
1082
9b2de807 1083 tracks.AddLast(part);
1084 }
1085
1086 QSortTracks( tracks, 0, tracks.GetEntriesFast() );
1087
1088 nTracks = tracks.GetEntriesFast();
1089
1090 if( !nTracks ) return 0;
1091 TObjArray *jets = new TObjArray(nTracks);
1092 TIter itrack(&tracks);
1093 while( nTracks ) {
1094 // 2- Start with the highest pT particle ...
1095 Float_t px,py,pz,pt;
1096 AliAODTrack* track = (AliAODTrack*)itrack.Next();
1097 if( !track ) continue;
1098 px = track->Px();
1099 py = track->Py();
1100 pz = track->Pz();
1101 pt = track->Pt();
1102 //jets->AddLast( new TLorentzVector(px, py, pz, pt) ); // Use the energy member to store Pt
1103 jets->AddLast( new AliAODJet(px,py,pz,pt) ); // Use the energy member to store Pt
1104 //TLorentzVector* jet = (TLorentzVector*)jets->Last();
1105 AliAODJet* jet = (AliAODJet*)jets->Last();
1106 jet->AddTrack(track);
1107 tracks.Remove( track );
1108
1109 TVector2 jetVect2(jet->Px(),jet->Py()); // OB
1110
1111 // 3- Go to the next highest pT particle not already included...
1112 AliAODTrack* track1;
1113 while ( (track1 = (AliAODTrack*)(itrack.Next())) ) {
1114 //Double_t dphi = TVector2::Phi_mpi_pi(jet->Phi()-track1->Phi()); // OB remove
1115 TVector2 trackVect2(track1->Px(),track1->Py());
1116 Double_t dphi = trackVect2.DeltaPhi(jetVect2);
1117
1118 Double_t r = TMath::Sqrt( (jet->Eta()-track1->Eta())*(jet->Eta()-track1->Eta()) +
1119 dphi*dphi );
1120
1121 if( r < fConeRadius ) {
1122 Double_t fPt = jet->E()+track1->Pt(); // Scalar sum of Pt
1123 // recalculating the centroid
1124 Double_t eta = jet->Eta()*jet->E()/fPt + track1->Eta()*track1->Pt()/fPt;
1125
1126 //Double_t phi = jet->Phi()*jet->E()/fPt + track1->Phi()*track1->Pt()/fPt; // OB - remove
1127 // OB - recalc phi via weighted 2vectors
1128 jetVect2.SetMagPhi(jet->E()/fPt,jetVect2.Phi());
1129 trackVect2.SetMagPhi(track1->Pt()/fPt,trackVect2.Phi());
1130
1131 TVector2 sumVect2 = jetVect2+trackVect2;
1132 Double_t phi = sumVect2.Phi();
1133
1134 //jet->SetPtEtaPhiE( 1., eta, phi, fPt );
1135 ((TLorentzVector*) jet->MomentumVector())->SetPtEtaPhiE(1,eta,phi,fPt);
1136
1137 jet->AddTrack(track1);
1138 tracks.Remove(track1);
1139 }
1140 }
1141
1142 tracks.Compress();
1143
1144 nTracks = tracks.GetEntries();
1145 // 4- Continue until all particles are in a jet.
1146 itrack.Reset();
1147 } // end while nTracks
1148
1149 // Convert to AODjets....
1150 Int_t njets = jets->GetEntriesFast();
1151
1152 TClonesArray* aodjets = new TClonesArray("AliAODJet",njets);
1153 aodjets->SetOwner(kTRUE);
1154
1155 Int_t count = 0;
1156
1157 for(Int_t ijet=0; ijet<njets; ++ijet) {
1158 //TLorentzVector* jet = (TLorentzVector*)jets->At(ijet);
1159 AliAODJet* jet = (AliAODJet*)jets->At(ijet);
1160
1161 if (jet->E() < fJetPtCut) continue;
1162 Float_t px, py,pz,en; // convert to 4-vector
1163 px = jet->E() * TMath::Cos(jet->Phi()); // Pt * cos(phi)
1164 py = jet->E() * TMath::Sin(jet->Phi()); // Pt * sin(phi)
1165 pz = jet->E() / TMath::Tan(2.0 * TMath::ATan(TMath::Exp(-jet->Eta())));
1166 en = TMath::Sqrt(px * px + py * py + pz * pz);
1167 jet->SetPxPyPzE(px,py,pz,en);
1168
1169 TClonesArray &tmpaodjets = *aodjets;
1170 new(tmpaodjets[count++]) AliAODJet(*jet);
1171 //aodjets->AddLast( new AliAODJet(*jet));
1172 //aodjets->AddLast( new AliAODJet(px, py, pz, en) );
1173 }
1174 // Order jets according to their pT .
1175 //QSortTracks( *aodjets, 0, aodjets->GetEntriesFast() ); // not for TClonesArray
1176
1177 if (fDebug>3) AliInfo(Form(" %d Charged jets found - after cuts %d \n",njets,count));
1178
1179 jets->Delete(); // OB: cleanup
1180 delete jets;
1181
1182
1183 return aodjets;
1184}
1185
1186//____________________________________________________________________
1187
1188TClonesArray* AliAnalysisTaskJetChem::GetPythiaJets(){
1189
1190 // return Pythia jets in jet acc
1191
1192 // note: present verion of AliMCEventHandler expects "AliAOD.root" (not "AliAODs.root"), otherwise galice.root will not be found in proper dir
1193 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
1194 if (!mcHandler) {
1195 Printf("ERROR: Could not retrieve MC event handler");
1196 return 0;
1197 }
1198
1199 AliMCEvent* mcEvent = mcHandler->MCEvent();
1200 if (!mcEvent) {
1201 Printf("ERROR: Could not retrieve MC event");
1202 return 0;
1203 }
1204
1205
1206 AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
1207 if(!pythiaGenHeader){
1208 Printf("ERROR: Could not retrieve pythiaGenHeader");
1209 return 0;
1210 }
1211
1212 //Get Jets from MC header
1213 Int_t nPythiaGenJets = pythiaGenHeader->NTriggerJets();
1214
1215 TClonesArray* aodjets = new TClonesArray("AliAODJet",nPythiaGenJets);
1216 aodjets->SetOwner(kTRUE);
1217
1218 int count = 0;
1219 for(int ip = 0; ip<nPythiaGenJets; ++ip){
1220 Float_t p[4];
1221 pythiaGenHeader->TriggerJet(ip,p);
1222 TVector3 tempVect(p[0],p[1],p[2]);
1223 if ( TMath::Abs(tempVect.Eta())>fJetEtaCut ) continue;
1224
1225 TClonesArray &tmpaodjets = *aodjets;
1226 new(tmpaodjets[count++]) AliAODJet(p[0],p[1],p[2],p[3]);
1227 //aodjets->AddLast(new AliAODJet(p[0],p[1],p[2],p[3]));
1228 }
1229
1230 return aodjets;
1231}
1232
1233
1234//____________________________________________________________________
1235void AliAnalysisTaskJetChem::QSortTracks(TObjArray &a, Int_t first, Int_t last)
1236{
1237 // Sort array of TObjArray of tracks by Pt using a quicksort algorithm.
1238
1239 static TObject *tmp;
1240 static int i; // "static" to save stack space
1241 int j;
1242
1243 while (last - first > 1) {
1244 i = first;
1245 j = last;
1246 for (;;) {
1247 while (++i < last && ((AliVParticle*)a[i])->Pt() > ((AliVParticle*)a[first])->Pt() )
1248 ;
1249 while (--j > first && ((AliVParticle*)a[j])->Pt() < ((AliVParticle*)a[first])->Pt() )
1250 ;
1251 if (i >= j)
1252 break;
1253
1254 tmp = a[i];
1255 a[i] = a[j];
1256 a[j] = tmp;
1257 }
1258 if (j == first) {
1259 ++first;
1260 continue;
1261 }
1262 tmp = a[first];
1263 a[first] = a[j];
1264 a[j] = tmp;
1265 if (j - first < last - (j + 1)) {
1266 QSortTracks(a, first, j);
1267 first = j + 1; // QSortTracks(j + 1, last);
1268 } else {
1269 QSortTracks(a, j + 1, last);
1270 last = j; // QSortTracks(first, j);
1271 }
1272 }
1273}
1274
1275//------------------------------------------------------------------
1276
1277TH1F* AliAnalysisTaskJetChem::CreatePIDhisto(const char* name){
1278
1279 // create histogram
1280
9280dfa6 1281 TH1F* result = new TH1F(name,"",60,0,60);
9b2de807 1282 result->SetOption("E");
1283
1284 // bin equal Geant ID
1285
1286 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(1),"photon");
1287 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(2),"e+");
1288 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(3),"e-");
1289 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(4),"e-neutrino");
1290 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(5),"mu+");
1291 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(6),"mu-");
1292 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(7),"pi0");
1293 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(8),"pi+");
1294 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(9),"pi-");
1295 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(10),"K long");
1296 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(11),"K+");
1297 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(12),"K-");
1298 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(13),"n");
1299 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(14),"p");
1300 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(15),"anti-proton");
1301 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(16),"K short");
1302 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(17),"eta");
1303 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(18),"Lambda");
1304 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(19),"Sigma+");
1305 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(20),"Sigma0");
1306 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(21),"Sigma-");
1307 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(22),"Xi0");
1308 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(23),"Xi-");
1309 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(24),"Omega-");
1310 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(25),"anti-neutron");
1311 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(26),"anti-Lambda");
1312 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(27),"anti-Sigma-");
1313 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(28),"anti-Sigma0");
1314 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(29),"anti-Sigma+");
1315 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(30),"anti-Xi0");
1316 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(31),"anti-Xi+");
1317 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(32),"anti-Omega+");
1318 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(33),"tau+");
1319 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(34),"tau-");
1320 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(35),"D+");
1321 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(36),"D-");
1322 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(37),"D0");
1323 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(38),"anti-D0");
1324 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(39),"Ds+");
1325 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(40),"anti Ds-");
1326 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(41),"Lamba_c+");
1327 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(42),"W+");
1328 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(43),"W-");
1329 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(44),"Z0");
1330 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(45),"d");
1331 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(46),"t");
1332 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(47),"alpha");
1333 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(48),"G_nu");
9280dfa6 1334 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGpm311Bin),"K0/#bar{K0}");
1335 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDG333Bin),"phi");
1336 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGpm313Bin),"K*(892)0");
1337 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGp323Bin),"K*(892)+");
1338 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGm323Bin),"K*(892)-");
1339 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGNeutrinoBin),"nu");
1340 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGCharmedBaryonBin),"charmed baryon");
1341 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGQuarkBin),"q/#bar{q}");
1342 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGDiQuarkBin),"q #bar{q}");
9b2de807 1343 result->GetXaxis()->LabelsOption("v"); // "u" ?
1344
1345 return result;
1346}
1347
1348//------------------------------------------------------------------
1349
9280dfa6 1350TH1F* AliAnalysisTaskJetChem::CreatePythiaIDhisto(const char* name){
9b2de807 1351
1352 // create histogram
1353
9280dfa6 1354 TH1F* result = new TH1F(name,"",22,0,22);
9b2de807 1355 result->SetOption("E");
1356
1357 result->GetXaxis()->SetBinLabel(kPythiaPIDP11Bin,"qq #rightarrow qq"); // ISUB = 11
1358 result->GetXaxis()->SetBinLabel(kPythiaPIDP12Bin,"q#bar{q} #rightarrow q#bar{q}"); // ISUB = 12
1359 result->GetXaxis()->SetBinLabel(kPythiaPIDP13Bin,"q#bar{q} #rightarrow gg"); // ISUB = 13
1360 result->GetXaxis()->SetBinLabel(kPythiaPIDP28Bin,"qg #rightarrow qg "); // ISUB = 28
1361 result->GetXaxis()->SetBinLabel(kPythiaPIDP53Bin,"gg #rightarrow q#bar{q}"); // ISUB = 53
1362 result->GetXaxis()->SetBinLabel(kPythiaPIDP68Bin,"gg #rightarrow gg"); // ISUB = 68
9280dfa6 1363 result->GetXaxis()->SetBinLabel(kPythiaPIDP92Bin,"SD"); // ISUB = 92
1364 result->GetXaxis()->SetBinLabel(kPythiaPIDP93Bin,"SD"); // ISUB = 93
1365 result->GetXaxis()->SetBinLabel(kPythiaPIDP94Bin,"DD"); // ISUB = 94
1366 result->GetXaxis()->SetBinLabel(kPythiaPIDP95Bin,"low pt (MPI)"); // ISUB = 95
1367 result->GetXaxis()->SetBinLabel(kPythiaPIDPOtherBin,"other"); // ISUB = XX
9b2de807 1368
1369 result->GetXaxis()->LabelsOption("v"); // "u" ?
1370
1371 return result;
1372}
1373
1374
1375//------------------------------------------------------------------
1376
1377void AliAnalysisTaskJetChem::FillPythiaIDhisto(TH1F* h, const Int_t PID){
1378
1379 // fill Pyhtia PID histogram
1380
1381 Int_t bin = -1;
1382
1383 if(PID == 11) bin = kPythiaPIDP11Bin;
1384 else if(PID == 12) bin = kPythiaPIDP12Bin;
1385 else if(PID == 13) bin = kPythiaPIDP13Bin;
1386 else if(PID == 28) bin = kPythiaPIDP28Bin;
1387 else if(PID == 53) bin = kPythiaPIDP53Bin;
1388 else if(PID == 68) bin = kPythiaPIDP68Bin;
9280dfa6 1389 else if(PID == 92) bin = kPythiaPIDP92Bin;
1390 else if(PID == 93) bin = kPythiaPIDP93Bin;
1391 else if(PID == 94) bin = kPythiaPIDP94Bin;
1392 else if(PID == 95) bin = kPythiaPIDP95Bin;
9b2de807 1393 else{
1394 if(PID != -1) AliInfo(Form("unknown PID %d",PID));
1395 return;
1396 }
1397
1398 h->Fill(h->GetBinCenter(bin),1);
1399}
1400
1401
1402//____________________________________________________________________
1403
1404
1405void AliAnalysisTaskJetChem::FillPIDhisto(TH1F* hist, Int_t pdg, Float_t weight){
1406
1407 // convert pdg code to Geant ID and fill corresponding bin
1408
1409 Int_t fGID = fpdgdb->ConvertPdgToGeant3(pdg);
1410
9280dfa6 1411 //cout<<" pdg "<<pdg<<" fGID "<<fGID<<endl;
1412
1413 if(TMath::Abs(pdg) == 311) fGID = kPDGpm311Bin;
1414 if(pdg == 333) fGID = kPDG333Bin;
1415 if(TMath::Abs(pdg) == 313) fGID = kPDGpm313Bin;
1416 if(pdg == 323) fGID = kPDGp323Bin;
1417 if(pdg == -323) fGID = kPDGm323Bin;
1418
1419 if(TMath::Abs(pdg)==12 || TMath::Abs(pdg)==14 || TMath::Abs(pdg)==16) fGID = kPDGNeutrinoBin;
1420
1421 if(TMath::Abs(pdg)==4122) fGID = kPDGCharmedBaryonBin;
1422
1423 if(1<=TMath::Abs(pdg) && TMath::Abs(pdg)<=6) fGID = kPDGQuarkBin;
1424
1425 if(TMath::Abs(pdg)==1103 || TMath::Abs(pdg)==2101 || TMath::Abs(pdg)==2103 || TMath::Abs(pdg)==2203 ||
1426 TMath::Abs(pdg)==3101 || TMath::Abs(pdg)==3103 || TMath::Abs(pdg)==3201 || TMath::Abs(pdg)==3203 ||
1427 TMath::Abs(pdg)==3303 || TMath::Abs(pdg)==4101 || TMath::Abs(pdg)==4103 || TMath::Abs(pdg)==4201 ||
1428 TMath::Abs(pdg)==4203 || TMath::Abs(pdg)==4301 || TMath::Abs(pdg)==4303 || TMath::Abs(pdg)==4403 ||
1429 TMath::Abs(pdg)==5101 || TMath::Abs(pdg)==5103 || TMath::Abs(pdg)==5201 || TMath::Abs(pdg)==5203 ||
1430 TMath::Abs(pdg)==5301 || TMath::Abs(pdg)==5303 || TMath::Abs(pdg)==5401 || TMath::Abs(pdg)==5403 ||
1431 TMath::Abs(pdg)==5503) fGID = kPDGDiQuarkBin;
1432
9b2de807 1433
1434 hist->Fill(fGID,weight);
1435
1436 if(fGID>hist->GetBinLowEdge(hist->GetNbinsX()+1) ||
1437 fGID<hist->GetBinLowEdge(1)){
1438
1439 AliError(Form("fGID %d for pdg %d exceeding histo limits ",fGID,pdg));
1440 }
1441
1442 if(fGID == 0){
9280dfa6 1443 AliInfo(Form("fGID 0 for pdg %d ",pdg));
9b2de807 1444 }
1445
1446}
1447
1448//____________________________________________________________________
1449void AliAnalysisTaskJetChem::CreateHistos(){
1450
1451 // make histos
1452
1453 fListOfHistos = new TList();
1454 fListOfHistos->SetOwner(kTRUE);
1455
9280dfa6 1456 fhPrimVertexNCont = new TH1F("hPrimVertexNCont","",52,-2,50);
1457 fhPrimVertexNCont->SetXTitle("");
1458 fhPrimVertexNCont->Sumw2();
1459 fListOfHistos->Add( fhPrimVertexNCont );
1460
9b2de807 1461 fhPrimVertexRho = new TH1F("hPrimVertexRho","",100,0,1);
1462 fhPrimVertexRho->SetXTitle("");
1463 fhPrimVertexRho->Sumw2();
1464 fListOfHistos->Add( fhPrimVertexRho );
1465
9280dfa6 1466 fhPrimVertexZ = new TH1F("hPrimVertexZ","",40,-20,20);
9b2de807 1467 fhPrimVertexZ->SetXTitle("");
1468 fhPrimVertexZ->Sumw2();
1469 fListOfHistos->Add( fhPrimVertexZ );
1470
1471 fhNJets = new TH1F("hNJets","",10, 0, 10);
1472 fhNJets->SetXTitle("# of jets");
1473 fhNJets->Sumw2();
1474 fListOfHistos->Add( fhNJets );
1475
1476 fhNJetsMC = new TH1F("hNJetsMC","",10,0,10);
1477 fhNJetsMC->SetXTitle("# of jets");
1478 fhNJetsMC->Sumw2();
1479 fListOfHistos->Add( fhNJetsMC );
1480
9280dfa6 1481 fhLeadingEta = new TH1F("hLeadingEta","Leading Jet eta",12,-0.6,0.6);
9b2de807 1482 fhLeadingEta->SetXTitle("eta");
1483 fhLeadingEta->SetYTitle("dN/deta");
1484 fhLeadingEta->Sumw2();
1485 fListOfHistos->Add(fhLeadingEta);
1486
9280dfa6 1487 fhLeadingNTracksVsEta = new TH2F("hLeadingNTracksVsEta","",20,-1.0,1.0,20,0,20);
1488 fhLeadingNTracksVsEta->SetXTitle("eta");
1489 fhLeadingNTracksVsEta->SetYTitle("# of tracks");
1490 fhLeadingNTracksVsEta->Sumw2();
1491 fListOfHistos->Add(fhLeadingNTracksVsEta);
1492
1493 fhLeadingPtVsEta = new TH2F("hLeadingPtVsEta","",20,-1.0,1.0,50,0,50);
1494 fhLeadingPtVsEta->SetXTitle("eta");
1495 fhLeadingPtVsEta->SetYTitle("# of tracks");
1496 fhLeadingPtVsEta->Sumw2();
1497 fListOfHistos->Add(fhLeadingPtVsEta);
1498
1499
1500 fhLeadingPhi = new TH1F("hLeadingPhi","Leading Jet phi",63,0,6.3);
9b2de807 1501 fhLeadingPhi->SetXTitle("phi");
1502 fhLeadingPhi->SetYTitle("dN/dphi");
1503 fhLeadingPhi->Sumw2();
1504 fListOfHistos->Add(fhLeadingPhi);
1505
9280dfa6 1506
9b2de807 1507 fhLeadingPt = new TH1F("hLeadingPt","leading Jet p_{T}",50,0,50);
1508 fhLeadingPt->SetXTitle("p_{T} (GeV/c)");
1509 fhLeadingPt->SetYTitle("dN/dp_{T} (1/GeV)");
1510 fhLeadingPt->Sumw2();
1511 fListOfHistos->Add( fhLeadingPt );
1512
1513 fhLeadingPtDiffr = new TH1F("hLeadingPtDiffr","leading Jet p_{T}",50,0,50);
9280dfa6 1514 fhLeadingPtDiffr->SetXTitle("P_{T} (GeV/c)");
9b2de807 1515 fhLeadingPtDiffr->SetYTitle("dN/dp_{T} (1/GeV)");
1516 fhLeadingPtDiffr->Sumw2();
1517 fListOfHistos->Add( fhLeadingPtDiffr );
1518
9280dfa6 1519 fhLeadingEtaMC = new TH1F("hLeadingEtaMC","Leading Jet eta",12,-0.6,0.6);
9b2de807 1520 fhLeadingEtaMC->SetXTitle("eta");
1521 fhLeadingEtaMC->SetYTitle("dN/deta");
1522 fhLeadingEtaMC->Sumw2();
1523 fListOfHistos->Add(fhLeadingEtaMC);
1524
9280dfa6 1525 fhLeadingPhiMC = new TH1F("hLeadingPhiMC","Leading Jet phi",63,0,6.3);
9b2de807 1526 fhLeadingPhiMC->SetXTitle("phi");
1527 fhLeadingPhiMC->SetYTitle("dN/dphi");
1528 fhLeadingPhiMC->Sumw2();
1529 fListOfHistos->Add(fhLeadingPhiMC);
1530
1531 fhLeadingPtMC = new TH1F("hLeadingPtMC","Leading Jet p_{T}",50,0,50);
1532 fhLeadingPtMC->SetXTitle("p_{T} (GeV/c)");
1533 fhLeadingPtMC->SetYTitle("dN/dp_{T} (1/GeV)");
1534 fhLeadingPtMC->Sumw2();
1535 fListOfHistos->Add( fhLeadingPtMC );
1536
1537 fhLeadingPtMCDiffr = new TH1F("hLeadingPtMCDiffr","Leading Jet p_{T}",50,0,50);
1538 fhLeadingPtMCDiffr->SetXTitle("p_{T} (GeV/c)");
1539 fhLeadingPtMCDiffr->SetYTitle("dN/dp_{T} (1/GeV)");
1540 fhLeadingPtMCDiffr->Sumw2();
1541 fListOfHistos->Add( fhLeadingPtMCDiffr );
1542
9280dfa6 1543 fhPhiEtaTracksNoCut = new TH2F("hPhiEtaTracksNoCut","phi vs eta tracks",20,-1.0,1.0,63,0,6.3);
1544 fhPhiEtaTracksNoCut->SetXTitle("eta");
1545 fhPhiEtaTracksNoCut->SetYTitle("phi");
1546 fhPhiEtaTracksNoCut->Sumw2();
1547 fListOfHistos->Add(fhPhiEtaTracksNoCut);
1548
1549 fhPtTracksNoCut = new TH1F("hPtTracksNoCut","p_{T} tracks",150,0,150);
1550 fhPtTracksNoCut->SetXTitle("p_{T} (GeV)");
1551 fhPtTracksNoCut->SetYTitle("dN/dp_{T} (1/GeV)");
1552 fhPtTracksNoCut->Sumw2();
1553 fListOfHistos->Add(fhPtTracksNoCut);
1554
1555 fhPhiEtaTracks = new TH2F("hPhiEtaTracks","phi vs eta tracks",20,-1.0,1.0,63,0,6.3);
1556 fhPhiEtaTracks->SetXTitle("eta");
1557 fhPhiEtaTracks->SetYTitle("phi");
1558 fhPhiEtaTracks->Sumw2();
1559 fListOfHistos->Add(fhPhiEtaTracks);
1560
1561 fhPtTracks = new TH1F("hPtTracks","p_{T} tracks",150,0,150);
9b2de807 1562 fhPtTracks->SetXTitle("P{T} (GeV)");
1563 fhPtTracks->SetYTitle("dN/dp_{T} (1/GeV)");
1564 fhPtTracks->Sumw2();
1565 fListOfHistos->Add(fhPtTracks);
1566
9280dfa6 1567 fhTrackMult = new TH1F("hTrackMult","",150,0,150);
1568 fhTrackMult->SetXTitle("n Tracks");
1569 fhTrackMult->SetYTitle("counts");
1570 fhTrackMult->Sumw2();
1571 fListOfHistos->Add(fhTrackMult);
1572
9b2de807 1573 fhEtaMCTracks = new TH1F("hEtaMCTracks","eta tracks",30,-1.5,1.5);
1574 fhEtaMCTracks->SetXTitle("eta");
1575 fhEtaMCTracks->SetYTitle("dN/deta");
1576 fhEtaMCTracks->Sumw2();
1577 fListOfHistos->Add(fhEtaMCTracks);
1578
9280dfa6 1579 fhPhiMCTracks = new TH1F("hPhiMCTracks","phi tracks",63,0,6.3);
9b2de807 1580 fhPhiMCTracks->SetXTitle("phi");
1581 fhPhiMCTracks->SetYTitle("dN/dphi");
1582 fhPhiMCTracks->Sumw2();
1583 fListOfHistos->Add(fhPhiMCTracks);
1584
1585 fhPtMCTracks = new TH1F("hPtMCTracks","p_{T} tracks",50,0,50);
1586 fhPtMCTracks->SetXTitle("p_{T} (GeV)");
1587 fhPtMCTracks->SetYTitle("dN/dp_{T} (1/GeV)");
1588 fhPtMCTracks->Sumw2();
1589 fListOfHistos->Add(fhPtMCTracks);
1590
1591 fhnTracksVsPtLeading = new TH2F("hnTracksVsPtLeading","",50,0.,50.,20,-0.5,19.5);
1592 fhnTracksVsPtLeading->SetXTitle("p_{T} (GeV/c)");
1593 fhnTracksVsPtLeading->SetYTitle("n tracks");
1594 fhnTracksVsPtLeading->Sumw2();
1595 fListOfHistos->Add( fhnTracksVsPtLeading );
1596
1597 fhdNdEtaPhiDist = new TH1F("hdNdEtaPhiDist","Charge particle density |#eta|< 1 vs #Delta#phi", 120, 0., 2.*TMath::Pi());
1598 fhdNdEtaPhiDist->SetXTitle("#Delta#phi");
1599 fhdNdEtaPhiDist->SetYTitle("dN{ch}/d#etad#phi");
1600 fhdNdEtaPhiDist->Sumw2();
1601 fListOfHistos->Add( fhdNdEtaPhiDist );
1602
1603 fhRegionSumPtMaxVsEt = new TH1F("hRegionSumPtMaxVsEt","P{T}^{90, max} vs Leading Jet P{T}",50,0.,50.);
1604 fhRegionSumPtMaxVsEt->SetXTitle("P{T} (GeV/c)");
1605 fhRegionSumPtMaxVsEt->Sumw2();
1606 fListOfHistos->Add( fhRegionSumPtMaxVsEt );
1607
1608 fhRegionMultMaxVsEt = new TH1F("hRegionMultMaxVsEt","N{ch}^{90, max} vs Leading Jet P{T}",50,0.,50.);
1609 fhRegionMultMaxVsEt->SetXTitle("E (GeV hRegionAveSumPtVsEt/c)");
1610 fhRegionMultMaxVsEt->Sumw2();
1611 fListOfHistos->Add( fhRegionMultMaxVsEt );
1612
1613 fhRegionSumPtMinVsEt = new TH1F("hRegionSumPtMinVsEt","P{T}^{90, min} vs Leading Jet P{T}",50,0.,50.);
1614 fhRegionSumPtMinVsEt->SetXTitle("P{T} (GeV/c)");
1615 fhRegionSumPtMinVsEt->Sumw2();
1616 fListOfHistos->Add( fhRegionSumPtMinVsEt );
1617
1618 fhRegionMultMinVsEt = new TH1F("hRegionMultMinVsEt","N{ch}^{90, min} vs Leading Jet P{T}",50,0.,50.);
1619 fhRegionMultMinVsEt->SetXTitle("E (GeV/c)");
1620 fhRegionMultMinVsEt->Sumw2();
1621 fListOfHistos->Add( fhRegionMultMinVsEt );
1622
1623 // V0s
1624
1625 fhNV0s = new TH1F("hNV0s","n V0s",50,0,50);
1626 fhNV0s->SetXTitle("n V0s");
1627 fhNV0s->Sumw2();
1628 fListOfHistos->Add(fhNV0s);
1629
1630 fhV0onFly = new TH1F("hV0onFly","on-the-fly V0",5,0,5);
1631 fhV0onFly->SetXTitle("is on-the-fly V0");
1632 fhV0onFly->Sumw2();
1633 fListOfHistos->Add(fhV0onFly);
1634
1635 fhV0DCADaughters = new TH1F("hV0DCADaughters","V0 DCA daughters",200,0,2.0);
1636 fhV0DCADaughters->SetXTitle("V0 DCA daughters");
1637 fhV0DCADaughters->Sumw2();
1638 fListOfHistos->Add(fhV0DCADaughters);
1639
1640 fhV0Radius = new TH1F("hV0Radius","V0 radius",2500,0,250);
1641 fhV0Radius->SetXTitle("V0 radius");
1642 fhV0Radius->Sumw2();
1643 fListOfHistos->Add(fhV0Radius);
1644
1645 fhV0DCAToVertex = new TH1F("hV0DCAToVertex","",100,0,10);
1646 fhV0DCAToVertex->SetXTitle("V0 DCA (cm)");
1647 fhV0DCAToVertex->Sumw2();
1648 fListOfHistos->Add(fhV0DCAToVertex);
1649
1650 fhV0DCAToVertexK0 = new TH1F("hV0DCAToVertexK0","",100,0,10);
1651 fhV0DCAToVertexK0->SetXTitle("V0 DCA (cm)");
1652 fhV0DCAToVertexK0->Sumw2();
1653 fListOfHistos->Add(fhV0DCAToVertexK0);
1654
1655 fhV0InvMassK0 = new TH1F("hV0InvMassK0","",2000,0,2);
1656 fhV0InvMassK0->SetXTitle("inv mass (GeV)");
1657 fhV0InvMassK0->Sumw2();
1658 fListOfHistos->Add(fhV0InvMassK0);
1659
1660 fhV0InvMassK0JetEvt = new TH1F("hV0InvMassK0JetEvt","",2000,0,2);
1661 fhV0InvMassK0JetEvt->SetXTitle("inv mass (GeV)");
1662 fhV0InvMassK0JetEvt->Sumw2();
1663 fListOfHistos->Add(fhV0InvMassK0JetEvt);
1664
1665 fhV0InvMassLambda = new TH1F("hV0InvMassLambda","",2000,0,2);
1666 fhV0InvMassLambda->SetXTitle("inv mass (GeV)");
1667 fhV0InvMassLambda->Sumw2();
1668 fListOfHistos->Add(fhV0InvMassLambda);
1669
1670 fhV0InvMassAntiLambda = new TH1F("hV0InvMassAntiLambda","",2000,0,2);
1671 fhV0InvMassAntiLambda->SetXTitle("inv mass (GeV)");
1672 fhV0InvMassAntiLambda->Sumw2();
1673 fListOfHistos->Add(fhV0InvMassAntiLambda);
1674
1675 fhV0InvMassLambdaJetEvt = new TH1F("hV0InvMassLambdaJetEvt","",2000,0,2);
1676 fhV0InvMassLambdaJetEvt->SetXTitle("inv mass (GeV)");
1677 fhV0InvMassLambdaJetEvt->Sumw2();
1678 fListOfHistos->Add(fhV0InvMassLambdaJetEvt);
1679
1680 fhV0InvMassAntiLambdaJetEvt = new TH1F("hV0InvMassAntiLambdaJetEvt","",2000,0,2);
1681 fhV0InvMassAntiLambdaJetEvt->SetXTitle("inv mass (GeV)");
1682 fhV0InvMassAntiLambdaJetEvt->Sumw2();
1683 fListOfHistos->Add(fhV0InvMassAntiLambdaJetEvt);
1684
1685 fhdROpanK0VsPt = new TH2F("hdROpanK0VsPt","V0 dR vs pt",100,0,10,100,0,10);
1686 fhdROpanK0VsPt->SetXTitle("opening angle R (rad)");
1687 fhdROpanK0VsPt->Sumw2();
1688 fListOfHistos->Add(fhdROpanK0VsPt);
1689
1690 fhdPhiJetV0 = new TH1F("hdPhiJetV0","",640,-3.2,3.2);
1691 fhdPhiJetV0->SetXTitle("#Delta #phi V0 - jet");
1692 fhdPhiJetV0->Sumw2();
1693 fListOfHistos->Add(fhdPhiJetV0);
1694
1695 fhdPhiJetK0 = new TH1F("hdPhiJetK0","",640,-3.2,3.2);
1696 fhdPhiJetK0->SetXTitle("#Delta #phi V0 - jet");
1697 fhdPhiJetK0->Sumw2();
1698 fListOfHistos->Add(fhdPhiJetK0);
1699
1700 fhdRJetK0 = new TH1F("hdRJetK0","dN/dR K0-jet",500,0,5);
1701 fhdRJetK0->SetXTitle("#Delta R K0 - jet");
9280dfa6 1702 fhdRJetK0->SetYTitle("1/N_{jet} dN/dR");
9b2de807 1703 fhdRJetK0->Sumw2();
1704 fListOfHistos->Add(fhdRJetK0);
1705
1706 fhdNdptV0 = new TH1F("hdNdptV0","dN/dpt V0",100,0,10);
1707 fhdNdptV0->SetXTitle("p_{T} (GeV/c)");
9280dfa6 1708 fhdNdptV0->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 1709 fhdNdptV0->Sumw2();
1710 fListOfHistos->Add(fhdNdptV0);
1711
1712 fhdNdptK0 = new TH1F("hdNdptK0","",100,0,10);
1713 fhdNdptK0->SetXTitle("p_{T} (GeV/c)");
9280dfa6 1714 fhdNdptK0->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 1715 fhdNdptK0->Sumw2();
1716 fListOfHistos->Add(fhdNdptK0);
1717
1718 fhPtVsEtaK0 = new TH2F("hPtVsEtaK0","",20,-1,1,100,0,10);
1719 fhPtVsEtaK0->SetXTitle("#eta");
1720 fhPtVsEtaK0->SetYTitle("p_{T} (GeV/c)");
1721 fhPtVsEtaK0->Sumw2();
1722 fListOfHistos->Add(fhPtVsEtaK0);
1723
1724 fhV0InvMassK0DCA = new TH1F("hV0InvMassK0DCA","",2000,0,2);
1725 fhV0InvMassK0DCA->SetXTitle("inv mass (GeV)");
1726 fhV0InvMassK0DCA->Sumw2();
1727 fListOfHistos->Add(fhV0InvMassK0DCA);
1728
1729 fhV0InvMassK0DCAdEdx = new TH1F("hV0InvMassK0DCAdEdx","",2000,0,2);
1730 fhV0InvMassK0DCAdEdx->SetXTitle("inv mass (GeV)");
1731 fhV0InvMassK0DCAdEdx->Sumw2();
1732 fListOfHistos->Add(fhV0InvMassK0DCAdEdx);
1733
1734 fhV0InvMassK0DCAPID = new TH1F("hV0InvMassK0DCAPID","",2000,0,2);
1735 fhV0InvMassK0DCAPID->SetXTitle("inv mass (GeV)");
1736 fhV0InvMassK0DCAPID->Sumw2();
1737 fListOfHistos->Add(fhV0InvMassK0DCAPID);
1738
1739 fhV0InvMassLambdaDCAdEdx = new TH1F("hV0InvMassLambdaDCAdEdx","",2000,0,2);
1740 fhV0InvMassLambdaDCAdEdx->SetXTitle("inv mass (GeV)");
1741 fhV0InvMassLambdaDCAdEdx->Sumw2();
1742 fListOfHistos->Add(fhV0InvMassLambdaDCAdEdx);
1743
1744 fhV0InvMassAntiLambdaDCAdEdx = new TH1F("hV0InvMassAntiLambdaDCAdEdx","",2000,0,2);
1745 fhV0InvMassAntiLambdaDCAdEdx->SetXTitle("inv mass (GeV)");
1746 fhV0InvMassAntiLambdaDCAdEdx->Sumw2();
1747 fListOfHistos->Add(fhV0InvMassAntiLambdaDCAdEdx);
1748
1749 fhdNdptK0DCA = new TH1F("hdNdptK0DCA","",100,0,10);
1750 fhdNdptK0DCA->SetXTitle("p_{T} (GeV/c)");
9280dfa6 1751 fhdNdptK0DCA->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 1752 fhdNdptK0DCA->Sumw2();
1753 fListOfHistos->Add(fhdNdptK0DCA);
1754
1755 fhdNdptK0DCAdEdx = new TH1F("hdNdptK0DCAdEdx","",100,0,10);
1756 fhdNdptK0DCAdEdx->SetXTitle("p_{T} (GeV/c)");
9280dfa6 1757 fhdNdptK0DCAdEdx->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 1758 fhdNdptK0DCAdEdx->Sumw2();
1759 fListOfHistos->Add(fhdNdptK0DCAdEdx);
1760
1761 fhV0InvMassK0Min = new TH1F("hV0InvMassK0Min","",2000,0,2);
1762 fhV0InvMassK0Min->SetXTitle("inv mass (GeV)");
1763 fhV0InvMassK0Min->Sumw2();
1764 fListOfHistos->Add(fhV0InvMassK0Min);
1765
1766 fhV0InvMassLambdaMin = new TH1F("hV0InvMassLambdaMin","",2000,0,2);
1767 fhV0InvMassLambdaMin->SetXTitle("inv mass (GeV)");
1768 fhV0InvMassLambdaMin->Sumw2();
1769 fListOfHistos->Add(fhV0InvMassLambdaMin);
1770
1771 fhV0InvMassAntiLambdaMin = new TH1F("hV0InvMassAntiLambdaMin","",2000,0,2);
1772 fhV0InvMassAntiLambdaMin->SetXTitle("inv mass (GeV)");
1773 fhV0InvMassAntiLambdaMin->Sumw2();
1774 fListOfHistos->Add(fhV0InvMassAntiLambdaMin);
1775
1776 fhV0InvMassK0Max = new TH1F("hV0InvMassK0Max","",2000,0,2);
1777 fhV0InvMassK0Max->SetXTitle("inv mass (GeV)");
1778 fhV0InvMassK0Max->Sumw2();
1779 fListOfHistos->Add(fhV0InvMassK0Max);
1780
1781 fhV0InvMassLambdaMax = new TH1F("hV0InvMassLambdaMax","",2000,0,2);
1782 fhV0InvMassLambdaMax->SetXTitle("inv mass (GeV)");
1783 fhV0InvMassLambdaMax->Sumw2();
1784 fListOfHistos->Add(fhV0InvMassLambdaMax);
1785
1786 fhV0InvMassAntiLambdaMax = new TH1F("hV0InvMassAntiLambdaMax","",2000,0,2);
1787 fhV0InvMassAntiLambdaMax->SetXTitle("inv mass (GeV)");
1788 fhV0InvMassAntiLambdaMax->Sumw2();
1789 fListOfHistos->Add(fhV0InvMassAntiLambdaMax);
1790
1791 fhV0InvMassK0Jet = new TH1F("hV0InvMassK0Jet","",2000,0,2);
1792 fhV0InvMassK0Jet->SetXTitle("inv mass (GeV)");
1793 fhV0InvMassK0Jet->Sumw2();
1794 fListOfHistos->Add(fhV0InvMassK0Jet);
1795
1796 fhV0InvMassLambdaJet = new TH1F("hV0InvMassLambdaJet","",2000,0,2);
1797 fhV0InvMassLambdaJet->SetXTitle("inv mass (GeV)");
1798 fhV0InvMassLambdaJet->Sumw2();
1799 fListOfHistos->Add(fhV0InvMassLambdaJet);
1800
1801 fhV0InvMassAntiLambdaJet = new TH1F("hV0InvMassAntiLambdaJet","",2000,0,2);
1802 fhV0InvMassAntiLambdaJet->SetXTitle("inv mass (GeV)");
1803 fhV0InvMassAntiLambdaJet->Sumw2();
1804 fListOfHistos->Add(fhV0InvMassAntiLambdaJet);
1805
1806 fhV0InvMassK0Lambda = new TH1F("hV0InvMassK0Lambda","",2000,0,2);
1807 fhV0InvMassK0Lambda->SetXTitle("inv mass (GeV)");
1808 fhV0InvMassK0Lambda->Sumw2();
1809 fListOfHistos->Add(fhV0InvMassK0Lambda);
1810
1811 fhdNdptK0JetEvt = new TH1F("hdNdptK0JetEvt","",100,0,10);
1812 fhdNdptK0JetEvt->SetXTitle("p_{T} (GeV/c)");
9280dfa6 1813 fhdNdptK0JetEvt->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 1814 fhdNdptK0JetEvt->Sumw2();
1815 fListOfHistos->Add(fhdNdptK0JetEvt);
1816
1817 fhdNdzK0 = new TH1F("hdNdzK0","",150,0,1.5);
1818 fhdNdzK0->SetXTitle("z");
9280dfa6 1819 fhdNdzK0->SetYTitle("1/N_{jet} dN/dz");
9b2de807 1820 fhdNdzK0->Sumw2();
1821 fListOfHistos->Add(fhdNdzK0);
1822
9280dfa6 1823 fhdNdzK05to10 = new TH1F("hdNdzK05to10","",150,0,1.5);
1824 fhdNdzK05to10->SetXTitle("z");
1825 fhdNdzK05to10->SetYTitle("1/N_{jet} dN/dz");
1826 fhdNdzK05to10->Sumw2();
1827 fListOfHistos->Add(fhdNdzK05to10);
1828
1829 fhdNdzK010to20 = new TH1F("hdNdzK010to20","",150,0,1.5);
1830 fhdNdzK010to20->SetXTitle("z");
1831 fhdNdzK010to20->SetYTitle("1/N_{jet} dN/dz");
1832 fhdNdzK010to20->Sumw2();
1833 fListOfHistos->Add(fhdNdzK010to20);
1834
1835 fhdNdzK020to30 = new TH1F("hdNdzK020to30","",150,0,1.5);
1836 fhdNdzK020to30->SetXTitle("z");
1837 fhdNdzK020to30->SetYTitle("1/N_{jet} dN/dz");
1838 fhdNdzK020to30->Sumw2();
1839 fListOfHistos->Add(fhdNdzK020to30);
1840
1841 fhdNdzK030to40 = new TH1F("hdNdzK030to40","",150,0,1.5);
1842 fhdNdzK030to40->SetXTitle("z");
1843 fhdNdzK030to40->SetYTitle("1/N_{jet} dN/dz");
1844 fhdNdzK030to40->Sumw2();
1845 fListOfHistos->Add(fhdNdzK030to40);
1846
1847 fhdNdzK040to60 = new TH1F("hdNdzK040to60","",150,0,1.5);
1848 fhdNdzK040to60->SetXTitle("z");
1849 fhdNdzK040to60->SetYTitle("1/N_{jet} dN/dz");
1850 fhdNdzK040to60->Sumw2();
1851 fListOfHistos->Add(fhdNdzK040to60);
1852
9b2de807 1853 fhdNdxiK0 = new TH1F("hdNdxiK0","",100,0,10);
1854 fhdNdxiK0->SetXTitle("xi");
9280dfa6 1855 fhdNdxiK0->SetYTitle("1/N_{jet} dN/dxi");
9b2de807 1856 fhdNdxiK0->Sumw2();
1857 fListOfHistos->Add(fhdNdxiK0);
1858
9280dfa6 1859 fhdNdzLambda = new TH1F("hdNdzLambda","",150,0,1.5);
1860 fhdNdzLambda->SetXTitle("z");
1861 fhdNdzLambda->SetYTitle("1/N_{jet} dN/dz");
1862 fhdNdzLambda->Sumw2();
1863 fListOfHistos->Add(fhdNdzLambda);
1864
1865 fhdNdzAntiLambda = new TH1F("hdNdzAntiLambda","",150,0,1.5);
1866 fhdNdzAntiLambda->SetXTitle("z");
1867 fhdNdzAntiLambda->SetYTitle("1/N_{jet} dN/dz");
1868 fhdNdzAntiLambda->Sumw2();
1869 fListOfHistos->Add(fhdNdzAntiLambda);
1870
9b2de807 1871 fhdNdzK0Max = new TH1F("hdNdzK0Max","",150,0,1.5);
1872 fhdNdzK0Max->SetXTitle("z");
9280dfa6 1873 fhdNdzK0Max->SetYTitle("1/N_{jet} dN/dz");
9b2de807 1874 fhdNdzK0Max->Sumw2();
1875 fListOfHistos->Add(fhdNdzK0Max);
1876
1877 fhdNdxiK0Max = new TH1F("hdNdxiK0Max","",100,0,10);
1878 fhdNdxiK0Max->SetXTitle("xi");
9280dfa6 1879 fhdNdxiK0Max->SetYTitle("1/N_{jet} dN/dxi");
9b2de807 1880 fhdNdxiK0Max->Sumw2();
1881 fListOfHistos->Add(fhdNdxiK0Max);
1882
1883 fhdNdzLambdaMax = new TH1F("hdNdzLambdaMax","",150,0,1.5);
1884 fhdNdzLambdaMax->SetXTitle("z");
9280dfa6 1885 fhdNdzLambdaMax->SetYTitle("1/N_{jet} dN/dz");
9b2de807 1886 fhdNdzLambdaMax->Sumw2();
1887 fListOfHistos->Add(fhdNdzLambdaMax);
1888
1889 fhdNdxiLambdaMax = new TH1F("hdNdxiLambdaMax","",700,0,7);
1890 fhdNdxiLambdaMax->SetXTitle("xi");
9280dfa6 1891 fhdNdxiLambdaMax->SetYTitle("1/N_{jet} dN/dxi");
9b2de807 1892 fhdNdxiLambdaMax->Sumw2();
1893 fListOfHistos->Add(fhdNdxiLambdaMax);
1894
1895 fhdNdptK0Max = new TH1F("hdNdptK0Max","",100,0,10);
1896 fhdNdptK0Max->SetXTitle("p_{T} (GeV/c)");
9280dfa6 1897 fhdNdptK0Max->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 1898 fhdNdptK0Max->Sumw2();
1899 fListOfHistos->Add(fhdNdptK0Max);
1900
1901 fhdNdptLambdaMax = new TH1F("hdNdptLambdaMax","",100,0,10);
1902 fhdNdptLambdaMax->SetXTitle("p_{T} (GeV/c)");
9280dfa6 1903 fhdNdptLambdaMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 1904 fhdNdptLambdaMax->Sumw2();
1905 fListOfHistos->Add(fhdNdptLambdaMax);
1906
1907 fhdNdzK0Min = new TH1F("hdNdzK0Min","",150,0,1.5);
1908 fhdNdzK0Min->SetXTitle("z");
9280dfa6 1909 fhdNdzK0Min->SetYTitle("1/N_{jet} dN/dz");
9b2de807 1910 fhdNdzK0Min->Sumw2();
1911 fListOfHistos->Add(fhdNdzK0Min);
1912
1913 fhdNdxiK0Min = new TH1F("hdNdxiK0Min","",100,0,10);
1914 fhdNdxiK0Min->SetXTitle("xi");
9280dfa6 1915 fhdNdxiK0Min->SetYTitle("1/N_{jet} dN/dxi");
9b2de807 1916 fhdNdxiK0Min->Sumw2();
1917 fListOfHistos->Add(fhdNdxiK0Min);
1918
1919 fhdNdzLambdaMin = new TH1F("hdNdzLambdaMin","",150,0,1.5);
1920 fhdNdzLambdaMin->SetXTitle("z");
9280dfa6 1921 fhdNdzLambdaMin->SetYTitle("1/N_{jet} dN/dz");
9b2de807 1922 fhdNdzLambdaMin->Sumw2();
1923 fListOfHistos->Add(fhdNdzLambdaMin);
1924
1925 fhdNdxiLambdaMin = new TH1F("hdNdxiLambdaMin","",700,0,7);
1926 fhdNdxiLambdaMin->SetXTitle("xi");
9280dfa6 1927 fhdNdxiLambdaMin->SetYTitle("1/N_{jet} dN/dxi");
9b2de807 1928 fhdNdxiLambdaMin->Sumw2();
1929 fListOfHistos->Add(fhdNdxiLambdaMin);
1930
1931 fhdNdptK0Min = new TH1F("hdNdptK0Min","",100,0,10);
1932 fhdNdptK0Min->SetXTitle("p_{T} (GeV/c)");
9280dfa6 1933 fhdNdptK0Min->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 1934 fhdNdptK0Min->Sumw2();
1935 fListOfHistos->Add(fhdNdptK0Min);
1936
1937 fhdNdptLambdaMin = new TH1F("hdNdptLambdaMin","",100,0,10);
1938 fhdNdptLambdaMin->SetXTitle("p_{T} (GeV/c)");
9280dfa6 1939 fhdNdptLambdaMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 1940 fhdNdptLambdaMin->Sumw2();
1941 fListOfHistos->Add(fhdNdptLambdaMin);
1942
1943 fhdNdzK0Jet = new TH1F("hdNdzK0Jet","",150,0,1.5);
1944 fhdNdzK0Jet->SetXTitle("z");
9280dfa6 1945 fhdNdzK0Jet->SetYTitle("1/N_{jet} dN/dz");
9b2de807 1946 fhdNdzK0Jet->Sumw2();
1947 fListOfHistos->Add(fhdNdzK0Jet);
1948
1949 fhdNdxiK0Jet = new TH1F("hdNdxiK0Jet","",100,0,10);
1950 fhdNdxiK0Jet->SetXTitle("xi");
9280dfa6 1951 fhdNdxiK0Jet->SetYTitle("1/N_{jet} dN/dxi");
9b2de807 1952 fhdNdxiK0Jet->Sumw2();
1953 fListOfHistos->Add(fhdNdxiK0Jet);
1954
1955 fhdNdzLambdaJet = new TH1F("hdNdzLambdaJet","",150,0,1.5);
1956 fhdNdzLambdaJet->SetXTitle("z");
9280dfa6 1957 fhdNdzLambdaJet->SetYTitle("1/N_{jet} dN/dz");
9b2de807 1958 fhdNdzLambdaJet->Sumw2();
1959 fListOfHistos->Add(fhdNdzLambdaJet);
1960
1961 fhdNdxiLambdaJet = new TH1F("hdNdxiLambdaJet","",700,0,7);
1962 fhdNdxiLambdaJet->SetXTitle("xi");
9280dfa6 1963 fhdNdxiLambdaJet->SetYTitle("1/N_{jet} dN/dxi");
9b2de807 1964 fhdNdxiLambdaJet->Sumw2();
1965 fListOfHistos->Add(fhdNdxiLambdaJet);
1966
1967 fhdNdptK0Jet = new TH1F("hdNdptK0Jet","",100,0,10);
1968 fhdNdptK0Jet->SetXTitle("p_{T} (GeV/c)");
9280dfa6 1969 fhdNdptK0Jet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 1970 fhdNdptK0Jet->Sumw2();
1971 fListOfHistos->Add(fhdNdptK0Jet);
1972
1973 fhdNdptLambdaJet = new TH1F("hdNdptLambdaJet","",100,0,10);
1974 fhdNdptLambdaJet->SetXTitle("p_{T} (GeV/c)");
9280dfa6 1975 fhdNdptLambdaJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 1976 fhdNdptLambdaJet->Sumw2();
1977 fListOfHistos->Add(fhdNdptLambdaJet);
1978
1979 fhdEdxVsMomV0 = new TH2F("hdEdxVsMomV0","",200,-10,10,200,0,2000);
1980 fhdEdxVsMomV0->SetXTitle("mom (GeV/c)");
1981 fhdEdxVsMomV0->SetYTitle("dE/dx");
1982 fhdEdxVsMomV0->Sumw2();
1983 fListOfHistos->Add(fhdEdxVsMomV0);
1984
1985 fhdEdxVsMomV0pidEdx = new TH2F("hdEdxVsMomV0pidEdx","",200,-10,10,200,0,2000);
1986 fhdEdxVsMomV0pidEdx->SetXTitle("mom (GeV/c)");
1987 fhdEdxVsMomV0pidEdx->SetYTitle("dE/dx");
1988 fhdEdxVsMomV0pidEdx->Sumw2();
1989 fListOfHistos->Add(fhdEdxVsMomV0pidEdx);
1990
1991 fhdEdxVsMomV0piPID = new TH2F("hdEdxVsMomV0piPID","",200,-10,10,200,0,2000);
1992 fhdEdxVsMomV0piPID->SetXTitle("mom (GeV/c)");
1993 fhdEdxVsMomV0piPID->SetYTitle("dE/dx");
1994 fhdEdxVsMomV0piPID->Sumw2();
1995 fListOfHistos->Add(fhdEdxVsMomV0piPID);
1996
1997 fhdPhiJetK0MC = new TH1F("hdPhiJetK0MC","",640,-3.2,3.2);
1998 fhdPhiJetK0MC->SetXTitle("#Delta #phi K0 - jet");
9280dfa6 1999 fhdPhiJetK0MC->SetYTitle("1/N_{jet} dN/dphi");
9b2de807 2000 fhdPhiJetK0MC->Sumw2();
2001 fListOfHistos->Add(fhdPhiJetK0MC);
2002
2003 fhdRJetK0MC = new TH1F("hdRJetK0MC","dN/R K0-jet",500,0,5);
2004 fhdRJetK0MC->SetXTitle("#Delta R K0 - jet");
9280dfa6 2005 fhdRJetK0MC->SetYTitle("1/N_{jet} dN/dR");
9b2de807 2006 fhdRJetK0MC->Sumw2();
2007 fListOfHistos->Add(fhdRJetK0MC);
2008
9280dfa6 2009 fhdRV0MC = new TH1F("hdRV0MC","",500,0.,1.);
2010 fhdRV0MC->SetXTitle("#Delta R");
2011 fhdRV0MC->SetYTitle("");
2012 fhdRV0MC->Sumw2();
2013 fListOfHistos->Add(fhdRV0MC);
2014
2015 fhdNdptchPiMCMax = new TH1F("hdNdptchPiMCMax","",100,0,10);
2016 fhdNdptchPiMCMax->SetXTitle("p_{T} (GeV/c)");
2017 fhdNdptchPiMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2018 fhdNdptchPiMCMax->Sumw2();
2019 fListOfHistos->Add(fhdNdptchPiMCMax);
2020
9b2de807 2021 fhdNdptK0MCMax = new TH1F("hdNdptK0MCMax","",100,0,10);
2022 fhdNdptK0MCMax->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2023 fhdNdptK0MCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 2024 fhdNdptK0MCMax->Sumw2();
2025 fListOfHistos->Add(fhdNdptK0MCMax);
2026
9280dfa6 2027 fhdNdptchKMCMax = new TH1F("hdNdptchKMCMax","",100,0,10);
2028 fhdNdptchKMCMax->SetXTitle("p_{T} (GeV/c)");
2029 fhdNdptchKMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2030 fhdNdptchKMCMax->Sumw2();
2031 fListOfHistos->Add(fhdNdptchKMCMax);
2032
2033 fhdNdptpMCMax = new TH1F("hdNdptpMCMax","",100,0,10);
2034 fhdNdptpMCMax->SetXTitle("p_{T} (GeV/c)");
2035 fhdNdptpMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2036 fhdNdptpMCMax->Sumw2();
2037 fListOfHistos->Add(fhdNdptpMCMax);
2038
2039 fhdNdptpBarMCMax = new TH1F("hdNdptpBarMCMax","",100,0,10);
2040 fhdNdptpBarMCMax->SetXTitle("p_{T} (GeV/c)");
2041 fhdNdptpBarMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2042 fhdNdptpBarMCMax->Sumw2();
2043 fListOfHistos->Add(fhdNdptpBarMCMax);
2044
2045 fhdNdptLambdaMCMax = new TH1F("hdNdptLambdaMCMax","",100,0,10);
2046 fhdNdptLambdaMCMax->SetXTitle("p_{T} (GeV/c)");
2047 fhdNdptLambdaMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2048 fhdNdptLambdaMCMax->Sumw2();
2049 fListOfHistos->Add(fhdNdptLambdaMCMax);
2050
2051 fhdNdptLambdaBarMCMax = new TH1F("hdNdptLambdaBarMCMax","",100,0,10);
2052 fhdNdptLambdaBarMCMax->SetXTitle("p_{T} (GeV/c)");
2053 fhdNdptLambdaBarMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2054 fhdNdptLambdaBarMCMax->Sumw2();
2055 fListOfHistos->Add(fhdNdptLambdaBarMCMax);
9b2de807 2056
9b2de807 2057
2058 fhdNdptchPiMCMin = new TH1F("hdNdptchPiMCMin","",100,0,10);
2059 fhdNdptchPiMCMin->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2060 fhdNdptchPiMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 2061 fhdNdptchPiMCMin->Sumw2();
2062 fListOfHistos->Add(fhdNdptchPiMCMin);
9280dfa6 2063
2064 fhdNdptK0MCMin = new TH1F("hdNdptK0MCMin","",100,0,10);
2065 fhdNdptK0MCMin->SetXTitle("p_{T} (GeV/c)");
2066 fhdNdptK0MCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2067 fhdNdptK0MCMin->Sumw2();
2068 fListOfHistos->Add(fhdNdptK0MCMin);
9b2de807 2069
2070 fhdNdptchKMCMin = new TH1F("hdNdptchKMCMin","",100,0,10);
2071 fhdNdptchKMCMin->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2072 fhdNdptchKMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 2073 fhdNdptchKMCMin->Sumw2();
2074 fListOfHistos->Add(fhdNdptchKMCMin);
2075
2076 fhdNdptpMCMin = new TH1F("hdNdptpMCMin","",100,0,10);
2077 fhdNdptpMCMin->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2078 fhdNdptpMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 2079 fhdNdptpMCMin->Sumw2();
2080 fListOfHistos->Add(fhdNdptpMCMin);
2081
2082 fhdNdptpBarMCMin = new TH1F("hdNdptpBarMCMin","",100,0,10);
2083 fhdNdptpBarMCMin->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2084 fhdNdptpBarMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 2085 fhdNdptpBarMCMin->Sumw2();
2086 fListOfHistos->Add(fhdNdptpBarMCMin);
2087
9280dfa6 2088 fhdNdptLambdaMCMin = new TH1F("hdNdptLambdaMCMin","",100,0,10);
2089 fhdNdptLambdaMCMin->SetXTitle("p_{T} (GeV/c)");
2090 fhdNdptLambdaMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2091 fhdNdptLambdaMCMin->Sumw2();
2092 fListOfHistos->Add(fhdNdptLambdaMCMin);
2093
2094 fhdNdptLambdaBarMCMin = new TH1F("hdNdptLambdaBarMCMin","",100,0,10);
2095 fhdNdptLambdaBarMCMin->SetXTitle("p_{T} (GeV/c)");
2096 fhdNdptLambdaBarMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2097 fhdNdptLambdaBarMCMin->Sumw2();
2098 fListOfHistos->Add(fhdNdptLambdaBarMCMin);
2099
2100 fhdNdptOmegaMCMin = new TH1F("hdNdptOmegaMCMin","",100,0,10);;
2101 fhdNdptOmegaMCMin->SetXTitle("p_{T} (GeV/c)");
2102 fhdNdptOmegaMCMin->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2103 fhdNdptOmegaMCMin->Sumw2();
2104 fListOfHistos->Add(fhdNdptOmegaMCMin);
2105
2106 fhdNdptOmegaBarMCMin = new TH1F("hdNdptOmegaBarMCMin","",100,0,10);;
2107 fhdNdptOmegaBarMCMin->SetXTitle("p_{T} (GeV/c)");
2108 fhdNdptOmegaBarMCMin->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2109 fhdNdptOmegaBarMCMin->Sumw2();
2110 fListOfHistos->Add(fhdNdptOmegaBarMCMin);
9b2de807 2111
9b2de807 2112 fhdNdptchPiMCJet = new TH1F("hdNdptchPiMCJet","",100,0,10);
2113 fhdNdptchPiMCJet->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2114 fhdNdptchPiMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 2115 fhdNdptchPiMCJet->Sumw2();
2116 fListOfHistos->Add(fhdNdptchPiMCJet);
2117
9280dfa6 2118 fhdNdptK0MCJet = new TH1F("hdNdptK0MCJet","",100,0,10);
2119 fhdNdptK0MCJet->SetXTitle("p_{T} (GeV/c)");
2120 fhdNdptK0MCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2121 fhdNdptK0MCJet->Sumw2();
2122 fListOfHistos->Add(fhdNdptK0MCJet);
2123
9b2de807 2124 fhdNdptchKMCJet = new TH1F("hdNdptchKMCJet","",100,0,10);
2125 fhdNdptchKMCJet->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2126 fhdNdptchKMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 2127 fhdNdptchKMCJet->Sumw2();
2128 fListOfHistos->Add(fhdNdptchKMCJet);
2129
2130 fhdNdptpMCJet = new TH1F("hdNdptpMCJet","",100,0,10);
2131 fhdNdptpMCJet->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2132 fhdNdptpMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 2133 fhdNdptpMCJet->Sumw2();
2134 fListOfHistos->Add(fhdNdptpMCJet);
2135
2136 fhdNdptpBarMCJet = new TH1F("hdNdptpBarMCJet","",100,0,10);
2137 fhdNdptpBarMCJet->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2138 fhdNdptpBarMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
9b2de807 2139 fhdNdptpBarMCJet->Sumw2();
2140 fListOfHistos->Add(fhdNdptpBarMCJet);
2141
9280dfa6 2142 fhdNdptLambdaMCJet = new TH1F("hdNdptLambdaMCJet","",100,0,10);
2143 fhdNdptLambdaMCJet->SetXTitle("p_{T} (GeV/c)");
2144 fhdNdptLambdaMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2145 fhdNdptLambdaMCJet->Sumw2();
2146 fListOfHistos->Add(fhdNdptLambdaMCJet);
9b2de807 2147
9280dfa6 2148 fhdNdptLambdaBarMCJet = new TH1F("hdNdptLambdaBarMCJet","",100,0,10);
2149 fhdNdptLambdaBarMCJet->SetXTitle("p_{T} (GeV/c)");
2150 fhdNdptLambdaBarMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2151 fhdNdptLambdaBarMCJet->Sumw2();
2152 fListOfHistos->Add(fhdNdptLambdaBarMCJet);
9b2de807 2153
2154 fhPIDMC = CreatePIDhisto("hPIDMC");
2155 fhPIDMC->Sumw2();
2156 fListOfHistos->Add(fhPIDMC);
2157
9280dfa6 2158 fhPIDMC_quarkEv = CreatePIDhisto("hPIDMC_quarkEv");
2159 fhPIDMC_quarkEv->Sumw2();
2160 fListOfHistos->Add(fhPIDMC_quarkEv);
2161
2162 fhPIDMC_gluonEv = CreatePIDhisto("hPIDMC_gluonEv");
2163 fhPIDMC_gluonEv->Sumw2();
2164 fListOfHistos->Add(fhPIDMC_gluonEv);
2165
9b2de807 2166 fhPIDMCAll = CreatePIDhisto("hPIDMCAll");
2167 fhPIDMCAll->Sumw2();
2168 fListOfHistos->Add(fhPIDMCAll);
2169
2170 fhPIDMCMin = CreatePIDhisto("hPIDMCMin");
2171 fhPIDMCMin->Sumw2();
2172 fListOfHistos->Add(fhPIDMCMin);
2173
2174 fhPIDMCJet = CreatePIDhisto("hPIDMCJet");
2175 fhPIDMCJet->Sumw2();
2176 fListOfHistos->Add(fhPIDMCJet);
9280dfa6 2177
2178 fhPIDMCMotherK0 = CreatePIDhisto("hPIDMCMotherK0");
2179 fhPIDMCMotherK0->Sumw2();
2180 fListOfHistos->Add(fhPIDMCMotherK0);
2181
2182 fhPIDMCGrandMotherK0 = CreatePIDhisto("hPIDMCGrandMotherK0");
2183 fhPIDMCGrandMotherK0->Sumw2();
2184 fListOfHistos->Add(fhPIDMCGrandMotherK0);
2185
2186 fhPIDMCMotherChK = CreatePIDhisto("hPIDMCMotherChK");
2187 fhPIDMCMotherChK->Sumw2();
2188 fListOfHistos->Add(fhPIDMCMotherChK);
2189
2190 fhPIDMCMotherK0Trans = CreatePIDhisto("hPIDMCMotherK0Trans");
2191 fhPIDMCMotherK0Trans->Sumw2();
2192 fListOfHistos->Add(fhPIDMCMotherK0Trans);
2193
2194 fhPIDMCGrandMotherK0Trans = CreatePIDhisto("hPIDMCGrandMotherK0Trans");
2195 fhPIDMCGrandMotherK0Trans->Sumw2();
2196 fListOfHistos->Add(fhPIDMCGrandMotherK0Trans);
2197
2198 fhPIDMCMotherChKTrans = CreatePIDhisto("hPIDMCMotherChKTrans");
2199 fhPIDMCMotherChKTrans->Sumw2();
2200 fListOfHistos->Add(fhPIDMCMotherChKTrans);
2201
9b2de807 2202 fhdNdptgammaMC = new TH1F("hdNdptgammaMC","",100,0,10);
2203 fhdNdptgammaMC->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2204 fhdNdptgammaMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
9b2de807 2205 fhdNdptgammaMC->Sumw2();
2206 fListOfHistos->Add(fhdNdptgammaMC);
2207
2208 fhdNdptchPiMC = new TH1F("hdNdptchPiMC","",100,0,10);;
2209 fhdNdptchPiMC->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2210 fhdNdptchPiMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
9b2de807 2211 fhdNdptchPiMC->Sumw2();
2212 fListOfHistos->Add(fhdNdptchPiMC);
2213
2214 fhdNdptpi0MC = new TH1F("hdNdptpi0MC","",100,0,10);;
2215 fhdNdptpi0MC->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2216 fhdNdptpi0MC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
9b2de807 2217 fhdNdptpi0MC->Sumw2();
2218 fListOfHistos->Add(fhdNdptpi0MC);
2219
2220 fhdNdptK0MC = new TH1F("hdNdptK0MC","",100,0,10);;
2221 fhdNdptK0MC->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2222 fhdNdptK0MC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
9b2de807 2223 fhdNdptK0MC->Sumw2();
2224 fListOfHistos->Add(fhdNdptK0MC);
2225
2226 fhdNdptchKMC = new TH1F("hdNdptchKMC","",100,0,10);;
2227 fhdNdptchKMC->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2228 fhdNdptchKMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
9b2de807 2229 fhdNdptchKMC->Sumw2();
2230 fListOfHistos->Add(fhdNdptchKMC);
2231
2232 fhdNdptpMC = new TH1F("hdNdptpMC","",100,0,10);;
2233 fhdNdptpMC->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2234 fhdNdptpMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
9b2de807 2235 fhdNdptpMC->Sumw2();
2236 fListOfHistos->Add(fhdNdptpMC);
2237
2238 fhdNdptpBarMC = new TH1F("hdNdptpBarMC","",100,0,10);;
2239 fhdNdptpBarMC->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2240 fhdNdptpBarMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
9b2de807 2241 fhdNdptpBarMC->Sumw2();
2242 fListOfHistos->Add(fhdNdptpBarMC);
2243
9280dfa6 2244 fhdNdptLambdaMC = new TH1F("hdNdptLambdaMC","",100,0,10);;
2245 fhdNdptLambdaMC->SetXTitle("p_{T} (GeV/c)");
2246 fhdNdptLambdaMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2247 fhdNdptLambdaMC->Sumw2();
2248 fListOfHistos->Add(fhdNdptLambdaMC);
2249
2250 fhdNdptLambdaBarMC = new TH1F("hdNdptLambdaBarMC","",100,0,10);;
2251 fhdNdptLambdaBarMC->SetXTitle("p_{T} (GeV/c)");
2252 fhdNdptLambdaBarMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2253 fhdNdptLambdaBarMC->Sumw2();
2254 fListOfHistos->Add(fhdNdptLambdaBarMC);
9b2de807 2255
2256 fhdNdptOmegaMC = new TH1F("hdNdptOmegaMC","",100,0,10);;
2257 fhdNdptOmegaMC->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2258 fhdNdptOmegaMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
9b2de807 2259 fhdNdptOmegaMC->Sumw2();
2260 fListOfHistos->Add(fhdNdptOmegaMC);
2261
2262 fhdNdptOmegaBarMC = new TH1F("hdNdptOmegaBarMC","",100,0,10);;
2263 fhdNdptOmegaBarMC->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2264 fhdNdptOmegaBarMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
9b2de807 2265 fhdNdptOmegaBarMC->Sumw2();
2266 fListOfHistos->Add(fhdNdptOmegaBarMC);
2267
2268 fhdNdxiMC = new TH1F("hdNdxiMC","",100,0,10);
2269 fhdNdxiMC->SetXTitle("#xi");
2270 fhdNdxiMC->Sumw2();
2271 fListOfHistos->Add(fhdNdxiMC);
2272
2273 fhdNdxiK0MC = new TH1F("hdNdxiK0MC","",100,0,10);
2274 fhdNdxiK0MC->SetXTitle("#xi");
2275 fhdNdxiK0MC->Sumw2();
2276 fListOfHistos->Add(fhdNdxiK0MC);
2277
2278 fhdNdxiK0MCJet = new TH1F("hdNdxiK0MCJet","",100,0,10);
2279 fhdNdxiK0MCJet->SetXTitle("#xi");
2280 fhdNdxiK0MCJet->Sumw2();
2281 fListOfHistos->Add(fhdNdxiK0MCJet);
2282
2283 fhdNdzK0MC = new TH1F("hdNdzK0MC","",150,0,1.5);
2284 fhdNdzK0MC->SetXTitle("z");
9280dfa6 2285 fhdNdzK0MC->SetYTitle("1/N_{jet} dN/dz");
9b2de807 2286 fhdNdzK0MC->Sumw2();
2287 fListOfHistos->Add(fhdNdzK0MC);
2288
2289 fhdNdzK0MCJet = new TH1F("hdNdzK0MCJet","",150,0,1.5);
2290 fhdNdzK0MCJet->SetXTitle("z");
9280dfa6 2291 fhdNdzK0MCJet->SetYTitle("1/N_{jet} dN/dz");
9b2de807 2292 fhdNdzK0MCJet->Sumw2();
2293 fListOfHistos->Add(fhdNdzK0MCJet);
2294
2295 fhdNdptK0MCJetEvt = new TH1F("hdNdptK0MCJetEvt","",100,0,10);;
2296 fhdNdptK0MCJetEvt->SetXTitle("p_{T} (GeV/c)");
9280dfa6 2297 fhdNdptK0MCJetEvt->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
9b2de807 2298 fhdNdptK0MCJetEvt->Sumw2();
2299 fListOfHistos->Add(fhdNdptK0MCJetEvt);
2300
2301 fhnJetsAODvsMC = new TH2F("hnJetsAODvsMC","",20,0,20,20,0.,20.);
2302 fhnJetsAODvsMC->SetXTitle("n jets MC");
2303 fhnJetsAODvsMC->SetXTitle("n jets AOD");
2304 fhnJetsAODvsMC->Sumw2();
2305 fListOfHistos->Add(fhnJetsAODvsMC);
2306
2307 fhLeadingPtAODvsMC = new TH2F("hLeadingPtAODvsMC","",20,0,20,20,0.,20.);
2308 fhLeadingPtAODvsMC->SetXTitle("p_{T} MC (GeV/c)");
2309 fhLeadingPtAODvsMC->SetYTitle("p_{T} AOD (GeV/c)");
2310 fhLeadingPtAODvsMC->Sumw2();
2311 fListOfHistos->Add(fhLeadingPtAODvsMC);
2312
2313 fhLeadingEtaAODvsMC = new TH2F("hLeadingEtaAODvsMC","",20,-1,1,20,-1.,1.);
2314 fhLeadingEtaAODvsMC->SetXTitle("#eta MC");
2315 fhLeadingEtaAODvsMC->SetYTitle("#eta AOD");
2316 fhLeadingEtaAODvsMC->Sumw2();
2317 fListOfHistos->Add(fhLeadingEtaAODvsMC);
2318
2319 fhLeadingPhiAODvsMC = new TH2F("hLeadingPhiAODvsMC","",63,0,6.3,63,0.,6.3);
2320 fhLeadingPhiAODvsMC->SetXTitle("#phi MC");
2321 fhLeadingPhiAODvsMC->SetYTitle("#phi AOD");
2322 fhLeadingPhiAODvsMC->Sumw2();
2323 fListOfHistos->Add(fhLeadingPhiAODvsMC);
2324
2325
2326 fhnTracksLeadingAODvsMC = new TH2F("hnTracksLeadingAODvsMC","",20,0.,20,20,-0.5,19.5);
2327 fhnTracksLeadingAODvsMC->SetXTitle("nTracks MC");
2328 fhnTracksLeadingAODvsMC->SetYTitle("nTracks AOD");
2329 fhnTracksLeadingAODvsMC->Sumw2();
2330 fListOfHistos->Add(fhnTracksLeadingAODvsMC);
2331
2332 fhLeadingdRAODMC = new TH1F("hLeadingdRAODMC","",40,0.,4);
2333 fhLeadingdRAODMC->SetXTitle("#Delta R");
2334 fhLeadingdRAODMC->Sumw2();
2335 fListOfHistos->Add(fhLeadingdRAODMC);
2336
2337
2338 fhLeadingPtAODvsMCdRcut = new TH2F("hLeadingPtAODvsMCdRcut","",40,0,20,40,0.,20.);
2339 fhLeadingPtAODvsMCdRcut->SetXTitle("p_{T} MC (GeV/c)");
2340 fhLeadingPtAODvsMCdRcut->SetYTitle("p_{T} AOD (GeV/c)");
2341 fhLeadingPtAODvsMCdRcut->Sumw2();
2342 fListOfHistos->Add(fhLeadingPtAODvsMCdRcut);
2343
2344
2345 fhdnTracksVsdPtLeadingAODMC = new TH2F("hdnTracksVsdPtLeadingAODMC","",40,-10.,10,20,-10.5,9.5);
2346 fhdnTracksVsdPtLeadingAODMC->SetXTitle("#Delta Pt AOD-MC (GeV/c)");
2347 fhdnTracksVsdPtLeadingAODMC->SetYTitle("#Delta N AOD-MC");
2348 fhdnTracksVsdPtLeadingAODMC->Sumw2();
2349 fListOfHistos->Add(fhdnTracksVsdPtLeadingAODMC);
2350
2351 fhnTracksJetVsPtAOD = new TH2F("hnTracksJetVsPtAOD","",50,0,50,20,-0.5,19.5);
2352 fhnTracksJetVsPtAOD->SetXTitle("p_{T} (GeV/c)");
2353 fhnTracksJetVsPtAOD->SetYTitle("nTracks");
2354 fhnTracksJetVsPtAOD->Sumw2();
2355 fListOfHistos->Add(fhnTracksJetVsPtAOD);
2356
2357 fhnTracksJetVsPtAODquarkEv = new TH2F("hnTracksJetVsPtAODquarkEv","",50,0,50,20,-0.5,19.5);
2358 fhnTracksJetVsPtAODquarkEv->SetXTitle("p_{T} (GeV/c)");
2359 fhnTracksJetVsPtAODquarkEv->SetYTitle("nTracks");
2360 fhnTracksJetVsPtAODquarkEv->Sumw2();
2361 fListOfHistos->Add(fhnTracksJetVsPtAODquarkEv);
2362
2363 fhRadiusJetVsPtAOD = new TH2F("hRadiusJetVsPtAOD","",50,0,50,10,0.,1);
2364 fhRadiusJetVsPtAOD->SetXTitle("p_{T} (GeV/c)");
2365 fhRadiusJetVsPtAOD->SetYTitle("radius");
2366 fhRadiusJetVsPtAOD->Sumw2();
2367 fListOfHistos->Add(fhRadiusJetVsPtAOD);
2368
2369 fhnTracksJetVsPtMC = new TH2F("hnTracksJetVsPtMC","",50,0,50,20,-0.5,19.5);
2370 fhnTracksJetVsPtMC->SetXTitle("p_{T} (GeV/c)");
2371 fhnTracksJetVsPtMC->SetYTitle("nTracks");
2372 fhnTracksJetVsPtMC->Sumw2();
2373 fListOfHistos->Add(fhnTracksJetVsPtMC);
2374
2375 fhnTracksJetVsPtMCquarkEv = new TH2F("hnTracksJetVsPtMCquarkEv","",50,0,50,20,-0.5,19.5);
2376 fhnTracksJetVsPtMCquarkEv->SetXTitle("p_{T} (GeV/c)");
2377 fhnTracksJetVsPtMCquarkEv->SetYTitle("nTracks");
2378 fhnTracksJetVsPtMCquarkEv->Sumw2();
2379 fListOfHistos->Add(fhnTracksJetVsPtMCquarkEv);
2380
2381 fhRadiusJetVsPtMC = new TH2F("hRadiusJetVsPtMC","",50,0,50,10,0.,1);
2382 fhRadiusJetVsPtMC->SetXTitle("p_{T} (GeV/c)");
2383 fhRadiusJetVsPtMC->SetYTitle("radius");
2384 fhRadiusJetVsPtMC->Sumw2();
2385 fListOfHistos->Add(fhRadiusJetVsPtMC);
2386
2387 fhnTracksJetVsPtMCK0 = new TH2F("hnTracksJetVsPtMCK0","",50,0,50,20,-0.5,19.5);
2388 fhnTracksJetVsPtMCK0->SetXTitle("p_{T} (GeV/c)");
2389 fhnTracksJetVsPtMCK0->SetYTitle("nTracks");
2390 fhnTracksJetVsPtMCK0->Sumw2();
2391 fListOfHistos->Add(fhnTracksJetVsPtMCK0);
2392
2393
2394 fhnTracksJetVsPtMCK0quarkEv = new TH2F("hnTracksJetVsPtMCK0quarkEv","",50,0,50,20,-0.5,19.5);
2395 fhnTracksJetVsPtMCK0quarkEv->SetXTitle("p_{T} (GeV/c)");
2396 fhnTracksJetVsPtMCK0quarkEv->SetYTitle("nTracks");
2397 fhnTracksJetVsPtMCK0quarkEv->Sumw2();
2398 fListOfHistos->Add(fhnTracksJetVsPtMCK0quarkEv);
2399
2400 fhRadiusJetVsPtMCK0 = new TH2F("hRadiusJetVsPtMCK0","",50,0,50,10,0.,1);
2401 fhRadiusJetVsPtMCK0->SetXTitle("p_{T} (GeV/c)");
2402 fhRadiusJetVsPtMCK0->SetYTitle("radius");
2403 fhRadiusJetVsPtMCK0->Sumw2();
2404 fListOfHistos->Add(fhRadiusJetVsPtMCK0);
2405
2406 fhnTracksJetVsPtAODK0 = new TH2F("hnTracksJetVsPtAODK0","",50,0,50,20,-0.5,19.5);
2407 fhnTracksJetVsPtAODK0->SetXTitle("p_{T} (GeV/c)");
2408 fhnTracksJetVsPtAODK0->SetYTitle("nTracks AODK0");
2409 fhnTracksJetVsPtAODK0->Sumw2();
2410 fListOfHistos->Add(fhnTracksJetVsPtAODK0);
2411
2412 fhnTracksJetVsPtAODK0quarkEv = new TH2F("hnTracksJetVsPtAODK0quarkEv","",50,0,50,20,-0.5,19.5);
2413 fhnTracksJetVsPtAODK0quarkEv->SetXTitle("p_{T} (GeV/c)");
2414 fhnTracksJetVsPtAODK0quarkEv->SetYTitle("nTracks AODK0quarkEv");
2415 fhnTracksJetVsPtAODK0quarkEv->Sumw2();
2416 fListOfHistos->Add(fhnTracksJetVsPtAODK0quarkEv);
2417
2418 fhRadiusJetVsPtAODK0 = new TH2F("hRadiusJetVsPtAODK0","",50,0,50,10,0.,1);
2419 fhRadiusJetVsPtAODK0->SetXTitle("p_{T} (GeV/c)");
2420 fhRadiusJetVsPtAODK0->SetYTitle("radius");
2421 fhRadiusJetVsPtAODK0->Sumw2();
2422 fListOfHistos->Add(fhRadiusJetVsPtAODK0);
2423
2424 fhnTracksJetVsPtAODpKch = new TH2F("hnTracksJetVsPtAODpKch","",50,0,50,20,-0.5,19.5);
2425 fhnTracksJetVsPtAODpKch->SetXTitle("p_{T} (GeV/c)");
2426 fhnTracksJetVsPtAODpKch->SetYTitle("nTracks AODpKch");
2427 fhnTracksJetVsPtAODpKch->Sumw2();
2428 fListOfHistos->Add(fhnTracksJetVsPtAODpKch);
2429
2430 fhRadiusJetVsPtAODpKch = new TH2F("hRadiusJetVsPtAODpKch","",50,0,50,20,-0.5,19.5);
2431 fhRadiusJetVsPtAODpKch->SetXTitle("p_{T} (GeV/c)");
2432 fhRadiusJetVsPtAODpKch->SetYTitle("Radius AODpKch");
2433 fhRadiusJetVsPtAODpKch->Sumw2();
2434 fListOfHistos->Add(fhRadiusJetVsPtAODpKch);
2435
2436 fhPythiaProcess = CreatePythiaIDhisto("hPythiaProcess");
2437 fListOfHistos->Add(fhPythiaProcess);
2438
2439 fhPythiaProcessK0 = CreatePythiaIDhisto("hPythiaProcessK0");
2440 fListOfHistos->Add(fhPythiaProcessK0);
2441
2442 fhPythiaProcessKch = CreatePythiaIDhisto("hPythiaProcessKch");
2443 fListOfHistos->Add(fhPythiaProcessKch);
2444
2445 fhPythiaProcessp = CreatePythiaIDhisto("hPythiaProcessp");
2446 fListOfHistos->Add(fhPythiaProcessp);
2447
2448 fhPythiaProcesspbar = CreatePythiaIDhisto("hPythiaProcesspbar");
2449 fListOfHistos->Add(fhPythiaProcesspbar);
9280dfa6 2450
2451 fhdNdzJets5to10 = new TH1F("hdNdzJets5to10","",25,0,1.25);
2452 fhdNdzJets5to10->SetXTitle("z");
2453 fhdNdzJets5to10->SetYTitle("dN/dz");
2454 fhdNdzJets5to10->Sumw2();
2455 fListOfHistos->Add(fhdNdzJets5to10);
2456
2457 fhdNdzJets10to20 = new TH1F("hdNdzJets10to20","",25,0,1.25);
2458 fhdNdzJets10to20->SetXTitle("z");
2459 fhdNdzJets10to20->SetYTitle("dN/dz");
2460 fhdNdzJets10to20->Sumw2();
2461 fListOfHistos->Add(fhdNdzJets10to20);
2462
2463 fhdNdzJets20to30 = new TH1F("hdNdzJets20to30","",25,0,1.25);
2464 fhdNdzJets20to30->SetXTitle("z");
2465 fhdNdzJets20to30->SetYTitle("dN/dz");
2466 fhdNdzJets20to30->Sumw2();
2467 fListOfHistos->Add(fhdNdzJets20to30);
2468
2469 fhdNdzJets30to40 = new TH1F("hdNdzJets30to40","",25,0,1.25);
2470 fhdNdzJets30to40->SetXTitle("z");
2471 fhdNdzJets30to40->SetYTitle("dN/dz");
2472 fhdNdzJets30to40->Sumw2();
2473 fListOfHistos->Add(fhdNdzJets30to40);
2474
2475 fhdNdzJets40to60 = new TH1F("hdNdzJets40to60","",25,0,1.25);
2476 fhdNdzJets40to60->SetXTitle("z");
2477 fhdNdzJets40to60->SetYTitle("dN/dz");
2478 fhdNdzJets40to60->Sumw2();
2479 fListOfHistos->Add(fhdNdzJets40to60);
2480
2481
2482 fhdNdxiJets5to10 = new TH1F("hdNdxiJets5to10","",70,0,7);
2483 fhdNdxiJets5to10->SetXTitle("z");
2484 fhdNdxiJets5to10->SetYTitle("dN/dz");
2485 fhdNdxiJets5to10->Sumw2();
2486 fListOfHistos->Add(fhdNdxiJets5to10);
2487
2488 fhdNdxiJets10to20 = new TH1F("hdNdxiJets10to20","",70,0,7);
2489 fhdNdxiJets10to20->SetXTitle("z");
2490 fhdNdxiJets10to20->SetYTitle("dN/dz");
2491 fhdNdxiJets10to20->Sumw2();
2492 fListOfHistos->Add(fhdNdxiJets10to20);
2493
2494 fhdNdxiJets20to30 = new TH1F("hdNdxiJets20to30","",70,0,7);
2495 fhdNdxiJets20to30->SetXTitle("z");
2496 fhdNdxiJets20to30->SetYTitle("dN/dz");
2497 fhdNdxiJets20to30->Sumw2();
2498 fListOfHistos->Add(fhdNdxiJets20to30);
2499
2500 fhdNdxiJets30to40 = new TH1F("hdNdxiJets30to40","",70,0,7);
2501 fhdNdxiJets30to40->SetXTitle("z");
2502 fhdNdxiJets30to40->SetYTitle("dN/dz");
2503 fhdNdxiJets30to40->Sumw2();
2504 fListOfHistos->Add(fhdNdxiJets30to40);
2505
2506 fhdNdxiJets40to60 = new TH1F("hdNdxiJets40to60","",70,0,7);
2507 fhdNdxiJets40to60->SetXTitle("z");
2508 fhdNdxiJets40to60->SetYTitle("dN/dz");
2509 fhdNdxiJets40to60->Sumw2();
2510 fListOfHistos->Add(fhdNdxiJets40to60);
2511
2512 fhdNdptTracksJetPt5to10 = new TH1F("hdNdptTracksJetPt5to10","",250,0,25);
2513 fhdNdptTracksJetPt5to10->SetXTitle("p_{T} (GeV)");
2514 fhdNdptTracksJetPt5to10->SetYTitle("dN/dp_{T} 1/GeV");
2515 fhdNdptTracksJetPt5to10->Sumw2();
2516 fListOfHistos->Add(fhdNdptTracksJetPt5to10);
2517
2518 fhdNdptTracksJetPt10to20 = new TH1F("hdNdptTracksJetPt10to20","",25,0,25);
2519 fhdNdptTracksJetPt10to20->SetXTitle("p_{T} (GeV)");
2520 fhdNdptTracksJetPt10to20->SetYTitle("dN/dp_{T} 1/GeV");
2521 fhdNdptTracksJetPt10to20->Sumw2();
2522 fListOfHistos->Add(fhdNdptTracksJetPt10to20);
2523
2524 fhdNdptTracksJetPt20to30 = new TH1F("hdNdptTracksJetPt20to30","",25,0,25);
2525 fhdNdptTracksJetPt20to30->SetXTitle("p_{T} (GeV)");
2526 fhdNdptTracksJetPt20to30->SetYTitle("dN/dp_{T} 1/GeV");
2527 fhdNdptTracksJetPt20to30->Sumw2();
2528 fListOfHistos->Add(fhdNdptTracksJetPt20to30);
2529
2530 fhdNdptTracksJetPt30to40 = new TH1F("hdNdptTracksJetPt30to40","",25,0,25);
2531 fhdNdptTracksJetPt30to40->SetXTitle("p_{T} (GeV)");
2532 fhdNdptTracksJetPt30to40->SetYTitle("dN/dp_{T} 1/GeV");
2533 fhdNdptTracksJetPt30to40->Sumw2();
2534 fListOfHistos->Add(fhdNdptTracksJetPt30to40);
2535
2536 fhdNdptTracksJetPt40to60 = new TH1F("hdNdptTracksJetPt40to60","",25,0,25);
2537 fhdNdptTracksJetPt40to60->SetXTitle("p_{T} (GeV)");
2538 fhdNdptTracksJetPt40to60->SetYTitle("dN/dp_{T} 1/GeV");
2539 fhdNdptTracksJetPt40to60->Sumw2();
2540 fListOfHistos->Add(fhdNdptTracksJetPt40to60);
2541
2542
9b2de807 2543 fh1Xsec = new TProfile("h1Xsec","xsec from pyxsec.root",1,0,1);
2544 fh1Xsec->SetXTitle("<#sigma>");
2545 fh1Xsec->Sumw2();
2546 fListOfHistos->Add( fh1Xsec );
2547
2548 fh1Trials = new TH1F("h1Trials","trials from pyxsec.root",1,0,1);
2549 fh1Trials->SetXTitle("#sum{ntrials}");
2550 fh1Trials->Sumw2();
2551 fListOfHistos->Add( fh1Trials );
2552
2553// fSettingsTree = new TTree("JetChemAnalysisSettings","Analysis Settings");
2554// fSettingsTree->Branch("fUseLOConeJets",&fUseLOConeJets,"UseLOConeJets/O");
2555// fSettingsTree->Branch("fUseLOConeMCJets",&fUseLOConeMCJets,"UseLOConeMCJets/O");
2556// fSettingsTree->Branch("fUsePythiaJets",&fUsePythiaJets,"UsePythiaJets/O");
2557// fSettingsTree->Branch("fConeRadius", &fConeRadius,"Rad/D");
2558// fSettingsTree->Branch("fTrackPtCutJF", &fTrackPtCutJF,"TrackPt/D");
2559// fSettingsTree->Branch("fFilterBitJF",&fFilterBitJF,"FilterBitJF/i");
2560// fSettingsTree->Branch("fRequireITSRefitJF",&fRequireITSRefitJF,"RequireITSRefitJF/O");
2561// fSettingsTree->Branch("fRejectK0TracksJF",&fRejectK0TracksJF,"RejectK0TracksJF/O");
2562// fSettingsTree->Branch("fJetPtCut",&fJetPtCut,"JetPtCut/D");
2563// fSettingsTree->Branch("fJetEtaCut",&fJetEtaCut,"JetEtaCut/D");
2564// fSettingsTree->Branch("fFilterBit",&fFilterBit,"FilterBit/i");
2565// fSettingsTree->Branch("fTrackPtCut",&fTrackPtCut,"TrackPtCut/D");
2566// fSettingsTree->Branch("fTrackEtaCut",&fTrackEtaCut,"TrackEtaCut/D");
2567// fSettingsTree->Branch("fUseOnFlyV0s",&fUseOnFlyV0s,"UseOnFlyV0s/O");
2568// fSettingsTree->Branch("fCutnSigdEdx",&fCutnSigdEdx,"CutnSigdEdx/D");
2569// fSettingsTree->Branch("fUseAODMCTracksForUE",&fUseAODMCTracksForUE,"UseAODMCTracksForUE/O");
2570// fSettingsTree->Branch("fAreaReg",&fAreaReg,"AreaReg/D");
2571// fSettingsTree->Branch("fAvgTrials",&fAvgTrials,"AvgTrials/D");
2572 // fListOfHistos->Add(fSettingsTree);
2573
2574}
2575
2576//____________________________________________________________________
2577void AliAnalysisTaskJetChem::Terminate(Option_t */*option*/){
2578
2579 // Terminate analysis
2580 //
2581
2582 WriteSettings();
2583
2584 // Update pointers reading them from the output slot
2585 fListOfHistos = dynamic_cast<TList*> (GetOutputData(1));
2586 if( !fListOfHistos ) {
2587 AliError("Histogram List is not available");
2588 return;
2589 }
2590
2591 fhLeadingPt = (TH1F*)fListOfHistos->FindObject("hLeadingPt");
2592 fhRegionSumPtMaxVsEt = (TH1F*)fListOfHistos->FindObject("hRegionSumPtMaxVsEt");
2593 fhRegionSumPtMinVsEt = (TH1F*)fListOfHistos->FindObject("hRegionSumPtMinVsEt");
2594 fhRegionMultMaxVsEt = (TH1F*)fListOfHistos->FindObject("hRegionMultMaxVsEt");
2595 fhRegionMultMinVsEt = (TH1F*)fListOfHistos->FindObject("hRegionMultMinVsEt");
2596
2597
2598 fhRegionSumPtMaxVsEt->Divide(fhLeadingPt);
2599 fhRegionSumPtMinVsEt->Divide(fhLeadingPt);
2600 fhRegionMultMaxVsEt->Divide(fhLeadingPt);
2601 fhRegionMultMinVsEt->Divide(fhLeadingPt);
2602
2603
2604 //Get Normalization
2605 fh1Xsec = (TProfile*) fListOfHistos->FindObject("fh1Xsec");
2606 fh1Trials = (TH1F*)fListOfHistos->FindObject("fh1Trials");
2607
2608 //std::cout<<" fh1Xsec "<<fh1Xsec<<" fh1Trials "<<fh1Trials<<std::endl;
2609
2610 if(fh1Xsec && fh1Trials){
2611
2612 Double_t xsec = fh1Xsec->GetBinContent(1);
2613 Double_t ntrials = fh1Trials->GetBinContent(1);
2614 Double_t normFactor = xsec/ntrials;
2615 Printf("xSec %f nTrials %f Norm %f \n",xsec,ntrials,normFactor);
2616 fhLeadingPt->Scale(normFactor);
2617 }
2618
2619
2620 if(fDebug > 1) AliInfo("End analysis");
2621}
2622
2623// ----------------------------------------------------------------------------
2624
2625void AliAnalysisTaskJetChem::WriteSettings(){
2626
2627 // fill settings tree
2628 // not on GRID !
2629
2630 //fSettingsTree->Fill();
2631
2632}
2633
2634// ---------------------------------------------------------------------------
2635
2636Bool_t AliAnalysisTaskJetChem::IsK0InvMass(const Double_t mass) const {
2637
2638 // K0 mass ? Use STAR params for bin counting & mass fit
2639
2640 const Double_t massK0 = 0.497; // from fits
2641 const Double_t sigmaK0 = 0.0046;
2642 const Double_t nSigmaSignal = 3.5; // STAR parameters for bin counting
2643
2644 if((massK0-nSigmaSignal*sigmaK0)<=mass &&
2645 mass<=(massK0 + nSigmaSignal*sigmaK0)) return kTRUE;
2646
2647 return kFALSE;
2648}
2649
2650// ---------------------------------------------------------------------------
2651
2652Bool_t AliAnalysisTaskJetChem::IsLambdaInvMass(const Double_t mass) const{
2653
2654 // Lambda mass ?
2655
2656
2657 if(1.1<mass && mass<1.13) return kTRUE; // FIXME - adjust range from fit
2658 return kFALSE;
2659}
2660
2661// ---------------------------------------------------------------------------
2662
9f77f0ce 2663Bool_t AliAnalysisTaskJetChem::IsAcceptedDCAK0(/*const Double_t dca*/) const{
9b2de807 2664
2665 // DCA cut
2666
2667 return kTRUE; // FIXME - adjust cut
2668}
2669
2670
2671// ---------------------------------------------------------------------------
2672
9f77f0ce 2673Bool_t AliAnalysisTaskJetChem::IsAcceptedDCALambda(/*const Double_t dca*/) const {
9b2de807 2674
2675 // DCA cut
2676
2677 return kTRUE; // FIXME - adjust cut
2678}
2679
2680// ---------------------------------------------------------------------------
2681
2682Bool_t AliAnalysisTaskJetChem::IsAccepteddEdx(const Double_t mom,const Double_t signal, AliPID::EParticleType n, const Double_t cutnSig) const{
2683
2684 // apply TPC dE/dx cut similar as in AliTPCpidESD
2685 // note: AliTPCidESD uses inner track param for momentum - not avaiable on AOD,
2686 // so we use global track momentum
2687 // should use separate parametrisation for MC and data, but probably ALEPH param & 7% resolution used here anyway not the last word
2688
2689
2690 const Double_t kBBMIP(50.);
2691 const Double_t kBBRes(0.07);
2692 //const Double_t kBBRange(5.);
2693 const Double_t kBBp1(0.76176e-1);
2694 const Double_t kBBp2(10.632);
2695 const Double_t kBBp3(0.13279e-4);
2696 const Double_t kBBp4(1.8631);
2697 const Double_t kBBp5(1.9479);
2698
2699 Double_t mass=AliPID::ParticleMass(n);
2700 Double_t betaGamma = mom/mass;
2701
2702 const Float_t kmeanCorrection =0.1;
2703 Double_t bb = AliExternalTrackParam::BetheBlochAleph(betaGamma,kBBp1,kBBp2,kBBp3,kBBp4,kBBp5);
2704 Double_t meanCorrection =(1+(bb-1)*kmeanCorrection);
2705 Double_t bethe = bb * meanCorrection; // expected
2706 Double_t sigma = bethe * kBBRes;
2707
2708
2709 Double_t dedx = signal/kBBMIP; // measured
2710
2711 Double_t nSig = (TMath::Abs(dedx - bethe))/sigma;
2712
2713 if(nSig > cutnSig) return kFALSE;
2714
2715 return kTRUE;
2716}
2717
2718// ----------------------------------------------------------------------------
2719
2720void AliAnalysisTaskJetChem::CheckV0s(AliAODJet* jetVect, Int_t maxPtRegionIndex,Bool_t& foundK0){
9280dfa6 2721
9b2de807 2722 // loop over AOD V0s, fill masses etc
2723
2724 Int_t nV0 = fAOD->GetNumberOfV0s();
2725 fhNV0s->Fill(nV0);
2726
2727 for(int i=0; i<fAOD->GetNumberOfV0s(); i++){ // loop over V0s
2728
2729 AliAODv0* v0 = fAOD->GetV0(i);
2730
2731 Bool_t isOnFly = v0->GetOnFlyStatus();
2732 isOnFly ? fhV0onFly->Fill(1) : fhV0onFly->Fill(0);
2733
2734 if((fUseOnFlyV0s && isOnFly) || (!fUseOnFlyV0s && !isOnFly) ){ // 'offline' V0s : using vertex tracks, on-fly: take decay vertex into account during mom fit
2735
2736 Double_t massK0 = v0->MassK0Short();
2737 Double_t massLambda = v0->MassLambda();
2738 Double_t massAntiLambda = v0->MassAntiLambda();
2739 Double_t radiusV0 = v0->RadiusV0();
9280dfa6 2740
9b2de807 2741 Double_t etaV0 = v0->Eta();
2742 Double_t fDCAV0ToVertex = v0->DcaV0ToPrimVertex();
2743 Double_t fDCADaughters = v0->DcaV0Daughters();
2744
2745
2746 Double_t v0Mom[3];
2747 v0->PxPyPz(v0Mom);
2748 TVector3 v0MomVect(v0Mom);
2749
2750 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0));
2751 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1));
2752
2753 Double_t posMom[3], negMom[3];
2754 trackPos->PxPyPz(posMom);
2755 trackNeg->PxPyPz(negMom);
2756
2757 TVector3 posMomVec(posMom);
2758 TVector3 negMomVec(negMom);
2759
2760 Double_t dROpanV0 = posMomVec.DeltaR(negMomVec);
2761
2762 AliAODPid* aodPidPos = trackPos->GetDetPid();
2763 AliAODPid* aodPidNeg = trackNeg->GetDetPid();
2764
2765 Double_t dEdxPos = aodPidPos->GetTPCsignal();
2766 Double_t dEdxNeg = aodPidNeg->GetTPCsignal();
2767
2768 Double_t momPos = trackPos->P();
2769 Double_t momNeg = trackNeg->P();
2770
2771 fhdEdxVsMomV0->Fill(momPos,dEdxPos);
2772 fhdEdxVsMomV0->Fill(momNeg,dEdxNeg);
2773
2774 Double_t pV0 = TMath::Sqrt(v0->Ptot2V0());
2775 Double_t ptV0 = TMath::Sqrt(v0->Pt2V0());
2776
2777 fhV0InvMassK0->Fill(massK0);
2778 fhV0InvMassLambda->Fill(massLambda);
2779 fhV0InvMassAntiLambda->Fill(massAntiLambda);
2780
2781 fhV0DCADaughters->Fill(fDCADaughters);
2782 fhV0Radius->Fill(radiusV0);
2783 fhV0DCAToVertex->Fill(fDCAV0ToVertex);
2784 fhdNdptV0->Fill(ptV0);
2785
2786
2787 // K0 signal before cuts
2788
2789 if(IsK0InvMass(massK0)){
2790
2791 fhdNdptK0->Fill(ptV0);
2792
2793 fhPtVsEtaK0->Fill(etaV0,ptV0);
2794
2795 Double_t dRV0MC = AssociateV0MC(&v0MomVect,310); // K0
2796 if(dRV0MC < 0) dRV0MC = 0.99;
2797 fhdRV0MC->Fill(dRV0MC);
2798 }
2799
9f77f0ce 2800 if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/)){
9b2de807 2801
2802 fhV0InvMassK0DCA->Fill(massK0);
2803 if(IsK0InvMass(massK0)) fhdNdptK0DCA->Fill(ptV0);
2804
2805
2806 if(IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)) fhdEdxVsMomV0pidEdx->Fill(momPos,dEdxPos);
2807 if(IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)) fhdEdxVsMomV0pidEdx->Fill(-1*momNeg,dEdxNeg);
2808
2809 if(IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) &&
2810 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
2811
2812 fhV0InvMassK0DCAdEdx->Fill(massK0);
2813
2814 if(IsK0InvMass(massK0)){
2815 fhdNdptK0DCAdEdx->Fill(ptV0);
2816 fhV0DCAToVertexK0->Fill(fDCAV0ToVertex);
2817 fhdROpanK0VsPt->Fill(ptV0,dROpanV0);
2818 foundK0 = kTRUE;
2819 }
2820 }
2821
2822
2823 // AOD pid - cuts strongly into signal
2824
2825 AliAODTrack::AODTrkPID_t mpPIDNeg = trackNeg->GetMostProbablePID();
2826 AliAODTrack::AODTrkPID_t mpPIDPos = trackPos->GetMostProbablePID();
2827
2828 if(mpPIDNeg == AliAODTrack::kPion) fhdEdxVsMomV0piPID->Fill(momPos,dEdxPos);
2829 if(mpPIDPos == AliAODTrack::kPion) fhdEdxVsMomV0piPID->Fill(-1*momNeg,dEdxNeg);
2830
2831
2832 if( (mpPIDNeg == AliAODTrack::kPion) && (mpPIDPos == AliAODTrack::kPion) ){
2833
2834 fhV0InvMassK0DCAPID->Fill(massK0);
2835 }
2836 }
2837
2838
2839 if(IsLambdaInvMass(massLambda) || (IsLambdaInvMass(massAntiLambda))){
2840
2841 fhV0InvMassK0Lambda->Fill(massK0);
2842 }
2843
2844
2845 // Lambda
2846
9f77f0ce 2847 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
9b2de807 2848 IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
2849 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
2850
2851 fhV0InvMassLambdaDCAdEdx->Fill(massLambda);
2852 }
2853
9f77f0ce 2854 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
9b2de807 2855 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
2856 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
2857
2858 fhV0InvMassAntiLambdaDCAdEdx->Fill(massAntiLambda);
2859 }
2860
2861
2862 // jet events
2863
2864 if(jetVect){
2865
2866 Int_t regionV0vect = IsTrackInsideRegion(jetVect,&v0MomVect);
2867
2868 // calc xi
2869 Double_t jetpt = jetVect->Pt();
2870 Double_t dPhiJetV0 = (jetVect->MomentumVector()->Vect()).DeltaPhi(v0MomVect);
2871 Double_t dRJetV0 = (jetVect->MomentumVector()->Vect()).DeltaR(v0MomVect);
2872
2873 Double_t z = pV0/jetpt;
2874 Double_t xi = TMath::Log(1/z);
2875
2876 fhdPhiJetV0->Fill(dPhiJetV0);
2877
2878 // K0
2879
9f77f0ce 2880 if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/) &&
9b2de807 2881 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) &&
2882 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
2883
2884 fhV0InvMassK0JetEvt->Fill(massK0);
2885
2886 if(IsK0InvMass(massK0)){
2887 fhdNdptK0JetEvt->Fill(ptV0);
9280dfa6 2888 fhdNdzK0->Fill(z);
9b2de807 2889 fhdNdxiK0->Fill(xi);
2890
2891 fhdPhiJetK0->Fill(dPhiJetV0);
2892 fhdRJetK0->Fill(dRJetV0);
9280dfa6 2893
2894 if(5<jetpt && jetpt<=10) fhdNdzK05to10->Fill(z);
2895 else if(10<jetpt && jetpt<=20) fhdNdzK010to20->Fill(z);
2896 else if(20<jetpt && jetpt<=30) fhdNdzK020to30->Fill(z);
2897 else if(30<jetpt && jetpt<=40) fhdNdzK030to40->Fill(z);
2898 else if(40<jetpt && jetpt<=60) fhdNdzK040to60->Fill(z);
9b2de807 2899 }
2900 }
9280dfa6 2901
9b2de807 2902
2903 // Lambda
2904
9f77f0ce 2905 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
9b2de807 2906 IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
2907 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
2908
2909 fhV0InvMassLambdaJetEvt->Fill(massLambda);
9280dfa6 2910
2911 if(IsLambdaInvMass(massLambda)){
2912 fhdNdzLambda->Fill(z);
2913 }
9b2de807 2914 }
2915
9f77f0ce 2916 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
9b2de807 2917 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
2918 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
2919
2920 fhV0InvMassAntiLambdaJetEvt->Fill(massAntiLambda);
9280dfa6 2921
2922 if(IsLambdaInvMass(massAntiLambda)){
2923 fhdNdzAntiLambda->Fill(z);
2924 }
9b2de807 2925 }
2926
2927 // fill histos max region
2928
2929 if(regionV0vect != 0 && regionV0vect == maxPtRegionIndex){ // max region
2930
2931 // K0
2932
9f77f0ce 2933 if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/) &&
9b2de807 2934 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) &&
2935 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){ // K0 cuts
2936
2937 fhV0InvMassK0Max->Fill(massK0);
2938
2939 if(IsK0InvMass(massK0)){
2940 fhdNdzK0Max->Fill(z);
2941 fhdNdxiK0Max->Fill(xi);
2942 fhdNdptK0Max->Fill(ptV0);
2943 }
2944 }
2945
2946 // Lambda
2947
9f77f0ce 2948 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
9b2de807 2949 IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
2950 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
2951
2952 fhV0InvMassLambdaMax->Fill(massLambda);
2953 }
2954
9f77f0ce 2955 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
9b2de807 2956 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
2957 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
2958
2959 fhV0InvMassAntiLambdaMax->Fill(massAntiLambda);
2960 }
2961
9f77f0ce 2962 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
24e00b4b 2963 ((IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
9280dfa6 2964 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)) ||
2965 (IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
2966 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)))){
2967
9b2de807 2968 if(IsLambdaInvMass(massLambda) || (IsLambdaInvMass(massAntiLambda))){
2969 fhdNdzLambdaMax->Fill(z);
2970 fhdNdxiLambdaMax->Fill(xi);
2971 fhdNdptLambdaMax->Fill(ptV0);
2972 }
2973 }
2974 }
2975
2976 // fill histos min region
2977
2978 if(regionV0vect != 0 && regionV0vect != maxPtRegionIndex){ // min region
2979
2980 // K0
2981
9f77f0ce 2982 if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/) &&
9b2de807 2983 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) &&
2984 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
2985
2986 fhV0InvMassK0Min->Fill(massK0);
2987
2988 if(IsK0InvMass(massK0)){
2989 fhdNdzK0Min->Fill(z);
2990 fhdNdxiK0Min->Fill(xi);
2991 fhdNdptK0Min->Fill(ptV0);
2992 }
2993 }
2994
2995
2996 // Lambda
2997
2998
9f77f0ce 2999 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
9b2de807 3000 IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
3001 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
3002
3003 fhV0InvMassLambdaMin->Fill(massLambda);
3004 }
3005
9f77f0ce 3006 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
9b2de807 3007 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
3008 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
3009
3010 fhV0InvMassAntiLambdaMin->Fill(massAntiLambda);
3011 }
3012
9f77f0ce 3013 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
24e00b4b 3014 ((IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
9280dfa6 3015 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)) ||
3016 (IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
3017 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)))){
9b2de807 3018
3019 if(IsLambdaInvMass(massLambda) || (IsLambdaInvMass(massAntiLambda))){
3020 fhdNdzLambdaMin->Fill(z);
3021 fhdNdxiLambdaMin->Fill(xi);
3022 fhdNdptLambdaMin->Fill(ptV0);
3023 }
3024 }
3025 }
3026
3027
3028 // jet region
3029
3030 if(regionV0vect == 0){ // jet region
3031
3032 //Double_t dRJetV0 = jetVect->DeltaR(v0MomVect);
3033
3034 if(dRJetV0 <= fConeRadius){
3035
3036 // fill histos jet region
3037
3038 // K0
3039
9f77f0ce 3040 if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/) &&
9b2de807 3041 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) &&
3042 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){ // K0 cuts
3043
3044 fhV0InvMassK0Jet->Fill(massK0);
3045
3046 if(IsK0InvMass(massK0)){
3047 fhdNdzK0Jet->Fill(z);
3048 fhdNdxiK0Jet->Fill(xi);
3049 fhdNdptK0Jet->Fill(ptV0);
3050 }
3051 }
3052
3053
3054 // Lambda
3055
9f77f0ce 3056 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
9b2de807 3057 IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
3058 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
3059
3060 fhV0InvMassLambdaJet->Fill(massLambda);
3061 }
9f77f0ce 3062 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
9b2de807 3063 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
3064 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
3065
3066 fhV0InvMassAntiLambdaJet->Fill(massAntiLambda);
3067 }
3068
9f77f0ce 3069 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
24e00b4b 3070 ((IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
9280dfa6 3071 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)) ||
3072 (IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
3073 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)))){
9b2de807 3074
3075 if(IsLambdaInvMass(massLambda) || (IsLambdaInvMass(massAntiLambda))){
3076 fhdNdzLambdaJet->Fill(z);
3077 fhdNdxiLambdaJet->Fill(xi);
3078 fhdNdptLambdaJet->Fill(ptV0);
3079 }
3080 }
3081 }
3082 }
3083 }
3084 }
3085 }
3086}
3087
3088// ----------------------------------------------------------------------------
3089
3090void AliAnalysisTaskJetChem::CheckMCParticles(AliAODJet* jetVect,Int_t maxPtRegionIndex,Bool_t& isK0Event){
3091
3092 // histos for generated MC
3093
3094 TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
3095
3096 if(!farray){
3097 AliInfo("no mcparticles branch");
3098 return;
3099 }
3100
3101 Int_t ntrks = farray->GetEntries();
3102 if (fDebug>1) AliInfo(Form("check MC particles, tracks %d \n",ntrks));
3103
3104 Int_t pythiaPID = GetPythiaProcessID();
3105
3106 isK0Event = kFALSE;
3107 Bool_t ispEvent = kFALSE;
3108 Bool_t ispBarEvent = kFALSE;
3109 Bool_t isKchEvent = kFALSE;
3110
9280dfa6 3111 Bool_t isQuarkHardScatteringEvent = IsQuarkHardScatteringEvent(pythiaPID);
3112 Bool_t isGluonHardScatteringEvent = IsGluonHardScatteringEvent(pythiaPID);
3113
9b2de807 3114 for(Int_t i =0 ; i<ntrks; i++){ // mc tracks loop
3115
3116 AliAODMCParticle* mctrk = (AliAODMCParticle*)farray->At(i);
3117
9280dfa6 3118 Double_t trackPt = mctrk->Pt();
3119
9b2de807 3120 //Cuts
3121 if (!(mctrk->IsPhysicalPrimary())) continue;
3122
9280dfa6 3123 if ((trackPt < fTrackPtCut) || (TMath::Abs(mctrk->Eta()) > fTrackEtaCut)) continue;
9b2de807 3124
3125 // OB PID histo
3126 Int_t pdg = mctrk->GetPdgCode();
3127
9280dfa6 3128
9b2de807 3129 FillPIDhisto(fhPIDMC,pdg);
9280dfa6 3130 if(isQuarkHardScatteringEvent) FillPIDhisto(fhPIDMC_quarkEv,pdg);
3131 if(isGluonHardScatteringEvent) FillPIDhisto(fhPIDMC_gluonEv,pdg);
3132
3133 if(pdg == 111) fhdNdptpi0MC->Fill(trackPt);
3134 if(pdg == 22) fhdNdptgammaMC->Fill(trackPt);
3135 if(TMath::Abs(pdg) == 211) fhdNdptchPiMC->Fill(trackPt);
3136 if(pdg == 310) fhdNdptK0MC->Fill(trackPt);
3137 if(TMath::Abs(pdg) == 321) fhdNdptchKMC->Fill(trackPt);
3138 if(pdg == 2212) fhdNdptpMC->Fill(trackPt);
3139 if(pdg == -2212) fhdNdptpBarMC->Fill(trackPt);
3140 if(pdg == 3122) fhdNdptLambdaMC->Fill(trackPt);
3141 if(pdg == -3122) fhdNdptLambdaBarMC->Fill(trackPt);
3142 if(pdg == 3332) fhdNdptOmegaMC->Fill(trackPt);
3143 if(pdg == -3332) fhdNdptOmegaBarMC->Fill(trackPt);
9b2de807 3144
9280dfa6 3145 if(pdg == 310) isK0Event = kTRUE;
3146 if(TMath::Abs(pdg) == 321) isKchEvent = kTRUE;
3147 if(pdg == 2212) ispEvent = kTRUE;
3148 if(pdg == -2212) ispBarEvent = kTRUE;
9b2de807 3149
9b2de807 3150
9280dfa6 3151 Int_t pdgMotherChK = -1;
3152 Int_t pdgMotherK0 = -1;
3153 Int_t pdgGrandMotherK0 = -1;
9b2de807 3154
9280dfa6 3155 if(TMath::Abs(pdg) == 321){ // chK
3156 Int_t labelMother = mctrk->GetMother();
3157 AliAODMCParticle* mctrkMother = NULL;
3158
3159 for(Int_t k=0 ; k<ntrks; k++){
3160
3161 AliAODMCParticle* mctrk2 = (AliAODMCParticle*)farray->At(k);
3162 if(mctrk2->GetLabel() == labelMother) mctrkMother = mctrk2;
3163 }
3164
3165 if(mctrkMother) pdgMotherChK = mctrkMother->GetPdgCode();
3166 FillPIDhisto(fhPIDMCMotherChK,pdgMotherChK);
3167
3168 //printf("pdgMotherChK %d \n",pdgMotherChK);
3169 }
3170
3171 if(pdg == 310){ // K0
3172 Int_t labelMother = mctrk->GetMother();
3173 AliAODMCParticle* mctrkMother = NULL;
3174
3175 for(Int_t k=0 ; k<ntrks; k++){
3176 AliAODMCParticle* mctrk2 = (AliAODMCParticle*)farray->At(k);
3177 if(mctrk2->GetLabel() == labelMother) mctrkMother = mctrk2;
3178 }
3179
3180 if(mctrkMother) pdgMotherK0 = mctrkMother->GetPdgCode();
3181 FillPIDhisto(fhPIDMCMotherK0,pdgMotherK0);
3182
3183 Int_t labelGrandMother = -1;
3184 if(mctrkMother) mctrkMother->GetMother();
3185 AliAODMCParticle* mctrkGrandMother = NULL;
3186
3187 for(Int_t k=0 ; k<ntrks; k++){
3188 AliAODMCParticle* mctrk2 = (AliAODMCParticle*)farray->At(k);
3189 if(mctrk2->GetLabel() == labelGrandMother) mctrkGrandMother = mctrk2;
3190 }
3191
3192 if(mctrkGrandMother) pdgGrandMotherK0 = mctrkGrandMother->GetPdgCode();
3193 FillPIDhisto(fhPIDMCGrandMotherK0,pdgGrandMotherK0);
3194 }
3195
9b2de807 3196 if(jetVect){ // jet event
3197
3198 FillPIDhisto(fhPIDMCAll,pdg);
3199
3200 TVector3 partVect(mctrk->Px(), mctrk->Py(), mctrk->Pz());
3201
3202 Int_t region = IsTrackInsideRegion(jetVect, &partVect );
3203
3204 //Double_t deltaPhi = jetVect->DeltaPhi(partVect); //+k270rad;
3205 Double_t deltaPhi = (jetVect->MomentumVector()->Vect()).DeltaPhi(partVect);
3206 if( deltaPhi > 2.*TMath::Pi() ) deltaPhi-= 2.*TMath::Pi();
3207 //Double_t deltaR = jetVect->DeltaR(partVect); //+k270rad;
3208 Double_t deltaR = (jetVect->MomentumVector()->Vect()).DeltaR(partVect);
3209
9b2de807 3210 // calc xi
3211 Double_t jetpt = jetVect->Pt();
9280dfa6 3212 //Double_t pV0 = partVect.Mag();
3213 //Double_t ptV0 = partVect.Pt();
3214
3215
3216 Double_t z = trackPt / jetpt;
9b2de807 3217 Double_t xi = TMath::Log(1/z);
3218
3219 if(!(mctrk->Charge() == 0 || mctrk->Charge()==-99)) fhdNdxiMC->Fill(xi);
3220
9280dfa6 3221 if(pdg == 310){ // K0
9b2de807 3222 fhdPhiJetK0MC->Fill(deltaPhi);
3223 fhdRJetK0MC->Fill(deltaR);
9b2de807 3224 fhdNdxiK0MC->Fill(xi);
3225 fhdNdzK0MC->Fill(z);
9280dfa6 3226 fhdNdptK0MCJetEvt->Fill(trackPt);
9b2de807 3227 }
3228
3229
3230 if(region != 0 && region == maxPtRegionIndex){ // max region
3231
9280dfa6 3232 if(TMath::Abs(pdg) == 211) fhdNdptchPiMCMax->Fill(trackPt);
3233 if(pdg == 310) fhdNdptK0MCMax->Fill(trackPt);
3234 if(TMath::Abs(pdg) == 321) fhdNdptchKMCMax->Fill(trackPt);
3235 if(pdg == 2212) fhdNdptpMCMax->Fill(trackPt);
3236 if(pdg == -2212) fhdNdptpBarMCMax->Fill(trackPt);
3237 if(pdg == 3122) fhdNdptLambdaMCMax->Fill(trackPt);
3238 if(pdg == -3122) fhdNdptLambdaBarMCMax->Fill(trackPt);
9b2de807 3239 }
3240
3241 if(region != 0 && region != maxPtRegionIndex){ // min region
3242
3243 FillPIDhisto(fhPIDMCMin,pdg);
3244
9280dfa6 3245 if(TMath::Abs(pdg) == 211) fhdNdptchPiMCMin->Fill(trackPt);
3246 if(pdg == 310) fhdNdptK0MCMin->Fill(trackPt);
3247 if(TMath::Abs(pdg) == 321) fhdNdptchKMCMin->Fill(trackPt);
3248 if(pdg == 2212) fhdNdptpMCMin->Fill(trackPt);
3249 if(pdg == -2212) fhdNdptpBarMCMin->Fill(trackPt);
3250 if(pdg == 3122) fhdNdptLambdaMCMin->Fill(trackPt);
3251 if(pdg == -3122) fhdNdptLambdaBarMCMin->Fill(trackPt);
3252
3253 if(pdg == 3332) fhdNdptOmegaMCMin->Fill(trackPt);
3254 if(pdg == -3332) fhdNdptOmegaBarMCMin->Fill(trackPt);
3255 }
9b2de807 3256
9280dfa6 3257 // trans region
3258 if(region != 0){
3259 if(TMath::Abs(pdg) == 321) FillPIDhisto(fhPIDMCMotherChKTrans,pdgMotherChK);
3260 if(pdg == 310){
3261 FillPIDhisto(fhPIDMCMotherK0Trans,pdgMotherK0);
3262 FillPIDhisto(fhPIDMCGrandMotherK0Trans,pdgGrandMotherK0);
9b2de807 3263 }
9b2de807 3264 }
9280dfa6 3265
9b2de807 3266 if(region == 0){ // jet region ?
3267
3268 FillPIDhisto(fhPIDMCJet,pdg);
3269
3270 //Double_t dRJetV0 = jetVect->DeltaR(partVect);
3271 Double_t dRJetV0 = (jetVect->MomentumVector()->Vect()).DeltaR(partVect);
3272
3273 if(dRJetV0 <= fConeRadius){
3274
9280dfa6 3275 if(pdg == 310){ // K0
9b2de807 3276
9280dfa6 3277 fhdNdptK0MCJet->Fill(trackPt);
9b2de807 3278 fhdNdxiK0MCJet->Fill(xi);
3279 fhdNdzK0MCJet->Fill(z);
3280 }
9280dfa6 3281
3282 if(TMath::Abs(pdg) == 211) fhdNdptchPiMCJet->Fill(trackPt);
3283 if(TMath::Abs(pdg) == 321) fhdNdptchKMCJet->Fill(trackPt);
3284 if(pdg == 2212) fhdNdptpMCJet->Fill(trackPt);
3285 if(pdg == -2212) fhdNdptpBarMCJet->Fill(trackPt);
3286 if(pdg == 3122) fhdNdptLambdaMCJet->Fill(trackPt);
3287 if(pdg == -3122) fhdNdptLambdaBarMCJet->Fill(trackPt);
9b2de807 3288 }
3289 }
3290 }
3291 }
3292
3293
3294 FillPythiaIDhisto(fhPythiaProcess,pythiaPID);
3295 if(isK0Event) FillPythiaIDhisto(fhPythiaProcessK0,pythiaPID);
3296 if(isKchEvent) FillPythiaIDhisto(fhPythiaProcessKch,pythiaPID);
3297 if(ispEvent) FillPythiaIDhisto(fhPythiaProcessp,pythiaPID);
3298 if(ispBarEvent) FillPythiaIDhisto(fhPythiaProcesspbar,pythiaPID);
3299
3300}
3301
3302// ----------------------------------------------------------------------------
3303
3304Double_t AliAnalysisTaskJetChem::AssociateV0MC(const TVector3* V0Mom,const Int_t pdg){
3305
3306 // find closest MC gen. particle for V0 vector
3307
3308 TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
3309
3310 if(!farray){
3311 AliInfo("no mcparticles branch");
3312 return -1;
3313 }
3314
3315 Double_t dRmin = -1;
3316
3317 Int_t ntrks = farray->GetEntries();
3318
3319 for(Int_t i =0 ; i<ntrks; i++){
3320
3321 AliAODMCParticle* mctrk = (AliAODMCParticle*)farray->At(i);
3322
3323 //Cuts
3324 if (!(mctrk->IsPhysicalPrimary())) continue;
3325
3326 Int_t pdgtrk = mctrk->GetPdgCode();
3327
3328 if(pdgtrk != pdg) continue;
3329
3330 TVector3 partVect(mctrk->Px(), mctrk->Py(), mctrk->Pz());
3331
3332 Double_t dR = V0Mom->DeltaR(partVect);
3333
3334 if(dRmin<0) dRmin = dR; // initialize
3335
3336 if(dR < dRmin) dRmin = dR;
3337
3338 }
3339
3340 return dRmin;
3341}
3342
3343// ----------------------------------------------------------------------------
3344
3345void AliAnalysisTaskJetChem::CompLeadingJets(AliAODJet* jetLeadingAOD,AliAODJet* jetLeadingMC,const Int_t pythiaPID,
3346 const Bool_t foundK0AOD,const Bool_t foundK0MC){
3347
3348 // leading jet properties
3349
3350 Double_t ptLeadingAOD = -1;
3351 Double_t etaLeadingAOD = -1;
3352 Double_t phiLeadingAOD = -1;
3353 Int_t nTracksLeadingAOD = -1;
3354
3355 Double_t ptLeadingMC = -1;
3356 Double_t etaLeadingMC = -1;
3357 Double_t phiLeadingMC = -1;
3358 Int_t nTracksLeadingMC = -1;
3359
3360 if(jetLeadingAOD){
3361 if(jetLeadingAOD->Pt()>fJetPtCut && TMath::Abs(jetLeadingAOD->Eta())<fJetEtaCut){
3362
3363 ptLeadingAOD = jetLeadingAOD->Pt();
3364 etaLeadingAOD = jetLeadingAOD->Eta();
3365 phiLeadingAOD = jetLeadingAOD->Phi();
3366 nTracksLeadingAOD = jetLeadingAOD->GetRefTracks()->GetEntriesFast();
3367
3368 Double_t radiusAOD = GetJetRadius(jetLeadingAOD,0.8);
3369 fhnTracksJetVsPtAOD->Fill(ptLeadingAOD,nTracksLeadingAOD);
3370 fhRadiusJetVsPtAOD->Fill(ptLeadingAOD,radiusAOD);
3371 if(IsQuarkHardScatteringEvent(pythiaPID)) fhnTracksJetVsPtAODquarkEv->Fill(ptLeadingAOD,nTracksLeadingAOD);
3372
3373 if(foundK0AOD){
3374
3375 fhnTracksJetVsPtAODK0->Fill(ptLeadingAOD,nTracksLeadingAOD);
3376 if(IsQuarkHardScatteringEvent(pythiaPID)) fhnTracksJetVsPtAODK0quarkEv->Fill(ptLeadingAOD,nTracksLeadingAOD);
3377 fhRadiusJetVsPtAODK0->Fill(ptLeadingAOD,radiusAOD);
3378 }
3379
3380
3381 // check if p/Kch in jet
3382
3383 Bool_t foundpKch = kFALSE;
3384 Int_t nTracksJet = jetLeadingAOD->GetRefTracks()->GetEntriesFast();
3385
3386 for(int i=0; i<nTracksJet; i++){
3387
3388 AliAODTrack* track = (AliAODTrack*) jetLeadingAOD->GetRefTracks()->At(i);
3389
3390 Double_t mom = track->P();
3391
3392 AliAODPid* aodPid = track->GetDetPid();
3393 Double_t dEdx = aodPid->GetTPCsignal();
3394
3395 if(IsAccepteddEdx(mom,dEdx,AliPID::kKaon,fCutnSigdEdx) ||
3396 IsAccepteddEdx(mom,dEdx,AliPID::kProton,fCutnSigdEdx)){
3397
3398 foundpKch = kTRUE;
3399 }
3400 } // track loop
3401
3402 if(foundpKch){
3403 fhnTracksJetVsPtAODpKch->Fill(ptLeadingAOD,nTracksLeadingAOD);
3404 fhRadiusJetVsPtAODpKch->Fill(ptLeadingAOD,radiusAOD);
3405 }
3406 }
3407 }
3408
3409
3410 if(jetLeadingMC){
3411 if(jetLeadingMC->Pt()>fJetPtCut && TMath::Abs(jetLeadingMC->Eta())<fJetEtaCut){
3412
3413 ptLeadingMC = jetLeadingMC->Pt();
3414 etaLeadingMC = jetLeadingMC->Eta();
3415 phiLeadingMC = jetLeadingMC->Phi();
3416 nTracksLeadingMC = jetLeadingMC->GetRefTracks()->GetEntriesFast();
3417
3418 Double_t radiusMC = GetJetRadius(jetLeadingMC,0.8);
3419 fhnTracksJetVsPtMC->Fill(ptLeadingMC,nTracksLeadingMC);
3420 fhRadiusJetVsPtMC->Fill(ptLeadingMC,radiusMC);
3421 if(IsQuarkHardScatteringEvent(pythiaPID)) fhnTracksJetVsPtMCquarkEv->Fill(ptLeadingMC,nTracksLeadingMC);
3422
3423 if(foundK0MC){
3424
3425 fhnTracksJetVsPtMCK0->Fill(ptLeadingMC,nTracksLeadingMC);
3426 if(IsQuarkHardScatteringEvent(pythiaPID)) fhnTracksJetVsPtMCK0quarkEv->Fill(ptLeadingMC,nTracksLeadingMC);
3427 fhRadiusJetVsPtMCK0->Fill(ptLeadingMC,radiusMC);
3428 }
3429 }
3430 }
3431
3432 if(jetLeadingAOD && jetLeadingMC){
3433
3434 //std::cout<<" comp: leading jetPt AOD "<<ptLeadingAOD<<" MC "<<ptLeadingMC<<std::endl;
3435 //if(jetLeadingAOD && jetLeadingMC)
3436 //std::cout<<" leading jet eta AOD "<<jetLeadingAOD->Eta()<<" MC "<<jetLeadingMC->Eta()<<std::endl;
3437
3438
3439 if(jetLeadingMC->Pt()>fJetPtCut && jetLeadingAOD->Pt()>fJetPtCut &&
3440 TMath::Abs(jetLeadingMC->Eta())<fJetEtaCut && TMath::Abs(jetLeadingAOD->Eta())<fJetEtaCut){
3441
3442 fhLeadingPtAODvsMC->Fill(ptLeadingMC,ptLeadingAOD);
3443 fhLeadingEtaAODvsMC->Fill(etaLeadingMC,etaLeadingAOD);
3444 fhLeadingPhiAODvsMC->Fill(phiLeadingMC,phiLeadingAOD);
3445 fhnTracksLeadingAODvsMC->Fill(nTracksLeadingMC,nTracksLeadingAOD);
3446
3447 TLorentzVector *mom4MC = jetLeadingMC->MomentumVector();
3448 TLorentzVector *mom4AOD = jetLeadingAOD->MomentumVector();
3449
3450 Double_t dR = mom4MC->DeltaR(*mom4AOD);
3451 fhLeadingdRAODMC->Fill(dR);
3452
3453 if(dR<0.4) fhLeadingPtAODvsMCdRcut->Fill(ptLeadingMC,ptLeadingAOD);
3454
3455 Double_t dPt = ptLeadingAOD - ptLeadingMC;
3456 Double_t dnTracks = nTracksLeadingAOD - nTracksLeadingMC;
3457
3458 fhdnTracksVsdPtLeadingAODMC->Fill(dPt,dnTracks);
3459 }
3460 }
3461
3462}
3463
3464// ---------------------------------------------------------------------------
3465
3466// void AliAnalysisTaskJetChem::CheckK0(){
3467
3468// TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
3469
3470// if(!farray){
3471// AliInfo("no mcparticles branch");
3472// return;
3473// }
3474
3475
3476// Int_t ntrks = farray->GetEntries();
3477
3478// for(Int_t i=0; i<ntrks; i++){ // trk loop
3479
3480// AliAODMCParticle* mctrk = (AliAODMCParticle*)farray->At(i);
3481
3482// Int_t pdg = mctrk->GetPdgCode();
3483
3484// if(pdg != 310) continue; // K0
3485
3486// cout<<" MC K0: physPrim "<<mctrk->IsPhysicalPrimary()<<endl;
3487
3488// Int_t labelMother = mctrk->GetMother();
3489
3490// cout<<" found K00, label mother "<<labelMother<<endl;
3491
3492// AliAODMCParticle* mctrkMother = NULL;
3493// Int_t pdgMother = -1;
3494
3495// for(Int_t k=0 ; k<ntrks; k++){
3496
3497// mctrkMother = (AliAODMCParticle*)farray->At(k);
3498// if(mctrkMother->GetLabel() == labelMother) break;
3499// }
3500
3501// pdgMother = mctrkMother->GetPdgCode();
3502// cout<<" K0 mother pdg "<<pdgMother<<" GID "<<fpdgdb->ConvertPdgToGeant3(pdgMother)<<" isPrimary "<<mctrkMother->IsPrimary()<<endl;
3503// //cout<<" mother name "<<mctrkMother->GetName()<<endl;
3504
3505// }
3506// }
3507
3508// ---------------------------------------------------------------------------
3509
3510// void AliAnalysisTaskJetChem::CheckK0Stack(){
3511
3512// AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
3513// if (!mcHandler) {
3514// Printf("ERROR: Could not retrieve MC event handler");
3515// return;
3516// }
3517
3518// AliMCEvent* mcEvent = mcHandler->MCEvent();
3519// if (!mcEvent) {
3520// Printf("ERROR: Could not retrieve MC event");
3521// return;
3522// }
3523
3524// AliStack* mcStack = mcEvent->Stack();//Load Stack
3525// AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
3526
3527// Int_t nTracksMC = mcStack->GetNtrack();
3528
3529// Bool_t foundK0 = kFALSE;
3530
3531// for (Int_t iTracks = 0; iTracks < nTracksMC; iTracks++) {
3532// //Cuts
3533// if(!(mcStack->IsPhysicalPrimary(iTracks))) continue;
3534
3535// TParticle* mctrk = mcStack->Particle(iTracks);
3536
3537// Int_t pdg = mctrk->GetPdgCode();
3538
3539// if ((mctrk->Pt() < 2*fTrackPtCut) || (TMath::Abs(mctrk->Eta()) > fTrackEtaCut )) continue;
3540
3541// if(pdg == 310) foundK0 = kTRUE;
3542
3543// if(pdg != 310) continue; // K0 short
3544
3545// cout<<" check K0 "<<endl;
3546
3547// Int_t indexMother = -999;
3548// TParticle* mctrkMother = mctrk;
3549// TParticle* mctrkFirstGenMother = NULL; // last mother which is not the primary proton
3550
3551// Int_t nAncestors = 0;
3552
3553// while(indexMother != -1){
3554
3555// indexMother = mctrkMother->GetFirstMother();
3556
3557// if(indexMother != -1){
3558// mctrkFirstGenMother = mctrkMother;
3559// mctrkMother = mcStack->Particle(indexMother);
3560// nAncestors++;
3561// }
3562
3563// cout<<" nAncestors "<<nAncestors<<" pdg mother "<<mctrkMother->GetPdgCode()<<" name "<<mctrkMother->GetName()<<endl;
3564// }
3565// cout<<" pdg firstGenMother "<<mctrkFirstGenMother->GetPdgCode()<<" name "<<mctrkFirstGenMother->GetName()<<endl;
3566
3567
3568
3569// cout<<" pythiaGenHeader "<<pythiaGenHeader<<endl;
3570// cout<<" Pythia Process type "<<pythiaGenHeader->ProcessType()<<" ptHard "<<pythiaGenHeader->GetPtHard()<<endl;
3571
3572
3573// fhPythiaProcess->Fill(pythiaGenHeader->ProcessType());
3574// if(foundK0) fhPythiaProcess_K0->Fill(pythiaGenHeader->ProcessType());
3575
3576// //Int_t indexGrandMother = mctrkMother->GetFirstMother();
3577// //cout<<" indexGrandMother "<<indexGrandMother<<endl;
3578
3579// //if(indexGrandMother>-1){
3580// // TParticle* mctrkGrandMother = mcStack->Particle(indexGrandMother);
3581// // cout<<" pdg grandMother "<<mctrkGrandMother->GetPdgCode()<<" name "<<mctrkGrandMother->GetName()<<endl;
3582// // }
3583
3584// }
3585// }
3586
3587// ---------------------------------------------------------------------------------
3588
3589Bool_t AliAnalysisTaskJetChem::IsQuarkHardScatteringEvent(const Int_t PID){
3590
3591 // Pythia Manual sec. 8.2.1 :
3592 // if Pythia PID = 92,93,94 event is diffractive
3593 // if Pythia PID = 13, 28, 68 hard scattering products are gg, qg
3594 // if Pythia PID = 11, 12, 53 hard scattering products are qq, q\bar{q}
3595
9280dfa6 3596 if(PID == 92 || PID == 93 || PID == 94 || PID==95) return kFALSE;
9b2de807 3597 else if(PID == 13 || PID == 28 || PID == 68) return kFALSE;
3598 else if(PID == 11 || PID == 12 || PID == 53) return kTRUE;
3599 else{
3600 AliInfo(Form("unknown Pythia PID %d",PID));
3601 }
3602
3603 return kFALSE;
3604}
3605
9280dfa6 3606
3607// ---------------------------------------------------------------------------------
3608
3609Bool_t AliAnalysisTaskJetChem::IsGluonHardScatteringEvent(const Int_t PID){
3610
3611 // Pythia Manual sec. 8.2.1 :
3612 // if Pythia PID = 92,93,94 event is diffractive
3613 // if Pythia PID = 95: low pt event (MPI)
3614 // if Pythia PID = 13, 28, 68 hard scattering products are gg, qg
3615 // if Pythia PID = 11, 12, 53 hard scattering products are qq, q\bar{q}
3616
3617
3618 if(PID == 92 || PID == 93 || PID == 94 || PID == 95) return kFALSE;
3619 else if(PID == 13 || PID == 68) return kTRUE;
3620 else if(PID == 28) return kFALSE; // mixed gq final state
3621 else if(PID == 11 || PID == 12 || PID == 53) return kFALSE;
3622 else{
3623 AliInfo(Form("unknown Pythia PID %d",PID));
3624 }
3625
3626 return kFALSE;
3627}
3628
9b2de807 3629// ---------------------------------------------------------------------------------
3630
3631Bool_t AliAnalysisTaskJetChem::IsDiffractiveEvent(const Int_t PID){
3632
3633 // Pythia Manual sec. 8.2.1 :
3634 // if Pythia PID = 92,93,94 event is diffractive
3635 // if Pythia PID = 13, 28, 68 hard scattering products are gg, qg
3636 // if Pythia PID = 11, 12, 53 hard scattering products are qq, q\bar{q}
3637
3638 if(PID == -1) return kFALSE;
3639
3640 if(PID == 13 || PID == 28 || PID == 68) return kFALSE;
3641 else if(PID == 11 || PID == 12 || PID == 53) return kFALSE;
9280dfa6 3642 else if(PID == 92 || PID == 93 || PID == 94 || PID==95) return kTRUE;
9b2de807 3643 else{
3644 AliInfo(Form("unknown Pythia PID %d",PID));
3645 }
3646
3647 return kFALSE;
3648
3649}
3650
3651
3652// ----------------------------------------------------------------------------------
3653
3654Int_t AliAnalysisTaskJetChem::GetPythiaProcessID(){
3655
3656 // Pythia PID for this event
3657
3658 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
3659 if (!mcHandler) {
3660 //Printf("ERROR: Could not retrieve MC event handler");
3661 return -1;
3662 }
3663
3664 AliMCEvent* mcEvent = mcHandler->MCEvent();
3665 if (!mcEvent) {
3568c3d6 3666 AliInfo("could not retrieve MC event");
9b2de807 3667 return -1;
3668 }
3669
3670 AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
3671
3568c3d6 3672 if (!pythiaGenHeader) {
3673 AliInfo("Could not retrieve pythiaEventHeader");
9b2de807 3674 return -1;
3675 }
3676
3677 Int_t pid = pythiaGenHeader->ProcessType();
3678
3679 return pid;
3680
3681}
9280dfa6 3682
3683// ----------------------------------------------------------------------------------
3684
3685void AliAnalysisTaskJetChem::GetJetTracksResum(TList* list, AliAODJet* jet, const Double_t radius){
3686
3687
3688 if(!jet) return; // no jet in acc in event
3689
3690 // list of AOD tracks in jet cone, using cone axis and distance axis-track (and not trackrefs)
3691
3692 Int_t nTracks = fAOD->GetNTracks();
3693
3694 if(!nTracks) return;
3695
3696 Double_t jetMom[3];
3697 jet->PxPyPz(jetMom);
3698 TVector3 jet3mom(jetMom);
3699
3700
3701 // standard cuts + ITS refit = stdrd PWG4 cut (compose them for productions with old ESDFilter task)
3702
3703 for (Int_t itrack=0; itrack<nTracks; itrack++) {
3704
3705 AliAODTrack* track = fAOD->GetTrack(itrack);
3706
3707 UInt_t status = track->GetStatus();
3708
3709 if(!track->TestFilterBit(fFilterBitJF)) continue; // track cut selection
3710 if(fRequireITSRefitJF && ((status&AliESDtrack::kITSrefit)==0)) continue;
3711
3712 if(TMath::Abs(track->Eta()) > fTrackEtaCut) continue;
3713 if( track->Pt() < fTrackPtCutJF ) continue;
3714
3715 Double_t trackMom[3];
3716 track->PxPyPz(trackMom);
3717 TVector3 track3mom(trackMom);
3718
3719 Double_t dR = jet3mom.DeltaR(track3mom);
3720
3721 if(dR<radius){
3722
3723 list->Add(track);
3724
3725 }
3726 }
3727}
3728
3729// ----------------------------------------------------------------------------------
3730
3731void AliAnalysisTaskJetChem::GetJetTracksTrackrefs(TList* list, AliAODJet* jet){
3732
3733 if(!jet) return; // no jet in acc in event
3734
3735 // list of AOD tracks in jet cone, using trackrefs
3736
3737 Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
3738
3739 if(!nTracks) return;
3740
3741 // standard cuts + ITS refit = stdrd PWG4 cut (compose them for productions with old ESDFilter task)
3742
3743 for (Int_t itrack=0; itrack<nTracks; itrack++) {
3744
3745 AliAODTrack* track = (AliAODTrack*) jet->GetRefTracks()->At(itrack);
3746
3747 UInt_t status = track->GetStatus();
3748
3749 if(!track->TestFilterBit(fFilterBitJF)) continue; // track cut selection
3750 if(fRequireITSRefitJF && ((status&AliESDtrack::kITSrefit)==0)) continue;
3751
3752 if(TMath::Abs(track->Eta()) > fTrackEtaCut) continue;
3753 if( track->Pt() < fTrackPtCutJF ) continue;
3754
3755 list->Add(track);
3756 }
3757
3758 //cout<<" trackrefs Size "<<nTracks<<" acc track list size "<<list->GetEntries()<<endl;
3759
3760}
3761
3762
3763// ----------------------------------------------------------------------------------
3764
3765void AliAnalysisTaskJetChem::FillReferenceFF(AliAODJet* jet){
3766
3767
3768 if(!jet) return;
3769
3770 TList* jetTracks = new TList(); // FIXME - avoid new/delete
3771 //GetJetTracksResum(jetTracks,jet,0.7);
3772 GetJetTracksTrackrefs(jetTracks,jet);
3773
3774 Double_t jetpt = jet->Pt();
3775
3776 TIter next(jetTracks);
3777 while(AliAODTrack* track = static_cast<AliAODTrack*>(next())){
3778
3779 Double_t trackpt = track->Pt();
3780 Double_t z = trackpt/jetpt;
3781 Double_t xi = TMath::Log(1/z);
3782
3783 //cout<<" trackpt "<<trackpt<<" jetpt "<<jetpt<<" z "<<z<<" xi "<<xi<<endl;
3784
3785 if(5<jetpt && jetpt<=10){
3786 fhdNdptTracksJetPt5to10->Fill(trackpt);
3787 fhdNdzJets5to10->Fill(z);
3788 fhdNdxiJets5to10->Fill(xi);
3789 }
3790 else if(10<jetpt && jetpt<=20){
3791 fhdNdptTracksJetPt10to20->Fill(trackpt);
3792 fhdNdzJets10to20->Fill(z);
3793 fhdNdxiJets10to20->Fill(xi);
3794 }
3795 else if(20<jetpt && jetpt<=30){
3796 fhdNdptTracksJetPt20to30->Fill(trackpt);
3797 fhdNdzJets20to30->Fill(z);
3798 fhdNdxiJets20to30->Fill(xi);
3799 }
3800 else if(30<jetpt && jetpt<=40){
3801 fhdNdptTracksJetPt30to40->Fill(trackpt);
3802 fhdNdzJets30to40->Fill(z);
3803 fhdNdxiJets30to40->Fill(xi);
3804 }
3805 else if(40<jetpt && jetpt<=60){
3806 fhdNdptTracksJetPt40to60->Fill(trackpt);
3807 fhdNdzJets40to60->Fill(z);
3808 fhdNdxiJets40to60->Fill(xi);
3809 }
3810 }
3811
3812
3813 delete jetTracks;
3814
3815}
3816
3817// -------------------------------------------------------------
3818
3819void AliAnalysisTaskJetChem::FillReferencePlotsTracks(){
3820
3821 // eta/phi & pt tracks before/after cuts
3822 // track multiplicity / evt
3823
3824 Int_t nTracks = fAOD->GetNTracks();
3825 Int_t countTracks = 0;
3826
3827 // standard cuts + ITS refit = stdrd PWG4 cut (compose them for productions with old ESDFilter task)
3828
3829 for(Int_t itrack=0; itrack<nTracks; itrack++){
3830
3831 AliAODTrack* track = fAOD->GetTrack(itrack);
3832
3833 Double_t trackPt = track->Pt();
3834 Double_t trackPhi = track->Phi();
3835 Double_t trackEta = track->Eta();
3836
3837 fhPhiEtaTracksNoCut->Fill(trackEta,trackPhi);
3838 fhPtTracksNoCut->Fill(trackPt);
3839
3840 UInt_t status = track->GetStatus();
3841
3842 if(!track->TestFilterBit(fFilterBitJF)) continue; // track cut selection
3843 if(fRequireITSRefitJF && ((status&AliESDtrack::kITSrefit)==0)) continue;
3844
3845 fhPhiEtaTracks->Fill(trackEta,trackPhi);
3846 fhPtTracks->Fill(trackPt);
3847
3848 if(TMath::Abs(track->Eta()) > fTrackEtaCut) continue;
3849 if( track->Pt() < fTrackPtCutJF ) continue;
3850
3851 countTracks++;
3852
3853 }
3854
3855 fhTrackMult->Fill(countTracks);
3856}
3857