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