1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
26 #include <TLorentzVector.h>
29 #include <TDatabasePDG.h>
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"
40 #include "AliGenPythiaEventHeader.h"
41 #include "AliAnalysisHelperJetTasks.h"
42 #include "AliESDtrack.h"
43 #include "AliExternalTrackParam.h"
45 //#include <iostream> // TEST!!!
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
53 // This class needs as input AOD with track and Jets
54 // the output is a list of histograms
56 // AOD can be either connected to the InputEventHandler
57 // for a chain of AOD files
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
67 ClassImp(AliAnalysisTaskJetChem)
69 ////////////////////////////////////////////////////////////////////////
72 //____________________________________________________________________
73 AliAnalysisTaskJetChem:: AliAnalysisTaskJetChem(const char* name): AliAnalysisTaskSE(name),
78 fDeltaAODBranchMC(""),
79 fAODBranchMC("jetsMC"),
86 fUseLOConeJets(kFALSE),
87 fUseLOConeMCJets(kFALSE),
88 fUsePythiaJets(kFALSE),
92 fRequireITSRefitJF(kFALSE),
93 fRejectK0TracksJF(kFALSE),
101 fUseAODMCTracksForUE(kFALSE),
104 fhPrimVertexNCont(0x0),
105 fhPrimVertexRho(0x0),
110 fhLeadingNTracksVsEta(0x0),
111 fhLeadingPtVsEta(0x0),
114 fhLeadingPtDiffr(0x0),
118 fhLeadingPtMCDiffr(0x0),
119 fhPhiEtaTracksNoCut(0x0),
120 fhPtTracksNoCut(0x0),
127 fhnTracksVsPtLeading(0x0),
128 fhdNdEtaPhiDist(0x0),
129 fhRegionSumPtMaxVsEt(0x0),
130 fhRegionMultMaxVsEt(0x0),
131 fhRegionSumPtMinVsEt(0x0),
132 fhRegionMultMinVsEt(0x0),
135 fhV0DCADaughters(0x0),
137 fhV0DCAToVertex(0x0),
138 fhV0DCAToVertexK0(0x0),
140 fhV0PtVsInvMassK0(0x0),
141 fhV0InvMassK0JetEvt(0x0),
142 fhV0InvMassLambda(0x0),
143 fhV0InvMassAntiLambda(0x0),
144 fhV0InvMassLambdaJetEvt(0x0),
145 fhV0InvMassAntiLambdaJetEvt(0x0),
153 fhV0InvMassK0DCA(0x0),
154 fhV0InvMassK0DCAdEdx(0x0),
155 fhV0PtVsInvMassK0DCAdEdx(0x0),
156 fhV0InvMassK0DCAPID(0x0),
157 fhV0InvMassLambdaDCAdEdx(0x0),
158 fhV0InvMassAntiLambdaDCAdEdx(0x0),
160 fhdNdptK0DCAdEdx(0x0),
161 fhV0InvMassK0Min(0x0),
162 fhV0InvMassLambdaMin(0x0),
163 fhV0InvMassAntiLambdaMin(0x0),
164 fhV0InvMassK0Max(0x0),
165 fhV0InvMassLambdaMax(0x0),
166 fhV0InvMassAntiLambdaMax(0x0),
167 fhV0InvMassK0Jet(0x0),
168 fhV0InvMassLambdaJet(0x0),
169 fhV0InvMassAntiLambdaJet(0x0),
170 fhV0InvMassK0Lambda(0x0),
171 fhdNdptK0JetEvt(0x0),
180 fhdNdzAntiLambda(0x0),
183 fhdNdzLambdaMax(0x0),
184 fhdNdxiLambdaMax(0x0),
186 fhdNdptLambdaMax(0x0),
189 fhdNdzLambdaMin(0x0),
190 fhdNdxiLambdaMin(0x0),
192 fhdNdptLambdaMin(0x0),
195 fhdNdzLambdaJet(0x0),
196 fhdNdxiLambdaJet(0x0),
198 fhdNdptLambdaJet(0x0),
200 fhdEdxVsMomV0pidEdx(0x0),
201 fhdEdxVsMomV0piPID(0x0),
205 fhdNdptchPiMCMax(0x0),
207 fhdNdptchKMCMax(0x0),
209 fhdNdptpBarMCMax(0x0),
210 fhdNdptLambdaMCMax(0x0),
211 fhdNdptLambdaBarMCMax(0x0),
212 fhdNdptchPiMCMin(0x0),
214 fhdNdptchKMCMin(0x0),
216 fhdNdptpBarMCMin(0x0),
217 fhdNdptLambdaMCMin(0x0),
218 fhdNdptLambdaBarMCMin(0x0),
219 fhdNdptOmegaMCMin(0x0),
220 fhdNdptOmegaBarMCMin(0x0),
221 fhdNdptchPiMCJet(0x0),
223 fhdNdptchKMCJet(0x0),
225 fhdNdptpBarMCJet(0x0),
226 fhdNdptLambdaMCJet(0x0),
227 fhdNdptLambdaBarMCJet(0x0),
229 fhPIDMC_quarkEv(0x0),
230 fhPIDMC_gluonEv(0x0),
234 fhPIDMCMotherK0(0x0),
235 fhPIDMCGrandMotherK0(0x0),
236 fhPIDMCMotherChK(0x0),
237 fhPIDMCMotherK0Trans(0x0),
238 fhPIDMCGrandMotherK0Trans(0x0),
239 fhPIDMCMotherChKTrans(0x0),
247 fhdNdptLambdaMC(0x0),
248 fhdNdptLambdaBarMC(0x0),
250 fhdNdptOmegaBarMC(0x0),
256 fhdNdptK0MCJetEvt(0x0),
258 fhLeadingPtAODvsMC(0x0),
259 fhLeadingEtaAODvsMC(0x0),
260 fhLeadingPhiAODvsMC(0x0),
261 fhnTracksLeadingAODvsMC(0x0),
262 fhLeadingdRAODMC(0x0),
263 fhLeadingPtAODvsMCdRcut(0x0),
264 fhdnTracksVsdPtLeadingAODMC(0x0),
265 fhnTracksJetVsPtAOD(0x0),
266 fhnTracksJetVsPtAODquarkEv(0x0),
267 fhRadiusJetVsPtAOD(0x0),
268 fhnTracksJetVsPtMC(0x0),
269 fhnTracksJetVsPtMCquarkEv(0x0),
270 fhRadiusJetVsPtMC(0x0),
271 fhnTracksJetVsPtMCK0(0x0),
272 fhnTracksJetVsPtMCK0quarkEv(0x0),
273 fhRadiusJetVsPtMCK0(0x0),
274 fhnTracksJetVsPtAODK0(0x0),
275 fhnTracksJetVsPtAODK0quarkEv(0x0),
276 fhRadiusJetVsPtAODK0(0x0),
277 fhnTracksJetVsPtAODpKch(0x0),
278 fhRadiusJetVsPtAODpKch(0x0),
279 fhPythiaProcess(0x0),
280 fhPythiaProcessK0(0x0),
281 fhPythiaProcessKch(0x0),
282 fhPythiaProcessp(0x0),
283 fhPythiaProcesspbar(0x0),
284 fhdNdzJets5to10(0x0),
285 fhdNdzJets10to20(0x0),
286 fhdNdzJets20to30(0x0),
287 fhdNdzJets30to40(0x0),
288 fhdNdzJets40to60(0x0),
289 fhdNdxiJets5to10(0x0),
290 fhdNdxiJets10to20(0x0),
291 fhdNdxiJets20to30(0x0),
292 fhdNdxiJets30to40(0x0),
293 fhdNdxiJets40to60(0x0),
294 fhdNdptTracksJetPt5to10(0x0),
295 fhdNdptTracksJetPt10to20(0x0),
296 fhdNdptTracksJetPt20to30(0x0),
297 fhdNdptTracksJetPt30to40(0x0),
298 fhdNdptTracksJetPt40to60(0x0),
302 // Default constructor
304 fAreaReg = 2*TMath::Pi()/6.0 * 2*fTrackEtaCut;
305 fpdgdb = TDatabasePDG::Instance();
307 // Output slot #1 writes into a TList container, 0 reserved for std AOD output
308 DefineOutput(1, TList::Class());
311 //______________________________________________________________
312 Bool_t AliAnalysisTaskJetChem::UserNotify()
315 // read the cross sections
316 // and number of trials from pyxsec.root
320 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
321 Float_t xsection = 0;
324 TFile *curfile = tree->GetCurrentFile();
326 Error("Notify","No current file");
329 if(!fh1Xsec||!fh1Trials){
330 Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
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;
342 //____________________________________________________________________
343 void AliAnalysisTaskJetChem::UserCreateOutputObjects()
345 // Create the output container
347 AliInfo("UserCreateOutPutObjects()");
353 PostData(1, fListOfHistos); // PostData at least once
357 //____________________________________________________________________
358 void AliAnalysisTaskJetChem::UserExec(Option_t */*option*/)
360 // Execute analysis for current event
362 if (fDebug > 3) AliInfo( " Processing event..." );
364 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
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
374 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
375 if( handler && handler->InheritsFrom("AliAODHandler") ) {
376 fAODjets = ((AliAODHandler*)handler)->GetAOD();
377 if (fDebug > 1) AliInfo(" ==== Jets from AliAODHandler");
381 if (fDebug > 1) AliInfo(" ==== Jets from AliAODInputHandler");
384 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
385 if( handler && handler->InheritsFrom("AliAODHandler") ) {
386 fAOD = ((AliAODHandler*)handler)->GetAOD();
388 if (fDebug > 1) AliInfo(" ==== Tracks and Jets from AliAODHandler");
390 AliFatal("I can't get any AOD Event Handler");
397 // fetch the pythia header info and get the trials
398 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
401 AliMCEvent* mcEvent = mcHandler->MCEvent();
403 AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
405 nTrials = pythiaGenHeader->Trials();
410 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
415 PostData(1, fListOfHistos);
418 //____________________________________________________________________
419 void AliAnalysisTaskJetChem::AnalyseEvent()
423 // TString firedTriggerClasses = fAOD->GetFiredTriggerClasses();
424 // AliInfo(Form("firedTriggerClasses %s",firedTriggerClasses.Data()));
425 // if(firedTriggerClasses.Length() > 0 && !firedTriggerClasses.Contains("CINT1B")) return;
427 AliAODVertex *primVertex = fAOD->GetPrimaryVertex(); // is this from tracks or from SPD ? SPD should only be used if no track vertex
429 AliInfo("no prim Vertex found - skip event");
430 fhPrimVertexNCont->Fill(-1);
434 TString primVtxName(primVertex->GetName());
436 if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
438 AliInfo("found TPC prim Vertex - skip event");
439 fhPrimVertexNCont->Fill(-1);
443 Int_t nVertContributors = primVertex->GetNContributors();
444 fhPrimVertexNCont->Fill(nVertContributors);
446 //cout<<" prim vertex name "<<primVertex->GetName()<<" nCont "<<nVertContributors<<endl;
448 if(nVertContributors<1){ // eventually check if not SPD vertex ???
449 AliInfo("prim vertex no contributors - skip event");
453 Double_t vertX = primVertex->GetX();
454 Double_t vertY = primVertex->GetY();
455 Double_t vertZ = primVertex->GetZ();
457 Double_t vertRho = TMath::Sqrt(vertX*vertX+vertY*vertY);
459 fhPrimVertexRho->Fill(vertRho);
460 fhPrimVertexZ->Fill(vertZ);
462 if(TMath::Abs(vertZ)>10){
463 AliInfo(Form("prim vertex z=%f - skip event",vertZ));
468 Int_t pythiaPID = GetPythiaProcessID();
470 // ------------------------------------------------
471 // Find Leading Jets 1,2,3
472 // (could be skipped if Jets are sort by Pt...)
475 //Int_t nTracksLeading = 0;
478 AliAODJet* leadingJetAOD = NULL; // non-zero if any jet in acc and leading jet survives pt cut
479 Int_t indexLeadingAOD = -1;
480 Double_t ptLeadingAOD = 0.;
481 Int_t indexMaxRegionAOD = 0; // initialize with '0', '+-1' transverse regions
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
489 Double_t ptLeadingAODAllEta = 0;
490 AliAODJet* leadingJetAODAllEta = NULL;
492 if(fUseLOConeJets) fArrayJetsAOD = FindChargedParticleJets();
493 else{ // Use jets on AOD/DeltaAOD
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());
500 AliFatal(" No jet-array! ");
505 if (fDebug > 1) AliInfo(" ==== AOD jets: Read Standard-AODs !");
506 if (fDebug > 1) AliInfo(Form(" ==== Reading Branch: %s ",fAODBranch.Data()));
508 fArrayJetsAOD = (TClonesArray*)fAODjets->FindListObject(fAODBranch.Data());
513 if(fUseLOConeMCJets) fArrayJetsMC = FindChargedParticleJetsMC();
514 else if(fUsePythiaJets) fArrayJetsMC = GetPythiaJets();
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());
522 AliFatal(" No jet-array! ");
527 if (fDebug > 1) AliInfo(" ==== MC jets: Read Standard-AODs !");
528 if (fDebug > 1) AliInfo(Form(" ==== Reading Branch: %s ",fAODBranchMC.Data()));
530 fArrayJetsMC = (TClonesArray*)fAODjets->FindListObject(fAODBranchMC.Data());
534 if(fArrayJetsAOD) nJetsAOD = fArrayJetsAOD->GetEntries();
535 if(fArrayJetsMC) nJetsMC = fArrayJetsMC->GetEntries();
537 fhNJets->Fill(nJetsAOD);
538 fhNJetsMC->Fill(nJetsMC);
539 fhnJetsAODvsMC->Fill(nJetsMC,nJetsAOD);
541 if(fDebug>1) AliInfo(Form("AOD %d jets",nJetsAOD));
544 // for xcheck: find leading jet in large eta range
545 for(Int_t i=0; i<nJetsAOD; i++){
547 AliAODJet* jet = (AliAODJet*) fArrayJetsAOD->At(i);
548 Double_t jetPt = jet->Pt();
550 if(jetPt > ptLeadingAODAllEta){
551 ptLeadingAODAllEta = jetPt;
552 leadingJetAODAllEta = jet;
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());
563 // find leading jet AOD
564 for(Int_t i=0; i<nJetsAOD; ++i){
566 AliAODJet* jet = (AliAODJet*) fArrayJetsAOD->At(i);
567 Double_t jetPt = jet->Pt();
568 Double_t jetEta = jet->Eta();
570 if((jetPt > ptLeadingAOD) && (TMath::Abs(jetEta)<fJetEtaCut)){
571 ptLeadingAOD = jetPt;
576 // find leading jet MC
577 for(Int_t i=0; i<nJetsMC; ++i){
579 AliAODJet* jet = (AliAODJet*) fArrayJetsMC->At(i);
580 Double_t jetPt = jet->Pt();
581 Double_t jetEta = jet->Eta();
583 if((jetPt > ptLeadingMC) && (TMath::Abs(jetEta)<fJetEtaCut)){
590 //cout<<" ptLeadingAOD "<<ptLeadingAOD<<" MC "<<ptLeadingMC<<endl;
592 if(indexLeadingAOD>=0){ // event with jet
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();
598 if(fDebug>1) AliInfo(Form("\n Pt Leading AOD Jet = %6.1f eta=%5.3f, nTracks %d",ptLeadingAOD,etaLeadingAOD,nTracksLeadingAOD));
600 fhLeadingEta->Fill(etaLeadingAOD);
602 if(TMath::Abs(etaLeadingAOD)<fJetEtaCut){ // leading jet eta cut
604 fhnTracksVsPtLeading->Fill(ptLeadingAOD,nTracksLeadingAOD);
605 fhLeadingPt->Fill(ptLeadingAOD);
606 if(IsDiffractiveEvent(pythiaPID)) fhLeadingPtDiffr->Fill(ptLeadingAOD);
608 if(ptLeadingAOD>fJetPtCut){ // leading jet pt cut
610 leadingJetAOD = (AliAODJet*) fArrayJetsAOD->At(indexLeadingAOD);
612 fhLeadingPhi->Fill(phiLeadingAOD);
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;
621 Int_t nTracks = fAOD->GetNTracks();
623 for (Int_t ipart=0; ipart<nTracks; ++ipart) {
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;
637 TVector3 partVect(part->Px(), part->Py(), part->Pz());
639 Int_t region = IsTrackInsideRegion(leadingJetAOD,&partVect );
642 sumPtRegionPosit += part->Pt();
646 sumPtRegionNegat += part->Pt();
651 // fill sumPt and mult density for transverse regions
653 if( sumPtRegionPosit > sumPtRegionNegat ) {
654 FillSumPtRegion( ptLeadingAOD, sumPtRegionPosit/fAreaReg, sumPtRegionNegat/fAreaReg );
657 FillSumPtRegion( ptLeadingAOD, sumPtRegionNegat/fAreaReg, sumPtRegionPosit/fAreaReg );
659 if (nTrackRegionPosit > nTrackRegionNegat ) {
660 FillMultRegion( ptLeadingAOD, nTrackRegionPosit/fAreaReg, nTrackRegionNegat/fAreaReg);
663 FillMultRegion( ptLeadingAOD, nTrackRegionNegat/fAreaReg, nTrackRegionPosit/fAreaReg);
666 indexMaxRegionAOD = (sumPtRegionPosit > sumPtRegionNegat) ? 1 : -1;
668 } // leading jet pt cut
669 } // leading jet eta cut
673 if(indexLeadingMC>=0){ // event with MC jet
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();
679 if(fDebug>1) AliInfo(Form("\n Pt Leading MC Jet = %6.1f eta=%5.3f, nTracks %d",ptLeadingMC,etaLeadingMC,nTracksLeadingMC));
681 fhLeadingEtaMC->Fill(etaLeadingMC);
683 if(TMath::Abs(etaLeadingMC)<fJetEtaCut){ // leading jet eta cut
685 fhLeadingPtMC->Fill(ptLeadingMC);
686 if(IsDiffractiveEvent(pythiaPID)) fhLeadingPtMCDiffr->Fill(ptLeadingMC);
688 if(ptLeadingMC>fJetPtCut){ // leading jet pt cut
690 leadingJetMC = (AliAODJet*) fArrayJetsMC->At(indexLeadingMC);
692 fhLeadingPhiMC->Fill(phiLeadingMC); // -pi to pi
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;
701 TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
703 Int_t ntrks = farray->GetEntries();
704 if (fDebug>1) AliInfo(Form("In UE MC analysis tracks %d \n",ntrks));
706 for(Int_t i =0 ; i < ntrks; i++){
708 AliAODMCParticle* mctrk = (AliAODMCParticle*)farray->At(i);
710 if (!(mctrk->IsPhysicalPrimary())) continue;
711 //if (!(mctrk->IsPrimary())) continue;
713 if (mctrk->Charge() == 0 || mctrk->Charge()==-99) continue;
715 if (mctrk->Pt() < fTrackPtCut ) continue;
716 if( TMath::Abs(mctrk->Eta()) > fTrackEtaCut ) continue;
718 Bool_t isHadron = TMath::Abs(mctrk->GetPdgCode())==211 ||
719 TMath::Abs(mctrk->GetPdgCode())==2212 ||
720 TMath::Abs(mctrk->GetPdgCode())==321;
722 if (!isHadron) continue;
724 TVector3 partVect(mctrk->Px(), mctrk->Py(), mctrk->Pz());
726 Int_t region = IsTrackInsideRegion(leadingJetMC,&partVect );
729 sumPtRegionPosit += mctrk->Pt();
733 sumPtRegionNegat += mctrk->Pt();
736 } // AliAODMCParticle loop
738 indexMaxRegionMC = (sumPtRegionPosit > sumPtRegionNegat) ? 1 : -1;
745 Bool_t foundK0AOD = kFALSE; // event with K0 ?
746 Bool_t foundK0MC = kFALSE; // event with K0 ?
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);
751 FillReferencePlotsTracks();
752 FillReferenceFF(leadingJetAOD);
755 if(fUseLOConeJets && fArrayJetsAOD){
756 fArrayJetsAOD->Delete(); // no 'Clear': AliAODjet contains TMomentum and TRefArray
757 delete fArrayJetsAOD;
759 if(fUseLOConeMCJets && fArrayJetsMC){
760 fArrayJetsMC->Delete(); // no 'Clear': AliAODjet contains TMomentum and TRefArray
765 // __________________________________________________________________
768 Double_t AliAnalysisTaskJetChem::GetJetRadius(const AliAODJet* jet, const Double_t energyFrac){
770 // calc jet radius containing fraction energyFrac of full jet pt
772 const Int_t kArraySize = 1000;
773 const Int_t kInitVal = -999;
775 Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
777 if(nTracks>kArraySize){
778 AliError(Form("nTracks in jet %d exceeds max array size",nTracks));
783 Double_t deltaR[kArraySize];
784 Double_t pt[kArraySize];
785 Int_t index[kArraySize];
786 for(int i=0; i<kArraySize; i++) index[i] = kInitVal;
791 for(int i=0; i<nTracks; i++){
793 AliAODTrack* track = (AliAODTrack*) jet->GetRefTracks()->At(i);
795 TLorentzVector *mom4Jet = jet->MomentumVector();
796 TVector3 mom3Jet = mom4Jet->Vect();
798 Double_t trackMom[3];
799 track->PxPyPz(trackMom);
800 TVector3 mom3Track(trackMom);
802 Double_t dR = mom3Jet.DeltaR(mom3Track);
810 //cout<<" ptTot "<<ptTot<<" jetPt "<<jet->Pt()<<endl; // Xcheck
812 TMath::Sort(nTracks,deltaR,index,kFALSE); // sort in decreasing order
816 for(int i=0; i<nTracks; i++){
818 Int_t ind = index[i];
820 Double_t ptTrack = pt[ind];
821 Double_t dR = deltaR[ind];
825 if(ptSum >= ptTot*energyFrac) return dR;
832 //____________________________________________________________________
833 void AliAnalysisTaskJetChem::FillSumPtRegion( Double_t leadingE, Double_t ptMax, Double_t ptMin )
835 // Fill sumPt of control regions
837 fhRegionSumPtMaxVsEt->Fill( leadingE, ptMax );
838 fhRegionSumPtMinVsEt->Fill( leadingE, ptMin );
841 //____________________________________________________________________
842 void AliAnalysisTaskJetChem::FillMultRegion(Double_t leadingE, Double_t nTrackPtmax, Double_t nTrackPtmin)
844 // Fill Nch multiplicity of control regions
846 fhRegionMultMaxVsEt->Fill( leadingE, nTrackPtmax );
847 fhRegionMultMinVsEt->Fill( leadingE, nTrackPtmin );
851 //____________________________________________________________________
852 Int_t AliAnalysisTaskJetChem::IsTrackInsideRegion(const AliAODJet* aodjetVect,const TVector3 *partVect)
854 // return de region in delta phi
855 // -1 negative delta phi
856 // 1 positive delta phi
859 TLorentzVector* jetVectLorentz = aodjetVect->MomentumVector();
860 TVector3 jetVect = jetVectLorentz->Vect();
862 static const Double_t k60rad = 60.*TMath::Pi()/180.;
863 static const Double_t k120rad = 120.*TMath::Pi()/180.;
866 if( TMath::Abs(partVect->Eta()) > fTrackEtaCut ) return 0;
867 // transverse regions
868 if (jetVect.DeltaPhi(*partVect) < -k60rad && jetVect.DeltaPhi(*partVect) > -k120rad ) region = -1;
869 if (jetVect.DeltaPhi(*partVect) > k60rad && jetVect.DeltaPhi(*partVect) < k120rad ) region = 1;
874 //____________________________________________________________________
875 TClonesArray* AliAnalysisTaskJetChem::FindChargedParticleJetsMC()
878 // loop over pt-ordered list of tracks, combine tracks within jet cone, recalc cone axis after each step
879 // based on implementation by Arian Abrahantes Quintana and Enesto Lopez
881 TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
883 AliInfo("no mcparticles branch");
887 Int_t nTracks = farray->GetEntries();
889 if( !nTracks ) return 0;
890 TObjArray tracks(nTracks);
892 for (Int_t ipart=0; ipart<nTracks; ++ipart) {
894 AliAODMCParticle* part = (AliAODMCParticle*)farray->At(ipart);
896 if(!part->IsPhysicalPrimary()) continue;
899 Int_t pdg = TMath::Abs(part->GetPdgCode());
900 if((pdg == 12 || pdg == 14 || pdg == 16)) return kFALSE;
902 if( !part->Charge() ) continue; // comment / uncomment here
903 fhEtaMCTracks->Fill(part->Eta());
905 if(TMath::Abs(part->Eta()) > fTrackEtaCut) continue;
906 fhPtMCTracks->Fill(part->Pt());
908 if( part->Pt() < fTrackPtCutJF ) continue;
909 fhPhiMCTracks->Fill(part->Phi());
911 tracks.AddLast(part);
915 QSortTracks( tracks, 0, tracks.GetEntriesFast() );
917 nTracks = tracks.GetEntriesFast();
919 if( !nTracks ) return 0;
920 TObjArray *jets = new TObjArray(nTracks);
921 TIter itrack(&tracks);
923 // 2- Start with the highest pT particle ...
926 AliAODMCParticle* track = (AliAODMCParticle*)itrack.Next();
928 if( !track ) continue;
933 jets->AddLast( new AliAODJet(px,py,pz,pt) ); // Use the energy member to store Pt
934 AliAODJet* jet = (AliAODJet*)jets->Last();
935 jet->AddTrack(track);
936 tracks.Remove(track);
939 TVector2 jetVect2(jet->Px(),jet->Py());
941 // 3- Go to the next highest pT particle not already included...
942 AliAODMCParticle* track1;
943 while ( (track1 = (AliAODMCParticle*)(itrack.Next())) ) {
944 TVector2 trackVect2(track1->Px(),track1->Py());
945 Double_t dphi = trackVect2.DeltaPhi(jetVect2);
947 Double_t r = TMath::Sqrt( (jet->Eta()-track1->Eta())*(jet->Eta()-track1->Eta()) +
950 if( r < fConeRadius ) {
951 Double_t fPt = jet->E()+track1->Pt(); // Scalar sum of Pt
952 // recalculating the centroid
953 Double_t eta = jet->Eta()*jet->E()/fPt + track1->Eta()*track1->Pt()/fPt;
955 jetVect2.SetMagPhi(jet->E()/fPt,jetVect2.Phi());
956 trackVect2.SetMagPhi(track1->Pt()/fPt,trackVect2.Phi());
958 TVector2 sumVect2 = jetVect2+trackVect2;
959 Double_t phi = sumVect2.Phi();
961 //jet->SetPtEtaPhiE( 1., eta, phi, fPt );
962 ((TLorentzVector*) jet->MomentumVector())->SetPtEtaPhiE(1,eta,phi,fPt);
964 jet->AddTrack(track1);
965 tracks.Remove(track1);
971 nTracks = tracks.GetEntries();
972 // 4- Continue until all particles are in a jet.
974 } // end while nTracks
976 // Convert to AODjets....
977 Int_t njets = jets->GetEntriesFast();
979 TClonesArray* aodjets = new TClonesArray("AliAODJet",njets);
980 aodjets->SetOwner(kTRUE);
983 for(Int_t ijet=0; ijet<njets; ++ijet) {
984 AliAODJet* jet = (AliAODJet*)jets->At(ijet);
986 if (jet->E() < fJetPtCut) continue;
987 Float_t px, py,pz,en; // convert to 4-vector
988 px = jet->E() * TMath::Cos(jet->Phi()); // Pt * cos(phi)
989 py = jet->E() * TMath::Sin(jet->Phi()); // Pt * sin(phi)
990 pz = jet->E() / TMath::Tan(2.0 * TMath::ATan(TMath::Exp(-jet->Eta())));
991 en = TMath::Sqrt(px * px + py * py + pz * pz);
992 jet->SetPxPyPzE(px,py,pz,en);
995 TClonesArray &tmpaodjets = *aodjets;
996 new(tmpaodjets[count++]) AliAODJet(*jet);
997 //aodjets->AddLast( new AliAODJet(*jet));
998 //aodjets->AddLast( new AliAODJet(px, py, pz, en) );
1000 // Order jets according to their pT .
1001 QSortTracks( *aodjets, 0, aodjets->GetEntriesFast() );
1004 //if (fDebug>3) AliInfo(Form(" %d Charged jets found\n",njets));
1006 jets->Delete(); // OB - should I cleanup or leave it to garbage collection ?
1012 // ___________________________________________________________________
1014 Bool_t AliAnalysisTaskJetChem::IsTrackFromK0(const Int_t indexTrack){
1016 // check wether track with index is from K0 (= V0 with proper inv mass)
1018 TClonesArray* tracks = fAOD->GetTracks();
1020 for(int i=0; i<fAOD->GetNumberOfV0s(); i++){ // loop over V0s
1022 AliAODv0* v0 = fAOD->GetV0(i);
1024 Bool_t isOnFly = v0->GetOnFlyStatus();
1026 if( (fUseOnFlyV0s && !isOnFly) || (!fUseOnFlyV0s && isOnFly) ) continue;
1028 Double_t massK0 = v0->MassK0Short();
1030 // FIXME: here should go cuts (at least V0 quality)
1032 if(IsK0InvMass(massK0)){
1034 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0));
1035 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1));
1037 Int_t indexPos = tracks->IndexOf(trackPos);
1038 Int_t indexNeg = tracks->IndexOf(trackNeg);
1041 if(indexPos == indexTrack){ return kTRUE;}
1042 if(indexNeg == indexTrack){ return kTRUE;}
1051 //____________________________________________________________________
1052 TClonesArray* AliAnalysisTaskJetChem::FindChargedParticleJets()
1056 // loop over pt-ordered list of tracks, combine tracks within jet cone, recalc cone axis after each step
1057 // based on implementation by Arian Abrahantes Quintana and Enesto Lopez
1058 // ref: PHYSICAL REVIEW D 65 092002, CDF Collaboration
1061 // 1 - Order all charged particles according to their pT .
1062 Int_t nTracks = fAOD->GetNTracks();
1064 if( !nTracks ) return 0;
1065 TObjArray tracks(nTracks);
1067 // standard cuts + ITS refit = stdrd PWG4 cut (compose them for productions with old ESDFilter task)
1069 for (Int_t ipart=0; ipart<nTracks; ++ipart) {
1070 AliAODTrack* part = fAOD->GetTrack( ipart );
1072 UInt_t status = part->GetStatus();
1074 if( !part->TestFilterBit(fFilterBitJF) ) continue; // track cut selection
1075 if(fRequireITSRefitJF && ((status&AliESDtrack::kITSrefit)==0)) continue;
1077 if(!part->Charge() ) continue;
1079 if(TMath::Abs(part->Eta()) > fTrackEtaCut) continue;
1081 if(fRejectK0TracksJF && IsTrackFromK0(ipart)) continue;
1083 if( part->Pt() < fTrackPtCutJF ) continue;
1085 tracks.AddLast(part);
1088 QSortTracks( tracks, 0, tracks.GetEntriesFast() );
1090 nTracks = tracks.GetEntriesFast();
1092 if( !nTracks ) return 0;
1093 TObjArray *jets = new TObjArray(nTracks);
1094 TIter itrack(&tracks);
1096 // 2- Start with the highest pT particle ...
1097 Float_t px,py,pz,pt;
1098 AliAODTrack* track = (AliAODTrack*)itrack.Next();
1099 if( !track ) continue;
1104 //jets->AddLast( new TLorentzVector(px, py, pz, pt) ); // Use the energy member to store Pt
1105 jets->AddLast( new AliAODJet(px,py,pz,pt) ); // Use the energy member to store Pt
1106 //TLorentzVector* jet = (TLorentzVector*)jets->Last();
1107 AliAODJet* jet = (AliAODJet*)jets->Last();
1108 jet->AddTrack(track);
1109 tracks.Remove( track );
1111 TVector2 jetVect2(jet->Px(),jet->Py()); // OB
1113 // 3- Go to the next highest pT particle not already included...
1114 AliAODTrack* track1;
1115 while ( (track1 = (AliAODTrack*)(itrack.Next())) ) {
1116 //Double_t dphi = TVector2::Phi_mpi_pi(jet->Phi()-track1->Phi()); // OB remove
1117 TVector2 trackVect2(track1->Px(),track1->Py());
1118 Double_t dphi = trackVect2.DeltaPhi(jetVect2);
1120 Double_t r = TMath::Sqrt( (jet->Eta()-track1->Eta())*(jet->Eta()-track1->Eta()) +
1123 if( r < fConeRadius ) {
1124 Double_t fPt = jet->E()+track1->Pt(); // Scalar sum of Pt
1125 // recalculating the centroid
1126 Double_t eta = jet->Eta()*jet->E()/fPt + track1->Eta()*track1->Pt()/fPt;
1128 //Double_t phi = jet->Phi()*jet->E()/fPt + track1->Phi()*track1->Pt()/fPt; // OB - remove
1129 // OB - recalc phi via weighted 2vectors
1130 jetVect2.SetMagPhi(jet->E()/fPt,jetVect2.Phi());
1131 trackVect2.SetMagPhi(track1->Pt()/fPt,trackVect2.Phi());
1133 TVector2 sumVect2 = jetVect2+trackVect2;
1134 Double_t phi = sumVect2.Phi();
1136 //jet->SetPtEtaPhiE( 1., eta, phi, fPt );
1137 ((TLorentzVector*) jet->MomentumVector())->SetPtEtaPhiE(1,eta,phi,fPt);
1139 jet->AddTrack(track1);
1140 tracks.Remove(track1);
1146 nTracks = tracks.GetEntries();
1147 // 4- Continue until all particles are in a jet.
1149 } // end while nTracks
1151 // Convert to AODjets....
1152 Int_t njets = jets->GetEntriesFast();
1154 TClonesArray* aodjets = new TClonesArray("AliAODJet",njets);
1155 aodjets->SetOwner(kTRUE);
1159 for(Int_t ijet=0; ijet<njets; ++ijet) {
1160 //TLorentzVector* jet = (TLorentzVector*)jets->At(ijet);
1161 AliAODJet* jet = (AliAODJet*)jets->At(ijet);
1163 if (jet->E() < fJetPtCut) continue;
1164 Float_t px, py,pz,en; // convert to 4-vector
1165 px = jet->E() * TMath::Cos(jet->Phi()); // Pt * cos(phi)
1166 py = jet->E() * TMath::Sin(jet->Phi()); // Pt * sin(phi)
1167 pz = jet->E() / TMath::Tan(2.0 * TMath::ATan(TMath::Exp(-jet->Eta())));
1168 en = TMath::Sqrt(px * px + py * py + pz * pz);
1169 jet->SetPxPyPzE(px,py,pz,en);
1171 TClonesArray &tmpaodjets = *aodjets;
1172 new(tmpaodjets[count++]) AliAODJet(*jet);
1173 //aodjets->AddLast( new AliAODJet(*jet));
1174 //aodjets->AddLast( new AliAODJet(px, py, pz, en) );
1176 // Order jets according to their pT .
1177 //QSortTracks( *aodjets, 0, aodjets->GetEntriesFast() ); // not for TClonesArray
1179 if (fDebug>3) AliInfo(Form(" %d Charged jets found - after cuts %d \n",njets,count));
1181 jets->Delete(); // OB: cleanup
1188 //____________________________________________________________________
1190 TClonesArray* AliAnalysisTaskJetChem::GetPythiaJets(){
1192 // return Pythia jets in jet acc
1194 // note: present verion of AliMCEventHandler expects "AliAOD.root" (not "AliAODs.root"), otherwise galice.root will not be found in proper dir
1195 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
1197 Printf("ERROR: Could not retrieve MC event handler");
1201 AliMCEvent* mcEvent = mcHandler->MCEvent();
1203 Printf("ERROR: Could not retrieve MC event");
1208 AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
1209 if(!pythiaGenHeader){
1210 Printf("ERROR: Could not retrieve pythiaGenHeader");
1214 //Get Jets from MC header
1215 Int_t nPythiaGenJets = pythiaGenHeader->NTriggerJets();
1217 TClonesArray* aodjets = new TClonesArray("AliAODJet",nPythiaGenJets);
1218 aodjets->SetOwner(kTRUE);
1221 for(int ip = 0; ip<nPythiaGenJets; ++ip){
1223 pythiaGenHeader->TriggerJet(ip,p);
1224 TVector3 tempVect(p[0],p[1],p[2]);
1225 if ( TMath::Abs(tempVect.Eta())>fJetEtaCut ) continue;
1227 TClonesArray &tmpaodjets = *aodjets;
1228 new(tmpaodjets[count++]) AliAODJet(p[0],p[1],p[2],p[3]);
1229 //aodjets->AddLast(new AliAODJet(p[0],p[1],p[2],p[3]));
1236 //____________________________________________________________________
1237 void AliAnalysisTaskJetChem::QSortTracks(TObjArray &a, Int_t first, Int_t last)
1239 // Sort array of TObjArray of tracks by Pt using a quicksort algorithm.
1241 static TObject *tmp;
1242 static int i; // "static" to save stack space
1245 while (last - first > 1) {
1249 while (++i < last && ((AliVParticle*)a[i])->Pt() > ((AliVParticle*)a[first])->Pt() )
1251 while (--j > first && ((AliVParticle*)a[j])->Pt() < ((AliVParticle*)a[first])->Pt() )
1267 if (j - first < last - (j + 1)) {
1268 QSortTracks(a, first, j);
1269 first = j + 1; // QSortTracks(j + 1, last);
1271 QSortTracks(a, j + 1, last);
1272 last = j; // QSortTracks(first, j);
1277 //------------------------------------------------------------------
1279 TH1F* AliAnalysisTaskJetChem::CreatePIDhisto(const char* name){
1283 TH1F* result = new TH1F(name,"",60,0,60);
1284 result->SetOption("E");
1286 // bin equal Geant ID
1288 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(1),"photon");
1289 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(2),"e+");
1290 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(3),"e-");
1291 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(4),"e-neutrino");
1292 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(5),"mu+");
1293 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(6),"mu-");
1294 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(7),"pi0");
1295 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(8),"pi+");
1296 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(9),"pi-");
1297 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(10),"K long");
1298 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(11),"K+");
1299 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(12),"K-");
1300 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(13),"n");
1301 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(14),"p");
1302 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(15),"anti-proton");
1303 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(16),"K short");
1304 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(17),"eta");
1305 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(18),"Lambda");
1306 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(19),"Sigma+");
1307 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(20),"Sigma0");
1308 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(21),"Sigma-");
1309 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(22),"Xi0");
1310 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(23),"Xi-");
1311 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(24),"Omega-");
1312 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(25),"anti-neutron");
1313 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(26),"anti-Lambda");
1314 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(27),"anti-Sigma-");
1315 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(28),"anti-Sigma0");
1316 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(29),"anti-Sigma+");
1317 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(30),"anti-Xi0");
1318 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(31),"anti-Xi+");
1319 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(32),"anti-Omega+");
1320 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(33),"tau+");
1321 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(34),"tau-");
1322 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(35),"D+");
1323 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(36),"D-");
1324 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(37),"D0");
1325 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(38),"anti-D0");
1326 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(39),"Ds+");
1327 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(40),"anti Ds-");
1328 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(41),"Lamba_c+");
1329 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(42),"W+");
1330 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(43),"W-");
1331 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(44),"Z0");
1332 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(45),"d");
1333 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(46),"t");
1334 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(47),"alpha");
1335 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(48),"G_nu");
1336 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGpm311Bin),"K0/#bar{K0}");
1337 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDG333Bin),"phi");
1338 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGpm313Bin),"K*(892)0");
1339 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGp323Bin),"K*(892)+");
1340 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGm323Bin),"K*(892)-");
1341 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGNeutrinoBin),"nu");
1342 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGCharmedBaryonBin),"charmed baryon");
1343 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGQuarkBin),"q/#bar{q}");
1344 result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGDiQuarkBin),"q #bar{q}");
1345 result->GetXaxis()->LabelsOption("v"); // "u" ?
1350 //------------------------------------------------------------------
1352 TH1F* AliAnalysisTaskJetChem::CreatePythiaIDhisto(const char* name){
1356 TH1F* result = new TH1F(name,"",22,0,22);
1357 result->SetOption("E");
1359 result->GetXaxis()->SetBinLabel(kPythiaPIDP11Bin,"qq #rightarrow qq"); // ISUB = 11
1360 result->GetXaxis()->SetBinLabel(kPythiaPIDP12Bin,"q#bar{q} #rightarrow q#bar{q}"); // ISUB = 12
1361 result->GetXaxis()->SetBinLabel(kPythiaPIDP13Bin,"q#bar{q} #rightarrow gg"); // ISUB = 13
1362 result->GetXaxis()->SetBinLabel(kPythiaPIDP28Bin,"qg #rightarrow qg "); // ISUB = 28
1363 result->GetXaxis()->SetBinLabel(kPythiaPIDP53Bin,"gg #rightarrow q#bar{q}"); // ISUB = 53
1364 result->GetXaxis()->SetBinLabel(kPythiaPIDP68Bin,"gg #rightarrow gg"); // ISUB = 68
1365 result->GetXaxis()->SetBinLabel(kPythiaPIDP92Bin,"SD"); // ISUB = 92
1366 result->GetXaxis()->SetBinLabel(kPythiaPIDP93Bin,"SD"); // ISUB = 93
1367 result->GetXaxis()->SetBinLabel(kPythiaPIDP94Bin,"DD"); // ISUB = 94
1368 result->GetXaxis()->SetBinLabel(kPythiaPIDP95Bin,"low pt (MPI)"); // ISUB = 95
1369 result->GetXaxis()->SetBinLabel(kPythiaPIDPOtherBin,"other"); // ISUB = XX
1371 result->GetXaxis()->LabelsOption("v"); // "u" ?
1377 //------------------------------------------------------------------
1379 void AliAnalysisTaskJetChem::FillPythiaIDhisto(TH1F* h, const Int_t PID){
1381 // fill Pyhtia PID histogram
1385 if(PID == 11) bin = kPythiaPIDP11Bin;
1386 else if(PID == 12) bin = kPythiaPIDP12Bin;
1387 else if(PID == 13) bin = kPythiaPIDP13Bin;
1388 else if(PID == 28) bin = kPythiaPIDP28Bin;
1389 else if(PID == 53) bin = kPythiaPIDP53Bin;
1390 else if(PID == 68) bin = kPythiaPIDP68Bin;
1391 else if(PID == 92) bin = kPythiaPIDP92Bin;
1392 else if(PID == 93) bin = kPythiaPIDP93Bin;
1393 else if(PID == 94) bin = kPythiaPIDP94Bin;
1394 else if(PID == 95) bin = kPythiaPIDP95Bin;
1396 if(PID != -1) AliInfo(Form("unknown PID %d",PID));
1400 h->Fill(h->GetBinCenter(bin),1);
1404 //____________________________________________________________________
1407 void AliAnalysisTaskJetChem::FillPIDhisto(TH1F* hist, Int_t pdg, Float_t weight){
1409 // convert pdg code to Geant ID and fill corresponding bin
1411 Int_t fGID = fpdgdb->ConvertPdgToGeant3(pdg);
1413 //cout<<" pdg "<<pdg<<" fGID "<<fGID<<endl;
1415 if(TMath::Abs(pdg) == 311) fGID = kPDGpm311Bin;
1416 if(pdg == 333) fGID = kPDG333Bin;
1417 if(TMath::Abs(pdg) == 313) fGID = kPDGpm313Bin;
1418 if(pdg == 323) fGID = kPDGp323Bin;
1419 if(pdg == -323) fGID = kPDGm323Bin;
1421 if(TMath::Abs(pdg)==12 || TMath::Abs(pdg)==14 || TMath::Abs(pdg)==16) fGID = kPDGNeutrinoBin;
1423 if(TMath::Abs(pdg)==4122) fGID = kPDGCharmedBaryonBin;
1425 if(1<=TMath::Abs(pdg) && TMath::Abs(pdg)<=6) fGID = kPDGQuarkBin;
1427 if(TMath::Abs(pdg)==1103 || TMath::Abs(pdg)==2101 || TMath::Abs(pdg)==2103 || TMath::Abs(pdg)==2203 ||
1428 TMath::Abs(pdg)==3101 || TMath::Abs(pdg)==3103 || TMath::Abs(pdg)==3201 || TMath::Abs(pdg)==3203 ||
1429 TMath::Abs(pdg)==3303 || TMath::Abs(pdg)==4101 || TMath::Abs(pdg)==4103 || TMath::Abs(pdg)==4201 ||
1430 TMath::Abs(pdg)==4203 || TMath::Abs(pdg)==4301 || TMath::Abs(pdg)==4303 || TMath::Abs(pdg)==4403 ||
1431 TMath::Abs(pdg)==5101 || TMath::Abs(pdg)==5103 || TMath::Abs(pdg)==5201 || TMath::Abs(pdg)==5203 ||
1432 TMath::Abs(pdg)==5301 || TMath::Abs(pdg)==5303 || TMath::Abs(pdg)==5401 || TMath::Abs(pdg)==5403 ||
1433 TMath::Abs(pdg)==5503) fGID = kPDGDiQuarkBin;
1436 hist->Fill(fGID,weight);
1438 if(fGID>hist->GetBinLowEdge(hist->GetNbinsX()+1) ||
1439 fGID<hist->GetBinLowEdge(1)){
1441 AliError(Form("fGID %d for pdg %d exceeding histo limits ",fGID,pdg));
1445 AliInfo(Form("fGID 0 for pdg %d ",pdg));
1450 //____________________________________________________________________
1451 void AliAnalysisTaskJetChem::CreateHistos(){
1455 fListOfHistos = new TList();
1456 fListOfHistos->SetOwner(kTRUE);
1458 fhPrimVertexNCont = new TH1F("hPrimVertexNCont","",52,-2,50);
1459 fhPrimVertexNCont->SetXTitle("");
1460 fhPrimVertexNCont->Sumw2();
1461 fListOfHistos->Add( fhPrimVertexNCont );
1463 fhPrimVertexRho = new TH1F("hPrimVertexRho","",100,0,1);
1464 fhPrimVertexRho->SetXTitle("");
1465 fhPrimVertexRho->Sumw2();
1466 fListOfHistos->Add( fhPrimVertexRho );
1468 fhPrimVertexZ = new TH1F("hPrimVertexZ","",40,-20,20);
1469 fhPrimVertexZ->SetXTitle("");
1470 fhPrimVertexZ->Sumw2();
1471 fListOfHistos->Add( fhPrimVertexZ );
1473 fhNJets = new TH1F("hNJets","",10, 0, 10);
1474 fhNJets->SetXTitle("# of jets");
1476 fListOfHistos->Add( fhNJets );
1478 fhNJetsMC = new TH1F("hNJetsMC","",10,0,10);
1479 fhNJetsMC->SetXTitle("# of jets");
1481 fListOfHistos->Add( fhNJetsMC );
1483 fhLeadingEta = new TH1F("hLeadingEta","Leading Jet eta",12,-0.6,0.6);
1484 fhLeadingEta->SetXTitle("eta");
1485 fhLeadingEta->SetYTitle("dN/deta");
1486 fhLeadingEta->Sumw2();
1487 fListOfHistos->Add(fhLeadingEta);
1489 fhLeadingNTracksVsEta = new TH2F("hLeadingNTracksVsEta","",20,-1.0,1.0,20,0,20);
1490 fhLeadingNTracksVsEta->SetXTitle("eta");
1491 fhLeadingNTracksVsEta->SetYTitle("# of tracks");
1492 fhLeadingNTracksVsEta->Sumw2();
1493 fListOfHistos->Add(fhLeadingNTracksVsEta);
1495 fhLeadingPtVsEta = new TH2F("hLeadingPtVsEta","",20,-1.0,1.0,50,0,50);
1496 fhLeadingPtVsEta->SetXTitle("eta");
1497 fhLeadingPtVsEta->SetYTitle("# of tracks");
1498 fhLeadingPtVsEta->Sumw2();
1499 fListOfHistos->Add(fhLeadingPtVsEta);
1502 fhLeadingPhi = new TH1F("hLeadingPhi","Leading Jet phi",63,0,6.3);
1503 fhLeadingPhi->SetXTitle("phi");
1504 fhLeadingPhi->SetYTitle("dN/dphi");
1505 fhLeadingPhi->Sumw2();
1506 fListOfHistos->Add(fhLeadingPhi);
1509 fhLeadingPt = new TH1F("hLeadingPt","leading Jet p_{T}",50,0,50);
1510 fhLeadingPt->SetXTitle("p_{T} (GeV/c)");
1511 fhLeadingPt->SetYTitle("dN/dp_{T} (1/GeV)");
1512 fhLeadingPt->Sumw2();
1513 fListOfHistos->Add( fhLeadingPt );
1515 fhLeadingPtDiffr = new TH1F("hLeadingPtDiffr","leading Jet p_{T}",50,0,50);
1516 fhLeadingPtDiffr->SetXTitle("P_{T} (GeV/c)");
1517 fhLeadingPtDiffr->SetYTitle("dN/dp_{T} (1/GeV)");
1518 fhLeadingPtDiffr->Sumw2();
1519 fListOfHistos->Add( fhLeadingPtDiffr );
1521 fhLeadingEtaMC = new TH1F("hLeadingEtaMC","Leading Jet eta",12,-0.6,0.6);
1522 fhLeadingEtaMC->SetXTitle("eta");
1523 fhLeadingEtaMC->SetYTitle("dN/deta");
1524 fhLeadingEtaMC->Sumw2();
1525 fListOfHistos->Add(fhLeadingEtaMC);
1527 fhLeadingPhiMC = new TH1F("hLeadingPhiMC","Leading Jet phi",63,0,6.3);
1528 fhLeadingPhiMC->SetXTitle("phi");
1529 fhLeadingPhiMC->SetYTitle("dN/dphi");
1530 fhLeadingPhiMC->Sumw2();
1531 fListOfHistos->Add(fhLeadingPhiMC);
1533 fhLeadingPtMC = new TH1F("hLeadingPtMC","Leading Jet p_{T}",50,0,50);
1534 fhLeadingPtMC->SetXTitle("p_{T} (GeV/c)");
1535 fhLeadingPtMC->SetYTitle("dN/dp_{T} (1/GeV)");
1536 fhLeadingPtMC->Sumw2();
1537 fListOfHistos->Add( fhLeadingPtMC );
1539 fhLeadingPtMCDiffr = new TH1F("hLeadingPtMCDiffr","Leading Jet p_{T}",50,0,50);
1540 fhLeadingPtMCDiffr->SetXTitle("p_{T} (GeV/c)");
1541 fhLeadingPtMCDiffr->SetYTitle("dN/dp_{T} (1/GeV)");
1542 fhLeadingPtMCDiffr->Sumw2();
1543 fListOfHistos->Add( fhLeadingPtMCDiffr );
1545 fhPhiEtaTracksNoCut = new TH2F("hPhiEtaTracksNoCut","phi vs eta tracks",20,-1.0,1.0,63,0,6.3);
1546 fhPhiEtaTracksNoCut->SetXTitle("eta");
1547 fhPhiEtaTracksNoCut->SetYTitle("phi");
1548 fhPhiEtaTracksNoCut->Sumw2();
1549 fListOfHistos->Add(fhPhiEtaTracksNoCut);
1551 fhPtTracksNoCut = new TH1F("hPtTracksNoCut","p_{T} tracks",150,0,150);
1552 fhPtTracksNoCut->SetXTitle("p_{T} (GeV)");
1553 fhPtTracksNoCut->SetYTitle("dN/dp_{T} (1/GeV)");
1554 fhPtTracksNoCut->Sumw2();
1555 fListOfHistos->Add(fhPtTracksNoCut);
1557 fhPhiEtaTracks = new TH2F("hPhiEtaTracks","phi vs eta tracks",20,-1.0,1.0,63,0,6.3);
1558 fhPhiEtaTracks->SetXTitle("eta");
1559 fhPhiEtaTracks->SetYTitle("phi");
1560 fhPhiEtaTracks->Sumw2();
1561 fListOfHistos->Add(fhPhiEtaTracks);
1563 fhPtTracks = new TH1F("hPtTracks","p_{T} tracks",150,0,150);
1564 fhPtTracks->SetXTitle("P{T} (GeV)");
1565 fhPtTracks->SetYTitle("dN/dp_{T} (1/GeV)");
1566 fhPtTracks->Sumw2();
1567 fListOfHistos->Add(fhPtTracks);
1569 fhTrackMult = new TH1F("hTrackMult","",150,0,150);
1570 fhTrackMult->SetXTitle("n Tracks");
1571 fhTrackMult->SetYTitle("counts");
1572 fhTrackMult->Sumw2();
1573 fListOfHistos->Add(fhTrackMult);
1575 fhEtaMCTracks = new TH1F("hEtaMCTracks","eta tracks",30,-1.5,1.5);
1576 fhEtaMCTracks->SetXTitle("eta");
1577 fhEtaMCTracks->SetYTitle("dN/deta");
1578 fhEtaMCTracks->Sumw2();
1579 fListOfHistos->Add(fhEtaMCTracks);
1581 fhPhiMCTracks = new TH1F("hPhiMCTracks","phi tracks",63,0,6.3);
1582 fhPhiMCTracks->SetXTitle("phi");
1583 fhPhiMCTracks->SetYTitle("dN/dphi");
1584 fhPhiMCTracks->Sumw2();
1585 fListOfHistos->Add(fhPhiMCTracks);
1587 fhPtMCTracks = new TH1F("hPtMCTracks","p_{T} tracks",50,0,50);
1588 fhPtMCTracks->SetXTitle("p_{T} (GeV)");
1589 fhPtMCTracks->SetYTitle("dN/dp_{T} (1/GeV)");
1590 fhPtMCTracks->Sumw2();
1591 fListOfHistos->Add(fhPtMCTracks);
1593 fhnTracksVsPtLeading = new TH2F("hnTracksVsPtLeading","",50,0.,50.,20,-0.5,19.5);
1594 fhnTracksVsPtLeading->SetXTitle("p_{T} (GeV/c)");
1595 fhnTracksVsPtLeading->SetYTitle("n tracks");
1596 fhnTracksVsPtLeading->Sumw2();
1597 fListOfHistos->Add( fhnTracksVsPtLeading );
1599 fhdNdEtaPhiDist = new TH1F("hdNdEtaPhiDist","Charge particle density |#eta|< 1 vs #Delta#phi", 120, 0., 2.*TMath::Pi());
1600 fhdNdEtaPhiDist->SetXTitle("#Delta#phi");
1601 fhdNdEtaPhiDist->SetYTitle("dN{ch}/d#etad#phi");
1602 fhdNdEtaPhiDist->Sumw2();
1603 fListOfHistos->Add( fhdNdEtaPhiDist );
1605 fhRegionSumPtMaxVsEt = new TH1F("hRegionSumPtMaxVsEt","P{T}^{90, max} vs Leading Jet P{T}",50,0.,50.);
1606 fhRegionSumPtMaxVsEt->SetXTitle("P{T} (GeV/c)");
1607 fhRegionSumPtMaxVsEt->Sumw2();
1608 fListOfHistos->Add( fhRegionSumPtMaxVsEt );
1610 fhRegionMultMaxVsEt = new TH1F("hRegionMultMaxVsEt","N{ch}^{90, max} vs Leading Jet P{T}",50,0.,50.);
1611 fhRegionMultMaxVsEt->SetXTitle("E (GeV hRegionAveSumPtVsEt/c)");
1612 fhRegionMultMaxVsEt->Sumw2();
1613 fListOfHistos->Add( fhRegionMultMaxVsEt );
1615 fhRegionSumPtMinVsEt = new TH1F("hRegionSumPtMinVsEt","P{T}^{90, min} vs Leading Jet P{T}",50,0.,50.);
1616 fhRegionSumPtMinVsEt->SetXTitle("P{T} (GeV/c)");
1617 fhRegionSumPtMinVsEt->Sumw2();
1618 fListOfHistos->Add( fhRegionSumPtMinVsEt );
1620 fhRegionMultMinVsEt = new TH1F("hRegionMultMinVsEt","N{ch}^{90, min} vs Leading Jet P{T}",50,0.,50.);
1621 fhRegionMultMinVsEt->SetXTitle("E (GeV/c)");
1622 fhRegionMultMinVsEt->Sumw2();
1623 fListOfHistos->Add( fhRegionMultMinVsEt );
1627 fhNV0s = new TH1F("hNV0s","n V0s",50,0,50);
1628 fhNV0s->SetXTitle("n V0s");
1630 fListOfHistos->Add(fhNV0s);
1632 fhV0onFly = new TH1F("hV0onFly","on-the-fly V0",5,0,5);
1633 fhV0onFly->SetXTitle("is on-the-fly V0");
1635 fListOfHistos->Add(fhV0onFly);
1637 fhV0DCADaughters = new TH1F("hV0DCADaughters","V0 DCA daughters",200,0,2.0);
1638 fhV0DCADaughters->SetXTitle("V0 DCA daughters");
1639 fhV0DCADaughters->Sumw2();
1640 fListOfHistos->Add(fhV0DCADaughters);
1642 fhV0Radius = new TH1F("hV0Radius","V0 radius",2500,0,250);
1643 fhV0Radius->SetXTitle("V0 radius");
1644 fhV0Radius->Sumw2();
1645 fListOfHistos->Add(fhV0Radius);
1647 fhV0DCAToVertex = new TH1F("hV0DCAToVertex","",100,0,10);
1648 fhV0DCAToVertex->SetXTitle("V0 DCA (cm)");
1649 fhV0DCAToVertex->Sumw2();
1650 fListOfHistos->Add(fhV0DCAToVertex);
1652 fhV0DCAToVertexK0 = new TH1F("hV0DCAToVertexK0","",100,0,10);
1653 fhV0DCAToVertexK0->SetXTitle("V0 DCA (cm)");
1654 fhV0DCAToVertexK0->Sumw2();
1655 fListOfHistos->Add(fhV0DCAToVertexK0);
1657 fhV0InvMassK0 = new TH1F("hV0InvMassK0","",2000,0,2);
1658 fhV0InvMassK0->SetXTitle("inv mass (GeV)");
1659 fhV0InvMassK0->Sumw2();
1660 fListOfHistos->Add(fhV0InvMassK0);
1662 fhV0PtVsInvMassK0 = new TH2F("hV0PtVsInvMassK0","",100,0.45,0.55,100,0.,50.);
1663 fhV0PtVsInvMassK0->SetXTitle("inv mass (GeV)");
1664 fhV0PtVsInvMassK0->SetYTitle("p_{T} (GeV)");
1665 fhV0PtVsInvMassK0->Sumw2();
1666 fListOfHistos->Add(fhV0PtVsInvMassK0);
1668 fhV0InvMassK0JetEvt = new TH1F("hV0InvMassK0JetEvt","",2000,0,2);
1669 fhV0InvMassK0JetEvt->SetXTitle("inv mass (GeV)");
1670 fhV0InvMassK0JetEvt->Sumw2();
1671 fListOfHistos->Add(fhV0InvMassK0JetEvt);
1673 fhV0InvMassLambda = new TH1F("hV0InvMassLambda","",2000,0,2);
1674 fhV0InvMassLambda->SetXTitle("inv mass (GeV)");
1675 fhV0InvMassLambda->Sumw2();
1676 fListOfHistos->Add(fhV0InvMassLambda);
1678 fhV0InvMassAntiLambda = new TH1F("hV0InvMassAntiLambda","",2000,0,2);
1679 fhV0InvMassAntiLambda->SetXTitle("inv mass (GeV)");
1680 fhV0InvMassAntiLambda->Sumw2();
1681 fListOfHistos->Add(fhV0InvMassAntiLambda);
1683 fhV0InvMassLambdaJetEvt = new TH1F("hV0InvMassLambdaJetEvt","",2000,0,2);
1684 fhV0InvMassLambdaJetEvt->SetXTitle("inv mass (GeV)");
1685 fhV0InvMassLambdaJetEvt->Sumw2();
1686 fListOfHistos->Add(fhV0InvMassLambdaJetEvt);
1688 fhV0InvMassAntiLambdaJetEvt = new TH1F("hV0InvMassAntiLambdaJetEvt","",2000,0,2);
1689 fhV0InvMassAntiLambdaJetEvt->SetXTitle("inv mass (GeV)");
1690 fhV0InvMassAntiLambdaJetEvt->Sumw2();
1691 fListOfHistos->Add(fhV0InvMassAntiLambdaJetEvt);
1693 fhdROpanK0VsPt = new TH2F("hdROpanK0VsPt","V0 dR vs pt",100,0,10,100,0,10);
1694 fhdROpanK0VsPt->SetXTitle("opening angle R (rad)");
1695 fhdROpanK0VsPt->Sumw2();
1696 fListOfHistos->Add(fhdROpanK0VsPt);
1698 fhdPhiJetV0 = new TH1F("hdPhiJetV0","",640,-3.2,3.2);
1699 fhdPhiJetV0->SetXTitle("#Delta #phi V0 - jet");
1700 fhdPhiJetV0->Sumw2();
1701 fListOfHistos->Add(fhdPhiJetV0);
1703 fhdPhiJetK0 = new TH1F("hdPhiJetK0","",640,-3.2,3.2);
1704 fhdPhiJetK0->SetXTitle("#Delta #phi V0 - jet");
1705 fhdPhiJetK0->Sumw2();
1706 fListOfHistos->Add(fhdPhiJetK0);
1708 fhdRJetK0 = new TH1F("hdRJetK0","dN/dR K0-jet",500,0,5);
1709 fhdRJetK0->SetXTitle("#Delta R K0 - jet");
1710 fhdRJetK0->SetYTitle("1/N_{jet} dN/dR");
1712 fListOfHistos->Add(fhdRJetK0);
1714 fhdNdptV0 = new TH1F("hdNdptV0","dN/dpt V0",100,0,10);
1715 fhdNdptV0->SetXTitle("p_{T} (GeV/c)");
1716 fhdNdptV0->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
1718 fListOfHistos->Add(fhdNdptV0);
1720 fhdNdptK0 = new TH1F("hdNdptK0","",100,0,10);
1721 fhdNdptK0->SetXTitle("p_{T} (GeV/c)");
1722 fhdNdptK0->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
1724 fListOfHistos->Add(fhdNdptK0);
1726 fhPtVsEtaK0 = new TH2F("hPtVsEtaK0","",20,-1,1,100,0,10);
1727 fhPtVsEtaK0->SetXTitle("#eta");
1728 fhPtVsEtaK0->SetYTitle("p_{T} (GeV/c)");
1729 fhPtVsEtaK0->Sumw2();
1730 fListOfHistos->Add(fhPtVsEtaK0);
1732 fhV0InvMassK0DCA = new TH1F("hV0InvMassK0DCA","",2000,0,2);
1733 fhV0InvMassK0DCA->SetXTitle("inv mass (GeV)");
1734 fhV0InvMassK0DCA->Sumw2();
1735 fListOfHistos->Add(fhV0InvMassK0DCA);
1737 fhV0InvMassK0DCAdEdx = new TH1F("hV0InvMassK0DCAdEdx","",2000,0,2);
1738 fhV0InvMassK0DCAdEdx->SetXTitle("inv mass (GeV)");
1739 fhV0InvMassK0DCAdEdx->Sumw2();
1740 fListOfHistos->Add(fhV0InvMassK0DCAdEdx);
1742 fhV0PtVsInvMassK0DCAdEdx = new TH2F("hV0PtVsInvMassK0DCAdEdx","",100,0.45,0.55,100,0.,50.);
1743 fhV0PtVsInvMassK0DCAdEdx->SetXTitle("inv mass (GeV)");
1744 fhV0PtVsInvMassK0DCAdEdx->SetYTitle("p_{T} (GeV)");
1745 fhV0PtVsInvMassK0DCAdEdx->Sumw2();
1746 fListOfHistos->Add(fhV0PtVsInvMassK0DCAdEdx);
1748 fhV0InvMassK0DCAPID = new TH1F("hV0InvMassK0DCAPID","",2000,0,2);
1749 fhV0InvMassK0DCAPID->SetXTitle("inv mass (GeV)");
1750 fhV0InvMassK0DCAPID->Sumw2();
1751 fListOfHistos->Add(fhV0InvMassK0DCAPID);
1753 fhV0InvMassLambdaDCAdEdx = new TH1F("hV0InvMassLambdaDCAdEdx","",2000,0,2);
1754 fhV0InvMassLambdaDCAdEdx->SetXTitle("inv mass (GeV)");
1755 fhV0InvMassLambdaDCAdEdx->Sumw2();
1756 fListOfHistos->Add(fhV0InvMassLambdaDCAdEdx);
1758 fhV0InvMassAntiLambdaDCAdEdx = new TH1F("hV0InvMassAntiLambdaDCAdEdx","",2000,0,2);
1759 fhV0InvMassAntiLambdaDCAdEdx->SetXTitle("inv mass (GeV)");
1760 fhV0InvMassAntiLambdaDCAdEdx->Sumw2();
1761 fListOfHistos->Add(fhV0InvMassAntiLambdaDCAdEdx);
1763 fhdNdptK0DCA = new TH1F("hdNdptK0DCA","",100,0,10);
1764 fhdNdptK0DCA->SetXTitle("p_{T} (GeV/c)");
1765 fhdNdptK0DCA->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
1766 fhdNdptK0DCA->Sumw2();
1767 fListOfHistos->Add(fhdNdptK0DCA);
1769 fhdNdptK0DCAdEdx = new TH1F("hdNdptK0DCAdEdx","",100,0,10);
1770 fhdNdptK0DCAdEdx->SetXTitle("p_{T} (GeV/c)");
1771 fhdNdptK0DCAdEdx->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
1772 fhdNdptK0DCAdEdx->Sumw2();
1773 fListOfHistos->Add(fhdNdptK0DCAdEdx);
1775 fhV0InvMassK0Min = new TH1F("hV0InvMassK0Min","",2000,0,2);
1776 fhV0InvMassK0Min->SetXTitle("inv mass (GeV)");
1777 fhV0InvMassK0Min->Sumw2();
1778 fListOfHistos->Add(fhV0InvMassK0Min);
1780 fhV0InvMassLambdaMin = new TH1F("hV0InvMassLambdaMin","",2000,0,2);
1781 fhV0InvMassLambdaMin->SetXTitle("inv mass (GeV)");
1782 fhV0InvMassLambdaMin->Sumw2();
1783 fListOfHistos->Add(fhV0InvMassLambdaMin);
1785 fhV0InvMassAntiLambdaMin = new TH1F("hV0InvMassAntiLambdaMin","",2000,0,2);
1786 fhV0InvMassAntiLambdaMin->SetXTitle("inv mass (GeV)");
1787 fhV0InvMassAntiLambdaMin->Sumw2();
1788 fListOfHistos->Add(fhV0InvMassAntiLambdaMin);
1790 fhV0InvMassK0Max = new TH1F("hV0InvMassK0Max","",2000,0,2);
1791 fhV0InvMassK0Max->SetXTitle("inv mass (GeV)");
1792 fhV0InvMassK0Max->Sumw2();
1793 fListOfHistos->Add(fhV0InvMassK0Max);
1795 fhV0InvMassLambdaMax = new TH1F("hV0InvMassLambdaMax","",2000,0,2);
1796 fhV0InvMassLambdaMax->SetXTitle("inv mass (GeV)");
1797 fhV0InvMassLambdaMax->Sumw2();
1798 fListOfHistos->Add(fhV0InvMassLambdaMax);
1800 fhV0InvMassAntiLambdaMax = new TH1F("hV0InvMassAntiLambdaMax","",2000,0,2);
1801 fhV0InvMassAntiLambdaMax->SetXTitle("inv mass (GeV)");
1802 fhV0InvMassAntiLambdaMax->Sumw2();
1803 fListOfHistos->Add(fhV0InvMassAntiLambdaMax);
1805 fhV0InvMassK0Jet = new TH1F("hV0InvMassK0Jet","",2000,0,2);
1806 fhV0InvMassK0Jet->SetXTitle("inv mass (GeV)");
1807 fhV0InvMassK0Jet->Sumw2();
1808 fListOfHistos->Add(fhV0InvMassK0Jet);
1810 fhV0InvMassLambdaJet = new TH1F("hV0InvMassLambdaJet","",2000,0,2);
1811 fhV0InvMassLambdaJet->SetXTitle("inv mass (GeV)");
1812 fhV0InvMassLambdaJet->Sumw2();
1813 fListOfHistos->Add(fhV0InvMassLambdaJet);
1815 fhV0InvMassAntiLambdaJet = new TH1F("hV0InvMassAntiLambdaJet","",2000,0,2);
1816 fhV0InvMassAntiLambdaJet->SetXTitle("inv mass (GeV)");
1817 fhV0InvMassAntiLambdaJet->Sumw2();
1818 fListOfHistos->Add(fhV0InvMassAntiLambdaJet);
1820 fhV0InvMassK0Lambda = new TH1F("hV0InvMassK0Lambda","",2000,0,2);
1821 fhV0InvMassK0Lambda->SetXTitle("inv mass (GeV)");
1822 fhV0InvMassK0Lambda->Sumw2();
1823 fListOfHistos->Add(fhV0InvMassK0Lambda);
1825 fhdNdptK0JetEvt = new TH1F("hdNdptK0JetEvt","",100,0,10);
1826 fhdNdptK0JetEvt->SetXTitle("p_{T} (GeV/c)");
1827 fhdNdptK0JetEvt->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
1828 fhdNdptK0JetEvt->Sumw2();
1829 fListOfHistos->Add(fhdNdptK0JetEvt);
1831 fhdNdzK0 = new TH1F("hdNdzK0","",150,0,1.5);
1832 fhdNdzK0->SetXTitle("z");
1833 fhdNdzK0->SetYTitle("1/N_{jet} dN/dz");
1835 fListOfHistos->Add(fhdNdzK0);
1837 fhdNdzK05to10 = new TH1F("hdNdzK05to10","",150,0,1.5);
1838 fhdNdzK05to10->SetXTitle("z");
1839 fhdNdzK05to10->SetYTitle("1/N_{jet} dN/dz");
1840 fhdNdzK05to10->Sumw2();
1841 fListOfHistos->Add(fhdNdzK05to10);
1843 fhdNdzK010to20 = new TH1F("hdNdzK010to20","",150,0,1.5);
1844 fhdNdzK010to20->SetXTitle("z");
1845 fhdNdzK010to20->SetYTitle("1/N_{jet} dN/dz");
1846 fhdNdzK010to20->Sumw2();
1847 fListOfHistos->Add(fhdNdzK010to20);
1849 fhdNdzK020to30 = new TH1F("hdNdzK020to30","",150,0,1.5);
1850 fhdNdzK020to30->SetXTitle("z");
1851 fhdNdzK020to30->SetYTitle("1/N_{jet} dN/dz");
1852 fhdNdzK020to30->Sumw2();
1853 fListOfHistos->Add(fhdNdzK020to30);
1855 fhdNdzK030to40 = new TH1F("hdNdzK030to40","",150,0,1.5);
1856 fhdNdzK030to40->SetXTitle("z");
1857 fhdNdzK030to40->SetYTitle("1/N_{jet} dN/dz");
1858 fhdNdzK030to40->Sumw2();
1859 fListOfHistos->Add(fhdNdzK030to40);
1861 fhdNdzK040to60 = new TH1F("hdNdzK040to60","",150,0,1.5);
1862 fhdNdzK040to60->SetXTitle("z");
1863 fhdNdzK040to60->SetYTitle("1/N_{jet} dN/dz");
1864 fhdNdzK040to60->Sumw2();
1865 fListOfHistos->Add(fhdNdzK040to60);
1867 fhdNdxiK0 = new TH1F("hdNdxiK0","",100,0,10);
1868 fhdNdxiK0->SetXTitle("xi");
1869 fhdNdxiK0->SetYTitle("1/N_{jet} dN/dxi");
1871 fListOfHistos->Add(fhdNdxiK0);
1873 fhdNdzLambda = new TH1F("hdNdzLambda","",150,0,1.5);
1874 fhdNdzLambda->SetXTitle("z");
1875 fhdNdzLambda->SetYTitle("1/N_{jet} dN/dz");
1876 fhdNdzLambda->Sumw2();
1877 fListOfHistos->Add(fhdNdzLambda);
1879 fhdNdzAntiLambda = new TH1F("hdNdzAntiLambda","",150,0,1.5);
1880 fhdNdzAntiLambda->SetXTitle("z");
1881 fhdNdzAntiLambda->SetYTitle("1/N_{jet} dN/dz");
1882 fhdNdzAntiLambda->Sumw2();
1883 fListOfHistos->Add(fhdNdzAntiLambda);
1885 fhdNdzK0Max = new TH1F("hdNdzK0Max","",150,0,1.5);
1886 fhdNdzK0Max->SetXTitle("z");
1887 fhdNdzK0Max->SetYTitle("1/N_{jet} dN/dz");
1888 fhdNdzK0Max->Sumw2();
1889 fListOfHistos->Add(fhdNdzK0Max);
1891 fhdNdxiK0Max = new TH1F("hdNdxiK0Max","",100,0,10);
1892 fhdNdxiK0Max->SetXTitle("xi");
1893 fhdNdxiK0Max->SetYTitle("1/N_{jet} dN/dxi");
1894 fhdNdxiK0Max->Sumw2();
1895 fListOfHistos->Add(fhdNdxiK0Max);
1897 fhdNdzLambdaMax = new TH1F("hdNdzLambdaMax","",150,0,1.5);
1898 fhdNdzLambdaMax->SetXTitle("z");
1899 fhdNdzLambdaMax->SetYTitle("1/N_{jet} dN/dz");
1900 fhdNdzLambdaMax->Sumw2();
1901 fListOfHistos->Add(fhdNdzLambdaMax);
1903 fhdNdxiLambdaMax = new TH1F("hdNdxiLambdaMax","",700,0,7);
1904 fhdNdxiLambdaMax->SetXTitle("xi");
1905 fhdNdxiLambdaMax->SetYTitle("1/N_{jet} dN/dxi");
1906 fhdNdxiLambdaMax->Sumw2();
1907 fListOfHistos->Add(fhdNdxiLambdaMax);
1909 fhdNdptK0Max = new TH1F("hdNdptK0Max","",100,0,10);
1910 fhdNdptK0Max->SetXTitle("p_{T} (GeV/c)");
1911 fhdNdptK0Max->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
1912 fhdNdptK0Max->Sumw2();
1913 fListOfHistos->Add(fhdNdptK0Max);
1915 fhdNdptLambdaMax = new TH1F("hdNdptLambdaMax","",100,0,10);
1916 fhdNdptLambdaMax->SetXTitle("p_{T} (GeV/c)");
1917 fhdNdptLambdaMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
1918 fhdNdptLambdaMax->Sumw2();
1919 fListOfHistos->Add(fhdNdptLambdaMax);
1921 fhdNdzK0Min = new TH1F("hdNdzK0Min","",150,0,1.5);
1922 fhdNdzK0Min->SetXTitle("z");
1923 fhdNdzK0Min->SetYTitle("1/N_{jet} dN/dz");
1924 fhdNdzK0Min->Sumw2();
1925 fListOfHistos->Add(fhdNdzK0Min);
1927 fhdNdxiK0Min = new TH1F("hdNdxiK0Min","",100,0,10);
1928 fhdNdxiK0Min->SetXTitle("xi");
1929 fhdNdxiK0Min->SetYTitle("1/N_{jet} dN/dxi");
1930 fhdNdxiK0Min->Sumw2();
1931 fListOfHistos->Add(fhdNdxiK0Min);
1933 fhdNdzLambdaMin = new TH1F("hdNdzLambdaMin","",150,0,1.5);
1934 fhdNdzLambdaMin->SetXTitle("z");
1935 fhdNdzLambdaMin->SetYTitle("1/N_{jet} dN/dz");
1936 fhdNdzLambdaMin->Sumw2();
1937 fListOfHistos->Add(fhdNdzLambdaMin);
1939 fhdNdxiLambdaMin = new TH1F("hdNdxiLambdaMin","",700,0,7);
1940 fhdNdxiLambdaMin->SetXTitle("xi");
1941 fhdNdxiLambdaMin->SetYTitle("1/N_{jet} dN/dxi");
1942 fhdNdxiLambdaMin->Sumw2();
1943 fListOfHistos->Add(fhdNdxiLambdaMin);
1945 fhdNdptK0Min = new TH1F("hdNdptK0Min","",100,0,10);
1946 fhdNdptK0Min->SetXTitle("p_{T} (GeV/c)");
1947 fhdNdptK0Min->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
1948 fhdNdptK0Min->Sumw2();
1949 fListOfHistos->Add(fhdNdptK0Min);
1951 fhdNdptLambdaMin = new TH1F("hdNdptLambdaMin","",100,0,10);
1952 fhdNdptLambdaMin->SetXTitle("p_{T} (GeV/c)");
1953 fhdNdptLambdaMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
1954 fhdNdptLambdaMin->Sumw2();
1955 fListOfHistos->Add(fhdNdptLambdaMin);
1957 fhdNdzK0Jet = new TH1F("hdNdzK0Jet","",150,0,1.5);
1958 fhdNdzK0Jet->SetXTitle("z");
1959 fhdNdzK0Jet->SetYTitle("1/N_{jet} dN/dz");
1960 fhdNdzK0Jet->Sumw2();
1961 fListOfHistos->Add(fhdNdzK0Jet);
1963 fhdNdxiK0Jet = new TH1F("hdNdxiK0Jet","",100,0,10);
1964 fhdNdxiK0Jet->SetXTitle("xi");
1965 fhdNdxiK0Jet->SetYTitle("1/N_{jet} dN/dxi");
1966 fhdNdxiK0Jet->Sumw2();
1967 fListOfHistos->Add(fhdNdxiK0Jet);
1969 fhdNdzLambdaJet = new TH1F("hdNdzLambdaJet","",150,0,1.5);
1970 fhdNdzLambdaJet->SetXTitle("z");
1971 fhdNdzLambdaJet->SetYTitle("1/N_{jet} dN/dz");
1972 fhdNdzLambdaJet->Sumw2();
1973 fListOfHistos->Add(fhdNdzLambdaJet);
1975 fhdNdxiLambdaJet = new TH1F("hdNdxiLambdaJet","",700,0,7);
1976 fhdNdxiLambdaJet->SetXTitle("xi");
1977 fhdNdxiLambdaJet->SetYTitle("1/N_{jet} dN/dxi");
1978 fhdNdxiLambdaJet->Sumw2();
1979 fListOfHistos->Add(fhdNdxiLambdaJet);
1981 fhdNdptK0Jet = new TH1F("hdNdptK0Jet","",100,0,10);
1982 fhdNdptK0Jet->SetXTitle("p_{T} (GeV/c)");
1983 fhdNdptK0Jet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
1984 fhdNdptK0Jet->Sumw2();
1985 fListOfHistos->Add(fhdNdptK0Jet);
1987 fhdNdptLambdaJet = new TH1F("hdNdptLambdaJet","",100,0,10);
1988 fhdNdptLambdaJet->SetXTitle("p_{T} (GeV/c)");
1989 fhdNdptLambdaJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
1990 fhdNdptLambdaJet->Sumw2();
1991 fListOfHistos->Add(fhdNdptLambdaJet);
1993 fhdEdxVsMomV0 = new TH2F("hdEdxVsMomV0","",200,-10,10,200,0,2000);
1994 fhdEdxVsMomV0->SetXTitle("mom (GeV/c)");
1995 fhdEdxVsMomV0->SetYTitle("dE/dx");
1996 fhdEdxVsMomV0->Sumw2();
1997 fListOfHistos->Add(fhdEdxVsMomV0);
1999 fhdEdxVsMomV0pidEdx = new TH2F("hdEdxVsMomV0pidEdx","",200,-10,10,200,0,2000);
2000 fhdEdxVsMomV0pidEdx->SetXTitle("mom (GeV/c)");
2001 fhdEdxVsMomV0pidEdx->SetYTitle("dE/dx");
2002 fhdEdxVsMomV0pidEdx->Sumw2();
2003 fListOfHistos->Add(fhdEdxVsMomV0pidEdx);
2005 fhdEdxVsMomV0piPID = new TH2F("hdEdxVsMomV0piPID","",200,-10,10,200,0,2000);
2006 fhdEdxVsMomV0piPID->SetXTitle("mom (GeV/c)");
2007 fhdEdxVsMomV0piPID->SetYTitle("dE/dx");
2008 fhdEdxVsMomV0piPID->Sumw2();
2009 fListOfHistos->Add(fhdEdxVsMomV0piPID);
2011 fhdPhiJetK0MC = new TH1F("hdPhiJetK0MC","",640,-3.2,3.2);
2012 fhdPhiJetK0MC->SetXTitle("#Delta #phi K0 - jet");
2013 fhdPhiJetK0MC->SetYTitle("1/N_{jet} dN/dphi");
2014 fhdPhiJetK0MC->Sumw2();
2015 fListOfHistos->Add(fhdPhiJetK0MC);
2017 fhdRJetK0MC = new TH1F("hdRJetK0MC","dN/R K0-jet",500,0,5);
2018 fhdRJetK0MC->SetXTitle("#Delta R K0 - jet");
2019 fhdRJetK0MC->SetYTitle("1/N_{jet} dN/dR");
2020 fhdRJetK0MC->Sumw2();
2021 fListOfHistos->Add(fhdRJetK0MC);
2023 fhdRV0MC = new TH1F("hdRV0MC","",500,0.,1.);
2024 fhdRV0MC->SetXTitle("#Delta R");
2025 fhdRV0MC->SetYTitle("");
2027 fListOfHistos->Add(fhdRV0MC);
2029 fhdNdptchPiMCMax = new TH1F("hdNdptchPiMCMax","",100,0,10);
2030 fhdNdptchPiMCMax->SetXTitle("p_{T} (GeV/c)");
2031 fhdNdptchPiMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2032 fhdNdptchPiMCMax->Sumw2();
2033 fListOfHistos->Add(fhdNdptchPiMCMax);
2035 fhdNdptK0MCMax = new TH1F("hdNdptK0MCMax","",100,0,10);
2036 fhdNdptK0MCMax->SetXTitle("p_{T} (GeV/c)");
2037 fhdNdptK0MCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2038 fhdNdptK0MCMax->Sumw2();
2039 fListOfHistos->Add(fhdNdptK0MCMax);
2041 fhdNdptchKMCMax = new TH1F("hdNdptchKMCMax","",100,0,10);
2042 fhdNdptchKMCMax->SetXTitle("p_{T} (GeV/c)");
2043 fhdNdptchKMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2044 fhdNdptchKMCMax->Sumw2();
2045 fListOfHistos->Add(fhdNdptchKMCMax);
2047 fhdNdptpMCMax = new TH1F("hdNdptpMCMax","",100,0,10);
2048 fhdNdptpMCMax->SetXTitle("p_{T} (GeV/c)");
2049 fhdNdptpMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2050 fhdNdptpMCMax->Sumw2();
2051 fListOfHistos->Add(fhdNdptpMCMax);
2053 fhdNdptpBarMCMax = new TH1F("hdNdptpBarMCMax","",100,0,10);
2054 fhdNdptpBarMCMax->SetXTitle("p_{T} (GeV/c)");
2055 fhdNdptpBarMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2056 fhdNdptpBarMCMax->Sumw2();
2057 fListOfHistos->Add(fhdNdptpBarMCMax);
2059 fhdNdptLambdaMCMax = new TH1F("hdNdptLambdaMCMax","",100,0,10);
2060 fhdNdptLambdaMCMax->SetXTitle("p_{T} (GeV/c)");
2061 fhdNdptLambdaMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2062 fhdNdptLambdaMCMax->Sumw2();
2063 fListOfHistos->Add(fhdNdptLambdaMCMax);
2065 fhdNdptLambdaBarMCMax = new TH1F("hdNdptLambdaBarMCMax","",100,0,10);
2066 fhdNdptLambdaBarMCMax->SetXTitle("p_{T} (GeV/c)");
2067 fhdNdptLambdaBarMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2068 fhdNdptLambdaBarMCMax->Sumw2();
2069 fListOfHistos->Add(fhdNdptLambdaBarMCMax);
2072 fhdNdptchPiMCMin = new TH1F("hdNdptchPiMCMin","",100,0,10);
2073 fhdNdptchPiMCMin->SetXTitle("p_{T} (GeV/c)");
2074 fhdNdptchPiMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2075 fhdNdptchPiMCMin->Sumw2();
2076 fListOfHistos->Add(fhdNdptchPiMCMin);
2078 fhdNdptK0MCMin = new TH1F("hdNdptK0MCMin","",100,0,10);
2079 fhdNdptK0MCMin->SetXTitle("p_{T} (GeV/c)");
2080 fhdNdptK0MCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2081 fhdNdptK0MCMin->Sumw2();
2082 fListOfHistos->Add(fhdNdptK0MCMin);
2084 fhdNdptchKMCMin = new TH1F("hdNdptchKMCMin","",100,0,10);
2085 fhdNdptchKMCMin->SetXTitle("p_{T} (GeV/c)");
2086 fhdNdptchKMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2087 fhdNdptchKMCMin->Sumw2();
2088 fListOfHistos->Add(fhdNdptchKMCMin);
2090 fhdNdptpMCMin = new TH1F("hdNdptpMCMin","",100,0,10);
2091 fhdNdptpMCMin->SetXTitle("p_{T} (GeV/c)");
2092 fhdNdptpMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2093 fhdNdptpMCMin->Sumw2();
2094 fListOfHistos->Add(fhdNdptpMCMin);
2096 fhdNdptpBarMCMin = new TH1F("hdNdptpBarMCMin","",100,0,10);
2097 fhdNdptpBarMCMin->SetXTitle("p_{T} (GeV/c)");
2098 fhdNdptpBarMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2099 fhdNdptpBarMCMin->Sumw2();
2100 fListOfHistos->Add(fhdNdptpBarMCMin);
2102 fhdNdptLambdaMCMin = new TH1F("hdNdptLambdaMCMin","",100,0,10);
2103 fhdNdptLambdaMCMin->SetXTitle("p_{T} (GeV/c)");
2104 fhdNdptLambdaMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2105 fhdNdptLambdaMCMin->Sumw2();
2106 fListOfHistos->Add(fhdNdptLambdaMCMin);
2108 fhdNdptLambdaBarMCMin = new TH1F("hdNdptLambdaBarMCMin","",100,0,10);
2109 fhdNdptLambdaBarMCMin->SetXTitle("p_{T} (GeV/c)");
2110 fhdNdptLambdaBarMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2111 fhdNdptLambdaBarMCMin->Sumw2();
2112 fListOfHistos->Add(fhdNdptLambdaBarMCMin);
2114 fhdNdptOmegaMCMin = new TH1F("hdNdptOmegaMCMin","",100,0,10);;
2115 fhdNdptOmegaMCMin->SetXTitle("p_{T} (GeV/c)");
2116 fhdNdptOmegaMCMin->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2117 fhdNdptOmegaMCMin->Sumw2();
2118 fListOfHistos->Add(fhdNdptOmegaMCMin);
2120 fhdNdptOmegaBarMCMin = new TH1F("hdNdptOmegaBarMCMin","",100,0,10);;
2121 fhdNdptOmegaBarMCMin->SetXTitle("p_{T} (GeV/c)");
2122 fhdNdptOmegaBarMCMin->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2123 fhdNdptOmegaBarMCMin->Sumw2();
2124 fListOfHistos->Add(fhdNdptOmegaBarMCMin);
2126 fhdNdptchPiMCJet = new TH1F("hdNdptchPiMCJet","",100,0,10);
2127 fhdNdptchPiMCJet->SetXTitle("p_{T} (GeV/c)");
2128 fhdNdptchPiMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2129 fhdNdptchPiMCJet->Sumw2();
2130 fListOfHistos->Add(fhdNdptchPiMCJet);
2132 fhdNdptK0MCJet = new TH1F("hdNdptK0MCJet","",100,0,10);
2133 fhdNdptK0MCJet->SetXTitle("p_{T} (GeV/c)");
2134 fhdNdptK0MCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2135 fhdNdptK0MCJet->Sumw2();
2136 fListOfHistos->Add(fhdNdptK0MCJet);
2138 fhdNdptchKMCJet = new TH1F("hdNdptchKMCJet","",100,0,10);
2139 fhdNdptchKMCJet->SetXTitle("p_{T} (GeV/c)");
2140 fhdNdptchKMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2141 fhdNdptchKMCJet->Sumw2();
2142 fListOfHistos->Add(fhdNdptchKMCJet);
2144 fhdNdptpMCJet = new TH1F("hdNdptpMCJet","",100,0,10);
2145 fhdNdptpMCJet->SetXTitle("p_{T} (GeV/c)");
2146 fhdNdptpMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2147 fhdNdptpMCJet->Sumw2();
2148 fListOfHistos->Add(fhdNdptpMCJet);
2150 fhdNdptpBarMCJet = new TH1F("hdNdptpBarMCJet","",100,0,10);
2151 fhdNdptpBarMCJet->SetXTitle("p_{T} (GeV/c)");
2152 fhdNdptpBarMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2153 fhdNdptpBarMCJet->Sumw2();
2154 fListOfHistos->Add(fhdNdptpBarMCJet);
2156 fhdNdptLambdaMCJet = new TH1F("hdNdptLambdaMCJet","",100,0,10);
2157 fhdNdptLambdaMCJet->SetXTitle("p_{T} (GeV/c)");
2158 fhdNdptLambdaMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2159 fhdNdptLambdaMCJet->Sumw2();
2160 fListOfHistos->Add(fhdNdptLambdaMCJet);
2162 fhdNdptLambdaBarMCJet = new TH1F("hdNdptLambdaBarMCJet","",100,0,10);
2163 fhdNdptLambdaBarMCJet->SetXTitle("p_{T} (GeV/c)");
2164 fhdNdptLambdaBarMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
2165 fhdNdptLambdaBarMCJet->Sumw2();
2166 fListOfHistos->Add(fhdNdptLambdaBarMCJet);
2168 fhPIDMC = CreatePIDhisto("hPIDMC");
2170 fListOfHistos->Add(fhPIDMC);
2172 fhPIDMC_quarkEv = CreatePIDhisto("hPIDMC_quarkEv");
2173 fhPIDMC_quarkEv->Sumw2();
2174 fListOfHistos->Add(fhPIDMC_quarkEv);
2176 fhPIDMC_gluonEv = CreatePIDhisto("hPIDMC_gluonEv");
2177 fhPIDMC_gluonEv->Sumw2();
2178 fListOfHistos->Add(fhPIDMC_gluonEv);
2180 fhPIDMCAll = CreatePIDhisto("hPIDMCAll");
2181 fhPIDMCAll->Sumw2();
2182 fListOfHistos->Add(fhPIDMCAll);
2184 fhPIDMCMin = CreatePIDhisto("hPIDMCMin");
2185 fhPIDMCMin->Sumw2();
2186 fListOfHistos->Add(fhPIDMCMin);
2188 fhPIDMCJet = CreatePIDhisto("hPIDMCJet");
2189 fhPIDMCJet->Sumw2();
2190 fListOfHistos->Add(fhPIDMCJet);
2192 fhPIDMCMotherK0 = CreatePIDhisto("hPIDMCMotherK0");
2193 fhPIDMCMotherK0->Sumw2();
2194 fListOfHistos->Add(fhPIDMCMotherK0);
2196 fhPIDMCGrandMotherK0 = CreatePIDhisto("hPIDMCGrandMotherK0");
2197 fhPIDMCGrandMotherK0->Sumw2();
2198 fListOfHistos->Add(fhPIDMCGrandMotherK0);
2200 fhPIDMCMotherChK = CreatePIDhisto("hPIDMCMotherChK");
2201 fhPIDMCMotherChK->Sumw2();
2202 fListOfHistos->Add(fhPIDMCMotherChK);
2204 fhPIDMCMotherK0Trans = CreatePIDhisto("hPIDMCMotherK0Trans");
2205 fhPIDMCMotherK0Trans->Sumw2();
2206 fListOfHistos->Add(fhPIDMCMotherK0Trans);
2208 fhPIDMCGrandMotherK0Trans = CreatePIDhisto("hPIDMCGrandMotherK0Trans");
2209 fhPIDMCGrandMotherK0Trans->Sumw2();
2210 fListOfHistos->Add(fhPIDMCGrandMotherK0Trans);
2212 fhPIDMCMotherChKTrans = CreatePIDhisto("hPIDMCMotherChKTrans");
2213 fhPIDMCMotherChKTrans->Sumw2();
2214 fListOfHistos->Add(fhPIDMCMotherChKTrans);
2216 fhdNdptgammaMC = new TH1F("hdNdptgammaMC","",100,0,10);
2217 fhdNdptgammaMC->SetXTitle("p_{T} (GeV/c)");
2218 fhdNdptgammaMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2219 fhdNdptgammaMC->Sumw2();
2220 fListOfHistos->Add(fhdNdptgammaMC);
2222 fhdNdptchPiMC = new TH1F("hdNdptchPiMC","",100,0,10);;
2223 fhdNdptchPiMC->SetXTitle("p_{T} (GeV/c)");
2224 fhdNdptchPiMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2225 fhdNdptchPiMC->Sumw2();
2226 fListOfHistos->Add(fhdNdptchPiMC);
2228 fhdNdptpi0MC = new TH1F("hdNdptpi0MC","",100,0,10);;
2229 fhdNdptpi0MC->SetXTitle("p_{T} (GeV/c)");
2230 fhdNdptpi0MC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2231 fhdNdptpi0MC->Sumw2();
2232 fListOfHistos->Add(fhdNdptpi0MC);
2234 fhdNdptK0MC = new TH1F("hdNdptK0MC","",100,0,10);;
2235 fhdNdptK0MC->SetXTitle("p_{T} (GeV/c)");
2236 fhdNdptK0MC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2237 fhdNdptK0MC->Sumw2();
2238 fListOfHistos->Add(fhdNdptK0MC);
2240 fhdNdptchKMC = new TH1F("hdNdptchKMC","",100,0,10);;
2241 fhdNdptchKMC->SetXTitle("p_{T} (GeV/c)");
2242 fhdNdptchKMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2243 fhdNdptchKMC->Sumw2();
2244 fListOfHistos->Add(fhdNdptchKMC);
2246 fhdNdptpMC = new TH1F("hdNdptpMC","",100,0,10);;
2247 fhdNdptpMC->SetXTitle("p_{T} (GeV/c)");
2248 fhdNdptpMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2249 fhdNdptpMC->Sumw2();
2250 fListOfHistos->Add(fhdNdptpMC);
2252 fhdNdptpBarMC = new TH1F("hdNdptpBarMC","",100,0,10);;
2253 fhdNdptpBarMC->SetXTitle("p_{T} (GeV/c)");
2254 fhdNdptpBarMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2255 fhdNdptpBarMC->Sumw2();
2256 fListOfHistos->Add(fhdNdptpBarMC);
2258 fhdNdptLambdaMC = new TH1F("hdNdptLambdaMC","",100,0,10);;
2259 fhdNdptLambdaMC->SetXTitle("p_{T} (GeV/c)");
2260 fhdNdptLambdaMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2261 fhdNdptLambdaMC->Sumw2();
2262 fListOfHistos->Add(fhdNdptLambdaMC);
2264 fhdNdptLambdaBarMC = new TH1F("hdNdptLambdaBarMC","",100,0,10);;
2265 fhdNdptLambdaBarMC->SetXTitle("p_{T} (GeV/c)");
2266 fhdNdptLambdaBarMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2267 fhdNdptLambdaBarMC->Sumw2();
2268 fListOfHistos->Add(fhdNdptLambdaBarMC);
2270 fhdNdptOmegaMC = new TH1F("hdNdptOmegaMC","",100,0,10);;
2271 fhdNdptOmegaMC->SetXTitle("p_{T} (GeV/c)");
2272 fhdNdptOmegaMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2273 fhdNdptOmegaMC->Sumw2();
2274 fListOfHistos->Add(fhdNdptOmegaMC);
2276 fhdNdptOmegaBarMC = new TH1F("hdNdptOmegaBarMC","",100,0,10);;
2277 fhdNdptOmegaBarMC->SetXTitle("p_{T} (GeV/c)");
2278 fhdNdptOmegaBarMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2279 fhdNdptOmegaBarMC->Sumw2();
2280 fListOfHistos->Add(fhdNdptOmegaBarMC);
2282 fhdNdxiMC = new TH1F("hdNdxiMC","",100,0,10);
2283 fhdNdxiMC->SetXTitle("#xi");
2285 fListOfHistos->Add(fhdNdxiMC);
2287 fhdNdxiK0MC = new TH1F("hdNdxiK0MC","",100,0,10);
2288 fhdNdxiK0MC->SetXTitle("#xi");
2289 fhdNdxiK0MC->Sumw2();
2290 fListOfHistos->Add(fhdNdxiK0MC);
2292 fhdNdxiK0MCJet = new TH1F("hdNdxiK0MCJet","",100,0,10);
2293 fhdNdxiK0MCJet->SetXTitle("#xi");
2294 fhdNdxiK0MCJet->Sumw2();
2295 fListOfHistos->Add(fhdNdxiK0MCJet);
2297 fhdNdzK0MC = new TH1F("hdNdzK0MC","",150,0,1.5);
2298 fhdNdzK0MC->SetXTitle("z");
2299 fhdNdzK0MC->SetYTitle("1/N_{jet} dN/dz");
2300 fhdNdzK0MC->Sumw2();
2301 fListOfHistos->Add(fhdNdzK0MC);
2303 fhdNdzK0MCJet = new TH1F("hdNdzK0MCJet","",150,0,1.5);
2304 fhdNdzK0MCJet->SetXTitle("z");
2305 fhdNdzK0MCJet->SetYTitle("1/N_{jet} dN/dz");
2306 fhdNdzK0MCJet->Sumw2();
2307 fListOfHistos->Add(fhdNdzK0MCJet);
2309 fhdNdptK0MCJetEvt = new TH1F("hdNdptK0MCJetEvt","",100,0,10);;
2310 fhdNdptK0MCJetEvt->SetXTitle("p_{T} (GeV/c)");
2311 fhdNdptK0MCJetEvt->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
2312 fhdNdptK0MCJetEvt->Sumw2();
2313 fListOfHistos->Add(fhdNdptK0MCJetEvt);
2315 fhnJetsAODvsMC = new TH2F("hnJetsAODvsMC","",20,0,20,20,0.,20.);
2316 fhnJetsAODvsMC->SetXTitle("n jets MC");
2317 fhnJetsAODvsMC->SetXTitle("n jets AOD");
2318 fhnJetsAODvsMC->Sumw2();
2319 fListOfHistos->Add(fhnJetsAODvsMC);
2321 fhLeadingPtAODvsMC = new TH2F("hLeadingPtAODvsMC","",20,0,20,20,0.,20.);
2322 fhLeadingPtAODvsMC->SetXTitle("p_{T} MC (GeV/c)");
2323 fhLeadingPtAODvsMC->SetYTitle("p_{T} AOD (GeV/c)");
2324 fhLeadingPtAODvsMC->Sumw2();
2325 fListOfHistos->Add(fhLeadingPtAODvsMC);
2327 fhLeadingEtaAODvsMC = new TH2F("hLeadingEtaAODvsMC","",20,-1,1,20,-1.,1.);
2328 fhLeadingEtaAODvsMC->SetXTitle("#eta MC");
2329 fhLeadingEtaAODvsMC->SetYTitle("#eta AOD");
2330 fhLeadingEtaAODvsMC->Sumw2();
2331 fListOfHistos->Add(fhLeadingEtaAODvsMC);
2333 fhLeadingPhiAODvsMC = new TH2F("hLeadingPhiAODvsMC","",63,0,6.3,63,0.,6.3);
2334 fhLeadingPhiAODvsMC->SetXTitle("#phi MC");
2335 fhLeadingPhiAODvsMC->SetYTitle("#phi AOD");
2336 fhLeadingPhiAODvsMC->Sumw2();
2337 fListOfHistos->Add(fhLeadingPhiAODvsMC);
2340 fhnTracksLeadingAODvsMC = new TH2F("hnTracksLeadingAODvsMC","",20,0.,20,20,-0.5,19.5);
2341 fhnTracksLeadingAODvsMC->SetXTitle("nTracks MC");
2342 fhnTracksLeadingAODvsMC->SetYTitle("nTracks AOD");
2343 fhnTracksLeadingAODvsMC->Sumw2();
2344 fListOfHistos->Add(fhnTracksLeadingAODvsMC);
2346 fhLeadingdRAODMC = new TH1F("hLeadingdRAODMC","",40,0.,4);
2347 fhLeadingdRAODMC->SetXTitle("#Delta R");
2348 fhLeadingdRAODMC->Sumw2();
2349 fListOfHistos->Add(fhLeadingdRAODMC);
2352 fhLeadingPtAODvsMCdRcut = new TH2F("hLeadingPtAODvsMCdRcut","",40,0,20,40,0.,20.);
2353 fhLeadingPtAODvsMCdRcut->SetXTitle("p_{T} MC (GeV/c)");
2354 fhLeadingPtAODvsMCdRcut->SetYTitle("p_{T} AOD (GeV/c)");
2355 fhLeadingPtAODvsMCdRcut->Sumw2();
2356 fListOfHistos->Add(fhLeadingPtAODvsMCdRcut);
2359 fhdnTracksVsdPtLeadingAODMC = new TH2F("hdnTracksVsdPtLeadingAODMC","",40,-10.,10,20,-10.5,9.5);
2360 fhdnTracksVsdPtLeadingAODMC->SetXTitle("#Delta Pt AOD-MC (GeV/c)");
2361 fhdnTracksVsdPtLeadingAODMC->SetYTitle("#Delta N AOD-MC");
2362 fhdnTracksVsdPtLeadingAODMC->Sumw2();
2363 fListOfHistos->Add(fhdnTracksVsdPtLeadingAODMC);
2365 fhnTracksJetVsPtAOD = new TH2F("hnTracksJetVsPtAOD","",50,0,50,20,-0.5,19.5);
2366 fhnTracksJetVsPtAOD->SetXTitle("p_{T} (GeV/c)");
2367 fhnTracksJetVsPtAOD->SetYTitle("nTracks");
2368 fhnTracksJetVsPtAOD->Sumw2();
2369 fListOfHistos->Add(fhnTracksJetVsPtAOD);
2371 fhnTracksJetVsPtAODquarkEv = new TH2F("hnTracksJetVsPtAODquarkEv","",50,0,50,20,-0.5,19.5);
2372 fhnTracksJetVsPtAODquarkEv->SetXTitle("p_{T} (GeV/c)");
2373 fhnTracksJetVsPtAODquarkEv->SetYTitle("nTracks");
2374 fhnTracksJetVsPtAODquarkEv->Sumw2();
2375 fListOfHistos->Add(fhnTracksJetVsPtAODquarkEv);
2377 fhRadiusJetVsPtAOD = new TH2F("hRadiusJetVsPtAOD","",50,0,50,10,0.,1);
2378 fhRadiusJetVsPtAOD->SetXTitle("p_{T} (GeV/c)");
2379 fhRadiusJetVsPtAOD->SetYTitle("radius");
2380 fhRadiusJetVsPtAOD->Sumw2();
2381 fListOfHistos->Add(fhRadiusJetVsPtAOD);
2383 fhnTracksJetVsPtMC = new TH2F("hnTracksJetVsPtMC","",50,0,50,20,-0.5,19.5);
2384 fhnTracksJetVsPtMC->SetXTitle("p_{T} (GeV/c)");
2385 fhnTracksJetVsPtMC->SetYTitle("nTracks");
2386 fhnTracksJetVsPtMC->Sumw2();
2387 fListOfHistos->Add(fhnTracksJetVsPtMC);
2389 fhnTracksJetVsPtMCquarkEv = new TH2F("hnTracksJetVsPtMCquarkEv","",50,0,50,20,-0.5,19.5);
2390 fhnTracksJetVsPtMCquarkEv->SetXTitle("p_{T} (GeV/c)");
2391 fhnTracksJetVsPtMCquarkEv->SetYTitle("nTracks");
2392 fhnTracksJetVsPtMCquarkEv->Sumw2();
2393 fListOfHistos->Add(fhnTracksJetVsPtMCquarkEv);
2395 fhRadiusJetVsPtMC = new TH2F("hRadiusJetVsPtMC","",50,0,50,10,0.,1);
2396 fhRadiusJetVsPtMC->SetXTitle("p_{T} (GeV/c)");
2397 fhRadiusJetVsPtMC->SetYTitle("radius");
2398 fhRadiusJetVsPtMC->Sumw2();
2399 fListOfHistos->Add(fhRadiusJetVsPtMC);
2401 fhnTracksJetVsPtMCK0 = new TH2F("hnTracksJetVsPtMCK0","",50,0,50,20,-0.5,19.5);
2402 fhnTracksJetVsPtMCK0->SetXTitle("p_{T} (GeV/c)");
2403 fhnTracksJetVsPtMCK0->SetYTitle("nTracks");
2404 fhnTracksJetVsPtMCK0->Sumw2();
2405 fListOfHistos->Add(fhnTracksJetVsPtMCK0);
2408 fhnTracksJetVsPtMCK0quarkEv = new TH2F("hnTracksJetVsPtMCK0quarkEv","",50,0,50,20,-0.5,19.5);
2409 fhnTracksJetVsPtMCK0quarkEv->SetXTitle("p_{T} (GeV/c)");
2410 fhnTracksJetVsPtMCK0quarkEv->SetYTitle("nTracks");
2411 fhnTracksJetVsPtMCK0quarkEv->Sumw2();
2412 fListOfHistos->Add(fhnTracksJetVsPtMCK0quarkEv);
2414 fhRadiusJetVsPtMCK0 = new TH2F("hRadiusJetVsPtMCK0","",50,0,50,10,0.,1);
2415 fhRadiusJetVsPtMCK0->SetXTitle("p_{T} (GeV/c)");
2416 fhRadiusJetVsPtMCK0->SetYTitle("radius");
2417 fhRadiusJetVsPtMCK0->Sumw2();
2418 fListOfHistos->Add(fhRadiusJetVsPtMCK0);
2420 fhnTracksJetVsPtAODK0 = new TH2F("hnTracksJetVsPtAODK0","",50,0,50,20,-0.5,19.5);
2421 fhnTracksJetVsPtAODK0->SetXTitle("p_{T} (GeV/c)");
2422 fhnTracksJetVsPtAODK0->SetYTitle("nTracks AODK0");
2423 fhnTracksJetVsPtAODK0->Sumw2();
2424 fListOfHistos->Add(fhnTracksJetVsPtAODK0);
2426 fhnTracksJetVsPtAODK0quarkEv = new TH2F("hnTracksJetVsPtAODK0quarkEv","",50,0,50,20,-0.5,19.5);
2427 fhnTracksJetVsPtAODK0quarkEv->SetXTitle("p_{T} (GeV/c)");
2428 fhnTracksJetVsPtAODK0quarkEv->SetYTitle("nTracks AODK0quarkEv");
2429 fhnTracksJetVsPtAODK0quarkEv->Sumw2();
2430 fListOfHistos->Add(fhnTracksJetVsPtAODK0quarkEv);
2432 fhRadiusJetVsPtAODK0 = new TH2F("hRadiusJetVsPtAODK0","",50,0,50,10,0.,1);
2433 fhRadiusJetVsPtAODK0->SetXTitle("p_{T} (GeV/c)");
2434 fhRadiusJetVsPtAODK0->SetYTitle("radius");
2435 fhRadiusJetVsPtAODK0->Sumw2();
2436 fListOfHistos->Add(fhRadiusJetVsPtAODK0);
2438 fhnTracksJetVsPtAODpKch = new TH2F("hnTracksJetVsPtAODpKch","",50,0,50,20,-0.5,19.5);
2439 fhnTracksJetVsPtAODpKch->SetXTitle("p_{T} (GeV/c)");
2440 fhnTracksJetVsPtAODpKch->SetYTitle("nTracks AODpKch");
2441 fhnTracksJetVsPtAODpKch->Sumw2();
2442 fListOfHistos->Add(fhnTracksJetVsPtAODpKch);
2444 fhRadiusJetVsPtAODpKch = new TH2F("hRadiusJetVsPtAODpKch","",50,0,50,20,-0.5,19.5);
2445 fhRadiusJetVsPtAODpKch->SetXTitle("p_{T} (GeV/c)");
2446 fhRadiusJetVsPtAODpKch->SetYTitle("Radius AODpKch");
2447 fhRadiusJetVsPtAODpKch->Sumw2();
2448 fListOfHistos->Add(fhRadiusJetVsPtAODpKch);
2450 fhPythiaProcess = CreatePythiaIDhisto("hPythiaProcess");
2451 fListOfHistos->Add(fhPythiaProcess);
2453 fhPythiaProcessK0 = CreatePythiaIDhisto("hPythiaProcessK0");
2454 fListOfHistos->Add(fhPythiaProcessK0);
2456 fhPythiaProcessKch = CreatePythiaIDhisto("hPythiaProcessKch");
2457 fListOfHistos->Add(fhPythiaProcessKch);
2459 fhPythiaProcessp = CreatePythiaIDhisto("hPythiaProcessp");
2460 fListOfHistos->Add(fhPythiaProcessp);
2462 fhPythiaProcesspbar = CreatePythiaIDhisto("hPythiaProcesspbar");
2463 fListOfHistos->Add(fhPythiaProcesspbar);
2465 fhdNdzJets5to10 = new TH1F("hdNdzJets5to10","",25,0,1.25);
2466 fhdNdzJets5to10->SetXTitle("z");
2467 fhdNdzJets5to10->SetYTitle("dN/dz");
2468 fhdNdzJets5to10->Sumw2();
2469 fListOfHistos->Add(fhdNdzJets5to10);
2471 fhdNdzJets10to20 = new TH1F("hdNdzJets10to20","",25,0,1.25);
2472 fhdNdzJets10to20->SetXTitle("z");
2473 fhdNdzJets10to20->SetYTitle("dN/dz");
2474 fhdNdzJets10to20->Sumw2();
2475 fListOfHistos->Add(fhdNdzJets10to20);
2477 fhdNdzJets20to30 = new TH1F("hdNdzJets20to30","",25,0,1.25);
2478 fhdNdzJets20to30->SetXTitle("z");
2479 fhdNdzJets20to30->SetYTitle("dN/dz");
2480 fhdNdzJets20to30->Sumw2();
2481 fListOfHistos->Add(fhdNdzJets20to30);
2483 fhdNdzJets30to40 = new TH1F("hdNdzJets30to40","",25,0,1.25);
2484 fhdNdzJets30to40->SetXTitle("z");
2485 fhdNdzJets30to40->SetYTitle("dN/dz");
2486 fhdNdzJets30to40->Sumw2();
2487 fListOfHistos->Add(fhdNdzJets30to40);
2489 fhdNdzJets40to60 = new TH1F("hdNdzJets40to60","",25,0,1.25);
2490 fhdNdzJets40to60->SetXTitle("z");
2491 fhdNdzJets40to60->SetYTitle("dN/dz");
2492 fhdNdzJets40to60->Sumw2();
2493 fListOfHistos->Add(fhdNdzJets40to60);
2496 fhdNdxiJets5to10 = new TH1F("hdNdxiJets5to10","",70,0,7);
2497 fhdNdxiJets5to10->SetXTitle("z");
2498 fhdNdxiJets5to10->SetYTitle("dN/dz");
2499 fhdNdxiJets5to10->Sumw2();
2500 fListOfHistos->Add(fhdNdxiJets5to10);
2502 fhdNdxiJets10to20 = new TH1F("hdNdxiJets10to20","",70,0,7);
2503 fhdNdxiJets10to20->SetXTitle("z");
2504 fhdNdxiJets10to20->SetYTitle("dN/dz");
2505 fhdNdxiJets10to20->Sumw2();
2506 fListOfHistos->Add(fhdNdxiJets10to20);
2508 fhdNdxiJets20to30 = new TH1F("hdNdxiJets20to30","",70,0,7);
2509 fhdNdxiJets20to30->SetXTitle("z");
2510 fhdNdxiJets20to30->SetYTitle("dN/dz");
2511 fhdNdxiJets20to30->Sumw2();
2512 fListOfHistos->Add(fhdNdxiJets20to30);
2514 fhdNdxiJets30to40 = new TH1F("hdNdxiJets30to40","",70,0,7);
2515 fhdNdxiJets30to40->SetXTitle("z");
2516 fhdNdxiJets30to40->SetYTitle("dN/dz");
2517 fhdNdxiJets30to40->Sumw2();
2518 fListOfHistos->Add(fhdNdxiJets30to40);
2520 fhdNdxiJets40to60 = new TH1F("hdNdxiJets40to60","",70,0,7);
2521 fhdNdxiJets40to60->SetXTitle("z");
2522 fhdNdxiJets40to60->SetYTitle("dN/dz");
2523 fhdNdxiJets40to60->Sumw2();
2524 fListOfHistos->Add(fhdNdxiJets40to60);
2526 fhdNdptTracksJetPt5to10 = new TH1F("hdNdptTracksJetPt5to10","",250,0,25);
2527 fhdNdptTracksJetPt5to10->SetXTitle("p_{T} (GeV)");
2528 fhdNdptTracksJetPt5to10->SetYTitle("dN/dp_{T} 1/GeV");
2529 fhdNdptTracksJetPt5to10->Sumw2();
2530 fListOfHistos->Add(fhdNdptTracksJetPt5to10);
2532 fhdNdptTracksJetPt10to20 = new TH1F("hdNdptTracksJetPt10to20","",25,0,25);
2533 fhdNdptTracksJetPt10to20->SetXTitle("p_{T} (GeV)");
2534 fhdNdptTracksJetPt10to20->SetYTitle("dN/dp_{T} 1/GeV");
2535 fhdNdptTracksJetPt10to20->Sumw2();
2536 fListOfHistos->Add(fhdNdptTracksJetPt10to20);
2538 fhdNdptTracksJetPt20to30 = new TH1F("hdNdptTracksJetPt20to30","",25,0,25);
2539 fhdNdptTracksJetPt20to30->SetXTitle("p_{T} (GeV)");
2540 fhdNdptTracksJetPt20to30->SetYTitle("dN/dp_{T} 1/GeV");
2541 fhdNdptTracksJetPt20to30->Sumw2();
2542 fListOfHistos->Add(fhdNdptTracksJetPt20to30);
2544 fhdNdptTracksJetPt30to40 = new TH1F("hdNdptTracksJetPt30to40","",25,0,25);
2545 fhdNdptTracksJetPt30to40->SetXTitle("p_{T} (GeV)");
2546 fhdNdptTracksJetPt30to40->SetYTitle("dN/dp_{T} 1/GeV");
2547 fhdNdptTracksJetPt30to40->Sumw2();
2548 fListOfHistos->Add(fhdNdptTracksJetPt30to40);
2550 fhdNdptTracksJetPt40to60 = new TH1F("hdNdptTracksJetPt40to60","",25,0,25);
2551 fhdNdptTracksJetPt40to60->SetXTitle("p_{T} (GeV)");
2552 fhdNdptTracksJetPt40to60->SetYTitle("dN/dp_{T} 1/GeV");
2553 fhdNdptTracksJetPt40to60->Sumw2();
2554 fListOfHistos->Add(fhdNdptTracksJetPt40to60);
2557 fh1Xsec = new TProfile("h1Xsec","xsec from pyxsec.root",1,0,1);
2558 fh1Xsec->SetXTitle("<#sigma>");
2560 fListOfHistos->Add( fh1Xsec );
2562 fh1Trials = new TH1F("h1Trials","trials from pyxsec.root",1,0,1);
2563 fh1Trials->SetXTitle("#sum{ntrials}");
2565 fListOfHistos->Add( fh1Trials );
2567 // fSettingsTree = new TTree("JetChemAnalysisSettings","Analysis Settings");
2568 // fSettingsTree->Branch("fUseLOConeJets",&fUseLOConeJets,"UseLOConeJets/O");
2569 // fSettingsTree->Branch("fUseLOConeMCJets",&fUseLOConeMCJets,"UseLOConeMCJets/O");
2570 // fSettingsTree->Branch("fUsePythiaJets",&fUsePythiaJets,"UsePythiaJets/O");
2571 // fSettingsTree->Branch("fConeRadius", &fConeRadius,"Rad/D");
2572 // fSettingsTree->Branch("fTrackPtCutJF", &fTrackPtCutJF,"TrackPt/D");
2573 // fSettingsTree->Branch("fFilterBitJF",&fFilterBitJF,"FilterBitJF/i");
2574 // fSettingsTree->Branch("fRequireITSRefitJF",&fRequireITSRefitJF,"RequireITSRefitJF/O");
2575 // fSettingsTree->Branch("fRejectK0TracksJF",&fRejectK0TracksJF,"RejectK0TracksJF/O");
2576 // fSettingsTree->Branch("fJetPtCut",&fJetPtCut,"JetPtCut/D");
2577 // fSettingsTree->Branch("fJetEtaCut",&fJetEtaCut,"JetEtaCut/D");
2578 // fSettingsTree->Branch("fFilterBit",&fFilterBit,"FilterBit/i");
2579 // fSettingsTree->Branch("fTrackPtCut",&fTrackPtCut,"TrackPtCut/D");
2580 // fSettingsTree->Branch("fTrackEtaCut",&fTrackEtaCut,"TrackEtaCut/D");
2581 // fSettingsTree->Branch("fUseOnFlyV0s",&fUseOnFlyV0s,"UseOnFlyV0s/O");
2582 // fSettingsTree->Branch("fCutnSigdEdx",&fCutnSigdEdx,"CutnSigdEdx/D");
2583 // fSettingsTree->Branch("fUseAODMCTracksForUE",&fUseAODMCTracksForUE,"UseAODMCTracksForUE/O");
2584 // fSettingsTree->Branch("fAreaReg",&fAreaReg,"AreaReg/D");
2585 // fSettingsTree->Branch("fAvgTrials",&fAvgTrials,"AvgTrials/D");
2586 // fListOfHistos->Add(fSettingsTree);
2590 //____________________________________________________________________
2591 void AliAnalysisTaskJetChem::Terminate(Option_t */*option*/){
2593 // Terminate analysis
2598 // Update pointers reading them from the output slot
2599 fListOfHistos = dynamic_cast<TList*> (GetOutputData(1));
2600 if( !fListOfHistos ) {
2601 AliError("Histogram List is not available");
2605 fhLeadingPt = (TH1F*)fListOfHistos->FindObject("hLeadingPt");
2606 fhRegionSumPtMaxVsEt = (TH1F*)fListOfHistos->FindObject("hRegionSumPtMaxVsEt");
2607 fhRegionSumPtMinVsEt = (TH1F*)fListOfHistos->FindObject("hRegionSumPtMinVsEt");
2608 fhRegionMultMaxVsEt = (TH1F*)fListOfHistos->FindObject("hRegionMultMaxVsEt");
2609 fhRegionMultMinVsEt = (TH1F*)fListOfHistos->FindObject("hRegionMultMinVsEt");
2612 fhRegionSumPtMaxVsEt->Divide(fhLeadingPt);
2613 fhRegionSumPtMinVsEt->Divide(fhLeadingPt);
2614 fhRegionMultMaxVsEt->Divide(fhLeadingPt);
2615 fhRegionMultMinVsEt->Divide(fhLeadingPt);
2619 fh1Xsec = (TProfile*) fListOfHistos->FindObject("fh1Xsec");
2620 fh1Trials = (TH1F*)fListOfHistos->FindObject("fh1Trials");
2622 //std::cout<<" fh1Xsec "<<fh1Xsec<<" fh1Trials "<<fh1Trials<<std::endl;
2624 if(fh1Xsec && fh1Trials){
2626 Double_t xsec = fh1Xsec->GetBinContent(1);
2627 Double_t ntrials = fh1Trials->GetBinContent(1);
2628 Double_t normFactor = xsec/ntrials;
2629 Printf("xSec %f nTrials %f Norm %f \n",xsec,ntrials,normFactor);
2630 fhLeadingPt->Scale(normFactor);
2634 if(fDebug > 1) AliInfo("End analysis");
2637 // ----------------------------------------------------------------------------
2639 void AliAnalysisTaskJetChem::WriteSettings(){
2641 // fill settings tree
2644 //fSettingsTree->Fill();
2648 // ---------------------------------------------------------------------------
2650 Bool_t AliAnalysisTaskJetChem::IsK0InvMass(const Double_t mass) const {
2652 // K0 mass ? Use STAR params for bin counting & mass fit
2654 const Double_t massK0 = 0.497; // from fits
2655 const Double_t sigmaK0 = 0.0046;
2656 const Double_t nSigmaSignal = 3.5; // STAR parameters for bin counting
2658 if((massK0-nSigmaSignal*sigmaK0)<=mass &&
2659 mass<=(massK0 + nSigmaSignal*sigmaK0)) return kTRUE;
2664 // ---------------------------------------------------------------------------
2666 Bool_t AliAnalysisTaskJetChem::IsLambdaInvMass(const Double_t mass) const{
2671 if(1.1<mass && mass<1.13) return kTRUE; // FIXME - adjust range from fit
2675 // ---------------------------------------------------------------------------
2677 Bool_t AliAnalysisTaskJetChem::IsAcceptedDCAK0(/*const Double_t dca*/) const{
2681 return kTRUE; // FIXME - adjust cut
2685 // ---------------------------------------------------------------------------
2687 Bool_t AliAnalysisTaskJetChem::IsAcceptedDCALambda(/*const Double_t dca*/) const {
2691 return kTRUE; // FIXME - adjust cut
2694 // ---------------------------------------------------------------------------
2696 Bool_t AliAnalysisTaskJetChem::IsAccepteddEdx(const Double_t mom,const Double_t signal, AliPID::EParticleType n, const Double_t cutnSig) const{
2698 // apply TPC dE/dx cut similar as in AliTPCpidESD
2699 // note: AliTPCidESD uses inner track param for momentum - not avaiable on AOD,
2700 // so we use global track momentum
2701 // should use separate parametrisation for MC and data, but probably ALEPH param & 7% resolution used here anyway not the last word
2704 const Double_t kBBMIP(50.);
2705 const Double_t kBBRes(0.07);
2706 //const Double_t kBBRange(5.);
2707 const Double_t kBBp1(0.76176e-1);
2708 const Double_t kBBp2(10.632);
2709 const Double_t kBBp3(0.13279e-4);
2710 const Double_t kBBp4(1.8631);
2711 const Double_t kBBp5(1.9479);
2713 Double_t mass=AliPID::ParticleMass(n);
2714 Double_t betaGamma = mom/mass;
2716 const Float_t kmeanCorrection =0.1;
2717 Double_t bb = AliExternalTrackParam::BetheBlochAleph(betaGamma,kBBp1,kBBp2,kBBp3,kBBp4,kBBp5);
2718 Double_t meanCorrection =(1+(bb-1)*kmeanCorrection);
2719 Double_t bethe = bb * meanCorrection; // expected
2720 Double_t sigma = bethe * kBBRes;
2723 Double_t dedx = signal/kBBMIP; // measured
2725 Double_t nSig = (TMath::Abs(dedx - bethe))/sigma;
2727 if(nSig > cutnSig) return kFALSE;
2732 // ----------------------------------------------------------------------------
2734 void AliAnalysisTaskJetChem::CheckV0s(AliAODJet* jetVect, Int_t maxPtRegionIndex,Bool_t& foundK0){
2736 // loop over AOD V0s, fill masses etc
2738 Int_t nV0 = fAOD->GetNumberOfV0s();
2741 for(int i=0; i<fAOD->GetNumberOfV0s(); i++){ // loop over V0s
2743 AliAODv0* v0 = fAOD->GetV0(i);
2745 Bool_t isOnFly = v0->GetOnFlyStatus();
2746 isOnFly ? fhV0onFly->Fill(1) : fhV0onFly->Fill(0);
2748 if((fUseOnFlyV0s && isOnFly) || (!fUseOnFlyV0s && !isOnFly) ){ // 'offline' V0s : using vertex tracks, on-fly: take decay vertex into account during mom fit
2750 Double_t massK0 = v0->MassK0Short();
2751 Double_t massLambda = v0->MassLambda();
2752 Double_t massAntiLambda = v0->MassAntiLambda();
2753 Double_t radiusV0 = v0->RadiusV0();
2755 Double_t etaV0 = v0->Eta();
2756 Double_t fDCAV0ToVertex = v0->DcaV0ToPrimVertex();
2757 Double_t fDCADaughters = v0->DcaV0Daughters();
2762 TVector3 v0MomVect(v0Mom);
2764 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0));
2765 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1));
2767 Double_t posMom[3], negMom[3];
2768 trackPos->PxPyPz(posMom);
2769 trackNeg->PxPyPz(negMom);
2771 TVector3 posMomVec(posMom);
2772 TVector3 negMomVec(negMom);
2774 Double_t dROpanV0 = posMomVec.DeltaR(negMomVec);
2776 AliAODPid* aodPidPos = trackPos->GetDetPid();
2777 AliAODPid* aodPidNeg = trackNeg->GetDetPid();
2779 Double_t dEdxPos = aodPidPos->GetTPCsignal();
2780 Double_t dEdxNeg = aodPidNeg->GetTPCsignal();
2782 Double_t momPos = trackPos->P();
2783 Double_t momNeg = trackNeg->P();
2785 fhdEdxVsMomV0->Fill(momPos,dEdxPos);
2786 fhdEdxVsMomV0->Fill(momNeg,dEdxNeg);
2788 Double_t pV0 = TMath::Sqrt(v0->Ptot2V0());
2789 Double_t ptV0 = TMath::Sqrt(v0->Pt2V0());
2791 fhV0InvMassK0->Fill(massK0);
2792 fhV0InvMassLambda->Fill(massLambda);
2793 fhV0InvMassAntiLambda->Fill(massAntiLambda);
2795 fhV0DCADaughters->Fill(fDCADaughters);
2796 fhV0Radius->Fill(radiusV0);
2797 fhV0DCAToVertex->Fill(fDCAV0ToVertex);
2798 fhdNdptV0->Fill(ptV0);
2800 fhV0PtVsInvMassK0->Fill(massK0,ptV0);
2803 // K0 signal before cuts
2805 if(IsK0InvMass(massK0)){
2807 fhdNdptK0->Fill(ptV0);
2809 fhPtVsEtaK0->Fill(etaV0,ptV0);
2811 Double_t dRV0MC = AssociateV0MC(&v0MomVect,310); // K0
2812 if(dRV0MC < 0) dRV0MC = 0.99;
2813 fhdRV0MC->Fill(dRV0MC);
2816 if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/)){
2818 fhV0InvMassK0DCA->Fill(massK0);
2819 if(IsK0InvMass(massK0)) fhdNdptK0DCA->Fill(ptV0);
2822 if(IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)) fhdEdxVsMomV0pidEdx->Fill(momPos,dEdxPos);
2823 if(IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)) fhdEdxVsMomV0pidEdx->Fill(-1*momNeg,dEdxNeg);
2825 if(IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) &&
2826 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
2828 fhV0InvMassK0DCAdEdx->Fill(massK0);
2829 fhV0PtVsInvMassK0DCAdEdx->Fill(massK0,ptV0);
2831 if(IsK0InvMass(massK0)){
2832 fhdNdptK0DCAdEdx->Fill(ptV0);
2833 fhV0DCAToVertexK0->Fill(fDCAV0ToVertex);
2834 fhdROpanK0VsPt->Fill(ptV0,dROpanV0);
2840 // AOD pid - cuts strongly into signal
2842 AliAODTrack::AODTrkPID_t mpPIDNeg = trackNeg->GetMostProbablePID();
2843 AliAODTrack::AODTrkPID_t mpPIDPos = trackPos->GetMostProbablePID();
2845 if(mpPIDNeg == AliAODTrack::kPion) fhdEdxVsMomV0piPID->Fill(momPos,dEdxPos);
2846 if(mpPIDPos == AliAODTrack::kPion) fhdEdxVsMomV0piPID->Fill(-1*momNeg,dEdxNeg);
2849 if( (mpPIDNeg == AliAODTrack::kPion) && (mpPIDPos == AliAODTrack::kPion) ){
2851 fhV0InvMassK0DCAPID->Fill(massK0);
2856 if(IsLambdaInvMass(massLambda) || (IsLambdaInvMass(massAntiLambda))){
2858 fhV0InvMassK0Lambda->Fill(massK0);
2864 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
2865 IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
2866 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
2868 fhV0InvMassLambdaDCAdEdx->Fill(massLambda);
2871 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
2872 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
2873 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
2875 fhV0InvMassAntiLambdaDCAdEdx->Fill(massAntiLambda);
2883 Int_t regionV0vect = IsTrackInsideRegion(jetVect,&v0MomVect);
2886 Double_t jetpt = jetVect->Pt();
2887 Double_t dPhiJetV0 = (jetVect->MomentumVector()->Vect()).DeltaPhi(v0MomVect);
2888 Double_t dRJetV0 = (jetVect->MomentumVector()->Vect()).DeltaR(v0MomVect);
2890 Double_t z = pV0/jetpt;
2891 Double_t xi = TMath::Log(1/z);
2893 fhdPhiJetV0->Fill(dPhiJetV0);
2897 if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/) &&
2898 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) &&
2899 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
2901 fhV0InvMassK0JetEvt->Fill(massK0);
2903 if(IsK0InvMass(massK0)){
2904 fhdNdptK0JetEvt->Fill(ptV0);
2906 fhdNdxiK0->Fill(xi);
2908 fhdPhiJetK0->Fill(dPhiJetV0);
2909 fhdRJetK0->Fill(dRJetV0);
2911 if(5<jetpt && jetpt<=10) fhdNdzK05to10->Fill(z);
2912 else if(10<jetpt && jetpt<=20) fhdNdzK010to20->Fill(z);
2913 else if(20<jetpt && jetpt<=30) fhdNdzK020to30->Fill(z);
2914 else if(30<jetpt && jetpt<=40) fhdNdzK030to40->Fill(z);
2915 else if(40<jetpt && jetpt<=60) fhdNdzK040to60->Fill(z);
2922 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
2923 IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
2924 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
2926 fhV0InvMassLambdaJetEvt->Fill(massLambda);
2928 if(IsLambdaInvMass(massLambda)){
2929 fhdNdzLambda->Fill(z);
2933 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
2934 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
2935 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
2937 fhV0InvMassAntiLambdaJetEvt->Fill(massAntiLambda);
2939 if(IsLambdaInvMass(massAntiLambda)){
2940 fhdNdzAntiLambda->Fill(z);
2944 // fill histos max region
2946 if(regionV0vect != 0 && regionV0vect == maxPtRegionIndex){ // max region
2950 if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/) &&
2951 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) &&
2952 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){ // K0 cuts
2954 fhV0InvMassK0Max->Fill(massK0);
2956 if(IsK0InvMass(massK0)){
2957 fhdNdzK0Max->Fill(z);
2958 fhdNdxiK0Max->Fill(xi);
2959 fhdNdptK0Max->Fill(ptV0);
2965 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
2966 IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
2967 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
2969 fhV0InvMassLambdaMax->Fill(massLambda);
2972 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
2973 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
2974 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
2976 fhV0InvMassAntiLambdaMax->Fill(massAntiLambda);
2979 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
2980 ((IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
2981 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)) ||
2982 (IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
2983 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)))){
2985 if(IsLambdaInvMass(massLambda) || (IsLambdaInvMass(massAntiLambda))){
2986 fhdNdzLambdaMax->Fill(z);
2987 fhdNdxiLambdaMax->Fill(xi);
2988 fhdNdptLambdaMax->Fill(ptV0);
2993 // fill histos min region
2995 if(regionV0vect != 0 && regionV0vect != maxPtRegionIndex){ // min region
2999 if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/) &&
3000 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) &&
3001 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
3003 fhV0InvMassK0Min->Fill(massK0);
3005 if(IsK0InvMass(massK0)){
3006 fhdNdzK0Min->Fill(z);
3007 fhdNdxiK0Min->Fill(xi);
3008 fhdNdptK0Min->Fill(ptV0);
3016 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
3017 IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
3018 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
3020 fhV0InvMassLambdaMin->Fill(massLambda);
3023 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
3024 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
3025 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
3027 fhV0InvMassAntiLambdaMin->Fill(massAntiLambda);
3030 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
3031 ((IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
3032 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)) ||
3033 (IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
3034 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)))){
3036 if(IsLambdaInvMass(massLambda) || (IsLambdaInvMass(massAntiLambda))){
3037 fhdNdzLambdaMin->Fill(z);
3038 fhdNdxiLambdaMin->Fill(xi);
3039 fhdNdptLambdaMin->Fill(ptV0);
3047 if(regionV0vect == 0){ // jet region
3049 //Double_t dRJetV0 = jetVect->DeltaR(v0MomVect);
3051 if(dRJetV0 <= fConeRadius){
3053 // fill histos jet region
3057 if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/) &&
3058 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) &&
3059 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){ // K0 cuts
3061 fhV0InvMassK0Jet->Fill(massK0);
3063 if(IsK0InvMass(massK0)){
3064 fhdNdzK0Jet->Fill(z);
3065 fhdNdxiK0Jet->Fill(xi);
3066 fhdNdptK0Jet->Fill(ptV0);
3073 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
3074 IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
3075 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
3077 fhV0InvMassLambdaJet->Fill(massLambda);
3079 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
3080 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
3081 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
3083 fhV0InvMassAntiLambdaJet->Fill(massAntiLambda);
3086 if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) &&
3087 ((IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) &&
3088 IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)) ||
3089 (IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) &&
3090 IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)))){
3092 if(IsLambdaInvMass(massLambda) || (IsLambdaInvMass(massAntiLambda))){
3093 fhdNdzLambdaJet->Fill(z);
3094 fhdNdxiLambdaJet->Fill(xi);
3095 fhdNdptLambdaJet->Fill(ptV0);
3105 // ----------------------------------------------------------------------------
3107 void AliAnalysisTaskJetChem::CheckMCParticles(AliAODJet* jetVect,Int_t maxPtRegionIndex,Bool_t& isK0Event){
3109 // histos for generated MC
3111 TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
3114 AliInfo("no mcparticles branch");
3118 Int_t ntrks = farray->GetEntries();
3119 if (fDebug>1) AliInfo(Form("check MC particles, tracks %d \n",ntrks));
3121 Int_t pythiaPID = GetPythiaProcessID();
3124 Bool_t ispEvent = kFALSE;
3125 Bool_t ispBarEvent = kFALSE;
3126 Bool_t isKchEvent = kFALSE;
3128 Bool_t isQuarkHardScatteringEvent = IsQuarkHardScatteringEvent(pythiaPID);
3129 Bool_t isGluonHardScatteringEvent = IsGluonHardScatteringEvent(pythiaPID);
3131 for(Int_t i =0 ; i<ntrks; i++){ // mc tracks loop
3133 AliAODMCParticle* mctrk = (AliAODMCParticle*)farray->At(i);
3135 Double_t trackPt = mctrk->Pt();
3138 if (!(mctrk->IsPhysicalPrimary())) continue;
3140 if ((trackPt < fTrackPtCut) || (TMath::Abs(mctrk->Eta()) > fTrackEtaCut)) continue;
3143 Int_t pdg = mctrk->GetPdgCode();
3146 FillPIDhisto(fhPIDMC,pdg);
3147 if(isQuarkHardScatteringEvent) FillPIDhisto(fhPIDMC_quarkEv,pdg);
3148 if(isGluonHardScatteringEvent) FillPIDhisto(fhPIDMC_gluonEv,pdg);
3150 if(pdg == 111) fhdNdptpi0MC->Fill(trackPt);
3151 if(pdg == 22) fhdNdptgammaMC->Fill(trackPt);
3152 if(TMath::Abs(pdg) == 211) fhdNdptchPiMC->Fill(trackPt);
3153 if(pdg == 310) fhdNdptK0MC->Fill(trackPt);
3154 if(TMath::Abs(pdg) == 321) fhdNdptchKMC->Fill(trackPt);
3155 if(pdg == 2212) fhdNdptpMC->Fill(trackPt);
3156 if(pdg == -2212) fhdNdptpBarMC->Fill(trackPt);
3157 if(pdg == 3122) fhdNdptLambdaMC->Fill(trackPt);
3158 if(pdg == -3122) fhdNdptLambdaBarMC->Fill(trackPt);
3159 if(pdg == 3332) fhdNdptOmegaMC->Fill(trackPt);
3160 if(pdg == -3332) fhdNdptOmegaBarMC->Fill(trackPt);
3162 if(pdg == 310) isK0Event = kTRUE;
3163 if(TMath::Abs(pdg) == 321) isKchEvent = kTRUE;
3164 if(pdg == 2212) ispEvent = kTRUE;
3165 if(pdg == -2212) ispBarEvent = kTRUE;
3168 Int_t pdgMotherChK = -1;
3169 Int_t pdgMotherK0 = -1;
3170 Int_t pdgGrandMotherK0 = -1;
3172 if(TMath::Abs(pdg) == 321){ // chK
3173 Int_t labelMother = mctrk->GetMother();
3174 AliAODMCParticle* mctrkMother = NULL;
3176 for(Int_t k=0 ; k<ntrks; k++){
3178 AliAODMCParticle* mctrk2 = (AliAODMCParticle*)farray->At(k);
3179 if(mctrk2->GetLabel() == labelMother) mctrkMother = mctrk2;
3182 if(mctrkMother) pdgMotherChK = mctrkMother->GetPdgCode();
3183 FillPIDhisto(fhPIDMCMotherChK,pdgMotherChK);
3185 //printf("pdgMotherChK %d \n",pdgMotherChK);
3188 if(pdg == 310){ // K0
3189 Int_t labelMother = mctrk->GetMother();
3190 AliAODMCParticle* mctrkMother = NULL;
3192 for(Int_t k=0 ; k<ntrks; k++){
3193 AliAODMCParticle* mctrk2 = (AliAODMCParticle*)farray->At(k);
3194 if(mctrk2->GetLabel() == labelMother) mctrkMother = mctrk2;
3197 if(mctrkMother) pdgMotherK0 = mctrkMother->GetPdgCode();
3198 FillPIDhisto(fhPIDMCMotherK0,pdgMotherK0);
3200 Int_t labelGrandMother = -1;
3201 if(mctrkMother) mctrkMother->GetMother();
3202 AliAODMCParticle* mctrkGrandMother = NULL;
3204 for(Int_t k=0 ; k<ntrks; k++){
3205 AliAODMCParticle* mctrk2 = (AliAODMCParticle*)farray->At(k);
3206 if(mctrk2->GetLabel() == labelGrandMother) mctrkGrandMother = mctrk2;
3209 if(mctrkGrandMother) pdgGrandMotherK0 = mctrkGrandMother->GetPdgCode();
3210 FillPIDhisto(fhPIDMCGrandMotherK0,pdgGrandMotherK0);
3213 if(jetVect){ // jet event
3215 FillPIDhisto(fhPIDMCAll,pdg);
3217 TVector3 partVect(mctrk->Px(), mctrk->Py(), mctrk->Pz());
3219 Int_t region = IsTrackInsideRegion(jetVect, &partVect );
3221 //Double_t deltaPhi = jetVect->DeltaPhi(partVect); //+k270rad;
3222 Double_t deltaPhi = (jetVect->MomentumVector()->Vect()).DeltaPhi(partVect);
3223 if( deltaPhi > 2.*TMath::Pi() ) deltaPhi-= 2.*TMath::Pi();
3224 //Double_t deltaR = jetVect->DeltaR(partVect); //+k270rad;
3225 Double_t deltaR = (jetVect->MomentumVector()->Vect()).DeltaR(partVect);
3228 Double_t jetpt = jetVect->Pt();
3229 //Double_t pV0 = partVect.Mag();
3230 //Double_t ptV0 = partVect.Pt();
3233 Double_t z = trackPt / jetpt;
3234 Double_t xi = TMath::Log(1/z);
3236 if(!(mctrk->Charge() == 0 || mctrk->Charge()==-99)) fhdNdxiMC->Fill(xi);
3238 if(pdg == 310){ // K0
3239 fhdPhiJetK0MC->Fill(deltaPhi);
3240 fhdRJetK0MC->Fill(deltaR);
3241 fhdNdxiK0MC->Fill(xi);
3242 fhdNdzK0MC->Fill(z);
3243 fhdNdptK0MCJetEvt->Fill(trackPt);
3247 if(region != 0 && region == maxPtRegionIndex){ // max region
3249 if(TMath::Abs(pdg) == 211) fhdNdptchPiMCMax->Fill(trackPt);
3250 if(pdg == 310) fhdNdptK0MCMax->Fill(trackPt);
3251 if(TMath::Abs(pdg) == 321) fhdNdptchKMCMax->Fill(trackPt);
3252 if(pdg == 2212) fhdNdptpMCMax->Fill(trackPt);
3253 if(pdg == -2212) fhdNdptpBarMCMax->Fill(trackPt);
3254 if(pdg == 3122) fhdNdptLambdaMCMax->Fill(trackPt);
3255 if(pdg == -3122) fhdNdptLambdaBarMCMax->Fill(trackPt);
3258 if(region != 0 && region != maxPtRegionIndex){ // min region
3260 FillPIDhisto(fhPIDMCMin,pdg);
3262 if(TMath::Abs(pdg) == 211) fhdNdptchPiMCMin->Fill(trackPt);
3263 if(pdg == 310) fhdNdptK0MCMin->Fill(trackPt);
3264 if(TMath::Abs(pdg) == 321) fhdNdptchKMCMin->Fill(trackPt);
3265 if(pdg == 2212) fhdNdptpMCMin->Fill(trackPt);
3266 if(pdg == -2212) fhdNdptpBarMCMin->Fill(trackPt);
3267 if(pdg == 3122) fhdNdptLambdaMCMin->Fill(trackPt);
3268 if(pdg == -3122) fhdNdptLambdaBarMCMin->Fill(trackPt);
3270 if(pdg == 3332) fhdNdptOmegaMCMin->Fill(trackPt);
3271 if(pdg == -3332) fhdNdptOmegaBarMCMin->Fill(trackPt);
3276 if(TMath::Abs(pdg) == 321) FillPIDhisto(fhPIDMCMotherChKTrans,pdgMotherChK);
3278 FillPIDhisto(fhPIDMCMotherK0Trans,pdgMotherK0);
3279 FillPIDhisto(fhPIDMCGrandMotherK0Trans,pdgGrandMotherK0);
3283 if(region == 0){ // jet region ?
3285 FillPIDhisto(fhPIDMCJet,pdg);
3287 //Double_t dRJetV0 = jetVect->DeltaR(partVect);
3288 Double_t dRJetV0 = (jetVect->MomentumVector()->Vect()).DeltaR(partVect);
3290 if(dRJetV0 <= fConeRadius){
3292 if(pdg == 310){ // K0
3294 fhdNdptK0MCJet->Fill(trackPt);
3295 fhdNdxiK0MCJet->Fill(xi);
3296 fhdNdzK0MCJet->Fill(z);
3299 if(TMath::Abs(pdg) == 211) fhdNdptchPiMCJet->Fill(trackPt);
3300 if(TMath::Abs(pdg) == 321) fhdNdptchKMCJet->Fill(trackPt);
3301 if(pdg == 2212) fhdNdptpMCJet->Fill(trackPt);
3302 if(pdg == -2212) fhdNdptpBarMCJet->Fill(trackPt);
3303 if(pdg == 3122) fhdNdptLambdaMCJet->Fill(trackPt);
3304 if(pdg == -3122) fhdNdptLambdaBarMCJet->Fill(trackPt);
3311 FillPythiaIDhisto(fhPythiaProcess,pythiaPID);
3312 if(isK0Event) FillPythiaIDhisto(fhPythiaProcessK0,pythiaPID);
3313 if(isKchEvent) FillPythiaIDhisto(fhPythiaProcessKch,pythiaPID);
3314 if(ispEvent) FillPythiaIDhisto(fhPythiaProcessp,pythiaPID);
3315 if(ispBarEvent) FillPythiaIDhisto(fhPythiaProcesspbar,pythiaPID);
3319 // ----------------------------------------------------------------------------
3321 Double_t AliAnalysisTaskJetChem::AssociateV0MC(const TVector3* V0Mom,const Int_t pdg){
3323 // find closest MC gen. particle for V0 vector
3325 TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
3328 AliInfo("no mcparticles branch");
3332 Double_t dRmin = -1;
3334 Int_t ntrks = farray->GetEntries();
3336 for(Int_t i =0 ; i<ntrks; i++){
3338 AliAODMCParticle* mctrk = (AliAODMCParticle*)farray->At(i);
3341 if (!(mctrk->IsPhysicalPrimary())) continue;
3343 Int_t pdgtrk = mctrk->GetPdgCode();
3345 if(pdgtrk != pdg) continue;
3347 TVector3 partVect(mctrk->Px(), mctrk->Py(), mctrk->Pz());
3349 Double_t dR = V0Mom->DeltaR(partVect);
3351 if(dRmin<0) dRmin = dR; // initialize
3353 if(dR < dRmin) dRmin = dR;
3360 // ----------------------------------------------------------------------------
3362 void AliAnalysisTaskJetChem::CompLeadingJets(AliAODJet* jetLeadingAOD,AliAODJet* jetLeadingMC,const Int_t pythiaPID,
3363 const Bool_t foundK0AOD,const Bool_t foundK0MC){
3365 // leading jet properties
3367 Double_t ptLeadingAOD = -1;
3368 Double_t etaLeadingAOD = -1;
3369 Double_t phiLeadingAOD = -1;
3370 Int_t nTracksLeadingAOD = -1;
3372 Double_t ptLeadingMC = -1;
3373 Double_t etaLeadingMC = -1;
3374 Double_t phiLeadingMC = -1;
3375 Int_t nTracksLeadingMC = -1;
3378 if(jetLeadingAOD->Pt()>fJetPtCut && TMath::Abs(jetLeadingAOD->Eta())<fJetEtaCut){
3380 ptLeadingAOD = jetLeadingAOD->Pt();
3381 etaLeadingAOD = jetLeadingAOD->Eta();
3382 phiLeadingAOD = jetLeadingAOD->Phi();
3383 nTracksLeadingAOD = jetLeadingAOD->GetRefTracks()->GetEntriesFast();
3385 Double_t radiusAOD = GetJetRadius(jetLeadingAOD,0.8);
3386 fhnTracksJetVsPtAOD->Fill(ptLeadingAOD,nTracksLeadingAOD);
3387 fhRadiusJetVsPtAOD->Fill(ptLeadingAOD,radiusAOD);
3388 if(IsQuarkHardScatteringEvent(pythiaPID)) fhnTracksJetVsPtAODquarkEv->Fill(ptLeadingAOD,nTracksLeadingAOD);
3392 fhnTracksJetVsPtAODK0->Fill(ptLeadingAOD,nTracksLeadingAOD);
3393 if(IsQuarkHardScatteringEvent(pythiaPID)) fhnTracksJetVsPtAODK0quarkEv->Fill(ptLeadingAOD,nTracksLeadingAOD);
3394 fhRadiusJetVsPtAODK0->Fill(ptLeadingAOD,radiusAOD);
3398 // check if p/Kch in jet
3400 Bool_t foundpKch = kFALSE;
3401 Int_t nTracksJet = jetLeadingAOD->GetRefTracks()->GetEntriesFast();
3403 for(int i=0; i<nTracksJet; i++){
3405 AliAODTrack* track = (AliAODTrack*) jetLeadingAOD->GetRefTracks()->At(i);
3407 Double_t mom = track->P();
3409 AliAODPid* aodPid = track->GetDetPid();
3410 Double_t dEdx = aodPid->GetTPCsignal();
3412 if(IsAccepteddEdx(mom,dEdx,AliPID::kKaon,fCutnSigdEdx) ||
3413 IsAccepteddEdx(mom,dEdx,AliPID::kProton,fCutnSigdEdx)){
3420 fhnTracksJetVsPtAODpKch->Fill(ptLeadingAOD,nTracksLeadingAOD);
3421 fhRadiusJetVsPtAODpKch->Fill(ptLeadingAOD,radiusAOD);
3428 if(jetLeadingMC->Pt()>fJetPtCut && TMath::Abs(jetLeadingMC->Eta())<fJetEtaCut){
3430 ptLeadingMC = jetLeadingMC->Pt();
3431 etaLeadingMC = jetLeadingMC->Eta();
3432 phiLeadingMC = jetLeadingMC->Phi();
3433 nTracksLeadingMC = jetLeadingMC->GetRefTracks()->GetEntriesFast();
3435 Double_t radiusMC = GetJetRadius(jetLeadingMC,0.8);
3436 fhnTracksJetVsPtMC->Fill(ptLeadingMC,nTracksLeadingMC);
3437 fhRadiusJetVsPtMC->Fill(ptLeadingMC,radiusMC);
3438 if(IsQuarkHardScatteringEvent(pythiaPID)) fhnTracksJetVsPtMCquarkEv->Fill(ptLeadingMC,nTracksLeadingMC);
3442 fhnTracksJetVsPtMCK0->Fill(ptLeadingMC,nTracksLeadingMC);
3443 if(IsQuarkHardScatteringEvent(pythiaPID)) fhnTracksJetVsPtMCK0quarkEv->Fill(ptLeadingMC,nTracksLeadingMC);
3444 fhRadiusJetVsPtMCK0->Fill(ptLeadingMC,radiusMC);
3449 if(jetLeadingAOD && jetLeadingMC){
3451 //std::cout<<" comp: leading jetPt AOD "<<ptLeadingAOD<<" MC "<<ptLeadingMC<<std::endl;
3452 //if(jetLeadingAOD && jetLeadingMC)
3453 //std::cout<<" leading jet eta AOD "<<jetLeadingAOD->Eta()<<" MC "<<jetLeadingMC->Eta()<<std::endl;
3456 if(jetLeadingMC->Pt()>fJetPtCut && jetLeadingAOD->Pt()>fJetPtCut &&
3457 TMath::Abs(jetLeadingMC->Eta())<fJetEtaCut && TMath::Abs(jetLeadingAOD->Eta())<fJetEtaCut){
3459 fhLeadingPtAODvsMC->Fill(ptLeadingMC,ptLeadingAOD);
3460 fhLeadingEtaAODvsMC->Fill(etaLeadingMC,etaLeadingAOD);
3461 fhLeadingPhiAODvsMC->Fill(phiLeadingMC,phiLeadingAOD);
3462 fhnTracksLeadingAODvsMC->Fill(nTracksLeadingMC,nTracksLeadingAOD);
3464 TLorentzVector *mom4MC = jetLeadingMC->MomentumVector();
3465 TLorentzVector *mom4AOD = jetLeadingAOD->MomentumVector();
3467 Double_t dR = mom4MC->DeltaR(*mom4AOD);
3468 fhLeadingdRAODMC->Fill(dR);
3470 if(dR<0.4) fhLeadingPtAODvsMCdRcut->Fill(ptLeadingMC,ptLeadingAOD);
3472 Double_t dPt = ptLeadingAOD - ptLeadingMC;
3473 Double_t dnTracks = nTracksLeadingAOD - nTracksLeadingMC;
3475 fhdnTracksVsdPtLeadingAODMC->Fill(dPt,dnTracks);
3481 // ---------------------------------------------------------------------------
3483 // void AliAnalysisTaskJetChem::CheckK0(){
3485 // TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
3488 // AliInfo("no mcparticles branch");
3493 // Int_t ntrks = farray->GetEntries();
3495 // for(Int_t i=0; i<ntrks; i++){ // trk loop
3497 // AliAODMCParticle* mctrk = (AliAODMCParticle*)farray->At(i);
3499 // Int_t pdg = mctrk->GetPdgCode();
3501 // if(pdg != 310) continue; // K0
3503 // cout<<" MC K0: physPrim "<<mctrk->IsPhysicalPrimary()<<endl;
3505 // Int_t labelMother = mctrk->GetMother();
3507 // cout<<" found K00, label mother "<<labelMother<<endl;
3509 // AliAODMCParticle* mctrkMother = NULL;
3510 // Int_t pdgMother = -1;
3512 // for(Int_t k=0 ; k<ntrks; k++){
3514 // mctrkMother = (AliAODMCParticle*)farray->At(k);
3515 // if(mctrkMother->GetLabel() == labelMother) break;
3518 // pdgMother = mctrkMother->GetPdgCode();
3519 // cout<<" K0 mother pdg "<<pdgMother<<" GID "<<fpdgdb->ConvertPdgToGeant3(pdgMother)<<" isPrimary "<<mctrkMother->IsPrimary()<<endl;
3520 // //cout<<" mother name "<<mctrkMother->GetName()<<endl;
3525 // ---------------------------------------------------------------------------
3527 // void AliAnalysisTaskJetChem::CheckK0Stack(){
3529 // AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
3530 // if (!mcHandler) {
3531 // Printf("ERROR: Could not retrieve MC event handler");
3535 // AliMCEvent* mcEvent = mcHandler->MCEvent();
3537 // Printf("ERROR: Could not retrieve MC event");
3541 // AliStack* mcStack = mcEvent->Stack();//Load Stack
3542 // AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
3544 // Int_t nTracksMC = mcStack->GetNtrack();
3546 // Bool_t foundK0 = kFALSE;
3548 // for (Int_t iTracks = 0; iTracks < nTracksMC; iTracks++) {
3550 // if(!(mcStack->IsPhysicalPrimary(iTracks))) continue;
3552 // TParticle* mctrk = mcStack->Particle(iTracks);
3554 // Int_t pdg = mctrk->GetPdgCode();
3556 // if ((mctrk->Pt() < 2*fTrackPtCut) || (TMath::Abs(mctrk->Eta()) > fTrackEtaCut )) continue;
3558 // if(pdg == 310) foundK0 = kTRUE;
3560 // if(pdg != 310) continue; // K0 short
3562 // cout<<" check K0 "<<endl;
3564 // Int_t indexMother = -999;
3565 // TParticle* mctrkMother = mctrk;
3566 // TParticle* mctrkFirstGenMother = NULL; // last mother which is not the primary proton
3568 // Int_t nAncestors = 0;
3570 // while(indexMother != -1){
3572 // indexMother = mctrkMother->GetFirstMother();
3574 // if(indexMother != -1){
3575 // mctrkFirstGenMother = mctrkMother;
3576 // mctrkMother = mcStack->Particle(indexMother);
3580 // cout<<" nAncestors "<<nAncestors<<" pdg mother "<<mctrkMother->GetPdgCode()<<" name "<<mctrkMother->GetName()<<endl;
3582 // cout<<" pdg firstGenMother "<<mctrkFirstGenMother->GetPdgCode()<<" name "<<mctrkFirstGenMother->GetName()<<endl;
3586 // cout<<" pythiaGenHeader "<<pythiaGenHeader<<endl;
3587 // cout<<" Pythia Process type "<<pythiaGenHeader->ProcessType()<<" ptHard "<<pythiaGenHeader->GetPtHard()<<endl;
3590 // fhPythiaProcess->Fill(pythiaGenHeader->ProcessType());
3591 // if(foundK0) fhPythiaProcess_K0->Fill(pythiaGenHeader->ProcessType());
3593 // //Int_t indexGrandMother = mctrkMother->GetFirstMother();
3594 // //cout<<" indexGrandMother "<<indexGrandMother<<endl;
3596 // //if(indexGrandMother>-1){
3597 // // TParticle* mctrkGrandMother = mcStack->Particle(indexGrandMother);
3598 // // cout<<" pdg grandMother "<<mctrkGrandMother->GetPdgCode()<<" name "<<mctrkGrandMother->GetName()<<endl;
3604 // ---------------------------------------------------------------------------------
3606 Bool_t AliAnalysisTaskJetChem::IsQuarkHardScatteringEvent(const Int_t PID){
3608 // Pythia Manual sec. 8.2.1 :
3609 // if Pythia PID = 92,93,94 event is diffractive
3610 // if Pythia PID = 13, 28, 68 hard scattering products are gg, qg
3611 // if Pythia PID = 11, 12, 53 hard scattering products are qq, q\bar{q}
3613 if(PID == 92 || PID == 93 || PID == 94 || PID==95) return kFALSE;
3614 else if(PID == 13 || PID == 28 || PID == 68) return kFALSE;
3615 else if(PID == 11 || PID == 12 || PID == 53) return kTRUE;
3617 AliInfo(Form("unknown Pythia PID %d",PID));
3624 // ---------------------------------------------------------------------------------
3626 Bool_t AliAnalysisTaskJetChem::IsGluonHardScatteringEvent(const Int_t PID){
3628 // Pythia Manual sec. 8.2.1 :
3629 // if Pythia PID = 92,93,94 event is diffractive
3630 // if Pythia PID = 95: low pt event (MPI)
3631 // if Pythia PID = 13, 28, 68 hard scattering products are gg, qg
3632 // if Pythia PID = 11, 12, 53 hard scattering products are qq, q\bar{q}
3635 if(PID == 92 || PID == 93 || PID == 94 || PID == 95) return kFALSE;
3636 else if(PID == 13 || PID == 68) return kTRUE;
3637 else if(PID == 28) return kFALSE; // mixed gq final state
3638 else if(PID == 11 || PID == 12 || PID == 53) return kFALSE;
3640 AliInfo(Form("unknown Pythia PID %d",PID));
3646 // ---------------------------------------------------------------------------------
3648 Bool_t AliAnalysisTaskJetChem::IsDiffractiveEvent(const Int_t PID){
3650 // Pythia Manual sec. 8.2.1 :
3651 // if Pythia PID = 92,93,94 event is diffractive
3652 // if Pythia PID = 13, 28, 68 hard scattering products are gg, qg
3653 // if Pythia PID = 11, 12, 53 hard scattering products are qq, q\bar{q}
3655 if(PID == -1) return kFALSE;
3657 if(PID == 13 || PID == 28 || PID == 68) return kFALSE;
3658 else if(PID == 11 || PID == 12 || PID == 53) return kFALSE;
3659 else if(PID == 92 || PID == 93 || PID == 94 || PID==95) return kTRUE;
3661 AliInfo(Form("unknown Pythia PID %d",PID));
3669 // ----------------------------------------------------------------------------------
3671 Int_t AliAnalysisTaskJetChem::GetPythiaProcessID(){
3673 // Pythia PID for this event
3675 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
3677 //Printf("ERROR: Could not retrieve MC event handler");
3681 AliMCEvent* mcEvent = mcHandler->MCEvent();
3683 AliInfo("could not retrieve MC event");
3687 AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
3689 if (!pythiaGenHeader) {
3690 AliInfo("Could not retrieve pythiaEventHeader");
3694 Int_t pid = pythiaGenHeader->ProcessType();
3700 // ----------------------------------------------------------------------------------
3702 void AliAnalysisTaskJetChem::GetJetTracksResum(TList* list, AliAODJet* jet, const Double_t radius){
3705 if(!jet) return; // no jet in acc in event
3707 // list of AOD tracks in jet cone, using cone axis and distance axis-track (and not trackrefs)
3709 Int_t nTracks = fAOD->GetNTracks();
3711 if(!nTracks) return;
3714 jet->PxPyPz(jetMom);
3715 TVector3 jet3mom(jetMom);
3718 // standard cuts + ITS refit = stdrd PWG4 cut (compose them for productions with old ESDFilter task)
3720 for (Int_t itrack=0; itrack<nTracks; itrack++) {
3722 AliAODTrack* track = fAOD->GetTrack(itrack);
3724 UInt_t status = track->GetStatus();
3726 if(!track->TestFilterBit(fFilterBitJF)) continue; // track cut selection
3727 if(fRequireITSRefitJF && ((status&AliESDtrack::kITSrefit)==0)) continue;
3729 if(TMath::Abs(track->Eta()) > fTrackEtaCut) continue;
3730 if( track->Pt() < fTrackPtCutJF ) continue;
3732 Double_t trackMom[3];
3733 track->PxPyPz(trackMom);
3734 TVector3 track3mom(trackMom);
3736 Double_t dR = jet3mom.DeltaR(track3mom);
3746 // ----------------------------------------------------------------------------------
3748 void AliAnalysisTaskJetChem::GetJetTracksTrackrefs(TList* list, AliAODJet* jet){
3750 if(!jet) return; // no jet in acc in event
3752 // list of AOD tracks in jet cone, using trackrefs
3754 Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
3756 if(!nTracks) return;
3758 // standard cuts + ITS refit = stdrd PWG4 cut (compose them for productions with old ESDFilter task)
3760 for (Int_t itrack=0; itrack<nTracks; itrack++) {
3762 AliAODTrack* track = (AliAODTrack*) jet->GetRefTracks()->At(itrack);
3764 UInt_t status = track->GetStatus();
3766 if(!track->TestFilterBit(fFilterBitJF)) continue; // track cut selection
3767 if(fRequireITSRefitJF && ((status&AliESDtrack::kITSrefit)==0)) continue;
3769 if(TMath::Abs(track->Eta()) > fTrackEtaCut) continue;
3770 if( track->Pt() < fTrackPtCutJF ) continue;
3775 //cout<<" trackrefs Size "<<nTracks<<" acc track list size "<<list->GetEntries()<<endl;
3780 // ----------------------------------------------------------------------------------
3782 void AliAnalysisTaskJetChem::FillReferenceFF(AliAODJet* jet){
3787 TList* jetTracks = new TList(); // FIXME - avoid new/delete
3788 //GetJetTracksResum(jetTracks,jet,0.7);
3789 GetJetTracksTrackrefs(jetTracks,jet);
3791 Double_t jetpt = jet->Pt();
3793 TIter next(jetTracks);
3794 while(AliAODTrack* track = static_cast<AliAODTrack*>(next())){
3796 Double_t trackpt = track->Pt();
3797 Double_t z = trackpt/jetpt;
3798 Double_t xi = TMath::Log(1/z);
3800 //cout<<" trackpt "<<trackpt<<" jetpt "<<jetpt<<" z "<<z<<" xi "<<xi<<endl;
3802 if(5<jetpt && jetpt<=10){
3803 fhdNdptTracksJetPt5to10->Fill(trackpt);
3804 fhdNdzJets5to10->Fill(z);
3805 fhdNdxiJets5to10->Fill(xi);
3807 else if(10<jetpt && jetpt<=20){
3808 fhdNdptTracksJetPt10to20->Fill(trackpt);
3809 fhdNdzJets10to20->Fill(z);
3810 fhdNdxiJets10to20->Fill(xi);
3812 else if(20<jetpt && jetpt<=30){
3813 fhdNdptTracksJetPt20to30->Fill(trackpt);
3814 fhdNdzJets20to30->Fill(z);
3815 fhdNdxiJets20to30->Fill(xi);
3817 else if(30<jetpt && jetpt<=40){
3818 fhdNdptTracksJetPt30to40->Fill(trackpt);
3819 fhdNdzJets30to40->Fill(z);
3820 fhdNdxiJets30to40->Fill(xi);
3822 else if(40<jetpt && jetpt<=60){
3823 fhdNdptTracksJetPt40to60->Fill(trackpt);
3824 fhdNdzJets40to60->Fill(z);
3825 fhdNdxiJets40to60->Fill(xi);
3834 // -------------------------------------------------------------
3836 void AliAnalysisTaskJetChem::FillReferencePlotsTracks(){
3838 // eta/phi & pt tracks before/after cuts
3839 // track multiplicity / evt
3841 Int_t nTracks = fAOD->GetNTracks();
3842 Int_t countTracks = 0;
3844 // standard cuts + ITS refit = stdrd PWG4 cut (compose them for productions with old ESDFilter task)
3846 for(Int_t itrack=0; itrack<nTracks; itrack++){
3848 AliAODTrack* track = fAOD->GetTrack(itrack);
3850 Double_t trackPt = track->Pt();
3851 Double_t trackPhi = track->Phi();
3852 Double_t trackEta = track->Eta();
3854 fhPhiEtaTracksNoCut->Fill(trackEta,trackPhi);
3855 fhPtTracksNoCut->Fill(trackPt);
3857 UInt_t status = track->GetStatus();
3859 if(!track->TestFilterBit(fFilterBitJF)) continue; // track cut selection
3860 if(fRequireITSRefitJF && ((status&AliESDtrack::kITSrefit)==0)) continue;
3862 fhPhiEtaTracks->Fill(trackEta,trackPhi);
3863 fhPtTracks->Fill(trackPt);
3865 if(TMath::Abs(track->Eta()) > fTrackEtaCut) continue;
3866 if( track->Pt() < fTrackPtCutJF ) continue;
3872 fhTrackMult->Fill(countTracks);