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