1 /*************************************************************************
4 * Task for Jet Chemistry Analysis in PWG4 Jet Task Force Train *
5 * Analysis of K0s, Lambda and Antilambda with and without Jetevents *
7 *************************************************************************/
9 /**************************************************************************
10 * Copyright(c) 1998-2012, ALICE Experiment at CERN, All rights reserved. *
12 * Author: The ALICE Off-line Project. *
13 * Contributors are mentioned in the code where appropriate. *
15 * Permission to use, copy, modify and distribute this software and its *
16 * documentation strictly for non-commercial purposes is hereby grante *
18 * without fee, provided that the above copyright notice appears in all *
19 * copies and that both the copyright notice and this permission notice *
20 * appear in the supporting documentation. The authors make no claims *
21 * about the suitability of this software for any purpose. It is *
22 * provided "as is" without express or implied warranty. *
23 **************************************************************************/
24 //Task for K0s, Lambda and Antilambda analysis in jets
25 //Author: Alice Zimmermann (zimmermann@physi.uni-heidelberg.de)
43 #include "THnSparse.h"
46 #include "AliAnalysisHelperJetTasks.h"
47 #include "TDatabasePDG.h"
49 #include "AliAnalysisManager.h"
50 #include "AliAODHandler.h"
51 #include "AliAODInputHandler.h"
52 #include "AliESDEvent.h"
53 #include "AliGenPythiaEventHeader.h"
54 #include "AliGenHijingEventHeader.h"
55 #include "AliGenEventHeader.h"
56 #include "TLorentzVector.h"
57 #include "AliAODEvent.h"
58 #include "AliAODJet.h"
60 #include "AliAODTrack.h"
61 #include "AliCentrality.h"
62 #include "AliAnalysisTaskSE.h"
63 #include "AliESDtrack.h"
64 #include "AliESDtrackCuts.h"
65 #include "AliESDEvent.h"
66 #include "AliESDInputHandler.h"
68 #include "AliPIDResponse.h"
69 #include "AliAODPid.h"
70 #include "AliExternalTrackParam.h"
71 #include "AliAnalysisTaskJetChem.h"
72 #include "AliPhysicsSelection.h"
73 #include "AliBackgroundSelection.h"
74 #include "AliInputEventHandler.h"
75 #include "AliAODMCHeader.h"
76 #include "AliAODPid.h"
77 #include "AliVEvent.h"
78 #include "AliAODMCParticle.h"
82 ClassImp(AliAnalysisTaskJetChem)
84 //____________________________________________________________________________
85 AliAnalysisTaskJetChem::AliAnalysisTaskJetChem()
86 : AliAnalysisTaskFragmentationFunction()
100 ,fCutV0cosPointAngle(0)
107 ,fCutDcaV0Daughters(0)
108 ,fCutDcaPosToPrimVertex(0)
109 ,fCutDcaNegToPrimVertex(0)
119 ,fFFHistosRecCutsK0Evt(0)
120 //,fFFHistosIMK0AllEvt(0)
121 //,fFFHistosIMK0Jet(0)
122 //,fFFHistosIMK0Cone(0)
126 // ,fFFHistosIMLaAllEvt(0)
127 // ,fFFHistosIMLaJet(0)
128 //,fFFHistosIMLaCone(0)
132 ,fListFeeddownLaCand(0)
133 ,fListFeeddownALaCand(0)
139 ,fListMCgenK0sCone(0)
141 ,fListMCgenALaCone(0)
142 ,IsArmenterosSelected(0)
143 // ,fFFHistosIMALaAllEvt(0)
144 // ,fFFHistosIMALaJet(0)
145 // ,fFFHistosIMALaCone(0)
161 ,fFFIMLaNBinsJetPt(0)
192 //,fh1trackPosNCls(0)
193 //,fh1trackNegNCls(0)
203 ,fh2ProperLifetimeK0sVsPtBeforeCut(0)
204 ,fh2ProperLifetimeK0sVsPtAfterCut(0)
206 ,fh1DcaV0Daughters(0)
207 ,fh1DcaPosToPrimVertex(0)
208 ,fh1DcaNegToPrimVertex(0)
209 ,fh2ArmenterosBeforeCuts(0)
210 ,fh2ArmenterosAfterCuts(0)
213 ,fh1PosDaughterCharge(0)
214 ,fh1NegDaughterCharge(0)
221 ,fhnInvMassEtaTrackPtK0s(0)
222 ,fhnInvMassEtaTrackPtLa(0)
223 ,fhnInvMassEtaTrackPtALa(0)
232 ,fh2MCEtagenK0Cone(0)
233 ,fh2MCEtagenLaCone(0)
234 ,fh2MCEtagenALaCone(0)
237 ,fh1IMALaConeSmear(0)
238 ,fhnrecMCHijingLaIncl(0)
239 ,fhnrecMCHijingLaCone(0)
240 ,fhnrecMCHijingALaIncl(0)
241 ,fhnrecMCHijingALaCone(0)
242 ,fhnrecMCInjectLaIncl(0)
243 ,fhnrecMCInjectLaCone(0)
244 ,fhnrecMCInjectALaIncl(0)
245 ,fhnrecMCInjectALaCone(0)
249 ,fhnMCrecK0ConeSmear(0)
250 ,fhnMCrecLaConeSmear(0)
251 ,fhnMCrecALaConeSmear(0)
252 ,fhnK0sSecContinCone(0)
253 ,fhnLaSecContinCone(0)
254 ,fhnALaSecContinCone(0)
276 ,fh1MCMultiplicityPrimary(0)
277 ,fh1MCMultiplicityTracks(0)
280 ,fhnFeedDownLaCone(0)
281 ,fhnFeedDownALaCone(0)
282 ,fh1MCProdRadiusK0s(0)
283 ,fh1MCProdRadiusLambda(0)
284 ,fh1MCProdRadiusAntiLambda(0)
288 ,fh1MCPtAntiLambda(0)
296 //,fh1MCRapAntiLambda(0)
300 ,fh1MCEtaAntiLambda(0)
303 // default constructor
306 //__________________________________________________________________________________________
307 AliAnalysisTaskJetChem::AliAnalysisTaskJetChem(const char *name)
308 : AliAnalysisTaskFragmentationFunction(name)
322 ,fCutV0cosPointAngle(0)
329 ,fCutDcaV0Daughters(0)
330 ,fCutDcaPosToPrimVertex(0)
331 ,fCutDcaNegToPrimVertex(0)
341 ,fFFHistosRecCutsK0Evt(0)
342 //,fFFHistosIMK0AllEvt(0)
343 //,fFFHistosIMK0Jet(0)
344 //,fFFHistosIMK0Cone(0)
348 //,fFFHistosIMLaAllEvt(0)
349 //,fFFHistosIMLaJet(0)
350 //,fFFHistosIMLaCone(0)
354 ,fListFeeddownLaCand(0)
355 ,fListFeeddownALaCand(0)
361 ,fListMCgenK0sCone(0)
363 ,fListMCgenALaCone(0)
364 ,IsArmenterosSelected(0)
365 //,fFFHistosIMALaAllEvt(0)
366 //,fFFHistosIMALaJet(0)
367 // ,fFFHistosIMALaCone(0)
383 ,fFFIMLaNBinsJetPt(0)
414 // ,fh1trackPosNCls(0)
415 // ,fh1trackNegNCls(0)
425 ,fh2ProperLifetimeK0sVsPtBeforeCut(0)
426 ,fh2ProperLifetimeK0sVsPtAfterCut(0)
428 ,fh1DcaV0Daughters(0)
429 ,fh1DcaPosToPrimVertex(0)
430 ,fh1DcaNegToPrimVertex(0)
431 ,fh2ArmenterosBeforeCuts(0)
432 ,fh2ArmenterosAfterCuts(0)
435 ,fh1PosDaughterCharge(0)
436 ,fh1NegDaughterCharge(0)
443 ,fhnInvMassEtaTrackPtK0s(0)
444 ,fhnInvMassEtaTrackPtLa(0)
445 ,fhnInvMassEtaTrackPtALa(0)
454 ,fh2MCEtagenK0Cone(0)
455 ,fh2MCEtagenLaCone(0)
456 ,fh2MCEtagenALaCone(0)
459 ,fh1IMALaConeSmear(0)
460 ,fhnrecMCHijingLaIncl(0)
461 ,fhnrecMCHijingLaCone(0)
462 ,fhnrecMCHijingALaIncl(0)
463 ,fhnrecMCHijingALaCone(0)
464 ,fhnrecMCInjectLaIncl(0)
465 ,fhnrecMCInjectLaCone(0)
466 ,fhnrecMCInjectALaIncl(0)
467 ,fhnrecMCInjectALaCone(0)
471 ,fhnMCrecK0ConeSmear(0)
472 ,fhnMCrecLaConeSmear(0)
473 ,fhnMCrecALaConeSmear(0)
474 ,fhnK0sSecContinCone(0)
475 ,fhnLaSecContinCone(0)
476 ,fhnALaSecContinCone(0)
498 ,fh1MCMultiplicityPrimary(0)
499 ,fh1MCMultiplicityTracks(0)
502 ,fhnFeedDownLaCone(0)
503 ,fhnFeedDownALaCone(0)
504 ,fh1MCProdRadiusK0s(0)
505 ,fh1MCProdRadiusLambda(0)
506 ,fh1MCProdRadiusAntiLambda(0)
510 ,fh1MCPtAntiLambda(0)
518 //,fh1MCRapAntiLambda(0)
522 ,fh1MCEtaAntiLambda(0)
528 DefineOutput(1,TList::Class());
531 //__________________________________________________________________________________________________________________________
532 AliAnalysisTaskJetChem::AliAnalysisTaskJetChem(const AliAnalysisTaskJetChem ©)
533 : AliAnalysisTaskFragmentationFunction()
535 ,fRandom(copy.fRandom)
536 ,fAnalysisMC(copy.fAnalysisMC)
537 ,fDeltaVertexZ(copy.fDeltaVertexZ)
538 ,fCutjetEta(copy.fCutjetEta)
539 ,fCuttrackNegNcls(copy.fCuttrackNegNcls)
540 ,fCuttrackPosNcls(copy.fCuttrackPosNcls)
541 ,fCutPostrackRap(copy.fCutPostrackRap)
542 ,fCutNegtrackRap(copy.fCutNegtrackRap)
543 ,fCutRap(copy.fCutRap)
544 ,fCutPostrackEta(copy.fCutPostrackEta)
545 ,fCutNegtrackEta(copy.fCutNegtrackEta)
546 ,fCutEta(copy.fCutEta)
547 ,fCutV0cosPointAngle(copy.fCutV0cosPointAngle)
548 ,fKinkDaughters(copy.fKinkDaughters)
549 ,fRequireTPCRefit(copy.fRequireTPCRefit)
550 ,fCutArmenteros(copy.fCutArmenteros)
551 ,fCutV0DecayMin(copy.fCutV0DecayMin)
552 ,fCutV0DecayMax(copy.fCutV0DecayMax)
553 ,fCutV0totMom(copy.fCutV0totMom)
554 ,fCutDcaV0Daughters(copy.fCutDcaV0Daughters)
555 ,fCutDcaPosToPrimVertex(copy.fCutDcaPosToPrimVertex)
556 ,fCutDcaNegToPrimVertex(copy.fCutDcaNegToPrimVertex)
557 ,fCutV0RadiusMin(copy.fCutV0RadiusMin)
558 ,fCutV0RadiusMax(copy.fCutV0RadiusMax)
559 ,fCutBetheBloch(copy.fCutBetheBloch)
560 ,fCutRatio(copy.fCutRatio)
561 ,fK0Type(copy.fK0Type)
562 ,fFilterMaskK0(copy.fFilterMaskK0)
563 ,fListK0s(copy.fListK0s)
564 ,fPIDResponse(copy.fPIDResponse)
565 ,fV0QAK0(copy.fV0QAK0)
566 ,fFFHistosRecCutsK0Evt(copy.fFFHistosRecCutsK0Evt)
567 //,fFFHistosIMK0AllEvt(copy.fFFHistosIMK0AllEvt)
568 //,fFFHistosIMK0Jet(copy.fFFHistosIMK0Jet)
569 //,fFFHistosIMK0Cone(copy.fFFHistosIMK0Cone)
570 ,fLaType(copy.fLaType)
571 ,fFilterMaskLa(copy.fFilterMaskLa)
572 ,fListLa(copy.fListLa)
573 //,fFFHistosIMLaAllEvt(copy.fFFHistosIMLaAllEvt)
574 //,fFFHistosIMLaJet(copy.fFFHistosIMLaJet)
575 //,fFFHistosIMLaCone(copy.fFFHistosIMLaCone)
576 ,fALaType(copy.fALaType)
577 ,fFilterMaskALa(copy.fFilterMaskALa)
578 ,fListALa(copy.fListALa)
579 ,fListFeeddownLaCand(copy.fListFeeddownLaCand)
580 ,fListFeeddownALaCand(copy.fListFeeddownALaCand)
581 ,jetConeFDLalist(copy.jetConeFDLalist)
582 ,jetConeFDALalist(copy.jetConeFDALalist)
583 ,fListMCgenK0s(copy.fListMCgenK0s)
584 ,fListMCgenLa(copy.fListMCgenLa)
585 ,fListMCgenALa(copy.fListMCgenALa)
586 ,fListMCgenK0sCone(copy.fListMCgenK0sCone)
587 ,fListMCgenLaCone(copy.fListMCgenLaCone)
588 ,fListMCgenALaCone(copy.fListMCgenALaCone)
589 ,IsArmenterosSelected(copy.IsArmenterosSelected)
590 //,fFFHistosIMALaAllEvt(copy.fFFHistosIMALaAllEvt)
591 //,fFFHistosIMALaJet(copy.fFFHistosIMALaJet)
592 //,fFFHistosIMALaCone(copy.fFFHistosIMALaCone)
593 ,fFFIMNBinsJetPt(copy.fFFIMNBinsJetPt)
594 ,fFFIMJetPtMin(copy.fFFIMJetPtMin)
595 ,fFFIMJetPtMax(copy.fFFIMJetPtMax)
596 ,fFFIMNBinsInvM(copy.fFFIMNBinsInvM)
597 ,fFFIMInvMMin(copy.fFFIMInvMMin)
598 ,fFFIMInvMMax(copy.fFFIMInvMMax)
599 ,fFFIMNBinsPt(copy.fFFIMNBinsPt)
600 ,fFFIMPtMin(copy.fFFIMPtMin)
601 ,fFFIMPtMax(copy.fFFIMPtMax)
602 ,fFFIMNBinsXi(copy.fFFIMNBinsXi)
603 ,fFFIMXiMin(copy.fFFIMXiMin)
604 ,fFFIMXiMax(copy.fFFIMXiMax)
605 ,fFFIMNBinsZ(copy.fFFIMNBinsZ)
606 ,fFFIMZMin(copy.fFFIMZMin)
607 ,fFFIMZMax(copy.fFFIMZMax)
608 ,fFFIMLaNBinsJetPt(copy.fFFIMLaNBinsJetPt)
609 ,fFFIMLaJetPtMin(copy.fFFIMLaJetPtMin)
610 ,fFFIMLaJetPtMax(copy.fFFIMLaJetPtMax)
611 ,fFFIMLaNBinsInvM(copy.fFFIMLaNBinsInvM)
612 ,fFFIMLaInvMMin(copy.fFFIMLaInvMMin)
613 ,fFFIMLaInvMMax(copy.fFFIMLaInvMMax)
614 ,fFFIMLaNBinsPt(copy.fFFIMLaNBinsPt)
615 ,fFFIMLaPtMin(copy.fFFIMLaPtMin)
616 ,fFFIMLaPtMax(copy.fFFIMLaPtMax)
617 ,fFFIMLaNBinsXi(copy.fFFIMLaNBinsXi)
618 ,fFFIMLaXiMin(copy.fFFIMLaXiMin)
619 ,fFFIMLaXiMax(copy.fFFIMLaXiMax)
620 ,fFFIMLaNBinsZ(copy.fFFIMLaNBinsZ)
621 ,fFFIMLaZMin(copy.fFFIMLaZMin)
622 ,fFFIMLaZMax(copy.fFFIMLaZMax)
623 ,fh1EvtAllCent(copy.fh1EvtAllCent)
625 ,fh1K0Mult(copy.fh1K0Mult)
626 ,fh1dPhiJetK0(copy.fh1dPhiJetK0)
627 ,fh1LaMult(copy.fh1LaMult)
628 ,fh1dPhiJetLa(copy.fh1dPhiJetLa)
629 ,fh1ALaMult(copy.fh1ALaMult)
630 ,fh1dPhiJetALa(copy.fh1dPhiJetALa)
631 ,fh1JetEta(copy.fh1JetEta)
632 ,fh1JetPhi(copy.fh1JetPhi)
633 ,fh2JetEtaPhi(copy.fh2JetEtaPhi)
634 //,fh1V0JetPt(copy.fh1V0JetPt)
635 ,fh1IMK0Cone(copy.fh1IMK0Cone)
636 ,fh1IMLaCone(copy.fh1IMLaCone)
637 ,fh1IMALaCone(copy.fh1IMALaCone)
638 ,fh2FFJetTrackEta(copy.fh2FFJetTrackEta)
639 //,fh1trackPosNCls(copy.fh1trackPosNCls)
640 //,fh1trackNegNCls(copy.fh1trackNegNCls)
641 ,fh1trackPosRap(copy.fh1trackPosRap)
642 ,fh1trackNegRap(copy.fh1trackNegRap)
643 //,fh1V0Rap(copy.fh1V0Rap)
644 ,fh1trackPosEta(copy.fh1trackPosEta)
645 ,fh1trackNegEta(copy.fh1trackNegEta)
646 ,fh1V0Eta(copy.fh1V0Eta)
647 //,fh1V0totMom(copy.fh1V0totMom)
648 ,fh1CosPointAngle(copy.fh1CosPointAngle)
649 ,fh1DecayLengthV0(copy.fh1DecayLengthV0)
650 ,fh2ProperLifetimeK0sVsPtBeforeCut(copy.fh2ProperLifetimeK0sVsPtBeforeCut)
651 ,fh2ProperLifetimeK0sVsPtAfterCut(copy.fh2ProperLifetimeK0sVsPtAfterCut)
652 ,fh1V0Radius(copy.fh1V0Radius)
653 ,fh1DcaV0Daughters(copy.fh1DcaV0Daughters)
654 ,fh1DcaPosToPrimVertex(copy.fh1DcaPosToPrimVertex)
655 ,fh1DcaNegToPrimVertex(copy.fh1DcaNegToPrimVertex)
656 ,fh2ArmenterosBeforeCuts(copy.fh2ArmenterosBeforeCuts)
657 ,fh2ArmenterosAfterCuts(copy.fh2ArmenterosAfterCuts)
658 ,fh2BBLaPos(copy.fh2BBLaPos)
659 ,fh2BBLaNeg(copy.fh2BBLaPos)
660 ,fh1PosDaughterCharge(copy.fh1PosDaughterCharge)
661 ,fh1NegDaughterCharge(copy.fh1NegDaughterCharge)
662 ,fh1PtMCK0s(copy.fh1PtMCK0s)
663 ,fh1PtMCLa(copy.fh1PtMCLa)
664 ,fh1PtMCALa(copy.fh1PtMCALa)
665 ,fh1EtaK0s(copy.fh1EtaK0s)
666 ,fh1EtaLa(copy.fh1EtaLa)
667 ,fh1EtaALa(copy.fh1EtaALa)
668 ,fhnInvMassEtaTrackPtK0s(copy.fhnInvMassEtaTrackPtK0s)
669 ,fhnInvMassEtaTrackPtLa(copy.fhnInvMassEtaTrackPtLa)
670 ,fhnInvMassEtaTrackPtALa(copy.fhnInvMassEtaTrackPtALa)
671 ,fh1TrackMultCone(copy.fh1TrackMultCone)
672 ,fh2TrackMultCone(copy.fh2TrackMultCone)
673 ,fhnNJK0(copy.fhnNJK0)
674 ,fhnNJLa(copy.fhnNJLa)
675 ,fhnNJALa(copy.fhnNJALa)
676 ,fh2MCgenK0Cone(copy.fh2MCgenK0Cone)
677 ,fh2MCgenLaCone(copy.fh2MCgenLaCone)
678 ,fh2MCgenALaCone(copy.fh2MCgenALaCone)
679 ,fh2MCEtagenK0Cone(copy.fh2MCEtagenK0Cone)
680 ,fh2MCEtagenLaCone(copy.fh2MCEtagenLaCone)
681 ,fh2MCEtagenALaCone(copy.fh2MCEtagenALaCone)
682 ,fh1IMK0ConeSmear(copy.fh1IMK0ConeSmear)
683 ,fh1IMLaConeSmear(copy.fh1IMLaConeSmear)
684 ,fh1IMALaConeSmear(copy.fh1IMALaConeSmear)
685 ,fhnrecMCHijingLaIncl(copy.fhnrecMCHijingLaIncl)
686 ,fhnrecMCHijingLaCone(copy.fhnrecMCHijingLaCone)
687 ,fhnrecMCHijingALaIncl(copy.fhnrecMCHijingALaIncl)
688 ,fhnrecMCHijingALaCone(copy.fhnrecMCHijingALaCone)
689 ,fhnrecMCInjectLaIncl(copy.fhnrecMCInjectLaIncl)
690 ,fhnrecMCInjectLaCone(copy.fhnrecMCInjectLaCone)
691 ,fhnrecMCInjectALaIncl(copy.fhnrecMCInjectALaIncl)
692 ,fhnrecMCInjectALaCone(copy.fhnrecMCInjectALaCone)
693 ,fhnMCrecK0Cone(copy.fhnMCrecK0Cone)
694 ,fhnMCrecLaCone(copy.fhnMCrecLaCone)
695 ,fhnMCrecALaCone(copy.fhnMCrecALaCone)
696 ,fhnMCrecK0ConeSmear(copy.fhnMCrecK0ConeSmear)
697 ,fhnMCrecLaConeSmear(copy.fhnMCrecLaConeSmear)
698 ,fhnMCrecALaConeSmear(copy.fhnMCrecALaConeSmear)
699 ,fhnK0sSecContinCone(copy.fhnK0sSecContinCone)
700 ,fhnLaSecContinCone(copy.fhnLaSecContinCone)
701 ,fhnALaSecContinCone(copy.fhnALaSecContinCone)
702 ,fhnK0sIncl(copy.fhnK0sIncl)
703 ,fhnK0sCone(copy.fhnK0sCone)
704 ,fhnLaIncl(copy.fhnLaIncl)
705 ,fhnLaCone(copy.fhnLaCone)
706 ,fhnALaIncl(copy.fhnALaIncl)
707 ,fhnALaCone(copy.fhnALaCone)
708 ,fhnK0sPC(copy.fhnK0sPC)
709 ,fhnLaPC(copy.fhnLaPC)
710 ,fhnALaPC(copy.fhnALaPC)
711 ,fhnK0sMCC(copy.fhnK0sMCC)
712 ,fhnLaMCC(copy.fhnLaMCC)
713 ,fhnALaMCC(copy.fhnALaMCC)
714 ,fhnK0sRC(copy.fhnK0sRC)
715 ,fhnLaRC(copy.fhnLaRC)
716 ,fhnALaRC(copy.fhnALaRC)
717 ,fhnK0sOC(copy.fhnK0sOC)
718 ,fhnLaOC(copy.fhnLaOC)
719 ,fhnALaOC(copy.fhnALaOC)
720 ,fh1AreaExcluded(copy.fh1AreaExcluded)
721 ,fh1MedianEta(copy.fh1MedianEta)
722 ,fh1JetPtMedian(copy.fh1JetPtMedian)
723 ,fh1MCMultiplicityPrimary(copy.fh1MCMultiplicityPrimary)
724 ,fh1MCMultiplicityTracks(copy.fh1MCMultiplicityTracks)
725 ,fhnFeedDownLa(copy.fhnFeedDownLa)
726 ,fhnFeedDownALa(copy.fhnFeedDownALa)
727 ,fhnFeedDownLaCone(copy.fhnFeedDownLaCone)
728 ,fhnFeedDownALaCone(copy.fhnFeedDownALaCone)
729 ,fh1MCProdRadiusK0s(copy.fh1MCProdRadiusK0s)
730 ,fh1MCProdRadiusLambda(copy.fh1MCProdRadiusLambda)
731 ,fh1MCProdRadiusAntiLambda(copy.fh1MCProdRadiusAntiLambda)
732 ,fh1MCPtV0s(copy.fh1MCPtV0s)
733 ,fh1MCPtK0s(copy.fh1MCPtK0s)
734 ,fh1MCPtLambda(copy.fh1MCPtLambda)
735 ,fh1MCPtAntiLambda(copy.fh1MCPtAntiLambda)
736 ,fh1MCXiPt(copy.fh1MCXiPt)
737 ,fh1MCXibarPt(copy.fh1MCXibarPt)
738 ,fh2MCEtaVsPtK0s(copy.fh2MCEtaVsPtK0s)
739 ,fh2MCEtaVsPtLa(copy.fh2MCEtaVsPtLa)
740 ,fh2MCEtaVsPtALa(copy.fh2MCEtaVsPtALa)
741 //,fh1MCRapK0s(copy.fh1MCRapK0s)
742 //,fh1MCRapLambda(copy.fh1MCRapLambda)
743 //,fh1MCRapAntiLambda(copy.fh1MCRapAntiLambda)
744 ,fh1MCEtaAllK0s(copy.fh1MCEtaAllK0s)
745 ,fh1MCEtaK0s(copy.fh1MCEtaK0s)
746 ,fh1MCEtaLambda(copy.fh1MCEtaLambda)
747 ,fh1MCEtaAntiLambda(copy.fh1MCEtaAntiLambda)
754 // _________________________________________________________________________________________________________________________________
755 AliAnalysisTaskJetChem& AliAnalysisTaskJetChem::operator=(const AliAnalysisTaskJetChem& o)
760 AliAnalysisTaskFragmentationFunction::operator=(o);
763 fAnalysisMC = o.fAnalysisMC;
764 fDeltaVertexZ = o.fDeltaVertexZ;
765 fCutjetEta = o.fCutjetEta;
766 fCuttrackNegNcls = o.fCuttrackNegNcls;
767 fCuttrackPosNcls = o.fCuttrackPosNcls;
768 fCutPostrackRap = o.fCutPostrackRap;
769 fCutNegtrackRap = o.fCutNegtrackRap;
771 fCutPostrackEta = o.fCutPostrackEta;
772 fCutNegtrackEta = o.fCutNegtrackEta;
774 fCutV0cosPointAngle = o.fCutV0cosPointAngle;
775 fKinkDaughters = o.fKinkDaughters;
776 fRequireTPCRefit = o.fRequireTPCRefit;
777 fCutArmenteros = o.fCutArmenteros;
778 fCutV0DecayMin = o.fCutV0DecayMin;
779 fCutV0DecayMax = o.fCutV0DecayMax;
780 fCutV0totMom = o.fCutV0totMom;
781 fCutDcaV0Daughters = o.fCutDcaV0Daughters;
782 fCutDcaPosToPrimVertex = o.fCutDcaPosToPrimVertex;
783 fCutDcaNegToPrimVertex = o.fCutDcaNegToPrimVertex;
784 fCutV0RadiusMin = o.fCutV0RadiusMin;
785 fCutV0RadiusMax = o.fCutV0RadiusMax;
786 fCutBetheBloch = o.fCutBetheBloch;
787 fCutRatio = o.fCutRatio;
789 fFilterMaskK0 = o.fFilterMaskK0;
790 fListK0s = o.fListK0s;
791 fPIDResponse = o.fPIDResponse;
793 fFFHistosRecCutsK0Evt = o.fFFHistosRecCutsK0Evt;
794 //fFFHistosIMK0AllEvt = o.fFFHistosIMK0AllEvt;
795 //fFFHistosIMK0Jet = o.fFFHistosIMK0Jet;
796 //fFFHistosIMK0Cone = o.fFFHistosIMK0Cone;
798 fFilterMaskLa = o.fFilterMaskLa;
800 //fFFHistosIMLaAllEvt = o.fFFHistosIMLaAllEvt;
801 //fFFHistosIMLaJet = o.fFFHistosIMLaJet;
802 //fFFHistosIMLaCone = o.fFFHistosIMLaCone;
803 fALaType = o.fALaType;
804 fFilterMaskALa = o.fFilterMaskALa;
805 fListFeeddownLaCand = o.fListFeeddownLaCand;
806 fListFeeddownALaCand = o.fListFeeddownALaCand;
807 jetConeFDLalist = o.jetConeFDLalist;
808 jetConeFDALalist = o.jetConeFDALalist;
809 fListMCgenK0s = o.fListMCgenK0s;
810 fListMCgenLa = o.fListMCgenLa;
811 fListMCgenALa = o.fListMCgenALa;
812 fListMCgenK0sCone = o.fListMCgenK0sCone;
813 fListMCgenLaCone = o.fListMCgenLaCone;
814 fListMCgenALaCone = o.fListMCgenALaCone;
815 IsArmenterosSelected = o.IsArmenterosSelected;
816 // fFFHistosIMALaAllEvt = o.fFFHistosIMALaAllEvt;
817 // fFFHistosIMALaJet = o.fFFHistosIMALaJet;
818 // fFFHistosIMALaCone = o.fFFHistosIMALaCone;
819 fFFIMNBinsJetPt = o.fFFIMNBinsJetPt;
820 fFFIMJetPtMin = o.fFFIMJetPtMin;
821 fFFIMJetPtMax = o.fFFIMJetPtMax;
822 fFFIMNBinsPt = o.fFFIMNBinsPt;
823 fFFIMPtMin = o.fFFIMPtMin;
824 fFFIMPtMax = o.fFFIMPtMax;
825 fFFIMNBinsXi = o.fFFIMNBinsXi;
826 fFFIMXiMin = o.fFFIMXiMin;
827 fFFIMXiMax = o.fFFIMXiMax;
828 fFFIMNBinsZ = o.fFFIMNBinsZ;
829 fFFIMZMin = o.fFFIMZMin;
830 fFFIMZMax = o.fFFIMZMax;
831 fFFIMLaNBinsJetPt = o.fFFIMLaNBinsJetPt;
832 fFFIMLaJetPtMin = o.fFFIMLaJetPtMin;
833 fFFIMLaJetPtMax = o.fFFIMLaJetPtMax;
834 fFFIMLaNBinsPt = o.fFFIMLaNBinsPt;
835 fFFIMLaPtMin = o.fFFIMLaPtMin;
836 fFFIMLaPtMax = o.fFFIMLaPtMax;
837 fFFIMLaNBinsXi = o.fFFIMLaNBinsXi;
838 fFFIMLaXiMin = o.fFFIMLaXiMin;
839 fFFIMLaXiMax = o.fFFIMLaXiMax;
840 fFFIMLaNBinsZ = o.fFFIMLaNBinsZ;
841 fFFIMLaZMin = o.fFFIMLaZMin;
842 fFFIMLaZMax = o.fFFIMLaZMax;
843 fh1EvtAllCent = o.fh1EvtAllCent;
845 fh1K0Mult = o.fh1K0Mult;
846 fh1dPhiJetK0 = o.fh1dPhiJetK0;
847 fh1LaMult = o.fh1LaMult;
848 fh1dPhiJetLa = o.fh1dPhiJetLa;
849 fh1ALaMult = o.fh1ALaMult;
850 fh1dPhiJetALa = o.fh1dPhiJetALa;
851 fh1JetEta = o.fh1JetEta;
852 fh1JetPhi = o.fh1JetPhi;
853 fh2JetEtaPhi = o.fh2JetEtaPhi;
854 //fh1V0JetPt = o.fh1V0JetPt;
855 fh1IMK0Cone = o.fh1IMK0Cone;
856 fh1IMLaCone = o.fh1IMLaCone;
857 fh1IMALaCone = o.fh1IMALaCone;
858 fh2FFJetTrackEta = o.fh2FFJetTrackEta;
859 //fh1trackPosNCls = o.fh1trackPosNCls;
860 //fh1trackNegNCls = o.fh1trackNegNCls;
861 fh1trackPosRap = o.fh1trackPosRap;
862 fh1trackNegRap = o.fh1trackNegRap;
863 //fh1V0Rap = o.fh1V0Rap;
864 fh1trackPosEta = o.fh1trackPosEta;
865 fh1trackNegEta = o.fh1trackNegEta;
866 fh1V0Eta = o.fh1V0Eta;
867 // fh1V0totMom = o.fh1V0totMom;
868 fh1CosPointAngle = o.fh1CosPointAngle;
869 fh1DecayLengthV0 = o.fh1DecayLengthV0;
870 fh2ProperLifetimeK0sVsPtBeforeCut = o.fh2ProperLifetimeK0sVsPtBeforeCut;
871 fh2ProperLifetimeK0sVsPtAfterCut= o.fh2ProperLifetimeK0sVsPtAfterCut;
872 fh1V0Radius = o.fh1V0Radius;
873 fh1DcaV0Daughters = o.fh1DcaV0Daughters;
874 fh1DcaPosToPrimVertex = o.fh1DcaPosToPrimVertex;
875 fh1DcaNegToPrimVertex = o.fh1DcaNegToPrimVertex;
876 fh2ArmenterosBeforeCuts = o.fh2ArmenterosBeforeCuts;
877 fh2ArmenterosAfterCuts = o.fh2ArmenterosAfterCuts;
878 fh2BBLaPos = o.fh2BBLaPos;
879 fh2BBLaNeg = o.fh2BBLaPos;
880 fh1PosDaughterCharge = o.fh1PosDaughterCharge;
881 fh1NegDaughterCharge = o.fh1NegDaughterCharge;
882 fh1PtMCK0s = o.fh1PtMCK0s;
883 fh1PtMCLa = o.fh1PtMCLa;
884 fh1PtMCALa = o.fh1PtMCALa;
885 fh1EtaK0s = o.fh1EtaK0s;
886 fh1EtaLa = o.fh1EtaLa;
887 fh1EtaALa = o.fh1EtaALa;
888 fhnInvMassEtaTrackPtK0s = o.fhnInvMassEtaTrackPtK0s;
889 fhnInvMassEtaTrackPtLa = o.fhnInvMassEtaTrackPtLa;
890 fhnInvMassEtaTrackPtALa = o.fhnInvMassEtaTrackPtALa;
891 fh1TrackMultCone = o.fh1TrackMultCone;
892 fh2TrackMultCone = o.fh2TrackMultCone;
895 fhnNJALa = o.fhnNJALa;
896 fh2MCgenK0Cone = o.fh2MCgenK0Cone;
897 fh2MCgenLaCone = o.fh2MCgenLaCone;
898 fh2MCgenALaCone = o.fh2MCgenALaCone;
899 fh2MCEtagenK0Cone = o.fh2MCEtagenK0Cone;
900 fh2MCEtagenLaCone = o.fh2MCEtagenLaCone;
901 fh2MCEtagenALaCone = o.fh2MCEtagenALaCone;
902 fh1IMK0ConeSmear = o.fh1IMK0ConeSmear;
903 fh1IMLaConeSmear = o.fh1IMLaConeSmear;
904 fh1IMALaConeSmear = o.fh1IMALaConeSmear;
905 fhnrecMCHijingLaIncl = o.fhnrecMCHijingLaIncl;
906 fhnrecMCHijingLaCone = o.fhnrecMCHijingLaCone;
907 fhnrecMCHijingALaIncl = o.fhnrecMCHijingALaIncl;
908 fhnrecMCHijingALaCone = o.fhnrecMCHijingALaCone;
909 fhnrecMCInjectLaIncl = o.fhnrecMCInjectLaIncl;
910 fhnrecMCInjectLaCone = o.fhnrecMCInjectLaCone;
911 fhnrecMCInjectALaIncl = o.fhnrecMCInjectALaIncl;
912 fhnrecMCInjectALaCone = o.fhnrecMCInjectALaCone;
913 fhnMCrecK0Cone = o.fhnMCrecK0Cone;
914 fhnMCrecLaCone = o.fhnMCrecLaCone;
915 fhnMCrecALaCone = o.fhnMCrecALaCone;
916 fhnMCrecK0ConeSmear = o.fhnMCrecK0ConeSmear;
917 fhnMCrecLaConeSmear = o.fhnMCrecLaConeSmear;
918 fhnMCrecALaConeSmear = o.fhnMCrecALaConeSmear;
919 fhnK0sSecContinCone = o.fhnK0sSecContinCone;
920 fhnLaSecContinCone = o.fhnLaSecContinCone;
921 fhnALaSecContinCone = o.fhnALaSecContinCone;
922 fhnK0sIncl = o.fhnK0sIncl;
923 fhnK0sCone = o.fhnK0sCone;
924 fhnLaIncl = o.fhnLaIncl;
925 fhnLaCone = o.fhnLaCone;
926 fhnALaIncl = o.fhnALaIncl;
927 fhnALaCone = o.fhnALaCone;
928 fhnK0sPC = o.fhnK0sPC;
930 fhnALaPC = o.fhnALaPC;
931 fhnK0sRC = o.fhnK0sRC;
933 fhnALaRC = o.fhnALaRC;
934 fhnK0sOC = o.fhnK0sOC;
936 fhnALaOC = o.fhnALaOC;
937 fh1AreaExcluded = o.fh1AreaExcluded;
938 fh1MedianEta = o.fh1MedianEta;
939 fh1JetPtMedian = o.fh1JetPtMedian;
940 fh1MCMultiplicityPrimary = o.fh1MCMultiplicityPrimary;
941 fh1MCMultiplicityTracks = o.fh1MCMultiplicityTracks;
942 fhnFeedDownLa = o.fhnFeedDownLa;
943 fhnFeedDownALa = o.fhnFeedDownALa;
944 fhnFeedDownLaCone = o.fhnFeedDownLaCone;
945 fhnFeedDownALaCone = o.fhnFeedDownALaCone;
946 fh1MCProdRadiusK0s = o.fh1MCProdRadiusK0s;
947 fh1MCProdRadiusLambda = o.fh1MCProdRadiusLambda;
948 fh1MCProdRadiusAntiLambda = o.fh1MCProdRadiusAntiLambda;
949 fh1MCPtV0s = o.fh1MCPtV0s;
950 fh1MCPtK0s = o.fh1MCPtK0s;
951 fh1MCPtLambda = o.fh1MCPtLambda;
952 fh1MCPtAntiLambda = o.fh1MCPtAntiLambda;
953 fh1MCXiPt = o.fh1MCXiPt;
954 fh1MCXibarPt = o.fh1MCXibarPt;
955 fh2MCEtaVsPtK0s = o.fh2MCEtaVsPtK0s;
956 fh2MCEtaVsPtLa = o.fh2MCEtaVsPtLa;
957 fh2MCEtaVsPtALa = o.fh2MCEtaVsPtALa;
958 //fh1MCRapK0s = o.fh1MCRapK0s;
959 //fh1MCRapLambda = o.fh1MCRapLambda;
960 //fh1MCRapAntiLambda = o.fh1MCRapAntiLambda;
961 fh1MCEtaAllK0s = o.fh1MCEtaAllK0s;
962 fh1MCEtaK0s = o.fh1MCEtaK0s;
963 fh1MCEtaLambda = o.fh1MCEtaLambda;
964 fh1MCEtaAntiLambda = o.fh1MCEtaAntiLambda;
970 //_______________________________________________
971 AliAnalysisTaskJetChem::~AliAnalysisTaskJetChem()
976 if(fListK0s) delete fListK0s;
977 if(fListLa) delete fListLa;
978 if(fListALa) delete fListALa;
979 if(fListFeeddownLaCand) delete fListFeeddownLaCand;
980 if(fListFeeddownALaCand) delete fListFeeddownALaCand;
981 if(jetConeFDLalist) delete jetConeFDLalist;
982 if(jetConeFDALalist) delete jetConeFDALalist;
983 if(fListMCgenK0s) delete fListMCgenK0s;
984 if(fListMCgenLa) delete fListMCgenLa;
985 if(fListMCgenALa) delete fListMCgenALa;
986 if(fRandom) delete fRandom;
989 //________________________________________________________________________________________________________________________________
990 AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::AliFragFuncHistosInvMass(const char* name,
991 Int_t nJetPt, Float_t jetPtMin, Float_t jetPtMax,
992 Int_t nInvMass, Float_t invMassMin, Float_t invMassMax,
993 Int_t nPt, Float_t ptMin, Float_t ptMax,
994 Int_t nXi, Float_t xiMin, Float_t xiMax,
995 Int_t nZ , Float_t zMin , Float_t zMax )
1000 ,fNBinsInvMass(nInvMass)
1001 ,fInvMassMin(invMassMin)
1002 ,fInvMassMax(invMassMax)
1018 // default constructor
1022 //______________________________________________________________________________________________________________
1023 AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::AliFragFuncHistosInvMass(const AliFragFuncHistosInvMass& copy)
1025 ,fNBinsJetPt(copy.fNBinsJetPt)
1026 ,fJetPtMin(copy.fJetPtMin)
1027 ,fJetPtMax(copy.fJetPtMax)
1028 ,fNBinsInvMass(copy.fNBinsInvMass)
1029 ,fInvMassMin(copy.fInvMassMin)
1030 ,fInvMassMax(copy.fInvMassMax)
1031 ,fNBinsPt(copy.fNBinsPt)
1032 ,fPtMin(copy.fPtMin)
1033 ,fPtMax(copy.fPtMax)
1034 ,fNBinsXi(copy.fNBinsXi)
1035 ,fXiMin(copy.fXiMin)
1036 ,fXiMax(copy.fXiMax)
1037 ,fNBinsZ(copy.fNBinsZ)
1040 ,fh3TrackPt(copy.fh3TrackPt)
1043 ,fh1JetPt(copy.fh1JetPt)
1044 ,fNameFF(copy.fNameFF)
1049 //______________________________________________________________________________________________________________________________________________________________________
1050 AliAnalysisTaskJetChem::AliFragFuncHistosInvMass& AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::operator=(const AliAnalysisTaskJetChem::AliFragFuncHistosInvMass& o)
1055 TObject::operator=(o);
1056 fNBinsJetPt = o.fNBinsJetPt;
1057 fJetPtMin = o.fJetPtMin;
1058 fJetPtMax = o.fJetPtMax;
1059 fNBinsInvMass = o.fNBinsInvMass;
1060 fInvMassMin = o.fInvMassMin;
1061 fInvMassMax = o.fInvMassMax;
1062 fNBinsPt = o.fNBinsPt;
1065 fNBinsXi = o.fNBinsXi;
1068 fNBinsZ = o.fNBinsZ;
1071 fh3TrackPt = o.fh3TrackPt;
1074 fh1JetPt = o.fh1JetPt;
1075 fNameFF = o.fNameFF;
1081 //___________________________________________________________________________
1082 AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::~AliFragFuncHistosInvMass()
1086 if(fh1JetPt) delete fh1JetPt;
1087 if(fh3TrackPt) delete fh3TrackPt;
1088 if(fh3Xi) delete fh3Xi;
1089 if(fh3Z) delete fh3Z;
1092 //_________________________________________________________________
1093 void AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::DefineHistos()
1097 fh1JetPt = new TH1F(Form("fh1FFJetPtIM%s", fNameFF.Data()),"",fNBinsJetPt,fJetPtMin,fJetPtMax);
1098 fh3TrackPt = new TH3F(Form("fh3FFTrackPtIM%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsInvMass, fInvMassMin, fInvMassMax, fNBinsPt, fPtMin, fPtMax);
1099 fh3Xi = new TH3F(Form("fh3FFXiIM%s", fNameFF.Data()),"", fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsInvMass, fInvMassMin, fInvMassMax, fNBinsXi, fXiMin, fXiMax);
1100 fh3Z = new TH3F(Form("fh3FFZIM%s", fNameFF.Data()),"", fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsInvMass, fInvMassMin, fInvMassMax, fNBinsZ, fZMin, fZMax);
1102 AliAnalysisTaskFragmentationFunction::SetProperties(fh1JetPt, "p_{t} (GeV/c)", "entries");
1103 AliAnalysisTaskJetChem::SetProperties(fh3TrackPt,"jet p_{t} (GeV/c)","inv Mass (GeV/c^2)","p_{t} (GeV/c)");
1104 AliAnalysisTaskJetChem::SetProperties(fh3Xi,"jet p_{t} (GeV/c)","inv Mass (GeV/c^2)","#xi");
1105 AliAnalysisTaskJetChem::SetProperties(fh3Z,"jet p_{t} (GeV/c)","inv Mass (GeV/c^2)","z");
1108 //________________________________________________________________________________________________________________________________
1109 void AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::FillFF(Float_t trackPt, Float_t invM, Float_t jetPt, Bool_t incrementJetPt)
1111 // fill FF, don't use TH3F anymore use THnSparse instead to save memory
1113 if(incrementJetPt) fh1JetPt->Fill(jetPt);
1114 //fh3TrackPt->Fill(jetPt,invM,trackPt);//Fill(x,y,z)
1117 if(jetPt>0) z = trackPt / jetPt;
1119 if(z>0) xi = TMath::Log(1/z);
1121 //fh3Xi->Fill(jetPt,invM,xi);
1122 //fh3Z->Fill(jetPt,invM,z);
1125 //___________________________________________________________________________________
1126 void AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::AddToOutput(TList* list) const
1128 // add histos to list
1130 list->Add(fh1JetPt);
1131 //list->Add(fh3TrackPt);
1137 //____________________________________________________
1138 void AliAnalysisTaskJetChem::UserCreateOutputObjects()
1140 // create output objects
1142 fRandom = new TRandom3(0);
1143 fRandom->SetSeed(0);
1145 if(fDebug > 1) Printf("AliAnalysisTaskJetChem::UserCreateOutputObjects()");
1147 // create list of tracks and jets
1149 fTracksRecCuts = new TList();
1150 fTracksRecCuts->SetOwner(kFALSE); //objects in TList wont be deleted when TList is deleted
1151 fJetsRecCuts = new TList();
1152 fJetsRecCuts->SetOwner(kFALSE);
1153 fBckgJetsRec = new TList();
1154 fBckgJetsRec->SetOwner(kFALSE);
1155 fListK0s = new TList();
1156 fListK0s->SetOwner(kFALSE);
1157 fListLa = new TList();
1158 fListLa->SetOwner(kFALSE);
1159 fListALa = new TList();
1160 fListALa->SetOwner(kFALSE);
1161 fListFeeddownLaCand = new TList(); //feeddown Lambda candidates
1162 fListFeeddownLaCand->SetOwner(kFALSE);
1163 fListFeeddownALaCand = new TList(); //feeddown Antilambda candidates
1164 fListFeeddownALaCand->SetOwner(kFALSE);
1165 jetConeFDLalist = new TList();
1166 jetConeFDLalist->SetOwner(kFALSE); //feeddown Lambda candidates in jet cone
1167 jetConeFDALalist = new TList();
1168 jetConeFDALalist->SetOwner(kFALSE); //feeddown Antilambda candidates in jet cone
1169 fListMCgenK0s = new TList(); //MC generated K0s
1170 fListMCgenK0s->SetOwner(kFALSE);
1171 fListMCgenLa = new TList(); //MC generated Lambdas
1172 fListMCgenLa->SetOwner(kFALSE);
1173 fListMCgenALa = new TList(); //MC generated Antilambdas
1174 fListMCgenALa->SetOwner(kFALSE);
1177 // Create histograms / output container
1179 fCommonHistList = new TList();
1180 fCommonHistList->SetOwner();
1182 Bool_t oldStatus = TH1::AddDirectoryStatus();
1183 TH1::AddDirectory(kFALSE);//By default (fAddDirectory = kTRUE), histograms are automatically added to the list of objects in memory
1185 // histograms inherited from AliAnalysisTaskFragmentationFunction
1187 fh1EvtSelection = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
1188 fh1EvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
1189 fh1EvtSelection->GetXaxis()->SetBinLabel(2,"event trigger selection: rejected");
1190 fh1EvtSelection->GetXaxis()->SetBinLabel(3,"event class: rejected");
1191 fh1EvtSelection->GetXaxis()->SetBinLabel(4,"vertex Ncontr: rejected");
1192 fh1EvtSelection->GetXaxis()->SetBinLabel(5,"vertex z: rejected");
1193 fh1EvtSelection->GetXaxis()->SetBinLabel(6,"vertex type: rejected");
1196 fh1EvtCent = new TH1F("fh1EvtCent","centrality",100,0.,100.);
1197 fh1VertexNContributors = new TH1F("fh1VertexNContributors", "Vertex N contributors", 11,-.5, 10.5);
1198 fh1VertexZ = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
1199 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
1200 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
1201 fh1Trials = new TH1F("fh1Trials","trials from pyxsec.root",1,0,1);
1202 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
1203 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
1204 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
1205 fh1nRecJetsCuts = new TH1F("fh1nRecJetsCuts","reconstructed jets per event",100,-0.5,99.5);
1207 // histograms JetChem task
1209 fh1EvtAllCent = new TH1F("fh1EvtAllCent","before centrality selection",100,0.,100.);
1210 fh1Evt = new TH1F("fh1Evt", "All events runned over", 3, 0.,1.);
1211 fh1EvtMult = new TH1F("fh1EvtMult","multiplicity",240,0.,240.);
1212 fh1K0Mult = new TH1F("fh1K0Mult","K0 multiplicity",100,0.,100.);//500. all
1213 fh1dPhiJetK0 = new TH1F("fh1dPhiJetK0","",64,-1,5.4);
1214 fh1LaMult = new TH1F("fh1LaMult","La multiplicity",100,0.,100.);
1215 fh1dPhiJetLa = new TH1F("fh1dPhiJetLa","",64,-1,5.4);
1216 fh1ALaMult = new TH1F("fh1ALaMult","ALa multiplicity",100,0.,100.);
1217 fh1dPhiJetALa = new TH1F("fh1dPhiJetALa","",64,-1,5.4);
1218 fh1JetEta = new TH1F("fh1JetEta","#eta distribution of all jets",40,-2.,2.);
1219 fh1JetPhi = new TH1F("fh1JetPhi","#phi distribution of all jets",63,0.,6.3);
1220 fh2JetEtaPhi = new TH2F("fh2JetEtaPhi","#eta and #phi distribution of all jets",400,-2.,2.,63,0.,6.3);
1223 //fh1V0JetPt = new TH1F("fh1V0JetPt","p_{T} distribution of all jets containing v0s",200,0.,200.);
1224 fh1IMK0Cone = new TH1F("fh1IMK0Cone","p_{T} distribution of all jets containing K0s candidates",19,5.,100.);
1225 fh1IMLaCone = new TH1F("fh1IMLaCone","p_{T} distribution of all jets containing #Lambda candidates",19,5.,100.);
1226 fh1IMALaCone = new TH1F("fh1IMALaCone","p_{T} distribution of all jets containing #bar{#Lambda} candidates",19,5.,100.);
1228 fh2FFJetTrackEta = new TH2F("fh2FFJetTrackEta","charged track eta distr. in jet cone",200,-1.,1.,40,0.,200.);
1229 //fh1trackPosNCls = new TH1F("fh1trackPosNCls","NTPC clusters positive daughters",10,0.,100.);
1230 //fh1trackNegNCls = new TH1F("fh1trackNegNCls","NTPC clusters negative daughters",10,0.,100.);
1231 fh1trackPosEta = new TH1F("fh1trackPosEta","eta positive daughters",100,-2.,2.);
1232 fh1trackNegEta = new TH1F("fh1trackNegEta","eta negative daughters",100,-2.,2.);
1233 fh1V0Eta = new TH1F("fh1V0Eta","V0 eta",60,-1.5,1.5);
1234 //fh1V0totMom = new TH1F("fh1V0totMom","V0 tot mom",100,0.,20.);
1235 fh1CosPointAngle = new TH1F("fh1CosPointAngle", "Cosine of V0's pointing angle",50,0.99,1.0);
1236 fh1DecayLengthV0 = new TH1F("fh1DecayLengthV0", "V0s decay Length;decay length(cm)",1200,0.,120.);
1237 fh2ProperLifetimeK0sVsPtBeforeCut = new TH2F("fh2ProperLifetimeK0sVsPtBeforeCut"," K0s ProperLifetime vs Pt; p_{T} (GeV/#it{c})",150,0.,15.,250,0.,250.);
1238 fh2ProperLifetimeK0sVsPtAfterCut = new TH2F("fh2ProperLifetimeK0sVsPtAfterCut"," K0s ProperLifetime vs Pt; p_{T} (GeV/#it{c})",150,0.,15.,250,0.,250.);
1239 fh1V0Radius = new TH1F("fh1V0Radius", "V0s Radius;Radius(cm)",200,0.,40.);
1240 fh1DcaV0Daughters = new TH1F("fh1DcaV0Daughters", "DCA between daughters;dca(cm)",200,0.,2.);
1241 fh1DcaPosToPrimVertex = new TH1F("fh1DcaPosToPrimVertex", "Positive V0 daughter;dca(cm)",100,0.,10.);
1242 fh1DcaNegToPrimVertex = new TH1F("fh1DcaNegToPrimVertex", "Negative V0 daughter;dca(cm)",100,0.,10.);
1243 fh2ArmenterosBeforeCuts = new TH2F("fh2ArmenterosBeforeCuts","Armenteros Podolanski Plot for K0s Candidates;#alpha;(p^{arm})_{T}/(GeV/#it{c})",200,-1.2,1.2,350,0.,0.35);
1244 fh2ArmenterosAfterCuts = new TH2F("fh2ArmenterosAfterCuts","Armenteros Podolanski Plot for K0s Candidates;#alpha;(p^{arm})_{T}/(GeV/#it{c});",200,-1.2,1.2,350,0.,0.35);
1245 fh2BBLaPos = new TH2F("fh2BBLaPos","PID of the positive daughter of La candidates; P (GeV); -dE/dx (keV/cm ?)",100,0,10,200,0,200);
1246 fh2BBLaNeg = new TH2F("fh2BBLaNeg","PID of the negative daughter of La candidates; P (GeV); -dE/dx (keV/cm ?)",100,0,10,200,0,200);
1247 fh1PosDaughterCharge = new TH1F("fh1PosDaughterCharge","charge of V0 positive daughters; V0 daughters",3,-2.,2.);
1248 fh1NegDaughterCharge = new TH1F("fh1NegDaughterCharge","charge of V0 negative daughters; V0 daughters",3,-2.,2.);
1249 fh1PtMCK0s = new TH1F("fh1PtMCK0s","Pt of MC rec K0s; #it{p}_{T} (GeV/#it{c})",200,0.,20.);
1250 fh1PtMCLa = new TH1F("fh1PtMCLa","Pt of MC rec La; #it{p}_{T} (GeV/#it{c})",200,0.,20.);
1251 fh1PtMCALa = new TH1F("fh1PtMCALa","Pt of MC rec ALa; #it{p}_{T} (GeV/#it{c})",200,0.,20.);
1252 fh1EtaK0s = new TH1F("fh1EtaK0s","K^{0}_{s} entries ;#eta",200,-1.,1.);
1253 fh1EtaLa = new TH1F("fh1EtaLa","#Lambda entries ;#eta",200,-1.,1.);
1254 fh1EtaALa = new TH1F("fh1EtaALa","#bar{#Lambda} entries ;#eta",200,-1.,1.);
1256 Int_t binsInvMassEtaTrackPtK0s[3] = {200, 200, 120};//eta,invM,trackPt
1257 Double_t xminInvMassEtaTrackPtK0s[3] = {-1.,0.3,0.};
1258 Double_t xmaxInvMassEtaTrackPtK0s[3] = {1.,0.7,12.};
1260 fhnInvMassEtaTrackPtK0s = new THnSparseF("fhnInvMassEtaTrackPtK0s","#eta; K0s invM (GeV/{#it{c}}^{2}); #it{p}_{T} (GeV/#it{c})",3,binsInvMassEtaTrackPtK0s,xminInvMassEtaTrackPtK0s,xmaxInvMassEtaTrackPtK0s);
1262 Int_t binsInvMassEtaTrackPtLa[3] = {200, 200, 120};//eta,invM,trackPt
1263 Double_t xminInvMassEtaTrackPtLa[3] = {-1.,1.05,0.};
1264 Double_t xmaxInvMassEtaTrackPtLa[3] = {1.,1.25,12.};
1266 fhnInvMassEtaTrackPtLa = new THnSparseF("fhnInvMassEtaTrackPtLa","#eta; #Lambda invM (GeV/{#it{c}}^{2}); #it{p}_{T} (GeV/#it{c})",3,binsInvMassEtaTrackPtLa,xminInvMassEtaTrackPtLa,xmaxInvMassEtaTrackPtLa);
1268 Int_t binsInvMassEtaTrackPtALa[3] = {200, 200, 120};//eta,invM,trackPt
1269 Double_t xminInvMassEtaTrackPtALa[3] = {-1.,1.05,0.};
1270 Double_t xmaxInvMassEtaTrackPtALa[3] = {1.,1.25,12.};
1272 fhnInvMassEtaTrackPtALa = new THnSparseF("fhnInvMassEtaTrackPtALa","#eta; #bar{#Lambda} invM (GeV/{#it{c}}^{2}); #it{p}_{T} (GeV/#it{c})",3,binsInvMassEtaTrackPtALa,xminInvMassEtaTrackPtALa,xmaxInvMassEtaTrackPtALa);
1274 Int_t binsK0sPC[4] = {19, 200, 200, 200};
1275 Double_t xminK0sPC[4] = {5.,0.3, 0., -1.};
1276 Double_t xmaxK0sPC[4] = {100.,0.7, 20., 1.};
1277 fhnK0sPC = new THnSparseF("fhnK0sPC","jet pT; K0s invM; particle pT; particle #eta",4,binsK0sPC,xminK0sPC,xmaxK0sPC);
1279 Int_t binsLaPC[4] = {19, 200, 200, 200};
1280 Double_t xminLaPC[4] = {5.,1.05, 0., -1.};
1281 Double_t xmaxLaPC[4] = {100.,1.25, 20., 1.};
1282 fhnLaPC = new THnSparseF("fhnLaPC","jet pT; #Lambda invM; particle pT; particle #eta",4,binsLaPC,xminLaPC,xmaxLaPC);
1284 Int_t binsALaPC[4] = {19, 200, 200, 200};
1285 Double_t xminALaPC[4] = {5.,1.05, 0., -1.};
1286 Double_t xmaxALaPC[4] = {100.,1.25, 20., 1.};
1287 fhnALaPC = new THnSparseF("fhnALaPC","jet pT; #bar#Lambda invM; particle pT; particle #eta",4,binsALaPC,xminALaPC,xmaxALaPC);
1289 Int_t binsK0sMCC[3] = {200, 200, 200};
1290 Double_t xminK0sMCC[3] = {0.3, 0., -1.};
1291 Double_t xmaxK0sMCC[3] = {0.7, 20., 1.};
1292 fhnK0sMCC = new THnSparseF("fhnK0sMCC","jet pT; K0s invM; particle pT; particle #eta",3,binsK0sMCC,xminK0sMCC,xmaxK0sMCC);
1294 Int_t binsLaMCC[3] = {200, 200, 200};
1295 Double_t xminLaMCC[3] = {1.05, 0., -1.};
1296 Double_t xmaxLaMCC[3] = {1.25, 20., 1.};
1297 fhnLaMCC = new THnSparseF("fhnLaMCC","jet pT; #Lambda invM; particle pT; particle #eta",3,binsLaMCC,xminLaMCC,xmaxLaMCC);
1299 Int_t binsALaMCC[3] = {200, 200, 200};
1300 Double_t xminALaMCC[3] = {1.05, 0., -1.};
1301 Double_t xmaxALaMCC[3] = {1.25, 20., 1.};
1302 fhnALaMCC = new THnSparseF("fhnALaMCC","jet pT; #bara#Lambda invM; particle pT; particle #eta",3,binsALaMCC,xminALaMCC,xmaxALaMCC);
1304 Int_t binsK0sRC[3] = {200, 200, 200};
1305 Double_t xminK0sRC[3] = {0.3, 0., -1.};
1306 Double_t xmaxK0sRC[3] = {0.7, 20., 1.};
1307 fhnK0sRC = new THnSparseF("fhnK0sRC","jet pT; K0s invM; particle pT; particle #eta",3,binsK0sRC,xminK0sRC,xmaxK0sRC);
1309 Int_t binsLaRC[3] = {200, 200, 200};
1310 Double_t xminLaRC[3] = {1.05, 0., -1.};
1311 Double_t xmaxLaRC[3] = {1.25, 20., 1.};
1312 fhnLaRC = new THnSparseF("fhnLaRC","jet pT; #Lambda invM; particle pT; particle #eta",3,binsLaRC,xminLaRC,xmaxLaRC);
1314 Int_t binsALaRC[3] = {200, 200, 200};
1315 Double_t xminALaRC[3] = {1.05, 0., -1.};
1316 Double_t xmaxALaRC[3] = {1.25, 20., 1.};
1317 fhnALaRC = new THnSparseF("fhnALaRC","jet pT; #bara#Lambda invM; particle pT; particle #eta",3,binsALaRC,xminALaRC,xmaxALaRC);
1320 Int_t binsK0sOC[3] = {200, 200, 200};
1321 Double_t xminK0sOC[3] = {0.3, 0., -1.};
1322 Double_t xmaxK0sOC[3] = {0.7, 20., 1.};
1323 fhnK0sOC = new THnSparseF("fhnK0sOC","jet pT; K0s invM; particle pT; particle #eta",3,binsK0sOC,xminK0sOC,xmaxK0sOC);
1325 Int_t binsLaOC[3] = {200, 200, 200};
1326 Double_t xminLaOC[3] = {1.05, 0., -1.};
1327 Double_t xmaxLaOC[3] = {1.25, 20., 1.};
1328 fhnLaOC = new THnSparseF("fhnLaOC","jet pT; #Lambda invM; particle pT; particle #eta",3,binsLaOC,xminLaOC,xmaxLaOC);
1330 Int_t binsALaOC[3] = {200, 200, 200};
1331 Double_t xminALaOC[3] = {1.05, 0., -1.};
1332 Double_t xmaxALaOC[3] = {1.25, 20., 1.};
1334 fhnALaOC = new THnSparseF("fhnALaOC","jet pT; #bara#Lambda invM; particle pT; particle #eta",3,binsALaOC,xminALaOC,xmaxALaOC);
1340 fh1AreaExcluded = new TH1F("fh1AreaExcluded","area excluded for selected jets in event acceptance",100,0.,5.);
1342 fh1MedianEta = new TH1F("fh1MedianEta","Median cluster axis ;#eta",200,-1.,1.);
1343 fh1JetPtMedian = new TH1F("fh1JetPtMedian"," (selected) jet it{p}_{T} distribution for MCC method; #GeV/it{c}",19,5.,100.);
1345 fh1TrackMultCone = new TH1F("fh1TrackMultCone","track multiplicity in jet cone; number of tracks",20,0.,50.);
1347 fh2TrackMultCone = new TH2F("fh2TrackMultCone","track multiplicity in jet cone vs. jet momentum; number of tracks; jet it{p}_{T} (GeV/it{c})",50,0.,50.,19,5.,100.);
1349 Int_t binsNJK0[3] = {200, 200, 200};
1350 Double_t xminNJK0[3] = {0.3, 0., -1.};
1351 Double_t xmaxNJK0[3] = {0.7, 20., 1.};
1352 fhnNJK0 = new THnSparseF("fhnNJK0","K0s candidates in events wo selected jets;",3,binsNJK0,xminNJK0,xmaxNJK0);
1354 Int_t binsNJLa[3] = {200, 200, 200};
1355 Double_t xminNJLa[3] = {1.05, 0., -1.};
1356 Double_t xmaxNJLa[3] = {1.25, 20., 1.};
1357 fhnNJLa = new THnSparseF("fhnNJLa","La candidates in events wo selected jets; ",3,binsNJLa,xminNJLa,xmaxNJLa);
1359 Int_t binsNJALa[3] = {200, 200, 200};
1360 Double_t xminNJALa[3] = {1.05, 0., -1.};
1361 Double_t xmaxNJALa[3] = {1.25, 20., 1.};
1362 fhnNJALa = new THnSparseF("fhnNJALa","ALa candidates in events wo selected jets; ",3,binsNJALa,xminNJALa,xmaxNJALa);
1364 fFFHistosRecCuts = new AliFragFuncHistos("RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1365 fFFNBinsPt, fFFPtMin, fFFPtMax,
1366 fFFNBinsXi, fFFXiMin, fFFXiMax,
1367 fFFNBinsZ , fFFZMin , fFFZMax);
1369 fV0QAK0 = new AliFragFuncQATrackHistos("V0QAK0",fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1370 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1371 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1372 fQATrackHighPtThreshold);
1374 fFFHistosRecCutsK0Evt = new AliFragFuncHistos("RecCutsK0Evt", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1375 fFFNBinsPt, fFFPtMin, fFFPtMax,
1376 fFFNBinsXi, fFFXiMin, fFFXiMax,
1377 fFFNBinsZ , fFFZMin , fFFZMax);
1380 fFFHistosIMK0AllEvt = new AliFragFuncHistosInvMass("K0AllEvt", fFFIMNBinsJetPt, fFFIMJetPtMin, fFFIMJetPtMax,
1381 fFFIMNBinsInvM,fFFIMInvMMin,fFFIMInvMMax,
1382 fFFIMNBinsPt, fFFIMPtMin, fFFIMPtMax,
1383 fFFIMNBinsXi, fFFIMXiMin, fFFIMXiMax,
1384 fFFIMNBinsZ , fFFIMZMin , fFFIMZMax);
1386 fFFHistosIMK0Jet = new AliFragFuncHistosInvMass("K0Jet", fFFIMNBinsJetPt, fFFIMJetPtMin, fFFIMJetPtMax,
1387 fFFIMNBinsInvM,fFFIMInvMMin,fFFIMInvMMax,
1388 fFFIMNBinsPt, fFFIMPtMin, fFFIMPtMax,
1389 fFFIMNBinsXi, fFFIMXiMin, fFFIMXiMax,
1390 fFFIMNBinsZ , fFFIMZMin , fFFIMZMax);
1392 fFFHistosIMK0Cone = new AliFragFuncHistosInvMass("K0Cone", fFFIMNBinsJetPt, fFFIMJetPtMin, fFFIMJetPtMax,
1393 fFFIMNBinsInvM,fFFIMInvMMin,fFFIMInvMMax,
1394 fFFIMNBinsPt, fFFIMPtMin, fFFIMPtMax,
1395 fFFIMNBinsXi, fFFIMXiMin, fFFIMXiMax,
1396 fFFIMNBinsZ , fFFIMZMin , fFFIMZMax);
1398 fFFHistosIMLaAllEvt = new AliFragFuncHistosInvMass("LaAllEvt", fFFIMLaNBinsJetPt, fFFIMLaJetPtMin, fFFIMLaJetPtMax,
1399 fFFIMLaNBinsInvM,fFFIMLaInvMMin,fFFIMLaInvMMax,
1400 fFFIMLaNBinsPt, fFFIMLaPtMin, fFFIMLaPtMax,
1401 fFFIMLaNBinsXi, fFFIMLaXiMin, fFFIMLaXiMax,
1402 fFFIMLaNBinsZ , fFFIMLaZMin , fFFIMLaZMax);
1404 fFFHistosIMLaJet = new AliFragFuncHistosInvMass("LaJet", fFFIMLaNBinsJetPt, fFFIMLaJetPtMin, fFFIMLaJetPtMax,
1405 fFFIMLaNBinsInvM,fFFIMLaInvMMin,fFFIMLaInvMMax,
1406 fFFIMLaNBinsPt, fFFIMLaPtMin, fFFIMLaPtMax,
1407 fFFIMLaNBinsXi, fFFIMLaXiMin, fFFIMLaXiMax,
1408 fFFIMLaNBinsZ , fFFIMLaZMin , fFFIMLaZMax);
1411 fFFHistosIMLaCone = new AliFragFuncHistosInvMass("LaCone", fFFIMLaNBinsJetPt, fFFIMLaJetPtMin, fFFIMLaJetPtMax,
1412 fFFIMLaNBinsInvM,fFFIMLaInvMMin,fFFIMLaInvMMax,
1413 fFFIMLaNBinsPt, fFFIMLaPtMin, fFFIMLaPtMax,
1414 fFFIMLaNBinsXi, fFFIMLaXiMin, fFFIMLaXiMax,
1415 fFFIMLaNBinsZ , fFFIMLaZMin , fFFIMLaZMax);
1418 fFFHistosIMALaAllEvt = new AliFragFuncHistosInvMass("ALaAllEvt", fFFIMLaNBinsJetPt, fFFIMLaJetPtMin, fFFIMLaJetPtMax,
1419 fFFIMLaNBinsInvM,fFFIMLaInvMMin,fFFIMLaInvMMax,
1420 fFFIMLaNBinsPt, fFFIMLaPtMin, fFFIMLaPtMax,
1421 fFFIMLaNBinsXi, fFFIMLaXiMin, fFFIMLaXiMax,
1422 fFFIMLaNBinsZ , fFFIMLaZMin , fFFIMLaZMax);
1424 fFFHistosIMALaJet = new AliFragFuncHistosInvMass("ALaJet", fFFIMLaNBinsJetPt, fFFIMLaJetPtMin, fFFIMLaJetPtMax,
1425 fFFIMLaNBinsInvM,fFFIMLaInvMMin,fFFIMLaInvMMax,
1426 fFFIMLaNBinsPt, fFFIMLaPtMin, fFFIMLaPtMax,
1427 fFFIMLaNBinsXi, fFFIMLaXiMin, fFFIMLaXiMax,
1428 fFFIMLaNBinsZ , fFFIMLaZMin , fFFIMLaZMax);
1430 fFFHistosIMALaCone = new AliFragFuncHistosInvMass("ALaCone", fFFIMLaNBinsJetPt, fFFIMLaJetPtMin, fFFIMLaJetPtMax,
1431 fFFIMLaNBinsInvM,fFFIMLaInvMMin,fFFIMLaInvMMax,
1432 fFFIMLaNBinsPt, fFFIMLaPtMin, fFFIMLaPtMax,
1433 fFFIMLaNBinsXi, fFFIMLaXiMin, fFFIMLaXiMax,
1434 fFFIMLaNBinsZ , fFFIMLaZMin , fFFIMLaZMax);
1441 fh2MCgenK0Cone = new TH2F("fh2MCgenK0Cone", "MC gen {K^{0}}^{s} #it{p}_{T} in cone around rec jet axis versus jet #it{p}_{T}; jet #it{p}_{T}",19,5.,100.,200,0.,20.);
1442 fh2MCgenLaCone = new TH2F("fh2MCgenLaCone", "MC gen #Lambda #it{p}_{T} in cone around rec jet axis versus jet #it{p}_{T} ; jet #it{p}_{T}",19,5.,100.,200,0.,20.);
1443 fh2MCgenALaCone = new TH2F("fh2MCgenALaCone", "MC gen #Antilambda #it{p}_{T} in cone around rec jet axis versus jet #it{p}_{T}; jet #it{p}_{T}",19,5.,100.,200,0.,20.);
1445 fh2MCgenK0Cone->GetYaxis()->SetTitle("MC gen K^{0}}^{s} #it{p}_{T}");
1446 fh2MCgenLaCone->GetYaxis()->SetTitle("MC gen #Lambda #it{p}_{T}");
1447 fh2MCgenALaCone->GetYaxis()->SetTitle("MC gen #Antilambda #it{p}_{T}");
1449 fh2MCEtagenK0Cone = new TH2F("fh2MCEtagenK0Cone","MC gen {K^{0}}^{s} #it{p}_{T} #eta distribution in jet cone;#eta",19,5.,100.,200,-1.,1.);
1450 fh2MCEtagenLaCone = new TH2F("fh2MCEtagenLaCone","MC gen #Lambda #it{p}_{T} #eta distribution in jet cone;#eta",19,5.,100.,200,-1.,1.);
1451 fh2MCEtagenALaCone = new TH2F("fh2MCEtagenALaCone","MC gen #Antilambda #it{p}_{T} #eta distribution in jet cone;#eta",19,5.,100.,200,-1.,1.);
1452 fh1IMK0ConeSmear = new TH1F("fh1IMK0ConeSmear","Smeared jet pt study for K0s-in-cone-jets; smeared jet #it{p}_{T}", 19,5.,100.);
1453 fh1IMLaConeSmear = new TH1F("fh1IMLaConeSmear","Smeared jet pt study for La-in-cone-jets; smeared jet #it{p}_{T}", 19,5.,100.);
1454 fh1IMALaConeSmear = new TH1F("fh1IMALaConeSmear","Smeared jet pt study for ALa-in-cone-jets; smeared jet #it{p}_{T}", 19,5.,100.);
1456 //8 new histograms: Cone, Incl, Lambda, Antilambda, Hijing, Injected:
1458 Int_t binsrecMCHijingLaIncl[3] = {200, 200, 200};
1459 Double_t xminrecMCHijingLaIncl[3] = {1.05, 0., -1.};
1460 Double_t xmaxrecMCHijingLaIncl[3] = {1.25, 20., 1.};
1461 fhnrecMCHijingLaIncl = new THnSparseF("fhnrecMCHijingLaIncl","La inv. mass; particle pT; particle #eta",3,binsrecMCHijingLaIncl,xminrecMCHijingLaIncl,xmaxrecMCHijingLaIncl);
1463 Int_t binsrecMCHijingLaCone[4] = {19, 200, 200, 200};
1464 Double_t xminrecMCHijingLaCone[4] = {5., 1.05, 0., -1.};
1465 Double_t xmaxrecMCHijingLaCone[4] = {100., 1.25, 20., 1.};
1466 fhnrecMCHijingLaCone = new THnSparseF("fhnrecMCHijingLaCone","La inv. mass; particle pT; particle #eta",4,binsrecMCHijingLaCone,xminrecMCHijingLaCone,xmaxrecMCHijingLaCone);
1468 Int_t binsrecMCHijingALaIncl[3] = {200, 200, 200};
1469 Double_t xminrecMCHijingALaIncl[3] = {1.05, 0., -1.};
1470 Double_t xmaxrecMCHijingALaIncl[3] = {1.25, 20., 1.};
1471 fhnrecMCHijingALaIncl = new THnSparseF("fhnrecMCHijingALaIncl","ALa inv. mass; particle pT; particle #eta",3,binsrecMCHijingALaIncl,xminrecMCHijingALaIncl,xmaxrecMCHijingALaIncl);
1473 Int_t binsrecMCHijingALaCone[4] = {19, 200, 200, 200};
1474 Double_t xminrecMCHijingALaCone[4] = {5., 1.05, 0., -1.};
1475 Double_t xmaxrecMCHijingALaCone[4] = {100., 1.25, 20., 1.};
1476 fhnrecMCHijingALaCone = new THnSparseF("fhnrecMCHijingALaCone","ALa inv. mass; particle pT; particle #eta",4,binsrecMCHijingALaCone,xminrecMCHijingALaCone,xmaxrecMCHijingALaCone);
1478 Int_t binsrecMCInjectLaIncl[3] = {200, 200, 200};
1479 Double_t xminrecMCInjectLaIncl[3] = {1.05, 0., -1.};
1480 Double_t xmaxrecMCInjectLaIncl[3] = {1.25, 20., 1.};
1481 fhnrecMCInjectLaIncl = new THnSparseF("fhnrecMCInjectLaIncl","La inv. mass; particle pT; particle #eta",3,binsrecMCInjectLaIncl,xminrecMCInjectLaIncl,xmaxrecMCInjectLaIncl);
1483 Int_t binsrecMCInjectLaCone[4] = {19, 200, 200, 200};
1484 Double_t xminrecMCInjectLaCone[4] = {5., 1.05, 0., -1.};
1485 Double_t xmaxrecMCInjectLaCone[4] = {100., 1.25, 20., 1.};
1486 fhnrecMCInjectLaCone = new THnSparseF("fhnrecMCInjectLaCone","La jet pT;inv. mass; particle pT; particle #eta",4,binsrecMCInjectLaCone,xminrecMCInjectLaCone,xmaxrecMCInjectLaCone);
1488 Int_t binsrecMCInjectALaIncl[3] = {200, 200, 200};
1489 Double_t xminrecMCInjectALaIncl[3] = {1.05, 0., -1.};
1490 Double_t xmaxrecMCInjectALaIncl[3] = {1.25, 20., 1.};
1491 fhnrecMCInjectALaIncl = new THnSparseF("fhnrecMCInjectALaIncl","ALa inv. mass; particle pT; particle #eta",3,binsrecMCInjectALaIncl,xminrecMCInjectALaIncl,xmaxrecMCInjectALaIncl);
1493 Int_t binsrecMCInjectALaCone[4] = {19, 200, 200, 200};
1494 Double_t xminrecMCInjectALaCone[4] = {5., 1.05, 0., -1.};
1495 Double_t xmaxrecMCInjectALaCone[4] = {100., 1.25, 20., 1.};
1496 fhnrecMCInjectALaCone = new THnSparseF("fhnrecMCInjectALaCone","ALa inv. mass; particle pT; particle #eta",4,binsrecMCInjectALaCone,xminrecMCInjectALaCone,xmaxrecMCInjectALaCone);
1499 Int_t binsMCrecK0Cone[4] = {19, 200, 200, 200};
1500 Double_t xminMCrecK0Cone[4] = {5.,0.3, 0., -1.};
1501 Double_t xmaxMCrecK0Cone[4] = {100.,0.7, 20., 1.};
1502 fhnMCrecK0Cone = new THnSparseF("fhnMCrecK0Cone", "MC rec {K^{0}}^{s} #it{p}_{T} in cone around jet axis matching MC gen particle; jet #it{p}_{T}; inv mass (GeV/#it{c}^{2};#it{p}_{T}",4,binsMCrecK0Cone,xminMCrecK0Cone,xmaxMCrecK0Cone);
1504 Int_t binsMCrecLaCone[4] = {19, 200, 200, 200};
1505 Double_t xminMCrecLaCone[4] = {5.,0.3, 0., -1.};
1506 Double_t xmaxMCrecLaCone[4] = {100.,0.7, 20., 1.};
1507 fhnMCrecLaCone = new THnSparseF("fhnMCrecLaCone", "MC rec {#Lambda #it{p}_{T} in cone around jet axis matching MC gen particle; jet #it{p}_{T}; inv mass (GeV/#it{c}^{2};#it{p}_{T}",4,binsMCrecLaCone,xminMCrecLaCone,xmaxMCrecLaCone);
1509 Int_t binsMCrecALaCone[4] = {19, 200, 200, 200};
1510 Double_t xminMCrecALaCone[4] = {5.,0.3, 0., -1.};
1511 Double_t xmaxMCrecALaCone[4] = {100.,0.7, 20., 1.};
1512 fhnMCrecALaCone = new THnSparseF("fhnMCrecALaCone", "MC rec {#bar{#Lambda} #it{p}_{T} in cone around jet axis matching MC gen particle; jet #it{p}_{T}; inv mass (GeV/#it{c}^{2};#it{p}_{T}",4,binsMCrecALaCone,xminMCrecALaCone,xmaxMCrecALaCone);
1514 Int_t binsMCrecK0ConeSmear[4] = {19, 200, 200, 200};
1515 Double_t xminMCrecK0ConeSmear[4] = {5.,0.3, 0., -1.};
1516 Double_t xmaxMCrecK0ConeSmear[4] = {100.,0.7, 20., 1.};
1517 fhnMCrecK0ConeSmear = new THnSparseF("fhnMCrecK0ConeSmear", "MC rec {K^{0}}^{s} #it{p}_{T} in cone around jet axis matching MC gen particle; jet #it{p}_{T}; inv mass (GeV/#it{c}^{2};#it{p}_{T}",4,binsMCrecK0ConeSmear,xminMCrecK0ConeSmear,xmaxMCrecK0ConeSmear);
1519 Int_t binsMCrecLaConeSmear[4] = {19, 200, 200, 200};
1520 Double_t xminMCrecLaConeSmear[4] = {5.,0.3, 0., -1.};
1521 Double_t xmaxMCrecLaConeSmear[4] = {100.,0.7, 20., 1.};
1522 fhnMCrecLaConeSmear = new THnSparseF("fhnMCrecLaConeSmear", "MC rec {#Lambda #it{p}_{T} in cone around jet axis matching MC gen particle; jet #it{p}_{T}; inv mass (GeV/#it{c}^{2};#it{p}_{T}",4,binsMCrecLaConeSmear,xminMCrecLaConeSmear,xmaxMCrecLaConeSmear);
1524 Int_t binsMCrecALaConeSmear[4] = {19, 200, 200, 200};
1525 Double_t xminMCrecALaConeSmear[4] = {5.,0.3, 0., -1.};
1526 Double_t xmaxMCrecALaConeSmear[4] = {100.,0.7, 20., 1.};
1527 fhnMCrecALaConeSmear = new THnSparseF("fhnMCrecALaConeSmear", "MC rec {#bar{#Lambda} #it{p}_{T} in cone around jet axis matching MC gen particle; jet #it{p}_{T}; inv mass (GeV/#it{c}^{2};#it{p}_{T}",4,binsMCrecALaConeSmear,xminMCrecALaConeSmear,xmaxMCrecALaConeSmear);
1529 Int_t binsK0sSecContinCone[3] = {19, 200, 200};
1530 Double_t xminK0sSecContinCone[3] = {5.,0., -1.};
1531 Double_t xmaxK0sSecContinCone[3] = {100.,20., 1.};
1532 fhnK0sSecContinCone = new THnSparseF("fhnK0sSecContinCone", "Secondary contamination {K^{0}}^{s} #it{p}_{T} in cone around jet axis; jet #it{p}_{T}; inv mass (GeV/#it{c}^{2};#it{p}_{T}",3,binsK0sSecContinCone,xminK0sSecContinCone,xmaxK0sSecContinCone);
1534 Int_t binsLaSecContinCone[3] = {19, 200, 200};
1535 Double_t xminLaSecContinCone[3] = {5.,0., -1.};
1536 Double_t xmaxLaSecContinCone[3] = {100.,20., 1.};
1537 fhnLaSecContinCone = new THnSparseF("fhnLaSecContinCone", "Secondary contamination {#Lambda #it{p}_{T} in cone around jet axis; jet #it{p}_{T}; inv mass (GeV/#it{c}^{2};#it{p}_{T}",3,binsLaSecContinCone,xminLaSecContinCone,xmaxLaSecContinCone);
1539 Int_t binsALaSecContinCone[3] = {19, 200, 200};
1540 Double_t xminALaSecContinCone[3] = {5.,0., -1.};
1541 Double_t xmaxALaSecContinCone[3] = {100.,20., 1.};
1542 fhnALaSecContinCone = new THnSparseF("fhnALaSecContinCone", "Secondary contamination {#bar{#Lambda} #it{p}_{T} in cone around jet axis; jet #it{p}_{T}; inv mass (GeV/#it{c}^{2};#it{p}_{T}",3,binsALaSecContinCone,xminALaSecContinCone,xmaxALaSecContinCone);
1544 Int_t binsK0sIncl[3] = {200, 200, 200};
1545 Double_t xminK0sIncl[3] = {0.3, 0., -1.};
1546 Double_t xmaxK0sIncl[3] = {0.7, 20., 1.};
1547 fhnK0sIncl = new THnSparseF("fhnK0sIncl","K0s inv. mass; particle pT; particle #eta",3,binsK0sIncl,xminK0sIncl,xmaxK0sIncl);
1549 Int_t binsK0sCone[4] = {19, 200, 200, 200};
1550 Double_t xminK0sCone[4] = {5.,0.3, 0., -1.};
1551 Double_t xmaxK0sCone[4] = {100.,0.7, 20., 1.};
1552 fhnK0sCone = new THnSparseF("fhnK0sCone","jet pT; K0s inv. mass; particle pT; particle #eta",4,binsK0sCone,xminK0sCone,xmaxK0sCone);
1554 Int_t binsLaIncl[3] = {200, 200, 200};
1555 Double_t xminLaIncl[3] = {1.05, 0., -1.};
1556 Double_t xmaxLaIncl[3] = {1.25, 20., 1.};
1557 fhnLaIncl = new THnSparseF("fhnLaIncl","La inv. mass; particle pT; particle #eta",3,binsLaIncl,xminLaIncl,xmaxLaIncl);
1559 Int_t binsLaCone[4] = {19, 200, 200, 200};
1560 Double_t xminLaCone[4] = {5.,1.05, 0., -1.};
1561 Double_t xmaxLaCone[4] = {100.,1.25, 20., 1.};
1562 fhnLaCone = new THnSparseF("fhnLaCone","jet pT; La inv. mass; particle pT; particle #eta",4,binsLaCone,xminLaCone,xmaxLaCone);
1564 Int_t binsALaIncl[3] = {200, 200, 200};
1565 Double_t xminALaIncl[3] = {1.05, 0., -1.};
1566 Double_t xmaxALaIncl[3] = {1.25, 20., 1.};
1567 fhnALaIncl = new THnSparseF("fhnALaIncl","ALa inv. mass; particle pT; particle #eta",3,binsALaIncl,xminALaIncl,xmaxALaIncl);
1569 Int_t binsALaCone[4] = {19, 200, 200, 200};
1570 Double_t xminALaCone[4] = {5.,1.05, 0., -1.};
1571 Double_t xmaxALaCone[4] = {100.,1.25, 20., 1.};
1572 fhnALaCone = new THnSparseF("fhnALaCone","jet pT; ALa inv. mass; particle pT; particle #eta",4,binsALaCone,xminALaCone,xmaxALaCone);
1574 fh1MCMultiplicityPrimary = new TH1F("fh1MCMultiplicityPrimary", "MC Primary Particles;NPrimary;Count", 201, -0.5, 200.5);
1575 fh1MCMultiplicityTracks = new TH1F("h1MCMultiplicityTracks", "MC Tracks;Ntracks;Count", 201, -0.5, 200.5);
1578 Int_t binsFeedDownLa[3] = {19, 200, 200};
1579 Double_t xminFeedDownLa[3] = {5.,1.05, 0.};
1580 Double_t xmaxFeedDownLa[3] = {100.,1.25, 20.};
1581 fhnFeedDownLa = new THnSparseF("fhnFeedDownLa","#Lambda stemming from feeddown from Xi(0/-)",3,binsFeedDownLa,xminFeedDownLa,xmaxFeedDownLa);
1583 Int_t binsFeedDownALa[3] = {19, 200, 200};
1584 Double_t xminFeedDownALa[3] = {5.,1.05, 0.};
1585 Double_t xmaxFeedDownALa[3] = {100.,1.25, 20.};
1586 fhnFeedDownALa = new THnSparseF("fhnFeedDownALa","#bar#Lambda stemming from feeddown from Xibar(0/+)",3,binsFeedDownALa,xminFeedDownALa,xmaxFeedDownALa);
1588 Int_t binsFeedDownLaCone[3] = {19, 200, 200};
1589 Double_t xminFeedDownLaCone[3] = {5.,1.05, 0.};
1590 Double_t xmaxFeedDownLaCone[3] = {100.,1.25, 20.};
1591 fhnFeedDownLaCone = new THnSparseF("fhnFeedDownLaCone","#Lambda stemming from feeddown from Xi(0/-) in jet cone",3,binsFeedDownLaCone,xminFeedDownLaCone,xmaxFeedDownLaCone);
1593 Int_t binsFeedDownALaCone[3] = {19, 200, 200};
1594 Double_t xminFeedDownALaCone[3] = {5.,1.05, 0.};
1595 Double_t xmaxFeedDownALaCone[3] = {100.,1.25, 20.};
1596 fhnFeedDownALaCone = new THnSparseF("fhnFeedDownALaCone","#bar#Lambda stemming from feeddown from Xibar(0/+) in jet cone",3,binsFeedDownALaCone,xminFeedDownALaCone,xmaxFeedDownALaCone);
1599 fh1MCProdRadiusK0s = new TH1F("fh1MCProdRadiusK0s","MC gen. MC K0s prod radius",200,0.,200.);
1600 fh1MCProdRadiusLambda = new TH1F("fh1MCProdRadiusLambda","MC gen. MC La prod radius",200,0.,200.);
1601 fh1MCProdRadiusAntiLambda = new TH1F("fh1MCProdRadiusAntiLambda","MC gen. MC ALa prod radius",200,0.,200.);
1603 // Pt and inv mass distributions
1605 fh1MCPtV0s = new TH1F("fh1MCPtV0s", "MC gen. V^{0} in rap range;#it{p}_{T} (GeV/#it{c})",200,0,20.);// 0.1 GeV/c steps
1606 fh1MCPtK0s = new TH1F("fh1MCPtK0s", "MC gen. K^{0}_{s} in eta range;#it{p}_{T} (GeV/#it{c})",200,0.,20.);
1607 fh1MCPtLambda = new TH1F("fh1MCPtLambda", "MC gen. #Lambda in rap range;#it{p}_{T} (GeV/#it{c})",200,0.,20.);
1608 fh1MCPtAntiLambda = new TH1F("fh1MCPtAntiLambda", "MC gen. #AntiLambda in rap range;#it{p}_{T} (GeV/#it{c})",200,0.,20.);
1609 fh1MCXiPt = new TH1F("fh1MCXiPt", "MC gen. #Xi^{-/o};#it{p}_{T} (GeV/#it{c})",200,0.,20.);
1610 fh1MCXibarPt = new TH1F("fh1MCXibarPt", "MC gen. #bar{#Xi}^{+/o};#it{p}_{T} (GeV/#it{c})",200,0.,20.);
1611 fh2MCEtaVsPtK0s = new TH2F("fh2MCEtaVsPtK0s","MC gen. K^{0}_{s} #eta; #it{p}_{T}",200,0.,20.,200,-1.,1.);
1612 fh2MCEtaVsPtLa = new TH2F("fh2MCEtaVsPtLa","MC gen. #Lambda #eta; #it{p}_{T}",200,0.,20.,200,-1.,1.);
1613 fh2MCEtaVsPtALa = new TH2F("fh2MCEtaVsPtALa","MC gen. #bar{#Lambda} #eta; #it{p}_{T}",200,0.,20.,200,-1.,1.);
1616 //fh1MCRapK0s = new TH1F("fh1MCRapK0s", "MC gen. K0s;rap with cut",200,-10,10);
1617 //fh1MCRapLambda = new TH1F("fh1MCRapLambda", "MC gen. #Lambda;rap",200,-10,10);
1618 //fh1MCRapAntiLambda = new TH1F("fh1MCRapAntiLambda", "MC gen. #bar{#Lambda};rap",200,-10,10);
1619 fh1MCEtaAllK0s = new TH1F("fh1MCEtaAllK0s", "MC gen. K0s;#eta",200,-1.,1.);
1620 fh1MCEtaK0s = new TH1F("fh1MCEtaK0s", "MC gen. K0s;#eta with cut",200,-1.,1.);
1621 fh1MCEtaLambda = new TH1F("fh1MCEtaLambda", "MC gen. #Lambda;#eta",200,-1.,1.);
1622 fh1MCEtaAntiLambda = new TH1F("fh1MCEtaAntiLambda", "MC gen. #bar{#Lambda};#eta",200,-1.,1.);
1624 fV0QAK0->DefineHistos();
1625 fFFHistosRecCuts->DefineHistos();
1626 fFFHistosRecCutsK0Evt->DefineHistos();
1627 /* fFFHistosIMK0AllEvt->DefineHistos();
1628 fFFHistosIMK0Jet->DefineHistos();
1629 fFFHistosIMK0Cone->DefineHistos();
1630 fFFHistosIMLaAllEvt->DefineHistos();
1631 fFFHistosIMLaJet->DefineHistos();
1632 fFFHistosIMLaCone->DefineHistos();
1633 fFFHistosIMALaAllEvt->DefineHistos();
1634 fFFHistosIMALaJet->DefineHistos();
1635 fFFHistosIMALaCone->DefineHistos();
1638 const Int_t saveLevel = 5;
1641 fCommonHistList->Add(fh1EvtAllCent);
1642 fCommonHistList->Add(fh1Evt);
1643 fCommonHistList->Add(fh1EvtSelection);
1644 fCommonHistList->Add(fh1EvtCent);
1645 fCommonHistList->Add(fh1VertexNContributors);
1646 fCommonHistList->Add(fh1VertexZ);
1647 fCommonHistList->Add(fh1Xsec);
1648 fCommonHistList->Add(fh1Trials);
1649 fCommonHistList->Add(fh1PtHard);
1650 fCommonHistList->Add(fh1PtHardTrials);
1651 fCommonHistList->Add(fh1nRecJetsCuts);
1652 fCommonHistList->Add(fh1EvtMult);
1653 fCommonHistList->Add(fh1K0Mult);
1654 fCommonHistList->Add(fh1dPhiJetK0);
1655 fCommonHistList->Add(fh1LaMult);
1656 fCommonHistList->Add(fh1dPhiJetLa);
1657 fCommonHistList->Add(fh1ALaMult);
1658 fCommonHistList->Add(fh1dPhiJetALa);
1659 fCommonHistList->Add(fh1JetEta);
1660 fCommonHistList->Add(fh1JetPhi);
1661 fCommonHistList->Add(fh2JetEtaPhi);
1662 //fCommonHistList->Add(fh1V0JetPt);
1663 fCommonHistList->Add(fh1IMK0Cone);
1664 fCommonHistList->Add(fh1IMLaCone);
1665 fCommonHistList->Add(fh1IMALaCone);
1666 fCommonHistList->Add(fh2FFJetTrackEta);
1667 // fCommonHistList->Add(fh1trackPosNCls);
1668 //fCommonHistList->Add(fh1trackNegNCls);
1669 fCommonHistList->Add(fh1trackPosEta);
1670 fCommonHistList->Add(fh1trackNegEta);
1671 fCommonHistList->Add(fh1V0Eta);
1672 // fCommonHistList->Add(fh1V0totMom);
1673 fCommonHistList->Add(fh1CosPointAngle);
1674 fCommonHistList->Add(fh1DecayLengthV0);
1675 fCommonHistList->Add(fh2ProperLifetimeK0sVsPtBeforeCut);
1676 fCommonHistList->Add(fh2ProperLifetimeK0sVsPtAfterCut);
1677 fCommonHistList->Add(fh1V0Radius);
1678 fCommonHistList->Add(fh1DcaV0Daughters);
1679 fCommonHistList->Add(fh1DcaPosToPrimVertex);
1680 fCommonHistList->Add(fh1DcaNegToPrimVertex);
1681 fCommonHistList->Add(fh2ArmenterosBeforeCuts);
1682 fCommonHistList->Add(fh2ArmenterosAfterCuts);
1683 fCommonHistList->Add(fh2BBLaPos);
1684 fCommonHistList->Add(fh2BBLaNeg);
1685 fCommonHistList->Add(fh1PosDaughterCharge);
1686 fCommonHistList->Add(fh1NegDaughterCharge);
1687 fCommonHistList->Add(fh1PtMCK0s);
1688 fCommonHistList->Add(fh1PtMCLa);
1689 fCommonHistList->Add(fh1PtMCALa);
1690 fCommonHistList->Add(fh1EtaK0s);
1691 fCommonHistList->Add(fh1EtaLa);
1692 fCommonHistList->Add(fh1EtaALa);
1693 fCommonHistList->Add(fhnInvMassEtaTrackPtK0s);
1694 fCommonHistList->Add(fhnInvMassEtaTrackPtLa);
1695 fCommonHistList->Add(fhnInvMassEtaTrackPtALa);
1696 fCommonHistList->Add(fh1TrackMultCone);
1697 fCommonHistList->Add(fh2TrackMultCone);
1698 fCommonHistList->Add(fhnNJK0);
1699 fCommonHistList->Add(fhnNJLa);
1700 fCommonHistList->Add(fhnNJALa);
1701 fCommonHistList->Add(fh2MCgenK0Cone);
1702 fCommonHistList->Add(fh2MCgenLaCone);
1703 fCommonHistList->Add(fh2MCgenALaCone);
1704 fCommonHistList->Add(fh2MCEtagenK0Cone);
1705 fCommonHistList->Add(fh2MCEtagenLaCone);
1706 fCommonHistList->Add(fh2MCEtagenALaCone);
1707 fCommonHistList->Add(fh1IMK0ConeSmear);
1708 fCommonHistList->Add(fh1IMLaConeSmear);
1709 fCommonHistList->Add(fh1IMALaConeSmear);
1710 fCommonHistList->Add(fhnrecMCHijingLaIncl);
1711 fCommonHistList->Add(fhnrecMCHijingLaCone);
1712 fCommonHistList->Add(fhnrecMCHijingALaIncl);
1713 fCommonHistList->Add(fhnrecMCHijingALaCone);
1714 fCommonHistList->Add(fhnrecMCInjectLaIncl);
1715 fCommonHistList->Add(fhnrecMCInjectLaCone);
1716 fCommonHistList->Add(fhnrecMCInjectALaIncl);
1717 fCommonHistList->Add(fhnrecMCInjectALaCone);
1718 fCommonHistList->Add(fhnrecMCHijingLaIncl);
1719 fCommonHistList->Add(fhnrecMCHijingLaCone);
1720 fCommonHistList->Add(fhnrecMCHijingALaIncl);
1721 fCommonHistList->Add(fhnrecMCHijingALaCone);
1722 fCommonHistList->Add(fhnrecMCInjectLaIncl);
1723 fCommonHistList->Add(fhnrecMCInjectLaCone);
1724 fCommonHistList->Add(fhnrecMCInjectALaIncl);
1725 fCommonHistList->Add(fhnrecMCInjectALaCone);
1726 fCommonHistList->Add(fhnMCrecK0Cone);
1727 fCommonHistList->Add(fhnMCrecLaCone);
1728 fCommonHistList->Add(fhnMCrecALaCone);
1729 fCommonHistList->Add(fhnMCrecK0ConeSmear);
1730 fCommonHistList->Add(fhnMCrecLaConeSmear);
1731 fCommonHistList->Add(fhnMCrecALaConeSmear);
1732 fCommonHistList->Add(fhnK0sSecContinCone);
1733 fCommonHistList->Add(fhnLaSecContinCone);
1734 fCommonHistList->Add(fhnALaSecContinCone);
1735 fCommonHistList->Add(fhnK0sIncl);
1736 fCommonHistList->Add(fhnK0sCone);
1737 fCommonHistList->Add(fhnLaIncl);
1738 fCommonHistList->Add(fhnLaCone);
1739 fCommonHistList->Add(fhnALaIncl);
1740 fCommonHistList->Add(fhnALaCone);
1741 fCommonHistList->Add(fhnK0sPC);
1742 fCommonHistList->Add(fhnLaPC);
1743 fCommonHistList->Add(fhnALaPC);
1744 fCommonHistList->Add(fhnK0sMCC);
1745 fCommonHistList->Add(fhnLaMCC);
1746 fCommonHistList->Add(fhnALaMCC);
1747 fCommonHistList->Add(fhnK0sRC);
1748 fCommonHistList->Add(fhnLaRC);
1749 fCommonHistList->Add(fhnALaRC);
1750 fCommonHistList->Add(fhnK0sOC);
1751 fCommonHistList->Add(fhnLaOC);
1752 fCommonHistList->Add(fhnALaOC);
1753 fCommonHistList->Add(fh1AreaExcluded);
1754 fCommonHistList->Add(fh1MedianEta);
1755 fCommonHistList->Add(fh1JetPtMedian);
1756 fCommonHistList->Add(fh1MCMultiplicityPrimary);
1757 fCommonHistList->Add(fh1MCMultiplicityTracks);
1758 fCommonHistList->Add(fhnFeedDownLa);
1759 fCommonHistList->Add(fhnFeedDownALa);
1760 fCommonHistList->Add(fhnFeedDownLaCone);
1761 fCommonHistList->Add(fhnFeedDownALaCone);
1762 fCommonHistList->Add(fh1MCProdRadiusK0s);
1763 fCommonHistList->Add(fh1MCProdRadiusLambda);
1764 fCommonHistList->Add(fh1MCProdRadiusAntiLambda);
1765 fCommonHistList->Add(fh1MCPtV0s);
1766 fCommonHistList->Add(fh1MCPtK0s);
1767 fCommonHistList->Add(fh1MCPtLambda);
1768 fCommonHistList->Add(fh1MCPtAntiLambda);
1769 fCommonHistList->Add(fh1MCXiPt);
1770 fCommonHistList->Add(fh1MCXibarPt);
1771 fCommonHistList->Add(fh2MCEtaVsPtK0s);
1772 fCommonHistList->Add(fh2MCEtaVsPtLa);
1773 fCommonHistList->Add(fh2MCEtaVsPtALa);
1774 //fCommonHistList->Add(fh1MCRapK0s);
1775 //fCommonHistList->Add(fh1MCRapLambda);
1776 //fCommonHistList->Add(fh1MCRapAntiLambda);
1777 fCommonHistList->Add(fh1MCEtaAllK0s);
1778 fCommonHistList->Add(fh1MCEtaK0s);
1779 fCommonHistList->Add(fh1MCEtaLambda);
1780 fCommonHistList->Add(fh1MCEtaAntiLambda);
1784 fV0QAK0->AddToOutput(fCommonHistList);
1785 fFFHistosRecCuts->AddToOutput(fCommonHistList);
1786 fFFHistosRecCutsK0Evt->AddToOutput(fCommonHistList);
1787 // fFFHistosIMK0AllEvt->AddToOutput(fCommonHistList);
1788 // fFFHistosIMK0Jet->AddToOutput(fCommonHistList);
1789 // fFFHistosIMK0Cone->AddToOutput(fCommonHistList);
1790 // fFFHistosIMLaAllEvt->AddToOutput(fCommonHistList);
1791 // fFFHistosIMLaJet->AddToOutput(fCommonHistList);
1792 // fFFHistosIMLaCone->AddToOutput(fCommonHistList);
1793 // fFFHistosIMALaAllEvt->AddToOutput(fCommonHistList);
1794 // fFFHistosIMALaJet->AddToOutput(fCommonHistList);
1795 // fFFHistosIMALaCone->AddToOutput(fCommonHistList);
1800 // =========== Switch on Sumw2 for all histos ===========
1801 for (Int_t i=0; i<fCommonHistList->GetEntries(); ++i){
1803 TH1 *h1 = dynamic_cast<TH1*>(fCommonHistList->At(i));
1805 if (h1) h1->Sumw2();//The error per bin will be computed as sqrt(sum of squares of weight) for each bin
1807 THnSparse *hnSparse = dynamic_cast<THnSparse*>(fCommonHistList->At(i));
1808 if(hnSparse) hnSparse->Sumw2();
1812 TH1::AddDirectory(oldStatus);
1813 PostData(1, fCommonHistList);
1816 //_______________________________________________
1817 void AliAnalysisTaskJetChem::UserExec(Option_t *)
1820 // Called for each event
1822 if(fDebug > 1) Printf("AliAnalysisTaskJetChem::UserExec()");
1824 if(fDebug > 1) Printf("Analysis event #%5d", (Int_t) fEntry);
1826 // Trigger selection
1827 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
1828 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
1831 //for AliPIDResponse:
1832 //AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
1833 //AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
1834 fPIDResponse = inputHandler->GetPIDResponse();
1836 if (!fPIDResponse){if(fDebug > 1) Printf("AliAnalysisTaskJetChem::UserExec(): fPIDResponse does not exist!"); return;}
1838 //std::cout<<"inputHandler->IsEventSelected(): "<<inputHandler->IsEventSelected()<<std::endl;
1839 //std::cout<<"fEvtSelectionMask: "<<fEvtSelectionMask<<std::endl;
1841 if(!(inputHandler->IsEventSelected() & fEvtSelectionMask)){
1842 //std::cout<<"########event rejected!!############"<<std::endl;
1843 fh1EvtSelection->Fill(1.);
1844 if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
1845 PostData(1, fCommonHistList);
1849 fESD = dynamic_cast<AliESDEvent*>(InputEvent());//casting of pointers for inherited class, only for ESDs
1851 if(fDebug>3) Printf("%s:%d ESDEvent not found in the input", (char*)__FILE__,__LINE__);
1854 fMCEvent = MCEvent();
1856 if(fDebug>3) Printf("%s:%d MCEvent not found in the input", (char*)__FILE__,__LINE__);
1859 // get AOD event from input/output
1860 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
1861 if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
1862 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
1863 if(fUseAODInputJets) fAODJets = fAOD;
1864 if (fDebug > 1) Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
1867 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
1868 if( handler && handler->InheritsFrom("AliAODHandler") ) {
1869 fAOD = ((AliAODHandler*)handler)->GetAOD();
1871 if (fDebug > 1) Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
1875 if(!fAODJets && !fUseAODInputJets){ // case we have AOD in input & output and want jets from output
1876 TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
1877 if( outHandler && outHandler->InheritsFrom("AliAODHandler") ){
1878 fAODJets = ((AliAODHandler*)outHandler)->GetAOD();
1879 if (fDebug > 1) Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
1883 if(fNonStdFile.Length()!=0){
1884 // case we have an AOD extension - fetch the jets from the extended output
1886 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
1887 fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
1889 if(fDebug>1)Printf("AODExtension not found for %s",fNonStdFile.Data());
1894 Printf("%s:%d AODEvent not found", (char*)__FILE__,__LINE__);
1898 Printf("%s:%d AODEvent with jet branch not found", (char*)__FILE__,__LINE__);
1902 //primary vertex position:
1903 AliAODVertex *myPrimaryVertex = NULL;
1904 myPrimaryVertex = (AliAODVertex*)fAOD->GetPrimaryVertex();
1905 if (!myPrimaryVertex) return;
1906 fh1Evt->Fill(1.);//fill in every event that was accessed with InputHandler
1908 // event selection *****************************************
1910 // *** vertex cut ***
1911 AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
1912 Int_t nTracksPrim = primVtx->GetNContributors();
1913 fh1VertexNContributors->Fill(nTracksPrim);
1915 if (fDebug > 1) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
1917 if(nTracksPrim <= 2){
1918 if (fDebug > 1) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__);
1919 fh1EvtSelection->Fill(3.);
1920 PostData(1, fCommonHistList);
1924 fh1VertexZ->Fill(primVtx->GetZ());
1926 if(TMath::Abs(primVtx->GetZ())>fMaxVertexZ){
1927 if (fDebug > 1) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
1928 fh1EvtSelection->Fill(4.);
1929 PostData(1, fCommonHistList);
1933 // accepts only events that have same "primary" and SPD vertex, special issue of LHC11h PbPb data
1935 //fAOD: pointer to global primary vertex
1937 const AliAODVertex* spdVtx = fAOD->GetPrimaryVertexSPD();
1939 if (TMath::Abs(spdVtx->GetZ() - primVtx->GetZ())>fDeltaVertexZ) { if (fDebug > 1) Printf("deltaZVertex: event REJECTED..."); return;}
1942 //check for vertex radius to be smaller than 1 cm, (that was first applied by Vit Kucera in his analysis)
1944 Double_t vtxX = primVtx->GetX();
1945 Double_t vtxY = primVtx->GetY();
1947 if(TMath::Sqrt(vtxX*vtxX + vtxY*vtxY)>=1){
1948 if (fDebug > 1) Printf("%s:%d primary vertex r = %f: event REJECTED...",(char*)__FILE__,__LINE__,TMath::Sqrt(vtxX*vtxX + vtxY*vtxY));
1953 TString primVtxName(primVtx->GetName());
1955 if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
1956 if (fDebug > 1) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
1957 fh1EvtSelection->Fill(5.);
1958 PostData(1, fCommonHistList);
1962 Bool_t selectedHelper = AliAnalysisHelperJetTasks::Selected();
1963 if(!selectedHelper){
1964 fh1EvtSelection->Fill(6.);
1965 PostData(1, fCommonHistList);
1969 // event selection *****************************************
1971 Double_t centPercent = -1;
1975 if(handler && handler->InheritsFrom("AliAODInputHandler")){
1977 centPercent = fAOD->GetHeader()->GetCentrality();
1979 //std::cout<<"centPercent: "<<centPercent<<std::endl;
1981 fh1EvtAllCent->Fill(centPercent);
1983 if(centPercent>10) cl = 2; //standard PWG-JE binning
1984 if(centPercent>30) cl = 3;
1985 if(centPercent>50) cl = 4;
1989 if(centPercent < 0) cl = -1;
1990 if(centPercent >= 0) cl = 1;
1991 if(centPercent > 10) cl = 2; //standard PWG-JE binning
1992 if(centPercent > 30) cl = 3;
1993 if(centPercent > 50) cl = 4;
1994 if(centPercent > 80) cl = 5; //takes centralities higher than my upper edge of 80%, not to be used
1999 cl = AliAnalysisHelperJetTasks::EventClass();
2001 if(fESD) centPercent = fESD->GetCentrality()->GetCentralityPercentile("V0M"); //ESD JetServices Task has the centrality binning 0-10,10-30,30-50,50-80
2002 fh1EvtAllCent->Fill(centPercent);
2005 if(cl!=fEventClass){ // event not in selected event class, reject event#########################################
2007 if (fDebug > 1) Printf("%s:%d event not in selected event class: event REJECTED ...",(char*)__FILE__,__LINE__);
2008 fh1EvtSelection->Fill(2.);
2009 PostData(1, fCommonHistList);
2012 }//end if fEventClass > 0
2015 if (fDebug > 1) Printf("%s:%d event ACCEPTED ...",(char*)__FILE__,__LINE__);
2018 //Printf("Analysis event #%5d", (Int_t) fEntry);
2020 fh1EvtSelection->Fill(0.);
2021 fh1EvtCent->Fill(centPercent);
2023 //___ get MC information __________________________________________________________________
2026 Double_t ptHard = 0.; //parton energy bins -> energy of particle
2027 Double_t nTrials = 1; // trials for MC trigger weight for real data
2030 AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
2031 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);//check usage of Pythia (pp) or Hijing (PbPb)
2032 AliGenHijingEventHeader* hijingGenHeader = 0x0;
2034 if(pythiaGenHeader){
2035 if(fDebug>3) Printf("%s:%d pythiaGenHeader found", (char*)__FILE__,__LINE__);
2036 nTrials = pythiaGenHeader->Trials();
2037 ptHard = pythiaGenHeader->GetPtHard();
2039 fh1PtHard->Fill(ptHard);
2040 fh1PtHardTrials->Fill(ptHard,nTrials);
2043 } else { // no pythia, hijing?
2045 if(fDebug>3) Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
2047 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
2048 if(!hijingGenHeader){
2049 if(fDebug>3) Printf("%s:%d no pythiaGenHeader or hjingGenHeader found", (char*)__FILE__,__LINE__);
2051 if(fDebug>3) Printf("%s:%d hijingGenHeader found", (char*)__FILE__,__LINE__);
2055 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
2058 //____ fetch jets _______________________________________________________________
2060 Int_t nJCuts = GetListOfJets(fJetsRecCuts, kJetsRecAcceptance);//fetch list with jets that survived all jet cuts: fJetsRecCuts
2062 Int_t nRecJetsCuts = 0; //number of reconstructed jets after jet cuts
2063 if(nJCuts>=0) nRecJetsCuts = fJetsRecCuts->GetEntries();
2064 if(fDebug>2)Printf("%s:%d Selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2065 if(nRecJetsCuts != nJCuts) Printf("%s:%d Mismatch selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2066 fh1nRecJetsCuts->Fill(nRecJetsCuts);
2069 //____ fetch background clusters ___________________________________________________
2070 if(fBranchRecBckgClusters.Length() != 0){
2072 Int_t nBJ = GetListOfBckgJets(fBckgJetsRec, kJetsRec);
2073 Int_t nRecBckgJets = 0;
2074 if(nBJ>=0) nRecBckgJets = fBckgJetsRec->GetEntries();
2075 if(fDebug>2)Printf("%s:%d Selected Rec background jets: %d %d",(char*)__FILE__,__LINE__,nBJ,nRecBckgJets);
2076 if(nBJ != nRecBckgJets) Printf("%s:%d Mismatch Selected Rec background jets: %d %d",(char*)__FILE__,__LINE__,nBJ,nRecBckgJets);
2080 //____ fetch reconstructed particles __________________________________________________________
2082 Int_t nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODCuts);//all tracks of event
2083 if(fDebug>2)Printf("%s:%d selected reconstructed tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,fTracksRecCuts->GetEntries());
2084 if(fTracksRecCuts->GetEntries() != nTCuts)
2085 Printf("%s:%d Mismatch selected reconstructed tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,fTracksRecCuts->GetEntries());
2086 fh1EvtMult->Fill(fTracksRecCuts->GetEntries());
2088 Int_t nK0s = GetListOfV0s(fListK0s,fK0Type,kK0,myPrimaryVertex,fAOD);//all V0s in event with K0s assumption
2090 if(fDebug>5){std::cout<<"fK0Type: "<<fK0Type<<" kK0: "<<kK0<<" myPrimaryVertex: "<<myPrimaryVertex<<" fAOD: "<<fAOD<<std::endl;}
2092 //std::cout<< "nK0s: "<<nK0s<<std::endl;
2094 if(fDebug>2)Printf("%s:%d Selected V0s after cuts: %d %d",(char*)__FILE__,__LINE__,nK0s,fListK0s->GetEntries());
2095 if(nK0s != fListK0s->GetEntries()) Printf("%s:%d Mismatch selected K0s: %d %d",(char*)__FILE__,__LINE__,nK0s,fListK0s->GetEntries());
2096 fh1K0Mult->Fill(fListK0s->GetEntries());
2099 Int_t nLa = GetListOfV0s(fListLa,fLaType,kLambda,myPrimaryVertex,fAOD);//all V0s in event with Lambda particle assumption
2100 if(fDebug>2)Printf("%s:%d Selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nLa,fListLa->GetEntries());
2101 if(nLa != fListLa->GetEntries()) Printf("%s:%d Mismatch selected La: %d %d",(char*)__FILE__,__LINE__,nLa,fListLa->GetEntries());
2102 fh1LaMult->Fill(fListLa->GetEntries());
2104 Int_t nALa = GetListOfV0s(fListALa,fALaType,kAntiLambda,myPrimaryVertex,fAOD);//all V0s in event with Antilambda particle assumption
2105 if(fDebug>2)Printf("%s:%d Selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nALa,fListALa->GetEntries());
2106 if(nALa != fListALa->GetEntries()) Printf("%s:%d Mismatch selected ALa: %d %d",(char*)__FILE__,__LINE__,nALa,fListALa->GetEntries());
2107 fh1ALaMult->Fill(fListALa->GetEntries());
2111 //fetch MC gen particles_______________________________________________________
2113 if(fAnalysisMC){ // here
2115 //fill feeddown histo for associated particles
2117 // Access MC generated particles, fill TLists and histograms :
2119 Int_t nMCgenK0s = GetListOfMCParticles(fListMCgenK0s,kK0,fAOD); //fill TList with MC generated primary true K0s (list to fill, particletype, mc aod event)
2120 if(nMCgenK0s != fListMCgenK0s->GetEntries()) Printf("%s:%d Mismatch selected MCgenK0s: %d %d",(char*)__FILE__,__LINE__,nMCgenK0s,fListMCgenK0s->GetEntries());
2123 for(Int_t it=0; it<fListMCgenK0s->GetSize(); ++it){ // loop MC generated K0s, filling histograms
2125 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenK0s->At(it));
2130 //Double_t fRapCurrentPart = MyRapidity(mcp0->E(),mcp0->Pz());
2131 Double_t fEtaCurrentPart = mcp0->Eta();
2132 Double_t fPtCurrentPart = mcp0->Pt();
2134 fh1MCEtaK0s->Fill(fEtaCurrentPart);
2135 //fh1MCRapK0s->Fill(fRapCurrentPart);
2136 fh1MCPtK0s->Fill(fPtCurrentPart);
2138 fh2MCEtaVsPtK0s->Fill(fPtCurrentPart,fEtaCurrentPart); //eta cut, physical primary selection and decay mode considered
2143 Int_t nMCgenLa = GetListOfMCParticles(fListMCgenLa,kLambda,fAOD); //fill TList with MC generated primary true Lambdas (list to fill, particletype, mc aod event)
2144 if(nMCgenLa != fListMCgenLa->GetEntries()) Printf("%s:%d Mismatch selected MCgenLa: %d %d",(char*)__FILE__,__LINE__,nMCgenLa,fListMCgenLa->GetEntries());
2147 for(Int_t it=0; it<fListMCgenLa->GetSize(); ++it){ // loop MC generated La, filling histograms
2149 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenLa->At(it));
2154 //Double_t fRapCurrentPart = MyRapidity(mcp0->E(),mcp0->Pz());
2155 Double_t fEtaCurrentPart = mcp0->Eta();
2156 Double_t fPtCurrentPart = mcp0->Pt();
2158 fh1MCEtaLambda->Fill(fEtaCurrentPart);
2159 //fh1MCRapLambda->Fill(fRapCurrentPart);
2160 fh1MCPtLambda->Fill(fPtCurrentPart);
2161 fh2MCEtaVsPtLa->Fill(fPtCurrentPart,fEtaCurrentPart); //eta cut, physical primary selection and decay mode considered
2166 Int_t nMCgenALa = GetListOfMCParticles(fListMCgenALa,kAntiLambda,fAOD); //fill TList with MC generated primary true Antilambdas (list to fill, particletype, mc aod event)
2167 if(nMCgenALa != fListMCgenALa->GetEntries()) Printf("%s:%d Mismatch selected MCgenALa: %d %d",(char*)__FILE__,__LINE__,nMCgenALa,fListMCgenALa->GetEntries());
2170 for(Int_t it=0; it<fListMCgenALa->GetSize(); ++it){ // loop MC generated ALa, filling histograms
2172 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenALa->At(it));
2175 //MC gen Antilambdas
2177 //Double_t fRapCurrentPart = MyRapidity(mcp0->E(),mcp0->Pz());
2178 Double_t fEtaCurrentPart = mcp0->Eta();
2179 Double_t fPtCurrentPart = mcp0->Pt();
2181 fh1MCEtaAntiLambda->Fill(fEtaCurrentPart);
2182 //fh1MCRapAntiLambda->Fill(fRapCurrentPart);
2183 fh1MCPtAntiLambda->Fill(fPtCurrentPart);
2184 fh2MCEtaVsPtALa->Fill(fPtCurrentPart,fEtaCurrentPart); //eta cut, physical primary selection and decay mode considered
2190 //loop over MC feeddown candidates in TList
2195 } //end MCAnalysis part for gen particles
2198 // ___ V0 QA + K0s + La + ALa pt spectra all events _______________________________________________
2200 Double_t lPrimaryVtxPosition[3];
2201 Double_t lV0Position[3];
2202 lPrimaryVtxPosition[0] = primVtx->GetX();
2203 lPrimaryVtxPosition[1] = primVtx->GetY();
2204 lPrimaryVtxPosition[2] = primVtx->GetZ();
2205 Double_t dRadiusExcludeCone = 2*GetFFRadius(); //2 times jet radius
2206 //------------------------------------------
2207 for(Int_t it=0; it<fListK0s->GetSize(); ++it){
2209 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0s->At(it));
2214 // VO's main characteristics to check the reconstruction cuts
2216 //Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2219 Double_t fV0Radius = -999;
2220 Double_t fDcaV0Daughters = v0->DcaV0Daughters();
2221 Double_t fDcaPosToPrimVertex = v0->DcaPosToPrimVertex();//IP of positive charged daughter
2222 Double_t fDcaNegToPrimVertex = v0->DcaNegToPrimVertex();//IP of negative charged daughter
2223 Int_t negDaughterpdg = 0;
2224 Int_t posDaughterpdg = 0;
2225 Int_t motherType = 0;
2228 Bool_t fPhysicalPrimary = kFALSE;//don't use IsPhysicalPrimary() anymore for MC analysis, use instead 2D distance from primary to secondary vertex
2229 Int_t MCv0PdgCode = 0;
2231 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0));
2232 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1));
2234 Double_t PosEta = trackPos->AliAODTrack::Eta();//daughter track charge is sometimes wrong here, account for that!!!
2235 Double_t NegEta = trackNeg->AliAODTrack::Eta();
2237 //Double_t trackPosNcls = trackPos->GetTPCNcls();//Get number of clusters for positive charged tracks
2238 //Double_t trackNegNcls = trackNeg->GetTPCNcls();//Get number of clusters for negative charged tracks
2240 CalculateInvMass(v0, kK0, invMK0s, trackPt); //function to calculate invMass with TLorentzVector class
2242 //Double_t jetPt = fFFIMJetPtMin; // assign pro forma jet energy
2243 //Double_t fRap = v0->RapK0Short();
2244 Double_t fEta = v0->PseudoRapV0();
2245 Bool_t bIsInCone = kFALSE;//init boolean, is not in any cone (OC)
2247 for(Int_t ij=0; ij<nRecJetsCuts; ++ij){ // loop over all jets in event
2249 AliAODJet* jet = (AliAODJet*) (fJetsRecCuts->At(ij));
2251 if(IsParticleInCone(jet, v0, dRadiusExcludeCone) == kTRUE) {bIsInCone = kTRUE;}
2255 if(bIsInCone==kFALSE){//K0s is not part of any selected jet in event
2256 Double_t vK0sOC[3] = {invMK0s,trackPt,fEta};
2257 fhnK0sOC->Fill(vK0sOC);
2260 Double_t fV0cosPointAngle = v0->CosPointingAngle(lPrimaryVtxPosition);
2262 lV0Position[0]= v0->DecayVertexV0X();
2263 lV0Position[1]= v0->DecayVertexV0Y();
2264 lV0Position[2]= v0->DecayVertexV0Z();
2268 fV0mom[0]=v0->MomV0X();
2269 fV0mom[1]=v0->MomV0Y();
2270 fV0mom[2]=v0->MomV0Z();
2271 // Double_t fV0TotalMomentum = TMath::Sqrt(fV0mom[0]*fV0mom[0]+fV0mom[1]*fV0mom[1]+fV0mom[2]*fV0mom[2]);
2272 Double_t fV0DecayLength = v0->DecayLengthV0(lPrimaryVtxPosition);
2273 fV0Radius = TMath::Sqrt(lV0Position[0]*lV0Position[0]+lV0Position[1]*lV0Position[1]);
2276 //fetch V0s outside of jet cones (outside of 2R):
2286 fV0QAK0->FillTrackQA(v0->Eta(), TVector2::Phi_0_2pi(v0->Phi()), v0->Pt());
2287 //fFFHistosIMK0AllEvt->FillFF(trackPt, invMK0s, jetPt, incrementJetPt);
2288 //fh1trackPosNCls->Fill(trackPosNcls);
2289 //fh1trackNegNCls->Fill(trackNegNcls);
2290 fh1EtaK0s->Fill(fEta);
2292 Double_t vK0sIncl[3] = {invMK0s,trackPt,fEta}; //fill all K0s in event into THnSparse of 3 dimensions
2293 fhnK0sIncl->Fill(vK0sIncl);
2297 TString generatorName;
2299 TList *listmc = fAOD->GetList();
2300 Bool_t mclabelcheck = MCLabelCheck(v0, kK0, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
2301 //if(fPhysicalPrimary == kFALSE)continue;
2302 //std::cout<<"mclabelcheck: "<<mclabelcheck<<std::endl;
2303 //std::cout<<"IsPhysicalPrimary: "<<fPhysicalPrimary<<std::endl;
2305 if(mclabelcheck == kFALSE)continue;
2307 Double_t vInvMassEtaTrackPtK0s[3] = {fEta,invMK0s,trackPt};
2308 fhnInvMassEtaTrackPtK0s->Fill(vInvMassEtaTrackPtK0s);//includes also feeddown particles, mainly phi particles whose decay products are considered here as primary
2311 fh1PtMCK0s->Fill(MCPt);
2315 fh1V0Eta->Fill(fEta);
2316 //fh1V0totMom->Fill(fV0TotalMomentum);
2317 fh1CosPointAngle->Fill(fV0cosPointAngle);
2318 fh1DecayLengthV0->Fill(fV0DecayLength);
2319 fh1V0Radius->Fill(fV0Radius);
2320 fh1DcaV0Daughters->Fill(fDcaV0Daughters);
2321 fh1DcaPosToPrimVertex->Fill(fDcaPosToPrimVertex);
2322 fh1DcaNegToPrimVertex->Fill(fDcaNegToPrimVertex);
2323 fh1trackPosEta->Fill(PosEta);
2324 fh1trackNegEta->Fill(NegEta);
2328 // __La pt spectra all events _______________________________________________
2331 for(Int_t it=0; it<fListLa->GetSize(); ++it){ // loop all La
2333 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListLa->At(it));
2336 // VO's main characteristics to check the reconstruction cuts
2337 // Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2340 Double_t fV0Radius = -999;
2341 Double_t fDcaV0Daughters = v0->DcaV0Daughters();
2342 Double_t fDcaPosToPrimVertex = v0->DcaPosToPrimVertex();//IP of positive charged daughter
2343 Double_t fDcaNegToPrimVertex = v0->DcaNegToPrimVertex();//IP of negative charged daughter
2344 Int_t negDaughterpdg = 0;
2345 Int_t posDaughterpdg = 0;
2346 Int_t motherType = 0;
2349 Bool_t fPhysicalPrimary = kFALSE;
2350 Int_t MCv0PdgCode = 0;
2351 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0));
2352 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1));
2354 //Double_t trackPosNcls = trackPos->GetTPCNcls();//Get number of clusters for positive charged tracks
2355 //Double_t trackNegNcls = trackNeg->GetTPCNcls();//Get number of clusters for negative charged tracks
2357 Double_t PosEta = trackPos->AliAODTrack::Eta();//daughter track charge is sometimes wrong here, account for that!!!
2358 Double_t NegEta = trackNeg->AliAODTrack::Eta();
2360 Double_t fEta = v0->PseudoRapV0();
2361 Bool_t bIsInCone = kFALSE;//init boolean, is not in any cone (OC)
2363 CalculateInvMass(v0, kLambda, invMLa, trackPt);//function to calculate invMass with TLorentzVector class
2366 for(Int_t ij=0; ij<nRecJetsCuts; ++ij){ // loop over all jets in event
2368 AliAODJet* jet = (AliAODJet*) (fJetsRecCuts->At(ij));
2370 if(IsParticleInCone(jet, v0, dRadiusExcludeCone) == kTRUE) {bIsInCone = kTRUE;
2374 if(bIsInCone == kFALSE){//success!
2375 Double_t vLaOC[3] = {invMLa, trackPt,fEta};
2376 fhnLaOC->Fill(vLaOC);
2379 // Double_t jetPt = fFFIMJetPtMin; // assign pro forma jet energy
2380 // Double_t fRap = v0->Y(3122);
2385 fV0mom[0]=v0->MomV0X();
2386 fV0mom[1]=v0->MomV0Y();
2387 fV0mom[2]=v0->MomV0Z();
2388 //Double_t fV0TotalMomentum = TMath::Sqrt(fV0mom[0]*fV0mom[0]+fV0mom[1]*fV0mom[1]+fV0mom[2]*fV0mom[2]);
2389 Double_t fV0DecayLength = v0->DecayLengthV0(lPrimaryVtxPosition);
2390 Double_t fV0cosPointAngle = v0->CosPointingAngle(lPrimaryVtxPosition);
2391 lV0Position[0]= v0->DecayVertexV0X();
2392 lV0Position[1]= v0->DecayVertexV0Y();
2393 lV0Position[2]= v0->DecayVertexV0Z();
2395 fV0Radius = TMath::Sqrt(lV0Position[0]*lV0Position[0]+lV0Position[1]*lV0Position[1]);
2397 //fFFHistosIMLaAllEvt->FillFF(trackPt, invMLa, jetPt, incrementJetPt);
2398 //fh1trackPosNCls->Fill(trackPosNcls);
2399 //fh1trackNegNCls->Fill(trackNegNcls);
2400 fh1EtaLa->Fill(fEta);
2402 Double_t vLaIncl[3] = {invMLa,trackPt,fEta};
2403 fhnLaIncl->Fill(vLaIncl);
2407 TString generatorName;
2409 TList* listmc = fAOD->GetList();
2410 Bool_t mclabelcheck = MCLabelCheck(v0, kLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
2411 if(mclabelcheck == kFALSE)continue;
2412 //if(fPhysicalPrimary == kFALSE)continue;
2414 if(generatorName == "Hijing"){
2415 Double_t vrecMCHijingLaIncl[3] = {invMLa,trackPt,fEta};
2416 fhnrecMCHijingLaIncl->Fill(vrecMCHijingLaIncl);
2419 if(isinjected == kTRUE){
2420 Double_t vrecMCInjectLaIncl[3] = {invMLa,trackPt,fEta};
2421 fhnrecMCInjectLaIncl->Fill(vrecMCInjectLaIncl);
2424 Double_t vInvMassEtaTrackPtLa[3] = {fEta,invMLa,trackPt};
2425 fhnInvMassEtaTrackPtLa->Fill(vInvMassEtaTrackPtLa);//includes also feed-down particles
2426 fh1PtMCLa->Fill(MCPt);
2429 fh1PtMCLa->Fill(MCPt);
2431 fh1V0Eta->Fill(fEta);
2432 //fh1V0totMom->Fill(fV0TotalMomentum);
2433 fh1CosPointAngle->Fill(fV0cosPointAngle);
2434 fh1DecayLengthV0->Fill(fV0DecayLength);
2435 fh1V0Radius->Fill(fV0Radius);
2436 fh1DcaV0Daughters->Fill(fDcaV0Daughters);
2437 fh1DcaPosToPrimVertex->Fill(fDcaPosToPrimVertex);
2438 fh1DcaNegToPrimVertex->Fill(fDcaNegToPrimVertex);
2439 fh1trackPosEta->Fill(PosEta);
2440 fh1trackNegEta->Fill(NegEta);
2443 // __ALa pt spectra all events _______________________________________________
2445 for(Int_t it=0; it<fListALa->GetSize(); ++it){ // loop all ALa
2447 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListALa->At(it));
2451 //VO's main characteristics to check the reconstruction cuts
2452 Double_t invMALa =0;
2454 Double_t fV0Radius = -999;
2455 Double_t fDcaV0Daughters = v0->DcaV0Daughters();
2456 Double_t fDcaPosToPrimVertex = v0->DcaPosToPrimVertex();//IP of positive charged daughter
2457 Double_t fDcaNegToPrimVertex = v0->DcaNegToPrimVertex();//IP of negative charged daughter
2458 Int_t negDaughterpdg = 0;
2459 Int_t posDaughterpdg = 0;
2460 Int_t motherType = 0;
2463 Bool_t fPhysicalPrimary = kFALSE;
2464 Int_t MCv0PdgCode = 0;
2466 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0));
2467 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1));
2469 Double_t PosEta = trackPos->AliAODTrack::Eta();//daughter track charge is sometimes wrong here, account for that!!!
2470 Double_t NegEta = trackNeg->AliAODTrack::Eta();
2472 Double_t fEta = v0->PseudoRapV0();
2473 Bool_t bIsInCone = kFALSE;//init boolean for OC
2476 CalculateInvMass(v0, kAntiLambda, invMALa, trackPt); //function to calculate invMass with TLorentzVector class
2478 for(Int_t ij=0; ij<nRecJetsCuts; ++ij){ // loop over all jets in event
2480 AliAODJet* jet = (AliAODJet*) (fJetsRecCuts->At(ij));
2482 if(IsParticleInCone(jet, v0, dRadiusExcludeCone) == kTRUE){
2487 if(bIsInCone == kFALSE){//success!
2488 Double_t vALaOC[3] = {invMALa, trackPt,fEta};
2489 fhnALaOC->Fill(vALaOC);
2492 //Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2493 //Double_t jetPt = fFFIMJetPtMin; // assign pro forma jet energy
2494 // Double_t fRap = v0->Y(-3122);
2499 fV0mom[0]=v0->MomV0X();
2500 fV0mom[1]=v0->MomV0Y();
2501 fV0mom[2]=v0->MomV0Z();
2502 //Double_t fV0TotalMomentum = TMath::Sqrt(fV0mom[0]*fV0mom[0]+fV0mom[1]*fV0mom[1]+fV0mom[2]*fV0mom[2]);
2504 Double_t fV0cosPointAngle = v0->CosPointingAngle(lPrimaryVtxPosition);
2505 lV0Position[0]= v0->DecayVertexV0X();
2506 lV0Position[1]= v0->DecayVertexV0Y();
2507 lV0Position[2]= v0->DecayVertexV0Z();
2508 Double_t fV0DecayLength = v0->DecayLengthV0(lPrimaryVtxPosition);
2509 fV0Radius = TMath::Sqrt(lV0Position[0]*lV0Position[0]+lV0Position[1]*lV0Position[1]);
2511 //fFFHistosIMALaAllEvt->FillFF(trackPt, invMALa, jetPt, incrementJetPt);
2512 //fh1trackPosNCls->Fill(trackPosNcls);
2513 //fh1trackNegNCls->Fill(trackNegNcls);
2514 fh1EtaALa->Fill(fEta);
2516 Double_t vALaIncl[3] = {invMALa,trackPt,fEta};
2517 fhnALaIncl->Fill(vALaIncl);
2520 TString generatorName;
2522 TList* listmc = fAOD->GetList();
2523 Bool_t mclabelcheck = MCLabelCheck(v0, kAntiLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
2524 if(mclabelcheck == kFALSE)continue;
2525 //if(fPhysicalPrimary == kFALSE)continue;
2527 if(generatorName == "Hijing"){
2528 Double_t vrecMCHijingALaIncl[3] = {invMALa,trackPt,fEta};
2529 fhnrecMCHijingALaIncl->Fill(vrecMCHijingALaIncl);
2532 if(isinjected == kTRUE){
2533 Double_t vrecMCInjectALaIncl[3] = {invMALa,trackPt,fEta};
2534 fhnrecMCInjectALaIncl->Fill(vrecMCInjectALaIncl);
2538 Double_t vInvMassEtaTrackPtALa[3] = {fEta,invMALa,trackPt};
2539 fhnInvMassEtaTrackPtALa->Fill(vInvMassEtaTrackPtALa);
2540 fh1PtMCALa->Fill(MCPt);
2543 fh1V0Eta->Fill(fEta);
2544 //fh1V0totMom->Fill(fV0TotalMomentum);
2545 fh1CosPointAngle->Fill(fV0cosPointAngle);
2546 fh1DecayLengthV0->Fill(fV0DecayLength);
2547 fh1V0Radius->Fill(fV0Radius);
2548 fh1DcaV0Daughters->Fill(fDcaV0Daughters);
2549 fh1DcaPosToPrimVertex->Fill(fDcaPosToPrimVertex);
2550 fh1DcaNegToPrimVertex->Fill(fDcaNegToPrimVertex);
2551 fh1trackPosEta->Fill(PosEta);
2552 fh1trackNegEta->Fill(NegEta);
2555 //_____no jets events______________________________________________________________________________________________________________________________________
2557 if(nRecJetsCuts == 0){//no jet events
2559 if(fDebug>6) { std::cout<<"################## nRecJetsCuts == 0 ###################"<<std::endl;
2560 std::cout<<"fListK0s->GetSize() in NJ event: "<<fListK0s->GetSize()<<std::endl;
2563 for(Int_t it=0; it<fListK0s->GetSize(); ++it){ // loop all K0s
2565 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0s->At(it));
2568 Double_t invMK0s =0;
2570 CalculateInvMass(v0, kK0, invMK0s, trackPt);
2571 Double_t fEta = v0->Eta();
2573 Double_t vNJK0[3] = {invMK0s,trackPt,fEta}; //fill all K0s in events wo selected jets
2574 fhnNJK0->Fill(vNJK0);
2578 for(Int_t it=0; it<fListLa->GetSize(); ++it){ // loop all La
2580 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListLa->At(it));
2585 CalculateInvMass(v0, kLambda, invMLa, trackPt);
2586 Double_t fEta = v0->Eta();
2588 Double_t vNJLa[3] = {invMLa,trackPt,fEta}; //fill all K0s in events wo selected jets
2589 fhnNJLa->Fill(vNJLa);
2593 for(Int_t it=0; it<fListALa->GetSize(); ++it){ // loop all ALa
2595 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListALa->At(it));
2598 Double_t invMALa =0;
2600 CalculateInvMass(v0, kAntiLambda, invMALa, trackPt);
2602 Double_t fEta = v0->Eta();
2604 Double_t vNJALa[3] = {invMALa,trackPt,fEta}; //fill all K0s in events wo selected jets
2605 fhnNJALa->Fill(vNJALa);
2612 //____ fill all jet related histos ________________________________________________________________________________________________________________________
2613 //##########################jet loop########################################################################################################################
2615 //fill jet histos in general
2616 for(Int_t ij=0; ij<nRecJetsCuts; ++ij){ // ij is an index running over the list of the reconstructed jets after cuts, all jets in event
2618 AliAODJet* jet = (AliAODJet*) (fJetsRecCuts->At(ij));
2620 Double_t jetPt = jet->Pt();
2621 Double_t jetEta = jet->Eta();
2622 Double_t jetPhi = jet->Phi();
2624 //if(ij==0){ // loop over leading jets for ij = 0, for ij>= 0 look into all jets
2626 if(ij>=0){//all jets in event
2628 TList* jettracklist = new TList();
2629 Double_t sumPt = 0.;
2630 Bool_t isBadJet = kFALSE;
2631 Int_t njetTracks = 0;
2633 if(GetFFRadius()<=0){
2634 GetJetTracksTrackrefs(jettracklist, jet, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);// list of jet tracks from trackrefs
2636 GetJetTracksPointing(fTracksRecCuts, jettracklist, jet, GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet); // fill list of tracks in cone around jet axis with cone Radius (= 0.4 standard)
2639 if(GetFFMinNTracks()>0 && jettracklist->GetSize() <= GetFFMinNTracks()) isBadJet = kTRUE; // reject jets with less tracks than fFFMinNTracks
2640 if(isBadJet) continue; // rejects jets in which no track has a track pt higher than 5 GeV/c (see AddTask macro)
2642 //Float_t fJetAreaMin = 0.6*TMath::Pi()*GetFFRadius()*GetFFRadius(); // minimum jet area cut
2644 //std::cout<<"GetFFRadius(): "<<GetFFRadius()<<std::endl;
2645 //std::cout<<"jet->EffectiveAreaCharged()"<<jet->EffectiveAreaCharged()<<std::endl;
2646 //std::cout<<"fJetAreaMin: "<<fJetAreaMin<<std::endl;
2648 //if(fDebug > 2) {if (jet->EffectiveAreaCharged() < fJetAreaMin) {std::cout<<" fCutjetArea cut removed a jet!!!!! Should not have to be done again!!"<<std::endl;}}// cut on jet area, already done by jet selection in FF task
2650 Double_t dAreaExcluded = TMath::Pi()*dRadiusExcludeCone*dRadiusExcludeCone; // area of the cone
2651 dAreaExcluded -= AreaCircSegment(dRadiusExcludeCone,fCutjetEta-jet->Eta()); // positive eta overhang
2652 dAreaExcluded -= AreaCircSegment(dRadiusExcludeCone,fCutjetEta+jet->Eta()); // negative eta overhang
2653 fh1AreaExcluded->Fill(dAreaExcluded);//histo contains all areas that are jet related and have to be excluded concerning OC UE pt spectrum normalisation by area
2655 fh1JetEta->Fill(jetEta);
2656 fh1JetPhi->Fill(jetPhi);
2657 fh2JetEtaPhi->Fill(jetEta,jetPhi);
2659 // printf("pT = %f, eta = %f, phi = %f, leadtr pt = %f\n, ",jetPt,jetEta,jetphi,leadtrack);
2661 for(Int_t it=0; it<jettracklist->GetSize(); ++it){//loop over all particles in jet
2663 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));//all tracks in jet cone
2664 if(!trackVP)continue;
2666 Float_t trackPt = trackVP->Pt();//transversal momentum of jet particle
2667 Float_t trackEta = trackVP->Eta();
2669 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2671 fFFHistosRecCuts->FillFF(trackPt, jetPt, incrementJetPt);//histo with tracks/jets after cut selection, for all events
2672 if(nK0s>0) fFFHistosRecCutsK0Evt->FillFF(trackPt, jetPt, incrementJetPt);//only for K0s events
2673 fh2FFJetTrackEta->Fill(trackEta,jetPt);
2678 njetTracks = jettracklist->GetSize();
2680 //____________________________________________________________________________________________________________________
2681 //strangeness constribution to jet cone
2685 TList *list = fAOD->GetList();
2686 AliAODMCHeader *mcHeadr=(AliAODMCHeader*)list->FindObject(AliAODMCHeader::StdBranchName());
2687 if(!mcHeadr)continue;
2689 Double_t mcXv=0., mcYv=0., mcZv=0.;//MC primary vertex position
2691 mcXv=mcHeadr->GetVtxX(); mcYv=mcHeadr->GetVtxY(); mcZv=mcHeadr->GetVtxZ(); // position of the MC primary vertex
2693 for(Int_t it=0; it<jettracklist->GetSize(); ++it){//loop over all tracks in the jet
2695 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));//track in jet cone
2696 if(!trackVP)continue;
2697 AliAODTrack *tr = dynamic_cast<AliAODTrack*> (trackVP); //fetch one jet track from the TList
2700 //get MC label information
2701 TList *mclist = fAOD->GetList();
2703 //fetch the MC stack
2704 TClonesArray* stackMC = (TClonesArray*)mclist->FindObject(AliAODMCParticle::StdBranchName()); //get MCAOD branch in data
2705 if (!stackMC) {Printf("ERROR: stack not available");}
2709 Int_t trackLabel = TMath::Abs(tr->GetLabel()); //fetch jet track label in MC stack
2711 AliAODMCParticle *part = dynamic_cast<AliAODMCParticle*>(stackMC->At(trackLabel)); //fetch MC gen. particle for rec. jet track
2713 if(!part)continue; //skip non-existing objects
2716 //Bool_t IsPhysicalPrimary = part->IsPhysicalPrimary();//not recommended to check, better use distance between primary vertex and secondary vertex
2718 Float_t fDistPrimaryMax = 0.01;
2719 // Get the distance between production point of the MC mother particle and the primary vertex
2721 Double_t dx = mcXv-part->Xv();//mc primary vertex - mc gen. v0 vertex
2722 Double_t dy = mcYv-part->Yv();
2723 Double_t dz = mcZv-part->Zv();
2725 Float_t fDistPrimary = TMath::Sqrt(dx*dx + dy*dy + dz*dz);
2726 Bool_t fPhysicalPrimary = (fDistPrimary < fDistPrimaryMax);
2728 // std::cout<<"fDistPrimary"<<fDistPrimary<<std::endl;
2729 // std::cout<<"fPhysicalPrimary"<<fPhysicalPrimary<<std::endl;
2731 if(!fPhysicalPrimary)continue;//rejects Kstar and other strong decaying particles from Secondary Contamination
2733 Bool_t isFromStrange = kFALSE;// flag to check whether particle has strange mother
2735 Int_t iMother = part->GetMother(); //get mother MC gen. particle label
2738 AliAODMCParticle *partM = dynamic_cast<AliAODMCParticle*>(stackMC->At(iMother)); //fetch mother of MC gen. particle
2739 if(!partM) continue;
2741 Int_t codeM = TMath::Abs(partM->GetPdgCode()); //mothers pdg code
2743 Int_t mfl = Int_t (codeM/ TMath::Power(10, Int_t(TMath::Log10(codeM)))); //asks for first number of mothers pdg code (strange particles always start with 3..)
2745 if (mfl == 3 && codeM != 3) isFromStrange = kTRUE;
2748 if(isFromStrange == kTRUE){
2750 Double_t trackPt = part->Pt();
2751 Double_t trackEta = part->Eta();
2752 //fh3StrContinCone->Fill(jetPt, trackPt, trackEta);//MC gen. particle parameters, but rec. jet pt
2754 }//isFromStrange is kTRUE
2756 }//end loop over jet tracks
2761 fh1TrackMultCone->Fill(njetTracks);
2762 fh2TrackMultCone->Fill(njetTracks,jetPt);
2766 // fQAJetHistosRecCutsLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
2768 for(Int_t it=0; it<fListK0s->GetSize(); ++it){ // loop all K0s
2770 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0s->At(it));
2771 if(!v0) continue;//rejection of events with no V0 vertex
2775 TVector3 v0MomVect(v0Mom);
2777 Double_t dPhiJetK0 = (jet->MomentumVector()->Vect()).DeltaPhi(v0MomVect);
2778 // Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2780 // if(incrementJetPt){fh1V0JetPt->Fill(jetPt);}
2782 Double_t invMK0s =0;
2784 CalculateInvMass(v0, kK0, invMK0s, trackPt); //function to calculate invMass with TLorentzVector class
2786 // fFFHistosIMK0Jet->FillFF(trackPt, invMK0s, jetPt, incrementJetPt);
2789 if(dPhiJetK0<fh1dPhiJetK0->GetXaxis()->GetXmin()) dPhiJetK0 += 2*TMath::Pi();
2790 fh1dPhiJetK0->Fill(dPhiJetK0);
2794 // if(fListK0s->GetSize() == 0){ // no K0: increment jet pt spectrum
2796 // Bool_t incrementJetPt = kTRUE;
2797 // fFFHistosIMK0Jet->FillFF(-1, -1, jetPt, incrementJetPt);
2800 //____fetch reconstructed K0s in cone around jet axis:_______________________________________________________________________________
2802 TList* jetConeK0list = new TList();
2804 Double_t sumPtK0 = 0.;
2806 Bool_t isBadJetK0 = kFALSE; // dummy, do not use
2809 GetTracksInCone(fListK0s, jetConeK0list, jet, GetFFRadius(), sumPtK0, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetK0); //reconstructed K0s in cone around jet axis
2811 if(fDebug>2)Printf("%s:%d nK0s total: %d, in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nK0s,jetConeK0list->GetEntries(),GetFFRadius());
2814 for(Int_t it=0; it<jetConeK0list->GetSize(); ++it){ // loop for K0s in jet cone
2816 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetConeK0list->At(it));
2819 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2820 Double_t invMK0s =0;
2825 CalculateInvMass(v0, kK0, invMK0s, trackPt); //function to calculate invMass with TLorentzVector class
2829 Double_t jetPtSmear = -1;
2830 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
2831 if(incrementJetPt == kTRUE){fh1IMK0ConeSmear->Fill(jetPtSmear);} //fill TH1F for normalization purposes
2834 if(incrementJetPt==kTRUE){
2835 fh1IMK0Cone->Fill(jetPt);}//normalisation by number of selected jets
2837 //fFFHistosIMK0Cone->FillFF(trackPt, invMK0s, jetPt, incrementJetPt);
2839 Double_t vK0sCone[4] = {jetPt, invMK0s,trackPt,fEta};
2840 fhnK0sCone->Fill(vK0sCone);
2844 if(jetConeK0list->GetSize() == 0){ // no K0: increment jet pt spectrum
2847 Bool_t incrementJetPt = kTRUE;//jets without K0s will be only filled in TH1F only once, so no increment needed
2848 //fFFHistosIMK0Cone->FillFF(-1, -1, jetPt, incrementJetPt);
2849 Double_t vK0sCone[4] = {jetPt, -1, -1, -1};
2850 fhnK0sCone->Fill(vK0sCone);
2852 if(incrementJetPt==kTRUE){
2853 fh1IMK0Cone->Fill(jetPt);}//normalisation by number of selected jets
2856 Double_t jetPtSmear = -1;
2857 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
2858 if(incrementJetPt == kTRUE){fh1IMK0ConeSmear->Fill(jetPtSmear);} //fill TH1F for normalization purposes
2862 //Random cones________________________________________________________________________
2864 if(ij==0){//fetch random cone V0s only once per event
2866 //______fetch random cones___________________________________________________________
2869 AliAODJet* jetRC = 0;
2870 jetRC = GetRandomCone(fJetsRecCuts, fCutjetEta, 2*GetFFRadius());//fetch one random cone for each event
2871 TList* fListK0sRC = new TList();//list for K0s in random cone (RC), one RC per event
2872 TList* fListLaRC = new TList();
2873 TList* fListALaRC = new TList();
2875 Double_t sumPtK0sRC = 0;
2876 Double_t sumPtLaRC = 0;
2877 Double_t sumPtALaRC = 0;
2878 Bool_t isBadJetK0sRC = kFALSE;
2879 Bool_t isBadJetLaRC = kFALSE;
2880 Bool_t isBadJetALaRC = kFALSE;
2885 GetTracksInCone(fListK0s, fListK0sRC, jetRC, GetFFRadius(), sumPtK0sRC, 0, 0, isBadJetK0sRC);
2887 for(Int_t it=0; it<fListK0sRC->GetSize(); ++it){ // loop for K0s in random cone
2889 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0sRC->At(it));
2892 Double_t invMK0s =0;
2897 CalculateInvMass(v0, kK0, invMK0s, trackPt); //function to calculate invMass with TLorentzVector class
2899 Double_t vK0sRC[3] = {invMK0s,trackPt,fEta};
2900 fhnK0sRC->Fill(vK0sRC);
2905 GetTracksInCone(fListLa, fListLaRC, jetRC, GetFFRadius(), sumPtLaRC, 0, 0, isBadJetLaRC);
2907 for(Int_t it=0; it<fListLaRC->GetSize(); ++it){ // loop for Lambdas in random cone
2909 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListLaRC->At(it));
2917 CalculateInvMass(v0, kLambda, invMLa, trackPt); //function to calculate invMass with TLorentzVector class
2919 Double_t vLaRC[3] = {invMLa,trackPt,fEta};
2920 fhnLaRC->Fill(vLaRC);
2925 GetTracksInCone(fListALa, fListALaRC, jetRC, GetFFRadius(), sumPtALaRC, 0, 0, isBadJetALaRC);
2927 for(Int_t it=0; it<fListALaRC->GetSize(); ++it){ // loop for Lambdas in random cone
2929 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListALaRC->At(it));
2932 Double_t invMALa =0;
2937 CalculateInvMass(v0, kAntiLambda, invMALa, trackPt); //function to calculate invMass with TLorentzVector class
2939 Double_t vALaRC[3] = {invMALa,trackPt,fEta};
2940 fhnALaRC->Fill(vALaRC);
2950 //fetch particles in perpendicular cone to estimate UE event contribution to particle spectrum
2951 //these perpendicular cone particle spectra serve to subtract the particles in jet cones, that are stemming from the Underlying event, on a statistical basis
2952 //for normalization the common jet pT spectrum is used: fh1IMK0Cone, fh1IMLaCone and fh1IMALaCone
2954 //____fetch reconstructed K0s in cone perpendicular to jet axis:_______________________________________________________________________________
2956 TList* jetPerpConeK0list = new TList();
2958 Double_t sumPerpPtK0 = 0.;
2960 GetTracksInPerpCone(fListK0s, jetPerpConeK0list, jet, GetFFRadius(), sumPerpPtK0); //reconstructed K0s in cone around jet axis
2962 if(fDebug>2)Printf("%s:%d nK0s total: %d, in perp jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nK0s,jetPerpConeK0list->GetEntries(),GetFFRadius());
2964 for(Int_t it=0; it<jetPerpConeK0list->GetSize(); ++it){ // loop for K0s in perpendicular cone
2966 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetPerpConeK0list->At(it));
2969 Double_t invMPerpK0s =0;
2974 CalculateInvMass(v0, kK0, invMPerpK0s, trackPt); //function to calculate invMass with TLorentzVector class
2975 Double_t vK0sPC[4] = {jetPt, invMPerpK0s,trackPt,fEta};
2977 fhnK0sPC->Fill(vK0sPC); //(x,y,z) //pay attention, this histogram contains the V0 content of both (+/- 90 degrees) perp. cones!!
2982 if(jetPerpConeK0list->GetSize() == 0){ // no K0s in jet cone
2984 Double_t vK0sPC[4] = {jetPt, -1, -1 , -999};//default values for case: no K0s is found in PC
2985 fhnK0sPC->Fill(vK0sPC);
2989 if(ij==0){//median cluster only once for event
2991 AliAODJet* medianCluster = GetMedianCluster();
2994 // ____ rec K0s in median cluster___________________________________________________________________________________________________________
2996 TList* jetMedianConeK0list = new TList();
2997 TList* jetMedianConeLalist = new TList();
2998 TList* jetMedianConeALalist = new TList();
3001 Double_t medianEta = medianCluster->Eta();
3003 if(TMath::Abs(medianEta)<=fCutjetEta){
3005 fh1MedianEta->Fill(medianEta);
3006 fh1JetPtMedian->Fill(jetPt); //for normalisation by total number of median cluster jets
3008 Double_t sumMedianPtK0 = 0.;
3010 Bool_t isBadJetK0Median = kFALSE; // dummy, do not use
3012 GetTracksInCone(fListK0s, jetMedianConeK0list, medianCluster, GetFFRadius(), sumMedianPtK0, 0., 0., isBadJetK0Median); //reconstructed K0s in median cone around jet axis
3013 //GetTracksInCone(fListK0s, jetConeK0list, jet, GetFFRadius(), sumPtK0, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetK0); //original use of function
3015 //cut parameters from Fragmentation Function task:
3016 //Float_t fFFMinLTrackPt; // reject jets with leading track with pt smaller than this value, use GetFFMinLTrackPt()
3017 //Float_t fFFMaxTrackPt; // reject jetscontaining any track with pt larger than this value, use GetFFMaxTrackPt()
3019 for(Int_t it=0; it<jetMedianConeK0list->GetSize(); ++it){ // loop for K0s in median cone
3021 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetMedianConeK0list->At(it));
3024 Double_t invMMedianK0s =0;
3029 CalculateInvMass(v0, kK0, invMMedianK0s, trackPt); //function to calculate invMass with TLorentzVector class
3030 Double_t vK0sMCC[3] = {invMMedianK0s,trackPt,fEta};
3031 fhnK0sMCC->Fill(vK0sMCC);
3035 if(jetMedianConeK0list->GetSize() == 0){ // no K0s in median cluster cone
3037 Double_t vK0sMCC[3] = {-1, -1, -999};
3038 fhnK0sMCC->Fill(vK0sMCC);
3042 //__________________________________________________________________________________________________________________________________________
3043 // ____ rec Lambdas in median cluster___________________________________________________________________________________________________________
3045 Double_t sumMedianPtLa = 0.;
3046 Bool_t isBadJetLaMedian = kFALSE; // dummy, do not use
3048 GetTracksInCone(fListLa, jetMedianConeLalist, medianCluster, GetFFRadius(), sumMedianPtLa, 0, 0, isBadJetLaMedian); //reconstructed Lambdas in median cone around jet axis
3050 //cut parameters from Fragmentation Function task:
3051 //Float_t fFFMinLTrackPt; // reject jets with leading track with pt smaller than this value, use GetFFMinLTrackPt()
3052 //Float_t fFFMaxTrackPt; // reject jets containing any track with pt larger than this value, use GetFFMaxTrackPt()
3054 for(Int_t it=0; it<jetMedianConeLalist->GetSize(); ++it){ // loop for Lambdas in perpendicular cone
3056 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetMedianConeLalist->At(it));
3059 Double_t invMMedianLa =0;
3064 CalculateInvMass(v0, kLambda, invMMedianLa, trackPt); //function to calculate invMass with TLorentzVector class
3066 Double_t vLaMCC[4] = {jetPt, invMMedianLa,trackPt,fEta};
3067 fhnLaMCC->Fill(vLaMCC);
3070 if(jetMedianConeLalist->GetSize() == 0){ // no Lambdas in median cluster cone
3072 Double_t vLaMCC[4] = {jetPt, -1, -1, -999};
3073 fhnLaMCC->Fill(vLaMCC);
3078 // ____ rec Antilambdas in median cluster___________________________________________________________________________________________________________
3081 Double_t sumMedianPtALa = 0.;
3083 Bool_t isBadJetALaMedian = kFALSE; // dummy, do not use
3085 GetTracksInCone(fListALa, jetMedianConeALalist, medianCluster, GetFFRadius(), sumMedianPtALa, 0, 0, isBadJetALaMedian); //reconstructed Antilambdas in median cone around jet axis
3088 //cut parameters from Fragmentation Function task:
3089 //Float_t fFFMinLTrackPt; // reject jets with leading track with pt smaller than this value, use GetFFMinLTrackPt()
3090 //Float_t fFFMaxTrackPt; // reject jets containing any track with pt larger than this value, use GetFFMaxTrackPt()
3092 for(Int_t it=0; it<jetMedianConeALalist->GetSize(); ++it){ // loop for Antilambdas in median cluster cone
3094 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetMedianConeALalist->At(it));
3097 Double_t invMMedianALa =0;
3103 CalculateInvMass(v0, kAntiLambda, invMMedianALa, trackPt); //function to calculate invMass with TLorentzVector class
3104 Double_t vALaMCC[4] = {jetPt, invMMedianALa,trackPt,fEta};
3105 fhnALaMCC->Fill(vALaMCC);
3109 if(jetMedianConeALalist->GetSize() == 0){ // no Antilambdas in median cluster cone
3111 Double_t vALaMCC[4] = {jetPt, -1, -1, -999};
3112 fhnALaMCC->Fill(vALaMCC);
3115 }//median cluster eta cut
3117 delete jetMedianConeK0list;
3118 delete jetMedianConeLalist;
3119 delete jetMedianConeALalist;
3121 }//if mediancluster is existing
3123 //_________________________________________________________________________________________________________________________________________
3125 //____fetch reconstructed Lambdas in cone perpendicular to jet axis:__________________________________________________________________________
3127 TList* jetPerpConeLalist = new TList();
3129 Double_t sumPerpPtLa = 0.;
3131 GetTracksInPerpCone(fListLa, jetPerpConeLalist, jet, GetFFRadius(), sumPerpPtLa); //reconstructed Lambdas in cone around jet axis //pay attention, this histogram contains the V0 content of both (+/- 90 degrees) perp. cones!!
3133 if(fDebug>2)Printf("%s:%d nLa total: %d, in perp jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nLa,jetPerpConeLalist->GetEntries(),GetFFRadius());
3135 for(Int_t it=0; it<jetPerpConeLalist->GetSize(); ++it){ // loop for Lambdas in perpendicular cone
3137 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetPerpConeLalist->At(it));
3140 Double_t invMPerpLa =0;
3145 CalculateInvMass(v0, kLambda, invMPerpLa, trackPt); //function to calculate invMass with TLorentzVector class
3146 Double_t vLaPC[4] = {jetPt, invMPerpLa,trackPt,fEta};
3147 fhnLaPC->Fill(vLaPC); //(x,y,z) //pay attention, this histogram contains the V0 content of both (+/- 90 degrees) perp. cones!!
3152 if(jetPerpConeLalist->GetSize() == 0){ // no Lambdas in jet
3154 Double_t vLaPC[4] = {jetPt, -1, -1 , -999};//default values for case: no K0s is found in PC
3155 fhnLaPC->Fill(vLaPC);
3161 //____fetch reconstructed Antilambdas in cone perpendicular to jet axis:___________________________________________________________________
3163 TList* jetPerpConeALalist = new TList();
3165 Double_t sumPerpPtALa = 0.;
3167 GetTracksInPerpCone(fListALa, jetPerpConeALalist, jet, GetFFRadius(), sumPerpPtALa); //reconstructed Antilambdas in cone around jet axis //pay attention, this histogram contains the V0 content of both (+/- 90 degrees) perp. cones!!
3169 if(fDebug>2)Printf("%s:%d nALa total: %d, in perp jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nALa,jetPerpConeALalist->GetEntries(),GetFFRadius());
3171 for(Int_t it=0; it<jetPerpConeALalist->GetSize(); ++it){ // loop for ALa in perpendicular cone
3173 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetPerpConeALalist->At(it));
3176 Double_t invMPerpALa =0;
3181 CalculateInvMass(v0, kAntiLambda, invMPerpALa, trackPt); //function to calculate invMass with TLorentzVector class
3182 Double_t vALaPC[4] = {jetPt, invMPerpALa,trackPt,fEta};
3183 fhnALaPC->Fill(vALaPC);
3188 if(jetPerpConeALalist->GetSize() == 0){ // no Antilambda
3190 Double_t vALaPC[4] = {jetPt, -1, -1, -999};
3191 fhnALaPC->Fill(vALaPC);
3211 //###########################################################################################################
3213 //__________________________________________________________________________________________________________________________________________
3217 //fill feeddown candidates from TList
3218 //std::cout<<"fListFeeddownLaCand entries: "<<fListFeeddownLaCand->GetSize()<<std::endl;
3220 Double_t sumPtFDLa = 0.;
3221 Bool_t isBadJetFDLa = kFALSE; // dummy, do not use
3223 GetTracksInCone(fListFeeddownLaCand, jetConeFDLalist, jet, GetFFRadius(), sumPtFDLa, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetFDLa);
3225 Double_t sumPtFDALa = 0.;
3226 Bool_t isBadJetFDALa = kFALSE; // dummy, do not use
3228 GetTracksInCone(fListFeeddownALaCand, jetConeFDALalist, jet, GetFFRadius(), sumPtFDALa, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetFDALa);
3230 //_________________________________________________________________
3231 for(Int_t it=0; it<fListFeeddownLaCand->GetSize(); ++it){
3233 AliAODv0* mcfd = dynamic_cast<AliAODv0*>(fListFeeddownLaCand->At(it));
3236 Double_t invMLaFDcand = 0;
3237 Double_t trackPt = 0;//pt of ass. particle, not used for the histos
3239 CalculateInvMass(mcfd, kLambda, invMLaFDcand, trackPt);
3241 //Get MC gen. Lambda transverse momentum
3242 TClonesArray *st = 0x0;
3245 TList *lt = fAOD->GetList();
3248 st = (TClonesArray*)lt->FindObject(AliAODMCParticle::StdBranchName()); //get MCAOD branch in data
3251 AliAODTrack *daughtertrack = (AliAODTrack *) (mcfd->GetSecondaryVtx()->GetDaughter(0));//fetch the first of the two daughter tracks
3252 Int_t AssLabel = TMath::Abs(daughtertrack->GetLabel());
3254 AliAODMCParticle *mcDaughterPart =(AliAODMCParticle*)st->UncheckedAt(AssLabel);
3256 Int_t v0lab = mcDaughterPart->GetMother();
3258 // Int_t v0lab= TMath::Abs(mcfd->GetLabel());//GetLabel doesn't work for AliAODv0 class!!! Only for AliAODtrack
3260 if((!v0lab) || (v0lab<0) || (v0lab > st->GetEntriesFast()))continue;//validity checks
3262 AliAODMCParticle *mcp=(AliAODMCParticle*)st->UncheckedAt(v0lab);
3264 Double_t genLaPt = mcp->Pt();
3266 //std::cout<<"Incl FD, genLaPt:"<<genLaPt<<std::endl;
3268 Double_t vFeedDownLa[3] = {5., invMLaFDcand, genLaPt};
3269 fhnFeedDownLa->Fill(vFeedDownLa);
3272 }//end loop over feeddown candidates for Lambda particles in jet cone
3273 //fetch MC truth in jet cones, denominator of rec. efficiency in jet cones
3274 //_________________________________________________________________
3275 for(Int_t it=0; it<jetConeFDLalist->GetSize(); ++it){
3277 AliAODv0* mcfd = dynamic_cast<AliAODv0*>(jetConeFDLalist->At(it));
3280 //std::cout<<"Cone, recLaPt:"<<mcfd->Pt()<<std::endl;
3282 Double_t invMLaFDcand = 0;
3283 Double_t trackPt = mcfd->Pt();//pt of ass. particle, not used for the histos
3285 CalculateInvMass(mcfd, kLambda, invMLaFDcand, trackPt);
3287 //Get MC gen. Lambda transverse momentum
3288 TClonesArray *st = 0x0;
3290 TList *lt = fAOD->GetList();
3293 st = (TClonesArray*)lt->FindObject(AliAODMCParticle::StdBranchName());
3295 AliAODTrack *daughtertrack = (AliAODTrack *) (mcfd->GetSecondaryVtx()->GetDaughter(0));//fetch the first of the two daughter tracks
3296 Int_t AssLabel = TMath::Abs(daughtertrack->GetLabel());
3298 AliAODMCParticle *mcDaughterPart =(AliAODMCParticle*)st->UncheckedAt(AssLabel);
3300 Int_t v0lab = mcDaughterPart->GetMother();
3302 //std::cout<<"v0lab: "<<v0lab<<std::endl;
3304 if((!v0lab) || (v0lab<0) || (v0lab > st->GetEntriesFast()))continue;//validity checks
3306 AliAODMCParticle *mcp=(AliAODMCParticle*)st->UncheckedAt(v0lab);
3308 Double_t genLaPt = mcp->Pt();
3311 //std::cout<<"Cone FD, genLaPt:"<<genLaPt<<std::endl;
3313 Double_t vFeedDownLaCone[3] = {jetPt, invMLaFDcand, genLaPt};
3314 fhnFeedDownLaCone->Fill(vFeedDownLaCone);
3317 }//end loop over feeddown candidates for Lambda particles in jet cone
3319 //_________________________________________________________________
3320 for(Int_t it=0; it<fListFeeddownALaCand->GetSize(); ++it){
3322 AliAODv0* mcfd = dynamic_cast<AliAODv0*>(fListFeeddownALaCand->At(it));
3325 Double_t invMALaFDcand = 0;
3326 Double_t trackPt = 0;//pt of ass. particle, not used for the histos
3328 CalculateInvMass(mcfd, kAntiLambda, invMALaFDcand, trackPt);
3330 //Get MC gen. Antilambda transverse momentum
3331 TClonesArray *st = 0x0;
3333 TList *lt = fAOD->GetList();
3336 st = (TClonesArray*)lt->FindObject(AliAODMCParticle::StdBranchName());
3338 AliAODTrack *daughtertrack = (AliAODTrack *) (mcfd->GetSecondaryVtx()->GetDaughter(0));//fetch the first of the two daughter tracks
3339 Int_t AssLabel = TMath::Abs(daughtertrack->GetLabel());
3341 AliAODMCParticle *mcDaughterPart =(AliAODMCParticle*)st->UncheckedAt(AssLabel);
3343 Int_t v0lab = mcDaughterPart->GetMother();
3346 if((!v0lab) || (v0lab<0) || (v0lab > st->GetEntriesFast()))continue;//validity checks
3348 AliAODMCParticle *mcp=(AliAODMCParticle*)st->UncheckedAt(v0lab);
3350 Double_t genALaPt = mcp->Pt();
3352 Double_t vFeedDownALa[3] = {5., invMALaFDcand, genALaPt};
3353 fhnFeedDownALa->Fill(vFeedDownALa);
3356 }//end loop over feeddown candidates for Antilambda particles
3358 //_________________________________________________________________
3359 //feeddown for Antilambdas from Xi(bar)+ and Xi(bar)0 in jet cone:
3361 for(Int_t it=0; it<jetConeFDALalist->GetSize(); ++it){
3363 AliAODv0* mcfd = dynamic_cast<AliAODv0*>(jetConeFDALalist->At(it));
3366 Double_t invMALaFDcand = 0;
3367 Double_t trackPt = 0;//pt of ass. particle, not used for the histos
3369 CalculateInvMass(mcfd, kAntiLambda, invMALaFDcand, trackPt);
3371 //Get MC gen. Antilambda transverse momentum
3372 TClonesArray *st = 0x0;
3374 TList *lt = fAOD->GetList();
3377 st = (TClonesArray*)lt->FindObject(AliAODMCParticle::StdBranchName());
3379 AliAODTrack *daughtertrack = (AliAODTrack *) (mcfd->GetSecondaryVtx()->GetDaughter(0));//fetch the first of the two daughter tracks
3380 Int_t AssLabel = TMath::Abs(daughtertrack->GetLabel());
3382 AliAODMCParticle *mcDaughterPart =(AliAODMCParticle*)st->UncheckedAt(AssLabel);
3384 Int_t v0lab = mcDaughterPart->GetMother();
3386 if((!v0lab) || (v0lab<0) || (v0lab > st->GetEntriesFast()))continue;//validity checks
3388 AliAODMCParticle *mcp=(AliAODMCParticle*)st->UncheckedAt(v0lab);
3390 Double_t genALaPt = mcp->Pt();
3392 Double_t vFeedDownALaCone[3] = {jetPt, invMALaFDcand, genALaPt};
3393 fhnFeedDownALaCone->Fill(vFeedDownALaCone);
3396 }//end loop over feeddown candidates for Antilambda particles in jet cone
3400 //____fetch MC generated K0s in cone around jet axis__(note: particles can stem from fragmentation but also from underlying event)________
3402 Double_t sumPtMCgenK0s = 0.;
3403 Bool_t isBadJetMCgenK0s = kFALSE; // dummy, do not use
3406 fListMCgenK0sCone = new TList(); //MC generated K0s in (only geometrical) jet cone (these are MC gen K0s falling geometrically into jet cone (R = 0.4) around jet axis, that was found by anti-kt jet finder, particles can stem from fragmentation but also from underlying event!!)
3408 //first: sampling MC gen K0s
3410 GetTracksInCone(fListMCgenK0s, fListMCgenK0sCone, jet, GetFFRadius(), sumPtMCgenK0s, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetMCgenK0s); //MC generated K0s in cone around jet axis
3412 if(fDebug>2)Printf("%s:%d nMCgenK0s in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,fListMCgenK0sCone->GetEntries(),GetFFRadius());
3415 for(Int_t it=0; it<fListMCgenK0sCone->GetSize(); ++it){ // loop MC generated K0s in cone around jet axis
3417 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenK0sCone->At(it));
3420 //Double_t fRapMCgenK0s = MyRapidity(mcp0->E(),mcp0->Pz());//get rec. particle in cone information
3421 Double_t fEtaMCgenK0s = mcp0->Eta();
3422 Double_t fPtMCgenK0s = mcp0->Pt();
3424 fh2MCgenK0Cone->Fill(jetPt,fPtMCgenK0s);
3425 fh2MCEtagenK0Cone->Fill(jetPt,fEtaMCgenK0s);
3429 //check whether the reconstructed K0s in jet cone are stemming from MC gen K0s (on MCgenK0s list):__________________________________________________
3431 for(Int_t ic=0; ic<jetConeK0list->GetSize(); ++ic){ //loop over all reconstructed K0s in jet cone
3433 //for(Int_t ic=0; ic<fListK0s->GetSize(); ++ic){ //loop over all reconstructed K0s -> previous definition of reconstruction efficiency, not sure what is the better one to choose
3435 Int_t negDaughterpdg;
3436 Int_t posDaughterpdg;
3439 Double_t fPtMCrecK0Match;
3440 Double_t invMK0Match;
3444 Bool_t fPhysicalPrimary = -1;
3445 Int_t MCv0PDGCode =0;
3446 Double_t jetPtSmear = -1;
3448 AliAODv0* v0c = dynamic_cast<AliAODv0*>(jetConeK0list->At(ic));//pointer to reconstructed K0s inside jet cone (cone is placed around reconstructed jet axis)
3450 //AliAODv0* v0c = dynamic_cast<AliAODv0*>(fListK0s->At(ic));//pointer to reconstructed K0s
3453 Bool_t daughtercheck = DaughterTrackCheck(v0c, nnum, pnum);//check daughter tracks have proper sign
3454 if(daughtercheck == kFALSE)continue;
3456 const AliAODTrack *trackMCNeg=(AliAODTrack *)(v0c->GetDaughter(nnum));
3457 const AliAODTrack *trackMCPos=(AliAODTrack *)(v0c->GetDaughter(pnum));
3459 TString generatorName;
3460 TList *listmc = fAOD->GetList();
3462 Bool_t mclabelcheck = MCLabelCheck(v0c, kK0, trackMCNeg, trackMCPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PDGCode, generatorName, isinjected);
3464 if(mclabelcheck == kFALSE)continue;
3465 if(fPhysicalPrimary == kFALSE)continue; //requirements for rec. V0 associated to MC true primary particle
3467 for(Int_t it=0; it<fListMCgenK0s->GetSize(); ++it){ // loop over MC generated K0s in event, check whether associated MC particle is part of it
3469 //for(Int_t it=0; it<fListMCgenK0sCone->GetSize(); ++it){//belongs to previous definition of rec. eff. of V0s within jet cone
3471 //Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3472 //AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenK0sCone->At(it));
3473 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenK0s->At(it));
3476 Bool_t particleMatching = IsParticleMatching(mcp0, v0Label);
3478 if(particleMatching == kFALSE)continue; //if reconstructed V0 particle doesn't match to the associated MC particle go to next stack entry
3479 CalculateInvMass(v0c, kK0, invMK0Match, fPtMCrecK0Match);
3480 Double_t fEta = v0c->Eta();
3481 Double_t fPtMCgenK0s = mcp0->Pt();//pt has to be always MC truth value!
3483 Double_t vMCrecK0Cone[4] = {jetPt, invMK0Match,fPtMCgenK0s,fEta};
3484 fhnMCrecK0Cone->Fill(vMCrecK0Cone); //fill matching rec. K0s in 3D histogram
3486 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear); //jetPt, cent, jetRadius, ptmintrack, &jetPtSmear
3488 Double_t vMCrecK0ConeSmear[4] = {jetPtSmear, invMK0Match,fPtMCgenK0s,fEta};
3489 fhnMCrecK0ConeSmear->Fill(vMCrecK0ConeSmear);
3491 //fill matching rec. K0s in 3D histogram, jet pT smeared according to deltaptjet distribution width
3494 } // end MCgenK0s / MCgenK0sCone loop
3497 //check the K0s daughters contamination of the jet tracks:
3499 TClonesArray *stackMC = 0x0;
3501 for(Int_t it=0; it<jettracklist->GetSize(); ++it){//loop over all particles in the jet
3503 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));//all tracks in jet cone
3504 if(!trackVP)continue;
3505 AliAODTrack *tr = dynamic_cast<AliAODTrack*> (trackVP); //fetch one jet track from the TList
3508 //get MC label information
3509 TList *mclist = fAOD->GetList(); //fetch the MC stack
3510 if(!mclist)continue;
3512 stackMC = (TClonesArray*)mclist->FindObject(AliAODMCParticle::StdBranchName()); //get MCAOD branch in data
3513 if (!stackMC) {Printf("ERROR: stack not available");}
3516 Int_t particleLabel = TMath::Abs(tr->GetLabel()); //fetch jet track label in MC stack
3518 //v0c is pointer to K0s candidate, is fetched already above, here it is just checked again whether daughters are properly ordered by their charge
3520 Bool_t daughterchecks = DaughterTrackCheck(v0c, nnum, pnum);
3522 if(daughterchecks == kFALSE)continue; //make sure that daughters are properly ordered
3524 const AliAODTrack *trackNeg=(AliAODTrack *)(v0c->GetDaughter(nnum)); //fetch v0 daughters of reconstructed K0s
3525 const AliAODTrack *trackPos=(AliAODTrack *)(v0c->GetDaughter(pnum));
3527 if(!trackNeg)continue;
3528 if(!trackPos)continue;
3530 Int_t negAssLabel = TMath::Abs(trackNeg->GetLabel()); //negative (reconstructed) charged track label in MC stack
3531 Int_t posAssLabel = TMath::Abs(trackPos->GetLabel()); //positive (reconstructed) charged track label in MC stack
3534 if(particleLabel == posAssLabel){ //check whether jet track and each of the rec. K0s daughters have same MC stack label -> are identical
3535 AliAODMCParticle* mctrackPos = dynamic_cast<AliAODMCParticle*>(stackMC->At(posAssLabel));
3536 if(!mctrackPos) continue;
3537 Double_t trackPosPt = mctrackPos->Pt();
3538 Double_t trackPosEta = mctrackPos->Eta();
3540 Double_t vK0sSecContinCone[3] = {jetPt, trackPosPt, trackPosEta};
3541 fhnK0sSecContinCone->Fill(vK0sSecContinCone);} //if it's the case, fill jet pt, daughter track pt and track eta in histo
3543 if(particleLabel == negAssLabel){
3544 AliAODMCParticle* mctrackNeg = dynamic_cast<AliAODMCParticle*>(stackMC->At(negAssLabel));
3545 if(!mctrackNeg) continue;
3546 Double_t trackNegPt = mctrackNeg->Pt();
3547 Double_t trackNegEta = mctrackNeg->Eta();
3549 Double_t vK0sSecContinCone[3] = {jetPt, trackNegPt, trackNegEta};
3550 fhnK0sSecContinCone->Fill(vK0sSecContinCone);} //if it's the case, fill jet pt, daughter track pt and track eta in histo
3558 } //end rec-K0-in-cone loop
3560 //________________________________________________________________________________________________________________________________________________________
3562 delete fListMCgenK0sCone;
3567 delete jetConeK0list;
3568 delete jetPerpConeK0list;
3569 delete jetPerpConeLalist;
3570 delete jetPerpConeALalist;
3573 //---------------La--------------------------------------------------------------------------------------------------------------------------------------------
3575 // fQAJetHistosRecCutsLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
3577 for(Int_t it=0; it<fListLa->GetSize(); ++it){ // loop all La
3579 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListLa->At(it));
3584 TVector3 v0MomVect(v0Mom);
3586 Double_t dPhiJetLa = (jet->MomentumVector()->Vect()).DeltaPhi(v0MomVect);
3591 CalculateInvMass(v0, kLambda, invMLa, trackPt); //function to calculate invMass with TLorentzVector class
3592 // Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3594 //if(incrementJetPt){fh1V0JetPt->Fill(jetPt);}
3596 //fFFHistosIMLaJet->FillFF(trackPt, invMLa, jetPt, incrementJetPt);
3598 if(dPhiJetLa<fh1dPhiJetLa->GetXaxis()->GetXmin()) dPhiJetLa += 2*TMath::Pi();
3599 fh1dPhiJetLa->Fill(dPhiJetLa);
3602 /* if(fListLa->GetSize() == 0){ // no La: increment jet pt spectrum
3604 Bool_t incrementJetPt = kTRUE;
3605 fFFHistosIMLaJet->FillFF(-1, -1, jetPt, incrementJetPt);
3609 // ____fetch rec. Lambdas in cone around jet axis_______________________________________________________________________________________
3611 TList* jetConeLalist = new TList();
3612 Double_t sumPtLa = 0.;
3613 Bool_t isBadJetLa = kFALSE; // dummy, do not use
3615 GetTracksInCone(fListLa, jetConeLalist, jet, GetFFRadius(), sumPtLa, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetLa);//method inherited from FF
3617 if(fDebug>2)Printf("%s:%d nLa total: %d, in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nLa,jetConeLalist->GetEntries(),GetFFRadius());
3619 for(Int_t it=0; it<jetConeLalist->GetSize(); ++it){ // loop La in jet cone
3621 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetConeLalist->At(it));
3627 Bool_t daughtercheck = DaughterTrackCheck(v0, nnum, pnum);
3628 if(daughtercheck == kFALSE)continue;
3635 CalculateInvMass(v0, kLambda, invMLa, trackPt); //function to calculate invMass with TLorentzVector class
3637 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;//needed for all histos, which serve for normalisation
3641 Int_t negDaughterpdg;
3642 Int_t posDaughterpdg;
3645 Double_t jetPtSmear = -1;
3647 Bool_t fPhysicalPrimary = -1;
3648 Int_t MCv0PDGCode =0;
3649 TString generatorName;
3651 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
3652 if(incrementJetPt == kTRUE){fh1IMLaConeSmear->Fill(jetPtSmear);
3654 const AliAODTrack *trackNeg=(AliAODTrack *)(v0->GetDaughter(nnum));
3655 const AliAODTrack *trackPos=(AliAODTrack *)(v0->GetDaughter(pnum));
3657 TList *listmc = fAOD->GetList();
3659 Bool_t mclabelcheck = MCLabelCheck(v0, kLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PDGCode, generatorName, isinjected);
3660 if(mclabelcheck == kFALSE)continue;
3662 //std::cout<<"generatorName: "<<generatorName<<std::endl;
3664 if(generatorName == "Hijing"){
3665 Double_t vrecMCHijingLaCone[4] = {jetPt, invMLa,trackPt,fEta};
3666 fhnrecMCHijingLaCone->Fill(vrecMCHijingLaCone);
3669 if(isinjected == kTRUE){
3670 Double_t vrecMCInjectLaCone[4] = {jetPt, invMLa,trackPt,fEta};
3671 fhnrecMCInjectLaCone->Fill(vrecMCInjectLaCone);
3674 }//fill TH1F for normalization purposes
3677 if(incrementJetPt==kTRUE){
3678 fh1IMLaCone->Fill(jetPt);}//normalisation by number of selected jets
3680 //fFFHistosIMLaCone->FillFF(trackPt, invMLa, jetPt, incrementJetPt);
3684 if(jetConeLalist->GetSize() == 0){ // no La: increment jet pt spectrum
3686 Bool_t incrementJetPt = kTRUE;
3687 // fFFHistosIMLaCone->FillFF(-1, -1, jetPt, incrementJetPt);
3688 Double_t vLaCone[4] = {jetPt, -1, -1, -1};
3689 fhnLaCone->Fill(vLaCone);
3691 if(incrementJetPt==kTRUE){
3692 fh1IMLaCone->Fill(jetPt);}//normalisation by number of selected jets
3695 Double_t jetPtSmear;
3696 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
3697 if(incrementJetPt == kTRUE){
3698 fh1IMLaConeSmear->Fill(jetPtSmear);
3707 //____fetch MC generated Lambdas in cone around jet axis__(particles can stem from fragmentation but also from underlying event)_____________
3709 Double_t sumPtMCgenLa = 0.;
3710 Bool_t isBadJetMCgenLa = kFALSE; // dummy, do not use
3712 //sampling MC gen. Lambdas in cone around reconstructed jet axis
3713 fListMCgenLaCone = new TList();
3715 GetTracksInCone(fListMCgenLa, fListMCgenLaCone, jet, GetFFRadius(), sumPtMCgenLa, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetMCgenLa);//fetch MC generated Lambdas in cone of resolution parameter R around jet axis
3717 if(fDebug>2)Printf("%s:%d nMCgenLa in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,fListMCgenLaCone->GetEntries(),GetFFRadius());
3719 for(Int_t it=0; it<fListMCgenLaCone->GetSize(); ++it){ // loop MC generated La in cone around jet axis
3721 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenLaCone->At(it));
3724 //Double_t fRapMCgenLa = MyRapidity(mcp0->E(),mcp0->Pz());
3725 Double_t fEtaMCgenLa = mcp0->Eta();
3726 Double_t fPtMCgenLa = mcp0->Pt();
3728 fh2MCgenLaCone->Fill(jetPt,fPtMCgenLa);
3729 fh2MCEtagenLaCone->Fill(jetPt,fEtaMCgenLa);
3733 //check whether the reconstructed La are stemming from MC gen La on fListMCgenLa List:__________________________________________________
3735 for(Int_t ic=0; ic<jetConeLalist->GetSize(); ++ic){//loop over all reconstructed La within jet cone, new definition
3737 //for(Int_t ic=0; ic<fListLa->GetSize(); ++ic){//old definition
3739 Int_t negDaughterpdg;
3740 Int_t posDaughterpdg;
3743 Double_t fPtMCrecLaMatch;
3744 Double_t invMLaMatch;
3748 Bool_t fPhysicalPrimary = -1;
3749 Int_t MCv0PDGCode =0;
3750 Double_t jetPtSmear = -1;
3751 TString generatorName;
3753 AliAODv0* v0c = dynamic_cast<AliAODv0*>(jetConeLalist->At(ic));//new definition
3758 Bool_t daughtercheck = DaughterTrackCheck(v0c, nnum, pnum);
3759 if(daughtercheck == kFALSE)continue;
3761 const AliAODTrack *trackMCNeg=(AliAODTrack *)(v0c->GetDaughter(nnum));
3762 const AliAODTrack *trackMCPos=(AliAODTrack *)(v0c->GetDaughter(pnum));
3764 TList *listmc = fAOD->GetList();
3766 Bool_t mclabelcheck = MCLabelCheck(v0c, kLambda, trackMCNeg, trackMCPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PDGCode, generatorName, isinjected);
3768 if(mclabelcheck == kFALSE)continue;
3769 if(fPhysicalPrimary == kFALSE)continue;
3771 for(Int_t it=0; it<fListMCgenLa->GetSize(); ++it){//new definition // loop over MC generated K0s in cone around jet axis
3773 // for(Int_t it=0; it<fListMCgenLaCone->GetSize(); ++it){//old definition // loop over MC generated La in cone around jet axis
3775 //Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3777 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenLa->At(it));//new definition
3778 //AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenLaCone->At(it));//old definition
3782 Bool_t particleMatching = IsParticleMatching(mcp0, v0Label);
3785 if(particleMatching == kFALSE)continue; //particle doesn't match on any associated MC gen particle in cone around rec jet axis
3787 CalculateInvMass(v0c, kLambda, invMLaMatch, fPtMCrecLaMatch);
3789 Double_t fPtMCgenLa = mcp0->Pt();
3790 Double_t fEta = v0c->Eta();//rec. MC particle
3791 Double_t vMCrecLaCone[4] = {jetPt, invMLaMatch,fPtMCgenLa,fEta};
3792 fhnMCrecLaCone->Fill(vMCrecLaCone);
3794 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
3796 Double_t vMCrecLaConeSmear[4] = {jetPtSmear, invMLaMatch,fPtMCgenLa,fEta};
3797 fhnMCrecLaConeSmear->Fill(vMCrecLaConeSmear); //fill matching rec. Lambdas in 3D histogram, jet pT smeared according to deltaptjet distribution width
3800 } // end MCgenLa loop
3802 //check the Lambda daughters contamination of the jet tracks://///////////////////////////////////////////////////////////////////////////////////////////
3804 TClonesArray *stackMC = 0x0;
3806 for(Int_t it=0; it<jettracklist->GetSize(); ++it){//loop over all particles in the jet
3808 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));//all tracks in jet cone
3809 if(!trackVP)continue;
3810 AliAODTrack *tr = dynamic_cast<AliAODTrack*> (trackVP); //fetch one jet track from the TList
3813 //get MC label information
3814 TList *mclist = fAOD->GetList(); //fetch the MC stack
3816 stackMC = (TClonesArray*)mclist->FindObject(AliAODMCParticle::StdBranchName()); //get MCAOD branch in data
3817 if (!stackMC) {Printf("ERROR: stack not available");}
3820 Int_t particleLabel = TMath::Abs(tr->GetLabel()); //fetch jet track label in MC stack
3822 Bool_t daughterchecks = DaughterTrackCheck(v0c, nnum, pnum);
3824 if(daughterchecks == kFALSE)continue; //make sure that daughters are properly ordered
3826 const AliAODTrack *trackNeg=(AliAODTrack *)(v0c->GetDaughter(nnum)); //fetch v0 daughters of reconstructed K0s
3827 const AliAODTrack *trackPos=(AliAODTrack *)(v0c->GetDaughter(pnum));
3829 Int_t negAssLabel = TMath::Abs(trackNeg->GetLabel()); //negative (reconstructed) charged track label in MC stack
3830 Int_t posAssLabel = TMath::Abs(trackPos->GetLabel()); //positive (reconstructed) charged track label in MC stack
3833 if(particleLabel == posAssLabel){ //check whether jet track and each of the rec. K0s daughters have same MC stack label -> are identical
3835 AliAODMCParticle* mctrackPos = dynamic_cast<AliAODMCParticle*>(stackMC->At(posAssLabel));
3836 if(!mctrackPos) continue;
3838 Double_t trackPosPt = trackPos->Pt();
3839 Double_t trackPosEta = trackPos->Eta();
3840 Double_t vLaSecContinCone[3] = {jetPt, trackPosPt, trackPosEta};
3841 fhnLaSecContinCone->Fill(vLaSecContinCone);
3843 } //if it's the case, fill jet pt, daughter track pt and track eta in histo
3846 if(particleLabel == negAssLabel){
3848 AliAODMCParticle* mctrackNeg = dynamic_cast<AliAODMCParticle*>(stackMC->At(negAssLabel));
3849 if(!mctrackNeg) continue;
3851 Double_t trackNegPt = trackNeg->Pt();
3852 Double_t trackNegEta = trackNeg->Eta();
3854 Double_t vLaSecContinCone[3] = {jetPt, trackNegPt, trackNegEta};
3855 fhnLaSecContinCone->Fill(vLaSecContinCone);
3858 } //if it's the case, fill jet pt, daughter track pt and track eta in histo
3863 } //end rec-La-in-cone loop
3864 //________________________________________________________________________________________________________________________________________________________
3866 delete fListMCgenLaCone;
3870 delete jetConeLalist;
3874 //---------------ALa-----------
3877 // fQAJetHistosRecCutsLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
3879 for(Int_t it=0; it<fListALa->GetSize(); ++it){ // loop all ALa
3881 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListALa->At(it));
3886 TVector3 v0MomVect(v0Mom);
3888 Double_t dPhiJetALa = (jet->MomentumVector()->Vect()).DeltaPhi(v0MomVect);
3890 Double_t invMALa =0;
3893 CalculateInvMass(v0, kAntiLambda, invMALa, trackPt); //function to calculate invMass with TLorentzVector class
3894 //Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3896 //if(incrementJetPt){fh1V0JetPt->Fill(jetPt);}
3898 //fFFHistosIMALaJet->FillFF(trackPt, invMALa, jetPt, incrementJetPt);
3900 if(dPhiJetALa<fh1dPhiJetALa->GetXaxis()->GetXmin()) dPhiJetALa += 2*TMath::Pi();
3901 fh1dPhiJetALa->Fill(dPhiJetALa);
3904 // if(fListALa->GetSize() == 0){ // no ALa: increment jet pt spectrum
3906 // Bool_t incrementJetPt = kTRUE;
3907 //fFFHistosIMALaJet->FillFF(-1, -1, jetPt, incrementJetPt);
3911 // ____fetch rec. Antilambdas in cone around jet axis_______________________________________________________________________________________
3913 TList* jetConeALalist = new TList();
3914 Double_t sumPtALa = 0.;
3915 Bool_t isBadJetALa = kFALSE; // dummy, do not use
3917 GetTracksInCone(fListALa, jetConeALalist, jet, GetFFRadius(), sumPtALa, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetALa);//method inherited from FF
3919 if(fDebug>2)Printf("%s:%d nALa total: %d, in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nALa,jetConeALalist->GetEntries(),GetFFRadius());
3921 for(Int_t it=0; it<jetConeALalist->GetSize(); ++it){ // loop ALa in jet cone
3923 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetConeALalist->At(it));
3930 Bool_t daughtercheck = DaughterTrackCheck(v0, nnum, pnum);
3931 if(daughtercheck == kFALSE)continue;
3934 Double_t invMALa =0;
3940 CalculateInvMass(v0, kAntiLambda, invMALa, trackPt); //function to calculate invMass with TLorentzVector class
3942 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3944 if(fAnalysisMC){ //jet pt smearing study for Antilambdas
3946 Int_t negDaughterpdg;
3947 Int_t posDaughterpdg;
3950 Double_t jetPtSmear = -1;
3952 Bool_t fPhysicalPrimary = -1;
3953 Int_t MCv0PDGCode =0;
3954 TString generatorName;
3956 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
3957 const AliAODTrack *trackNeg=(AliAODTrack *)(v0->GetDaughter(nnum));
3958 const AliAODTrack *trackPos=(AliAODTrack *)(v0->GetDaughter(pnum));
3960 TList *listmc = fAOD->GetList();
3962 Bool_t mclabelcheck = MCLabelCheck(v0, kAntiLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PDGCode, generatorName, isinjected);
3963 if(mclabelcheck == kFALSE)continue;
3965 //std::cout<<"generatorName: "<<generatorName<<std::endl;
3967 if(generatorName == "Hijing"){
3968 Double_t vrecMCHijingALaCone[4] = {jetPt, invMALa,trackPt,fEta};
3969 fhnrecMCHijingALaCone->Fill(vrecMCHijingALaCone);
3972 if(isinjected == kTRUE){
3973 Double_t vrecMCInjectALaCone[4] = {jetPt, invMALa,trackPt,fEta};
3974 fhnrecMCInjectALaCone->Fill(vrecMCInjectALaCone);
3977 if(incrementJetPt == kTRUE){fh1IMALaConeSmear->Fill(jetPtSmear);} //fill TH1F for normalization purposes
3980 if(incrementJetPt==kTRUE){
3981 fh1IMALaCone->Fill(jetPt);}//normalisation by number of selected jets
3983 //fFFHistosIMALaCone->FillFF(trackPt, invMALa, jetPt, incrementJetPt);
3984 Double_t vALaCone[4] = {jetPt, invMALa,trackPt,fEta};
3985 fhnALaCone->Fill(vALaCone);
3988 if(jetConeALalist->GetSize() == 0){ // no ALa: increment jet pt spectrum
3990 Bool_t incrementJetPt = kTRUE;
3992 if(incrementJetPt==kTRUE){
3993 fh1IMALaCone->Fill(jetPt);}//normalisation by number of selected jets
3995 //fFFHistosIMALaCone->FillFF(-1, -1, jetPt, incrementJetPt);
3996 Double_t vALaCone[4] = {jetPt, -1, -1, -1};
3997 fhnALaCone->Fill(vALaCone);
4000 Double_t jetPtSmear;
4001 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
4002 if(incrementJetPt == kTRUE)fh1IMALaConeSmear->Fill(jetPtSmear);}
4008 //____fetch MC generated Antilambdas in cone around jet axis__(particles can stem from fragmentation but also from underlying event)_____________
4010 Double_t sumPtMCgenALa = 0.;
4011 Bool_t isBadJetMCgenALa = kFALSE; // dummy, do not use
4013 //sampling MC gen Antilambdas in cone around reconstructed jet axis
4014 fListMCgenALaCone = new TList();
4016 GetTracksInCone(fListMCgenALa, fListMCgenALaCone, jet, GetFFRadius(), sumPtMCgenALa, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetMCgenALa);//MC generated K0s in cone around jet axis
4018 if(fDebug>2)Printf("%s:%d nMCgenALa in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,fListMCgenALaCone->GetEntries(),GetFFRadius());
4020 for(Int_t it=0; it<fListMCgenALaCone->GetSize(); ++it){ // loop MC generated La in cone around jet axis
4022 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenALaCone->At(it));
4025 //Double_t fRapMCgenALa = MyRapidity(mcp0->E(),mcp0->Pz());
4026 Double_t fEtaMCgenALa = mcp0->Eta();
4027 Double_t fPtMCgenALa = mcp0->Pt();
4029 fh2MCgenALaCone->Fill(jetPt,fPtMCgenALa);
4030 fh2MCEtagenALaCone->Fill(jetPt,fEtaMCgenALa);
4034 //check whether the reconstructed ALa are stemming from MC gen ALa on MCgenALa List:__________________________________________________
4036 for(Int_t ic=0; ic<jetConeALalist->GetSize(); ++ic){//loop over all reconstructed ALa
4038 Int_t negDaughterpdg;
4039 Int_t posDaughterpdg;
4042 Double_t fPtMCrecALaMatch;
4043 Double_t invMALaMatch;
4047 Bool_t fPhysicalPrimary = -1;
4048 Int_t MCv0PDGCode =0;
4049 Double_t jetPtSmear = -1;
4050 TString generatorName;
4052 AliAODv0* v0c = dynamic_cast<AliAODv0*>(jetConeALalist->At(ic));
4055 Bool_t daughtercheck = DaughterTrackCheck(v0c, nnum, pnum);
4056 if(daughtercheck == kFALSE)continue;
4058 const AliAODTrack *trackMCNeg=(AliAODTrack *)(v0c->GetDaughter(nnum));
4059 const AliAODTrack *trackMCPos=(AliAODTrack *)(v0c->GetDaughter(pnum));
4061 TList *listmc = fAOD->GetList();
4062 if(!listmc)continue;
4064 Bool_t mclabelcheck = MCLabelCheck(v0c, kAntiLambda, trackMCNeg, trackMCPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PDGCode, generatorName, isinjected);
4066 if(mclabelcheck == kFALSE)continue;
4067 if(fPhysicalPrimary == kFALSE)continue;
4069 for(Int_t it=0; it<fListMCgenALa->GetSize(); ++it){ // loop over MC generated Antilambdas in cone around jet axis
4071 //Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4073 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenALa->At(it));
4076 Bool_t particleMatching = IsParticleMatching(mcp0, v0Label);
4078 if(particleMatching == kFALSE)continue; //particle doesn't match on any associated MC gen particle in cone around rec jet axis
4080 CalculateInvMass(v0c, kAntiLambda, invMALaMatch, fPtMCrecALaMatch);
4082 Double_t fPtMCgenALa = mcp0->Pt();
4083 Double_t fEta = v0c->Eta();
4084 Double_t vMCrecALaCone[4] = {jetPt, invMALaMatch,fPtMCgenALa,fEta};
4085 fhnMCrecALaCone->Fill(vMCrecALaCone); //fill matching rec. Antilambda in 3D histogram
4087 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
4089 Double_t vMCrecALaConeSmear[4] = {jetPtSmear, invMALaMatch,fPtMCgenALa,fEta};
4090 fhnMCrecALaConeSmear->Fill(vMCrecALaConeSmear); //fill matching rec. Antilambda in 3D histogram
4092 } // end MCgenALa loop
4096 //check the Antilambda daughters contamination of the jet tracks:
4098 TClonesArray *stackMC = 0x0;
4100 for(Int_t it=0; it<jettracklist->GetSize(); ++it){//loop over all particles in the jet
4102 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));//all tracks in jet cone
4103 if(!trackVP)continue;
4104 AliAODTrack *tr = dynamic_cast<AliAODTrack*> (trackVP); //fetch one jet track from the TList
4107 //get MC label information
4108 TList *mclist = fAOD->GetList(); //fetch the MC stack
4109 if(!mclist)continue;
4111 stackMC = (TClonesArray*)mclist->FindObject(AliAODMCParticle::StdBranchName()); //get MCAOD branch in data
4112 if (!stackMC) {Printf("ERROR: stack not available");}
4115 Int_t particleLabel = TMath::Abs(tr->GetLabel()); //fetch jet track label in MC stack
4117 Bool_t daughterchecks = DaughterTrackCheck(v0c, nnum, pnum);
4119 if(daughterchecks == kFALSE)continue; //make sure that daughters are properly ordered
4121 const AliAODTrack *trackNeg=(AliAODTrack *)(v0c->GetDaughter(nnum)); //fetch v0 daughters of reconstructed K0s
4122 const AliAODTrack *trackPos=(AliAODTrack *)(v0c->GetDaughter(pnum));
4123 if(!trackPos)continue;
4124 if(!trackNeg)continue;
4126 Int_t negAssLabel = TMath::Abs(trackNeg->GetLabel()); //negative (reconstructed) charged track label in MC stack
4127 Int_t posAssLabel = TMath::Abs(trackPos->GetLabel()); //positive (reconstructed) charged track label in MC stack
4129 if(!negAssLabel)continue;
4130 if(!posAssLabel)continue;
4132 if(particleLabel == posAssLabel){ //check whether jet track and each of the rec. K0s daughters have same MC stack label -> are identical
4133 AliAODMCParticle* mctrackPos = dynamic_cast<AliAODMCParticle*>(stackMC->At(posAssLabel));
4134 if(!mctrackPos) continue;
4136 Double_t trackPosPt = trackPos->Pt();
4137 Double_t trackPosEta = trackPos->Eta();
4138 if(!trackPosPt)continue;
4139 if(!trackPosEta)continue;
4141 Double_t vLaSecContinCone[3] = {jetPt, trackPosPt, trackPosEta};
4142 fhnLaSecContinCone->Fill(vLaSecContinCone);
4146 //fh3SecContinCone->Fill(jetPt, trackPosPt, trackPosEta);
4147 } //if it's the case, fill jet pt, daughter track pt and track eta in histo
4149 if(particleLabel == negAssLabel){
4151 AliAODMCParticle* mctrackNeg = dynamic_cast<AliAODMCParticle*>(stackMC->At(negAssLabel));
4152 if(!mctrackNeg) continue;
4154 Double_t trackNegPt = trackNeg->Pt();
4155 Double_t trackNegEta = trackNeg->Eta();
4157 if(!trackNegPt)continue;
4158 if(!trackNegEta)continue;
4160 Double_t vLaSecContinCone[3] = {jetPt, trackNegPt, trackNegEta};
4161 fhnLaSecContinCone->Fill(vLaSecContinCone);
4163 //fh3SecContinCone->Fill(jetPt, trackNegPt, trackNegEta);
4164 } //if it's the case, fill jet pt, daughter track pt and track eta in histo
4168 } //end rec-ALa-in-cone loop
4169 //________________________________________________________________________________________________________________________________________________________
4171 delete fListMCgenALaCone;
4175 delete jetConeALalist;
4176 delete jettracklist; //had been initialised at jet loop beginning
4178 }//end of if 'leading' or 'all jet' requirement
4184 fTracksRecCuts->Clear();
4185 fJetsRecCuts->Clear();
4186 fBckgJetsRec->Clear();
4190 fListFeeddownLaCand->Clear();
4191 fListFeeddownALaCand->Clear();
4192 jetConeFDLalist->Clear();
4193 jetConeFDALalist->Clear();
4194 fListMCgenK0s->Clear();
4195 fListMCgenLa->Clear();
4196 fListMCgenALa->Clear();
4200 PostData(1, fCommonHistList);
4205 // ____________________________________________________________________________________________
4206 void AliAnalysisTaskJetChem::SetProperties(TH3F* h,const char* x, const char* y, const char* z)
4208 //Set properties of histos (x,y and z title)
4213 h->GetXaxis()->SetTitleColor(1);
4214 h->GetYaxis()->SetTitleColor(1);
4215 h->GetZaxis()->SetTitleColor(1);
4219 //________________________________________________________________________________________________________________________________________
4220 Bool_t AliAnalysisTaskJetChem::AcceptBetheBloch(AliAODv0 *v0, AliPIDResponse *PIDResponse, const Int_t particletype) //dont use for MC Analysis
4226 const AliAODTrack *ntracktest=(AliAODTrack *)v0->GetDaughter(nnum);
4227 if(ntracktest->Charge() > 0){nnum = 0; pnum = 1;}
4229 const AliAODTrack *trackNeg=(AliAODTrack *)(v0->GetDaughter(nnum));
4230 const AliAODTrack *trackPos=(AliAODTrack *)(v0->GetDaughter(pnum));
4232 //Check if both tracks are available
4233 if (!trackPos || !trackNeg) {
4234 Printf("strange analysis::UserExec:: Error:Could not retrieve one of the daughter tracks\n");
4238 //remove like sign V0s
4239 if ( trackPos->Charge() == trackNeg->Charge() ){
4240 //if(fDebug>1) Printf("%s:%d found like-sign V0", (char*)__FILE__,__LINE__);
4245 Double_t nsig_p = 0; //number of sigmas that positive daughter track has got in TPC pid information
4246 Double_t nsig_n = 0;
4248 const AliAODPid *pid_p=trackPos->GetDetPid(); // returns fDetPID, more detailed or detector specific pid information
4249 const AliAODPid *pid_n=trackNeg->GetDetPid();
4251 if(!pid_p)return kFALSE;
4252 if(!pid_n)return kFALSE;
4256 if(particletype == 1) //PID cut on positive charged Lambda daughters (only those with pt < 1 GeV/c)
4259 nsig_p=PIDResponse->NumberOfSigmasTPC(trackPos,AliPID::kProton);
4260 Double_t protonPt = trackPos->Pt();
4261 if ((TMath::Abs(nsig_p) >= fCutBetheBloch) && (fCutBetheBloch >0) && (protonPt < 1)) return kFALSE;
4270 if(particletype == 2)
4272 nsig_n=PIDResponse->NumberOfSigmasTPC(trackNeg,AliPID::kProton);
4273 Double_t antiprotonPt = trackNeg->Pt();
4274 if ((TMath::Abs(nsig_n) >= fCutBetheBloch) && (fCutBetheBloch >0) && (antiprotonPt < 1)) return kFALSE;
4282 //___________________________________________________________________
4283 Bool_t AliAnalysisTaskJetChem::IsK0InvMass(const Double_t mass) const
4285 // K0 mass ? Use FF histo limits
4287 if(fFFIMInvMMin <= mass && mass < fFFIMInvMMax) return kTRUE;
4291 //___________________________________________________________________
4292 Bool_t AliAnalysisTaskJetChem::IsLaInvMass(const Double_t mass) const
4294 // La mass ? Use FF histo limits
4297 if(fFFIMLaInvMMin <= mass && mass < fFFIMLaInvMMax) return kTRUE;
4302 //_____________________________________________________________________________________
4303 Int_t AliAnalysisTaskJetChem::GetListOfV0s(TList *list, const Int_t type, const Int_t particletype, AliAODVertex* primVertex, AliAODEvent* aod)
4305 // fill list of V0s selected according to type
4308 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
4313 if(fDebug>5){std::cout<<"AliAnalysisTaskJetChem::GetListOfV0s(): type: "<<type<<" particletype: "<<particletype<<"aod: "<<aod<<std::endl;
4314 if(type==kTrackUndef){std::cout<<"AliAnalysisTaskJetChem::GetListOfV0s(): kTrackUndef!! "<<std::endl;}
4318 if(type==kTrackUndef) return 0;
4320 if(!primVertex) return 0;
4322 Double_t lPrimaryVtxPosition[3];
4323 Double_t lV0Position[3];
4324 lPrimaryVtxPosition[0] = primVertex->GetX();
4325 lPrimaryVtxPosition[1] = primVertex->GetY();
4326 lPrimaryVtxPosition[2] = primVertex->GetZ();
4328 if(fDebug>5){ std::cout<<"AliAnalysisTaskJetChem::GetListOfV0s(): aod->GetNumberOfV0s: "<<aod->GetNumberOfV0s()<<std::endl; }
4331 for(int i=0; i<aod->GetNumberOfV0s(); i++){ // loop over V0s
4334 AliAODv0* v0 = aod->GetV0(i);
4338 std::cout << std::endl
4339 << "Warning in AliAnalysisTaskJetChem::GetListOfV0s:" << std::endl
4340 << "v0 = " << v0 << std::endl;
4344 Bool_t isOnFly = v0->GetOnFlyStatus();
4346 if(!isOnFly && (type == kOnFly || type == kOnFlyPID || type == kOnFlydEdx || type == kOnFlyPrim)) continue;
4347 if( isOnFly && (type == kOffl || type == kOfflPID || type == kOffldEdx || type == kOfflPrim)) continue;
4349 Int_t motherType = -1;
4350 //Double_t v0CalcMass = 0; //mass of MC v0
4351 Double_t MCPt = 0; //pt of MC v0
4353 Double_t pp[3]={0,0,0}; //3-momentum positive charged track
4354 Double_t pm[3]={0,0,0}; //3-momentum negative charged track
4355 Double_t v0mom[3]={0,0,0};
4366 Bool_t daughtercheck = DaughterTrackCheck(v0, nnum, pnum);
4368 if(daughtercheck == kFALSE)continue;
4370 const AliAODTrack *trackNeg=(AliAODTrack *)(v0->GetDaughter(nnum));
4371 const AliAODTrack *trackPos=(AliAODTrack *)(v0->GetDaughter(pnum));
4374 ///////////////////////////////////////////////////////////////////////////////////
4376 //calculate InvMass for every V0 particle assumption (Kaon=1,Lambda=2,Antilambda=3)
4377 switch(particletype){
4379 CalculateInvMass(v0, kK0, invM, trackPt); //function to calculate invMass with TLorentzVector class
4383 CalculateInvMass(v0, kLambda, invM, trackPt);
4387 CalculateInvMass(v0, kAntiLambda, invM, trackPt);
4391 std::cout<<"***NO VALID PARTICLETYPE***"<<std::endl;
4396 /////////////////////////////////////////////////////////////
4397 //V0 and track Cuts:
4400 if(fDebug>7){if(!(IsK0InvMass(invMK0s)) && !(IsLaInvMass(invMLa)) && !(IsLaInvMass(invMALa))){std::cout<<"AliAnalysisTaskJetChem::GetListOfV0s: invM not in selected mass window "<<std::endl;}}
4402 if(!(IsK0InvMass(invMK0s)) && !(IsLaInvMass(invMLa)) && !(IsLaInvMass(invMALa)))continue;
4404 // Double_t PosEta = trackPos->AliAODTrack::Eta();//daughter track charge is sometimes wrong here, account for that!!!
4405 // Double_t NegEta = trackNeg->AliAODTrack::Eta();
4407 Double_t PosEta = trackPos->Eta();//daughter track charge is sometimes wrong here, account for that!!!
4408 Double_t NegEta = trackNeg->Eta();
4410 Double_t PosCharge = trackPos->Charge();
4411 Double_t NegCharge = trackNeg->Charge();
4413 if((trackPos->Charge() == 1) && (trackNeg->Charge() == -1)) //Fill daughters charge into histo to check if they are symmetric distributed
4414 { fh1PosDaughterCharge->Fill(PosCharge);
4415 fh1NegDaughterCharge->Fill(NegCharge);
4418 //DistOverTotMom_in_2D___________
4420 Float_t fMassK0s = TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass();
4421 Float_t fMassLambda = TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass();
4424 AliAODVertex* primVtx = fAOD->GetPrimaryVertex(); // get the primary vertex
4425 Double_t dPrimVtxPos[3]; // primary vertex position {x,y,z}
4426 primVtx->GetXYZ(dPrimVtxPos);
4428 Float_t fPtV0 = TMath::Sqrt(v0->Pt2V0()); // transverse momentum of V0
4429 Double_t dSecVtxPos[3]; // V0 vertex position {x,y,z}
4430 v0->GetSecondaryVtx(dSecVtxPos);
4431 Double_t dDecayPath[3];
4432 for (Int_t iPos = 0; iPos < 3; iPos++)
4433 dDecayPath[iPos] = dSecVtxPos[iPos]-dPrimVtxPos[iPos]; // vector of the V0 path
4434 Float_t fDecLen2D = TMath::Sqrt(dDecayPath[0]*dDecayPath[0]+dDecayPath[1]*dDecayPath[1]); //transverse path length R
4435 Float_t fROverPt = fDecLen2D/fPtV0; // R/pT
4437 Float_t fMROverPtK0s = fMassK0s*fROverPt; // m*R/pT
4438 Float_t fMROverPtLambda = fMassLambda*fROverPt; // m*R/pT
4440 //___________________
4441 //Double_t fRap = -999;//init values
4442 Double_t fEta = -999;
4443 Double_t fV0cosPointAngle = -999;
4444 Double_t fV0DecayLength = v0->DecayLengthV0(lPrimaryVtxPosition);
4448 fV0mom[0]=v0->MomV0X();
4449 fV0mom[1]=v0->MomV0Y();
4450 fV0mom[2]=v0->MomV0Z();
4452 Double_t fV0TotalMomentum = TMath::Sqrt(fV0mom[0]*fV0mom[0]+fV0mom[1]*fV0mom[1]+fV0mom[2]*fV0mom[2]);
4453 // const Double_t K0sPDGmass = 0.497614;
4454 // const Double_t LambdaPDGmass = 1.115683;
4456 const Double_t K0sPDGmass = TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass();
4457 const Double_t LambdaPDGmass = TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass();
4459 Double_t fDistOverTotMomK0s = 0;
4460 Double_t fDistOverTotMomLa = 0;
4462 //calculate proper lifetime of particles in 3D (not recommended anymore)
4464 if(particletype == kK0){
4466 fDistOverTotMomK0s = fV0DecayLength * K0sPDGmass;
4467 fDistOverTotMomK0s /= (fV0TotalMomentum+1e-10);
4470 if((particletype == kLambda)||(particletype == kAntiLambda)){
4472 fDistOverTotMomLa = fV0DecayLength * LambdaPDGmass;
4473 fDistOverTotMomLa /= (fV0TotalMomentum+1e-10);
4476 //TPC cluster (not used anymore) and TPCRefit cuts
4478 //Double_t trackPosNcls = trackPos->GetTPCNcls();//Get number of clusters for positive charged tracks
4479 //Double_t trackNegNcls = trackNeg->GetTPCNcls();//Get number of clusters for negative charged tracks
4481 if(fRequireTPCRefit==kTRUE){//if kTRUE: accept only if daughter track is refitted in TPC!!
4482 Bool_t isPosTPCRefit = (trackPos->AliAODTrack::IsOn(AliESDtrack::kTPCrefit));
4483 Bool_t isNegTPCRefit = (trackNeg->AliAODTrack::IsOn(AliESDtrack::kTPCrefit));
4484 if (!isPosTPCRefit)continue;
4485 if (!isNegTPCRefit)continue;
4488 if(fKinkDaughters==kFALSE){//if kFALSE: no acceptance of kink daughters
4489 AliAODVertex* ProdVtxDaughtersPos = (AliAODVertex*) (trackPos->AliAODTrack::GetProdVertex());
4490 Char_t isAcceptKinkDaughtersPos = ProdVtxDaughtersPos->GetType();
4491 if(isAcceptKinkDaughtersPos==AliAODVertex::kKink)continue;
4493 AliAODVertex* ProdVtxDaughtersNeg = (AliAODVertex*) (trackNeg->AliAODTrack::GetProdVertex());
4494 Char_t isAcceptKinkDaughtersNeg = ProdVtxDaughtersNeg->GetType();
4495 if(isAcceptKinkDaughtersNeg==AliAODVertex::kKink)continue;
4499 Double_t fV0Radius = -999;
4500 Double_t fDcaV0Daughters = v0->DcaV0Daughters();
4501 Double_t fDcaPosToPrimVertex = v0->DcaPosToPrimVertex();//IP of positive charged daughter
4502 Double_t fDcaNegToPrimVertex = v0->DcaNegToPrimVertex();//IP of negative charged daughter
4503 Double_t avDecayLengthK0s = 2.6844;
4504 Double_t avDecayLengthLa = 7.89;
4506 //Float_t fCTauK0s = 2.6844; // [cm] c tau of K0S
4507 //Float_t fCTauLambda = 7.89; // [cm] c tau of Lambda and Antilambda
4509 fV0cosPointAngle = v0->CosPointingAngle(lPrimaryVtxPosition);
4510 lV0Position[0]= v0->DecayVertexV0X();
4511 lV0Position[1]= v0->DecayVertexV0Y();
4512 lV0Position[2]= v0->DecayVertexV0Z();
4514 fV0Radius = TMath::Sqrt(lV0Position[0]*lV0Position[0]+lV0Position[1]*lV0Position[1]);
4516 if(particletype == kK0) {//fRap = v0->RapK0Short();
4517 fEta = v0->PseudoRapV0();}
4518 if(particletype == kLambda) {//fRap = v0->RapLambda();
4519 fEta = v0->PseudoRapV0();}
4520 if(particletype == kAntiLambda) {//fRap = v0->Y(-3122);
4521 fEta = v0->PseudoRapV0();}
4524 //cut on 3D DistOverTotMom: (not used anymore)
4525 //if((particletype == kLambda)||(particletype == kAntiLambda)){if(fDistOverTotMomLa >= (fCutV0DecayMax * avDecayLengthLa)) continue;}
4527 //cut on K0s applied below after all other cuts for histo fill purposes..
4529 //cut on 2D DistOverTransMom: (recommended from Iouri)
4530 if((particletype == kLambda)||(particletype == kAntiLambda)){if(fMROverPtLambda > (fCutV0DecayMax * avDecayLengthLa))continue;}//fCutV0DecayMax set to 5 in AddTask macro
4532 //Armenteros Podolanski Plot for K0s:////////////////////////////
4534 Double_t ArmenterosAlpha=-999;
4535 Double_t ArmenterosPt=-999;
4541 if(particletype == kK0){
4543 pp[0]=v0->MomPosX();
4544 pp[1]=v0->MomPosY();
4545 pp[2]=v0->MomPosZ();
4547 pm[0]=v0->MomNegX();
4548 pm[1]=v0->MomNegY();
4549 pm[2]=v0->MomNegZ();
4552 v0mom[0]=v0->MomV0X();
4553 v0mom[1]=v0->MomV0Y();
4554 v0mom[2]=v0->MomV0Z();
4556 TVector3 v0Pos(pp[0],pp[1],pp[2]);
4557 TVector3 v0Neg(pm[0],pm[1],pm[2]);
4558 TVector3 v0totMom(v0mom[0], v0mom[1], v0mom[2]); //vector for tot v0 momentum
4560 //PosPt = v0Pos.Perp(v0totMom); //longitudinal momentum of positive charged daughter track
4561 PosPl = v0Pos.Dot(v0totMom)/v0totMom.Mag(); //transversal momentum of positive charged daughter track
4563 //NegPt = v0Neg.Perp(v0totMom); //longitudinal momentum of negative charged daughter track
4564 NegPl = v0Neg.Dot(v0totMom)/v0totMom.Mag(); //transversal momentum of nergative charged daughter track
4566 ArmenterosAlpha = 1.-2./(1+(PosPl/NegPl));
4567 ArmenterosPt= v0->PtArmV0();
4571 if(particletype == kK0){//only cut on K0s histos
4572 if(IsArmenterosSelected == 1){// Armenteros Cut to reject Lambdas contamination in K0s inv. massspectrum
4573 fh2ArmenterosBeforeCuts->Fill(ArmenterosAlpha,ArmenterosPt);
4577 //some more cuts on v0s and daughter tracks:
4580 if((TMath::Abs(PosEta)>fCutPostrackEta) || (TMath::Abs(NegEta)>fCutNegtrackEta))continue; //Daughters pseudorapidity cut
4581 if (fV0cosPointAngle < fCutV0cosPointAngle) continue; //cospointangle cut
4583 //if(TMath::Abs(fRap) > fCutRap)continue; //V0 Rapidity Cut
4584 if(TMath::Abs(fEta) > fCutEta) continue; //V0 Eta Cut
4585 if (fDcaV0Daughters > fCutDcaV0Daughters)continue;
4586 if ((fDcaPosToPrimVertex < fCutDcaPosToPrimVertex) || (fDcaNegToPrimVertex < fCutDcaNegToPrimVertex))continue;
4587 if ((fV0Radius < fCutV0RadiusMin) || (fV0Radius > fCutV0RadiusMax))continue;
4589 const AliAODPid *pid_p1=trackPos->GetDetPid();
4590 const AliAODPid *pid_n1=trackNeg->GetDetPid();
4593 if(particletype == kLambda){
4594 // if(AcceptBetheBloch(v0, fPIDResponse, 1) == kFALSE){std::cout<<"******PID cut rejects Lambda!!!************"<<std::endl;}
4595 if(AcceptBetheBloch(v0, fPIDResponse, 1) == kFALSE)continue;
4596 fh2BBLaPos->Fill(pid_p1->GetTPCmomentum(),pid_p1->GetTPCsignal());//positive lambda daughter
4597 fh2BBLaNeg->Fill(pid_n1->GetTPCmomentum(),pid_n1->GetTPCsignal());//negative lambda daughter
4599 //Double_t phi = v0->Phi();
4600 //Double_t massLa = v0->MassLambda();
4602 //printf("La: i = %d, m = %f, pT = %f, eta = %f, phi = %f\n, ",i,massLa,trackPt,fEta,phi);
4606 if(particletype == kAntiLambda){
4608 if(AcceptBetheBloch(v0, fPIDResponse, 2) == kFALSE)continue;
4609 fh2BBLaPos->Fill(pid_p1->GetTPCmomentum(),pid_p1->GetTPCsignal());//positive antilambda daughter
4610 fh2BBLaNeg->Fill(pid_n1->GetTPCmomentum(),pid_n1->GetTPCsignal());//negative antilambda daughter
4615 //Armenteros cut on K0s:
4616 if(particletype == kK0){
4617 if(IsArmenterosSelected == 1){// Armenteros Cut to reject Lambdas contamination in K0s inv. massspectrum
4619 if((ArmenterosPt<=(TMath::Abs(fCutArmenteros*ArmenterosAlpha))) && (fCutArmenteros!=-999))continue; //Cuts out Lambda contamination in K0s histos
4620 fh2ArmenterosAfterCuts->Fill(ArmenterosAlpha,ArmenterosPt);
4624 //not used anymore in 3D, z component of total momentum has bad resolution, cut instead in 2D and use pT
4625 //Proper Lifetime Cut: DecayLength3D * PDGmass / |p_tot| < 3*2.68cm (ctau(betagamma=1)) ; |p|/mass = beta*gamma
4626 //////////////////////////////////////////////
4629 //cut on 2D DistOverTransMom
4630 if(particletype == kK0){//the cut on Lambdas you can find above
4632 fh2ProperLifetimeK0sVsPtBeforeCut->Fill(trackPt,fMROverPtK0s); //fill these histos after all other cuts
4633 if(fMROverPtK0s > (fCutV0DecayMax * avDecayLengthK0s))continue;
4634 fh2ProperLifetimeK0sVsPtAfterCut->Fill(trackPt,fMROverPtK0s);
4636 //Double_t phi = v0->Phi();
4637 // Double_t massK0s = v0->MassK0Short();
4638 //printf("K0S: i = %d, m = %f, pT = %f, eta = %f, phi = %f\n",i,invMK0s,trackPt,fEta,phi);
4640 //test std::cout<<" Index accepted K0s candidate in list of V0s in event: "<<i<<" m: "<<invMK0s<<" pT: "<<trackPt<<" eta: "<<fEta<<" phi: "<<v0->Phi()<<std::endl;
4643 //MC Associated V0 particles: (reconstructed particles associated with MC truth (MC truth: true primary MC generated particle))
4646 if(fAnalysisMC){// begin MC part
4648 Int_t negDaughterpdg = 0;
4649 Int_t posDaughterpdg = 0;
4651 Bool_t fPhysicalPrimary = -1; //v0 physical primary check
4652 Int_t MCv0PdgCode = 0;
4653 Bool_t mclabelcheck = kFALSE;
4655 TList *listmc = aod->GetList(); //AliAODEvent* is inherited from AliVEvent*, listmc is pointer to reconstructed event in MC list, member of AliAODEvent
4657 if(!listmc)continue;
4659 if((particletype == kLambda) || (particletype == kAntiLambda)){// at this point the v0 candidates already survived all V0 cuts, for the MC analysis they still have to survive the association checks in the following block
4661 //feeddown-correction for Lambda/Antilambda particles
4662 //feedddown comes mainly from charged and neutral Xi particles
4663 //feeddown from Sigma decays so quickly that it's not possible to distinguish from primary Lambdas with detector
4664 //feeddown for K0s from phi decays is neglectible
4665 //TH2F* fh2FeedDownMatrix = 0x0; //histo for feeddown already decleared above
4668 //first for all Lambda and Antilambda candidates____________________________________________________________________
4669 TString generatorName;
4671 if(particletype == kLambda){
4673 mclabelcheck = MCLabelCheck(v0, kLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
4676 if((motherType == 3312)||(motherType == 3322)){//mother of v0 is neutral or negative Xi
4678 fListFeeddownLaCand->Add(v0); //fill TList with ass. particles, stemming from feeddown from Xi(bar) decays
4682 if(particletype == kAntiLambda){
4684 mclabelcheck = MCLabelCheck(v0, kAntiLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
4686 if((motherType == -3312)||(motherType == -3322)){
4687 fListFeeddownALaCand->Add(v0); //fill TList with ass. particles, stemming from feeddown from Xi(bar) decays
4693 //_only true primary particles survive the following checks_______________________________________________________________________________________________
4694 TString generatorName;
4696 if(particletype == kK0){
4698 mclabelcheck = MCLabelCheck(v0, kK0, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
4699 if(mclabelcheck == kFALSE)continue;
4701 if(particletype == kLambda){
4703 mclabelcheck = MCLabelCheck(v0, kLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
4704 if(mclabelcheck == kFALSE)continue;
4706 if(particletype == kAntiLambda){
4708 mclabelcheck = MCLabelCheck(v0, kAntiLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
4709 if(mclabelcheck == kFALSE)continue;
4712 if(fPhysicalPrimary != 1)continue; //V0 candidate (K0s, Lambda or Antilambda) must be physical primary, this means there is no mother particle existing
4721 Int_t nPart=list->GetSize();
4724 } // end GetListOfV0s()
4726 // -------------------------------------------------------------------------------------------------------
4728 void AliAnalysisTaskJetChem::CalculateInvMass(AliAODv0* v0vtx, const Int_t particletype, Double_t& invM, Double_t& trackPt){
4738 Double_t pp[3]={0,0,0}; //3-momentum positive charged track
4739 Double_t pm[3]={0,0,0}; //3-momentum negative charged track
4741 const Double_t massPi = 0.13957018; //better use PDG code at this point
4742 const Double_t massP = 0.93827203;
4747 TLorentzVector vector; //lorentzvector V0 particle
4748 TLorentzVector fourmom1;//lorentzvector positive daughter
4749 TLorentzVector fourmom2;//lorentzvector negative daughter
4751 //--------------------------------------------------------------
4753 AliAODTrack *trackPos = (AliAODTrack *) (v0vtx->GetSecondaryVtx()->GetDaughter(0));//index 0 defined as positive charged track in AliESDFilter
4755 if( trackPos->Charge() == 1 ){
4757 pp[0]=v0vtx->MomPosX();
4758 pp[1]=v0vtx->MomPosY();
4759 pp[2]=v0vtx->MomPosZ();
4761 pm[0]=v0vtx->MomNegX();
4762 pm[1]=v0vtx->MomNegY();
4763 pm[2]=v0vtx->MomNegZ();
4766 if( trackPos->Charge() == -1 ){
4768 pm[0]=v0vtx->MomPosX();
4769 pm[1]=v0vtx->MomPosY();
4770 pm[2]=v0vtx->MomPosZ();
4772 pp[0]=v0vtx->MomNegX();
4773 pp[1]=v0vtx->MomNegY();
4774 pp[2]=v0vtx->MomNegZ();
4777 if (particletype == kK0){ // case K0s
4778 mass1 = massPi;//positive particle
4779 mass2 = massPi;//negative particle
4780 } else if (particletype == kLambda){ // case Lambda
4781 mass1 = massP;//positive particle
4782 mass2 = massPi;//negative particle
4783 } else if (particletype == kAntiLambda){ //case AntiLambda
4784 mass1 = massPi;//positive particle
4785 mass2 = massP; //negative particle
4788 fourmom1.SetXYZM(pp[0],pp[1],pp[2],mass1);//positive track
4789 fourmom2.SetXYZM(pm[0],pm[1],pm[2],mass2);//negative track
4790 vector=fourmom1 + fourmom2;
4793 trackPt = vector.Pt();
4795 /*// don't apply AliAODv0 methods to get the inv. mass for the OnFly finder, since the daughter labels are sometimes switched!!!! For Offline V0 finder no problem
4797 if(particletype == kK0){
4798 std::cout << "invMK0s: " << invM <<std::endl;
4799 std::cout << "v0vtx->MassK0Short(): " << v0vtx->MassK0Short() << std::endl;
4800 std::cout << " " <<std::endl;
4801 //invM = v0vtx->MassK0Short();
4804 if(particletype == kLambda){
4805 std::cout << "invMLambda: " << invM <<std::endl;
4806 std::cout << "v0vtx->MassMassLambda(): " << v0vtx->MassLambda() << std::endl;
4807 std::cout << " " <<std::endl;
4808 //invM = v0vtx->MassLambda();
4811 if(particletype == kAntiLambda){
4812 std::cout << "invMAntiLambda: " << invM <<std::endl;
4813 std::cout << "v0vtx->MassAntiLambda(): " << v0vtx->MassAntiLambda() << std::endl;
4814 std::cout << " " <<std::endl;
4815 //invM = v0vtx->MassAntiLambda();
4823 //_____________________________________________________________________________________
4824 Int_t AliAnalysisTaskJetChem::GetListOfMCParticles(TList *outputlist, const Int_t particletype, AliAODEvent *mcaodevent) //(list to fill here e.g. fListMCgenK0s, particle species to search for)
4827 outputlist->Clear();
4829 TClonesArray *stack = 0x0;
4830 Double_t mcXv=0., mcYv=0., mcZv=0.;//MC vertex position
4833 // get MC generated particles
4835 Int_t fPdgcodeCurrentPart = 0; //pdg code current particle
4836 //Double_t fRapCurrentPart = 0; //get rapidity
4837 //Double_t fPtCurrentPart = 0; //get transverse momentum
4838 Double_t fEtaCurrentPart = 0; //get pseudorapidity
4840 //variable for check: physical primary particle
4841 //Bool_t IsPhysicalPrimary = -1;
4842 //Int_t index = 0; //check number of injected particles
4843 //****************************
4844 // Start loop over MC particles
4846 TList *lst = mcaodevent->GetList();
4849 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
4853 stack = (TClonesArray*)lst->FindObject(AliAODMCParticle::StdBranchName());
4855 Printf("ERROR: stack not available");
4859 AliAODMCHeader *mcHdr=(AliAODMCHeader*)lst->FindObject(AliAODMCHeader::StdBranchName());
4860 if(!mcHdr)return -1;
4862 mcXv=mcHdr->GetVtxX(); mcYv=mcHdr->GetVtxY(); mcZv=mcHdr->GetVtxZ(); // position of the MC primary vertex
4865 ntrk=stack->GetEntriesFast();
4867 //if(TMath::Abs(mcZv)>10)return; //i also cut at the reconstructed particles - here i also want to cut for a second time on z vertex (?) -> could be possible bias because of resolution effects on edges of acceptance, also the case for pseudorapidity...
4870 for (Int_t iMc = 0; iMc < ntrk; iMc++) { //loop over mc generated particles
4873 AliAODMCParticle *p0=(AliAODMCParticle*)stack->UncheckedAt(iMc);
4875 //Printf("ERROR: particle with label %d not found in stack (mc loop)", iMc);
4878 fPdgcodeCurrentPart = p0->GetPdgCode();
4880 // Keep only K0s, Lambda and AntiLambda, Xi and Phi:
4881 //if ( (fPdgcodeCurrentPart != 310 ) && (fPdgcodeCurrentPart != 3122 ) && (fPdgcodeCurrentPart != -3122 ) && (fPdgcodeCurrentPart != 3312 ) && (fPdgcodeCurrentPart != -3312) && (fPdgcodeCurrentPart != -333) ) continue;
4885 //Rejection of Pythia injected particles with David Chinellatos method - not the latest method, better Method with TString from MC generator in IsInjected() function below!
4887 /* if( (p0->GetStatus()==21) ||
4888 ((p0->GetPdgCode() == 443) &&
4889 (p0->GetMother() == -1) &&
4890 (p0->GetDaughter(0) == (iMc))) ){ index++; }
4892 if(p0->GetStatus()==21){std::cout<< "hello !!!!" <<std::endl;}
4894 std::cout<< "MC particle status: " << p0->GetStatus() <<std::endl;
4898 //if(index>=1){std::cout<< "MC particle status: " << p0->GetStatus() <<std::endl;}//if first injected MC particle was found, the Status is printed out for this and every following MC particle
4901 //injected particles could be from GenBox (single high pt tracks) or jet related tracks, both generated from PYTHIA MC generator
4903 //Check: MC particle mother
4905 //for feed-down checks
4906 /* //MC gen particles
4907 Int_t iMother = p0->GetMother(); //Motherparticle of V0 candidate (e.g. phi particle,..)
4909 AliAODMCParticle *partM = (AliAODMCParticle*)stack->UncheckedAt(iMother);
4911 if(partM) codeM = TMath::Abs(partM->GetPdgCode());
4914 3312 Xi- -3312 Xibar+
4915 3322 Xi0 -3322 Xibar0
4918 if((codeM == 3312)||(codeM == 3322))// feeddown for Lambda coming from Xi- and Xi0
4924 /* //Check: MC gen. particle decays via 2-pion decay? -> only to be done for the rec. particles !! (-> branching ratio ~ 70 % for K0s -> pi+ pi-)
4926 Int_t daughter0Label = p0->GetDaughter(0);
4927 AliAODMCParticle *mcDaughter0 = (AliAODMCParticle *)stack->UncheckedAt(daughter0Label);
4928 if(daughter0Label >= 0)
4929 {daughter0Type = mcDaughter0->GetPdgCode();}
4931 Int_t daughter1Label = p0->GetDaughter(1);
4932 AliAODMCParticle *mcDaughter1 = (AliAODMCParticle *)stack->UncheckedAt(daughter1Label);
4934 if(daughter1Label >= 1)
4935 {daughter1Type = mcDaughter1->GetPdgCode();} //requirement that daughters are pions is only done for the reconstructed V0s in GetListofV0s() below
4940 // Keep only K0s, Lambda and AntiLambda:
4941 if ( (fPdgcodeCurrentPart != 310 ) && (fPdgcodeCurrentPart != 3122 ) && (fPdgcodeCurrentPart != -3122 )) continue;
4942 // Check: Is physical primary
4944 //do not use anymore: //IsPhysicalPrimary = p0->IsPhysicalPrimary();
4945 //if(!IsPhysicalPrimary)continue;
4947 Float_t fDistPrimaryMax = 0.01; // [cm] max distance of production point to the primary vertex (criterion for choice of MC particles considered as primary)
4949 // Get the distance between production point of the MC mother particle and the primary vertex
4951 Double_t dx = mcXv-p0->Xv();//mc primary vertex - mc gen. v0 vertex
4952 Double_t dy = mcYv-p0->Yv();
4953 Double_t dz = mcZv-p0->Zv();
4955 Double_t fDistPrimary = TMath::Sqrt(dx*dx + dy*dy + dz*dz);
4956 Bool_t fPhysicalPrimary = (fDistPrimary < fDistPrimaryMax);
4958 if(!fPhysicalPrimary)continue;
4960 //if(fPhysicalPrimary){std::cout<<"hello**********************"<<std::endl;}
4962 /* std::cout<<"dx: "<<dx<<std::endl;
4963 std::cout<<"dy: "<<dy<<std::endl;
4964 std::cout<<"dz: "<<dz<<std::endl;
4966 std::cout<<"start: "<<std::endl;
4967 std::cout<<"mcXv: "<<mcXv<<std::endl;
4968 std::cout<<"mcYv: "<<mcYv<<std::endl;
4969 std::cout<<"mcZv: "<<mcZv<<std::endl;
4971 std::cout<<"p0->Xv(): "<<p0->Xv()<<std::endl;
4972 std::cout<<"p0->Yv(): "<<p0->Yv()<<std::endl;
4973 std::cout<<"p0->Zv(): "<<p0->Zv()<<std::endl;
4974 std::cout<<" "<<std::endl;
4975 std::cout<<"fDistPrimary"<<fDistPrimary<<std::endl;
4976 std::cout<<"fPhysicalPrimary"<<fPhysicalPrimary<<std::endl;
4979 //Is close enough to primary vertex to be considered as primary-like?
4981 //fRapCurrentPart = MyRapidity(p0->E(),p0->Pz());
4982 fEtaCurrentPart = p0->Eta();
4983 //fPtCurrentPart = p0->Pt();
4985 if (TMath::Abs(fEtaCurrentPart) < fCutEta){
4986 // if (TMath::Abs(fRapCurrentPart) > fCutRap)continue; //rap cut for crosschecks
4988 if(particletype == kK0){ //MC gen. K0s
4989 if (fPdgcodeCurrentPart==310){
4990 outputlist->Add(p0);
4994 if(particletype == kLambda){ //MC gen. Lambdas
4995 if (fPdgcodeCurrentPart==3122) {
4996 outputlist->Add(p0);
5000 if(particletype == kAntiLambda){
5001 if (fPdgcodeCurrentPart==-3122) { //MC gen. Antilambdas
5002 outputlist->Add(p0);
5007 }//end loop over MC generated particle
5009 Int_t nMCPart=outputlist->GetSize();
5016 //---------------------------------------------------------------------------
5018 Bool_t AliAnalysisTaskJetChem::FillFeeddownMatrix(TList* fListFeeddownCand, Int_t particletype)
5021 // Define Feeddown matrix
5022 Double_t lFeedDownMatrix [100][100];
5023 // FeedDownMatrix [Lambda Bin][Xi Bin];
5025 //Initialize entries of matrix:
5026 for(Int_t ilb = 0; ilb<100; ilb++){
5027 for(Int_t ixb = 0; ixb<100; ixb++){
5028 lFeedDownMatrix[ilb][ixb]=0; //first lambda bins, xi bins
5033 //----------------------------------------------------------------------------
5035 Double_t AliAnalysisTaskJetChem::MyRapidity(Double_t rE, Double_t rPz) const
5037 // Local calculation for rapidity
5038 return 0.5*TMath::Log((rE+rPz)/(rE-rPz+1.e-13));
5040 //----------------------------------------------------------------------------
5042 // ________________________________________________________________________________________________________________________//function to get the MC gen. jet particles
5044 void AliAnalysisTaskJetChem::GetTracksInCone(TList* inputlist, TList* outputlist, const AliAODJet* jet,
5045 const Double_t radius, Double_t& sumPt, const Double_t minPt, const Double_t maxPt, Bool_t& isBadPt)
5047 // fill list of tracks in cone around jet axis
5050 Bool_t isBadMaxPt = kFALSE;
5051 Bool_t isBadMinPt = kTRUE;
5055 jet->PxPyPz(jetMom);
5056 TVector3 jet3mom(jetMom);
5058 //if(jetets < jetetscutr)continue;
5060 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
5062 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
5064 Double_t trackMom[3];
5065 track->PxPyPz(trackMom);
5066 TVector3 track3mom(trackMom);
5068 Double_t dR = jet3mom.DeltaR(track3mom);
5072 outputlist->Add(track);
5074 sumPt += track->Pt();
5076 if(maxPt>0 && track->Pt()>maxPt) isBadMaxPt = kTRUE; // reject jets containing any track with pt larger than this value, use GetFFMaxTrackPt()
5077 if(minPt>0 && track->Pt()>minPt) isBadMinPt = kFALSE; // reject jets with leading track with pt smaller than this value, use GetFFMinLTrackPt()
5083 if(minPt>0 && isBadMinPt) isBadPt = kTRUE; //either the jet is bad because of too small leading track pt.. (probability to be purely combinatorial jet is too high to accept it)
5084 if(maxPt>0 && isBadMaxPt) isBadPt = kTRUE; //..or because of leading track with too high pt (could be fake track)
5090 //____________________________________________________________________________________________________________________
5093 void AliAnalysisTaskJetChem::GetTracksInPerpCone(TList* inputlist, TList* outputlist, const AliAODJet* jet,
5094 const Double_t radius, Double_t& sumPerpPt)
5096 // fill list of tracks in two cones around jet axis rotated in phi +/- 90 degrees
5098 Double_t jetMom[3]; //array for entries in TVector3
5099 Double_t perpjetplusMom[3]; //array for entries in TVector3
5100 Double_t perpjetnegMom[3];
5104 jet->PxPyPz(jetMom); //get 3D jet momentum
5105 Double_t jetPerpPt = jet->Pt(); //original jet pt, invariant under rotations
5106 Double_t jetPhi = jet->Phi(); //original jet phi
5108 Double_t jetPerpposPhi = jetPhi + ((TMath::Pi())*0.5);//get new perp. jet axis phi clockwise
5109 Double_t jetPerpnegPhi = jetPhi - ((TMath::Pi())*0.5);//get new perp. jet axis phi counterclockwise
5111 TVector3 jet3mom(jetMom); //3-Vector for original rec. jet axis
5113 //Double_t phitest = jet3mom.Phi();
5115 perpjetplusMom[0]=(TMath::Cos(jetPerpposPhi)*jetPerpPt); //x coordinate (sidewards - when looking in beam direction)
5116 perpjetplusMom[1]=(TMath::Sin(jetPerpposPhi)*jetPerpPt); //y coordinate (upwards - when looking in beam direction)
5117 perpjetplusMom[2]=jetMom[2]; //z coordinate (along beam axis), invariant under azimuthal rotation
5119 perpjetnegMom[0]=(TMath::Cos(jetPerpnegPhi)*jetPerpPt); //x coordinate (sidewards - when looking in beam direction)
5120 perpjetnegMom[1]=(TMath::Sin(jetPerpnegPhi)*jetPerpPt); //y coordinate (upwards - when looking in beam direction)
5121 perpjetnegMom[2]=jetMom[2]; //z coordinate (along beam axis), invariant under azimuthal rotation
5124 TVector3 perpjetplus3mom(perpjetplusMom); //3-Vector for new perp. jet axis, clockwise rotated
5125 TVector3 perpjetneg3mom(perpjetnegMom); //3-Vector for new perp. jet axis, counterclockwise rotated
5127 //for crosscheck TVector3 rotation method
5129 //Double_t jetMomplusTest[3];
5130 //Double_t jetMomminusTest[3];
5132 //jet3mom.RotateZ(TMath::Pi()*0.5);//rotate original jet axis around +90 degrees in phi
5134 //perpjetminus3momTest = jet3mom.RotateZ((-1)*TMath::Pi()*0.5);
5136 // jet3mom.RotateZ(TMath::Pi()*0.5);
5137 // jet3mom.RotateZ((-1)*TMath::Pi()*0.5);
5139 //jetMomplusTest[0] = jet3mom.X(); //fetching perp. axis coordinates
5140 //jetMomplusTest[1] = jet3mom.Y();
5141 //jetMomplusTest[2] = jet3mom.Z();
5143 //TVector3 perpjetplus3momTest(jetMomplusTest); //new TVector3 for +90deg rotated jet axis with rotation method from ROOT
5144 //TVector3 perpjetminus3momTest(jetMomminusTest); //new TVector3 for -90deg rotated jet axis with rotation method from ROOT
5147 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){ //collect V0 content in perp cone, rotated clockwise
5149 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack)); //inputlist is fListK0s, all reconstructed K0s in event
5150 if(!track){std::cout<<"K0s track not found!!!"<<std::endl; continue;}
5152 Double_t trackMom[3];//3-mom of V0 particle
5153 track->PxPyPz(trackMom);
5154 TVector3 track3mom(trackMom);
5156 Double_t dR = perpjetplus3mom.DeltaR(track3mom);
5160 outputlist->Add(track); // output list is jetPerpConeK0list
5162 sumPerpPt += track->Pt();
5169 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){//collect V0 content in perp cone, rotated counterclockwise
5171 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack)); //inputlist is fListK0s, all reconstructed K0s in event
5172 if(!track){std::cout<<"K0s track not found!!!"<<std::endl; continue;}
5174 Double_t trackMom[3];//3-mom of V0 particle
5175 track->PxPyPz(trackMom);
5176 TVector3 track3mom(trackMom);
5178 Double_t dR = perpjetneg3mom.DeltaR(track3mom);
5182 outputlist->Add(track); // output list is jetPerpConeK0list
5184 sumPerpPt += track->Pt();
5190 // pay attention: this list contains the double amount of V0s, found in both cones
5191 // before using it, devide spectra by 2!!!
5192 sumPerpPt = sumPerpPt*0.5; //correct to do this?
5200 // _______________________________________________________________________________________________________________________________________________________
5202 Bool_t AliAnalysisTaskJetChem::MCLabelCheck(AliAODv0* v0, Int_t particletype,const AliAODTrack* trackNeg, const AliAODTrack* trackPos, TList *listmc, Int_t& negDaughterpdg, Int_t& posDaughterpdg, Int_t& motherType, Int_t& v0Label, Double_t& MCPt, Bool_t& fPhysicalPrimary, Int_t& MCv0PDGCode, TString& generatorName, Bool_t& isinjected){
5204 if(!v0)return kFALSE;
5206 TClonesArray *stackmc = 0x0;
5207 stackmc = (TClonesArray*)listmc->FindObject(AliAODMCParticle::StdBranchName()); //get MCAOD branch in data
5210 Printf("ERROR: stack not available");
5215 Int_t negAssLabel = TMath::Abs(trackNeg->GetLabel()); //negative (reconstructed) charged track label in MC stack
5216 Int_t posAssLabel = TMath::Abs(trackPos->GetLabel()); //positive (reconstructed) charged track label in MC stack
5218 //injected particle checks
5224 AliAODMCHeader *header=(AliAODMCHeader*)listmc->FindObject(AliAODMCHeader::StdBranchName());
5225 if(!header)return kFALSE;
5227 mcXv=header->GetVtxX(); mcYv=header->GetVtxY(); mcZv=header->GetVtxZ();
5233 if(negAssLabel>=0 && negAssLabel < stackmc->GetEntriesFast() && posAssLabel>=0 && posAssLabel < stackmc->GetEntriesFast()){//safety check if label has valid value of stack
5235 AliAODMCParticle *mcNegPart =(AliAODMCParticle*)stackmc->UncheckedAt(negAssLabel);//fetch the, with one MC truth track associated (reconstructed), negative charged track
5236 v0Label = mcNegPart->GetMother();// mother of negative charged particle is v0, get v0 label here
5237 negDaughterpdg = mcNegPart->GetPdgCode();
5238 AliAODMCParticle *mcPosPart =(AliAODMCParticle*)stackmc->UncheckedAt(posAssLabel);//fetch the, with one MC truth track associated (reconstructed), positive charged track
5239 Int_t v0PosLabel = mcPosPart->GetMother(); //get mother label of positive charged track label
5240 posDaughterpdg = mcPosPart->GetPdgCode();
5242 if(v0Label >= 0 && v0Label < stackmc->GetEntriesFast() && v0Label == v0PosLabel){//first v0 mc label check, then: check if both daughters are stemming from same particle
5244 AliAODMCParticle *mcv0 = (AliAODMCParticle *)stackmc->UncheckedAt(v0Label); //fetch MC ass. particle to v0 (mother of the both charged daughter tracks)
5246 //do not use anymore:
5247 //fPhysicalPrimary = mcv0->IsPhysicalPrimary();
5249 Float_t fDistPrimaryMax = 0.01; // [cm] max distance of production point to the primary vertex (criterion for choice of MC particles considered as primary)
5251 // Get the distance between production point of the MC mother particle and the primary vertex
5253 Double_t dx = mcXv-mcv0->Xv();//mc primary vertex - mc particle production vertex
5254 Double_t dy = mcYv-mcv0->Yv();
5255 Double_t dz = mcZv-mcv0->Zv();
5257 Float_t fDistPrimary = TMath::Sqrt(dx*dx + dy*dy + dz*dz);
5258 fPhysicalPrimary = kFALSE;//init
5260 fPhysicalPrimary = (fDistPrimary < fDistPrimaryMax);
5261 MCv0PDGCode = mcv0->GetPdgCode();
5263 //if(fPhysicalPrimary == kTRUE){//look only at physical primary particles
5265 isinjected = IsTrackInjected(v0Label, header, stackmc, generatorName); //requires AliAODv0 instead of AliVTrack
5267 //trackinjected is kFALSE if it is either Hijing or has no generator name
5268 // std::cout<<" "<<std::endl;
5269 // std::cout<<"#### next particle: ####"<<std::endl;
5270 //std::cout<<"Is track injected: "<< trackinjected <<std::endl;
5271 // std::cout<<"pdg code: "<<MCv0PDGCode<<std::endl;
5272 // std::cout<<"v0Label: "<<v0Label<<std::endl;
5274 MCPt = mcv0->Pt();//for MC data, always use MC gen. pt for any pt distributions, also for the spectra, used for normalisation
5275 //for feed-down checks later
5277 Int_t motherLabel = mcv0->GetMother(); //get mother particle label of v0 particle
5278 // std::cout<<"motherLabel: "<<motherLabel<<std::endl;
5280 if(motherLabel >= 0 && v0Label < stackmc->GetEntriesFast()) //do safety check for mother label
5282 AliAODMCParticle *mcMother = (AliAODMCParticle *)stackmc->UncheckedAt(motherLabel); //get mother particle
5283 motherType = mcMother->GetPdgCode(); //get PDG code of mother
5286 Double_t XibarPt = 0.;
5288 if(particletype == kLambda){
5289 if((motherType == 3312)||(motherType == 3322)){ //if v0 mother is Xi0 or Xi- fill MC gen. pt in FD La histogram
5290 XiPt = mcMother->Pt();
5291 fh1MCXiPt->Fill(XiPt);
5294 if(particletype == kAntiLambda){
5295 if((motherType == -3312)||(motherType == -3322)){ //if v0 mother is Xibar0 or Xibar+ fill MC gen. pt in FD ALa histogram
5296 XibarPt = mcMother->Pt();
5297 fh1MCXibarPt->Fill(XibarPt);
5303 //pdg code checks etc..
5305 if(particletype == kK0){
5307 if(TMath::Abs(posDaughterpdg) != 211){return kFALSE;}//one or both of the daughters are not a pion
5308 if(TMath::Abs(negDaughterpdg) != 211){return kFALSE;}
5310 if(MCv0PDGCode != 310) {return kFALSE;}
5313 if(particletype == kLambda){
5314 if(MCv0PDGCode != 3122)return kFALSE;//if particle is not Antilambda, v0 is rejected
5315 if(posDaughterpdg != 2212)return kFALSE;
5316 if(negDaughterpdg != -211)return kFALSE; //pdg code check for Lambda daughters
5318 //{if((motherType == 3312)||(motherType == 3322)){continue;}//if Xi0 and Xi- is motherparticle of Lambda, particle is rejected, pay attention, most possible Xi-, Xi0 and Omega- are not distributed physically and are much more abundant than expected by physics //}
5321 if(particletype == kAntiLambda){
5322 if(MCv0PDGCode != -3122)return kFALSE;
5323 if(posDaughterpdg != 211)return kFALSE;
5324 if(negDaughterpdg !=-2212)return kFALSE; //pdg code check for Antilambda daughters
5327 //{if((motherType == -3312)||(motherType == -3322)){continue;}//if bar{Xi0} and Xi+ is motherparticle of Antilambda, particle is rejected
5331 return kTRUE; //check was successful
5332 }//end mc v0 label check
5333 }// end of stack label check
5338 return kFALSE; //check wasn't successful
5340 //________________________________________________________________________________________________________________________________________________________
5343 Bool_t AliAnalysisTaskJetChem::IsParticleMatching(const AliAODMCParticle* mcp0, const Int_t v0Label){
5345 const Int_t mcp0label = mcp0->GetLabel();
5347 if(v0Label == mcp0label)return kTRUE;
5352 //_______________________________________________________________________________________________________________________________________________________
5354 Bool_t AliAnalysisTaskJetChem::DaughterTrackCheck(AliAODv0* v0, Int_t& nnum, Int_t& pnum){
5357 if(v0->GetNDaughters() != 2) return kFALSE;//case v0 has more or less than 2 daughters, avoids seg. break at some AOD files //reason?
5360 // safety check of input parameters
5363 if(fDebug > 1){std::cout << std::endl
5364 << "Warning in AliAnalysisTaskJetChem::DaughterTrackCheck:" << std::endl
5365 << "v0 = " << v0 << std::endl;}
5371 //Daughters track check: its Luke Hanrattys method to check daughters charge
5377 AliAODTrack *ntracktest =(AliAODTrack*)(v0->GetDaughter(nnum));
5379 if(ntracktest == NULL)
5381 if(fDebug > 1){std::cout << std::endl
5382 << "Warning in AliAnalysisTaskJetChem::DaughterTrackCheck:" << std::endl
5383 << "ntracktest = " << ntracktest << std::endl;}
5388 if(ntracktest->Charge() > 0)
5394 const AliAODTrack *trackNeg=(AliAODTrack *)(v0->GetDaughter(nnum));
5395 const AliAODTrack *trackPos=(AliAODTrack *)(v0->GetDaughter(pnum));
5397 //Check if both tracks are available
5398 if (!trackPos || !trackNeg) {
5399 if(fDebug > 1) Printf("strange analysis::UserExec:: Error:Could not retrieve one of the daughter tracks\n");
5404 //remove like sign V0s
5405 if ( trackPos->Charge() == trackNeg->Charge() ){
5406 //if(fDebug>1) Printf("%s:%d found like-sign V0", (char*)__FILE__,__LINE__);
5414 //______________________________________________________________________
5415 TString AliAnalysisTaskJetChem::GetGenerator(Int_t label, AliAODMCHeader* header){
5416 Int_t nsumpart=0;//number of particles
5417 TList *lh=header->GetCocktailHeaders();//TList with all generator headers
5418 Int_t nh=lh->GetEntries();//number of entries in TList with all headers
5420 for(Int_t i=0;i<nh;i++){
5421 AliGenEventHeader* gh=(AliGenEventHeader*)lh->At(i);
5422 TString genname=gh->GetName();//name of particle generator
5423 Int_t npart=gh->NProduced();//number of stable or undecayed particles in MC stack block (?)
5424 if(label>=nsumpart && label<(nsumpart+npart)) return genname;
5431 //____________________________________________________________________________________________________________-
5433 Int_t AliAnalysisTaskJetChem::SplitCocktail(AliAODMCParticle *mcv0, Int_t v0Label, AliAODMCHeader *header, TClonesArray *arrayMC){//info in TString should be available from 2011 data productions on..
5435 if(!mcv0){std::cout << " !mcv0 " << std::endl;return 1;}
5436 if(!header){std::cout << " !header " << std::endl;return 1;}
5437 if(!arrayMC){std::cout << " !arrayMC " << std::endl;return 1;}
5439 //comment: all MC truth particles are sorted in the MC stack, first comes a block with all hijing produced particles, then several blocks for particletype specific injected particles, after this comes a (not-ordered) block with all decay-products
5440 //the complete amount of MC truth particles produced by several sources is named 'cocktail'
5442 //std::cout<<"v0 label: "<<v0Label<<std::endl;
5444 if(v0Label < 0) {std::cout<<"v0Label is negative!"<<std::endl;return 1;} //if label is negative, this particle is in the first part of the MC stack, it can be either Hijing or injected, but it is a primary particle
5446 TString generatorName = GetGenerator(v0Label,header);//this function returns a string with the generator name, used to produce this certain particle
5447 TString empty="";//no header was found
5449 std::cout << " TString generatorName: " << generatorName << std::endl;
5451 std::cout << " FIRST CALL " << generatorName << std::endl;
5453 while(generatorName.IsWhitespace()){
5454 AliAODMCParticle *mcpart= (AliAODMCParticle*)arrayMC->At(v0Label);
5455 if(!mcpart){return 1;}
5456 Int_t mother = mcpart->GetMother();
5458 generatorName = GetGenerator(mother,header);//see function directly below..
5459 std::cout << "Testing " << generatorName << " " << std::endl;
5462 std::cout << " FINAL CALL " << generatorName << std::endl;
5464 //std::transform(generatorName.begin(), generatorName.end(), generatorName.begin(), ::tolower); //convert TString bbb into lower case, to avoid that TString could be written in lower or upper case
5466 if(generatorName.Contains("ijing")){std::cout << " particle is Hijing!! " << std::endl; return 0;}//if TString returns something with "ijing" return this method with 0 -> select out all HIJING particles, all others return with "1"
5471 //_____________________________________________________________________
5472 void AliAnalysisTaskJetChem::GetTrackPrimaryGenerator(Int_t lab, AliAODMCHeader *header,TClonesArray *arrayMC,TString &nameGen){
5474 // method to check if a v0 comes from a given generator
5476 nameGen=GetGenerator(lab,header);
5478 // Int_t countControl=0;
5480 while(nameGen.IsWhitespace()){
5481 AliAODMCParticle *mcpart= (AliAODMCParticle*)arrayMC->At(lab);//get MC generated particle for AliAODv0 particle
5483 printf("AliAnalysisTaskJetChem::IsTrackInjected - BREAK: No valid AliAODMCParticle at label %i\n",lab);
5486 Int_t mother = mcpart->GetMother();
5489 printf("AliAnalysisTaskJetChem::IsTrackInjected - BREAK: Reached primary particle without valid mother\n");
5493 nameGen=GetGenerator(mother,header);
5496 // if(countControl>=10){ // 10 = arbitrary number; protection from infinite loops
5497 // printf("AliVertexingHFUtils::IsTrackInjected - BREAK: Protection from infinite loop active\n");
5507 //---------------------------------------------------------------------------------------------------------------------
5508 Bool_t AliAnalysisTaskJetChem::IsTrackInjected(Int_t lab, AliAODMCHeader *header,TClonesArray *arrayMC, TString& nameGen){
5509 // method to check if a v0 particle comes from the signal event or from the underlying Hijing event
5512 GetTrackPrimaryGenerator(lab, header, arrayMC, nameGen);
5514 if(nameGen.IsWhitespace() || nameGen.Contains("ijing")) return kFALSE;//particle has either no info about generator or is Hijing particle, so it is not injected
5516 //std::cout<<"generator name: "<<nameGen<<std::endl;
5521 //_________________________________________________________________________________________________________________________________________
5522 Double_t AliAnalysisTaskJetChem::SmearJetPt(Double_t jetPt, Int_t /*cent*/, Double_t /*jetRadius*/, Double_t /*ptmintrack*/, Double_t& jetPtSmear){
5524 static TF1 fsmear("f1","[0]*exp(-1*(x-[1])*(x-[1])/(2*[2]*[2]))",-100.,100.); //smearing according to gaussian function in between +/- 10 GeV/c
5528 /* if(cent>10) cl = 2;
5533 fsmear.SetParameters(1,0,11.19);//for 2010 PbPb jets, R=0.4, ptmintrack = 0.15 GeV/c, cent 00-10%, delta-pt width estimated via single track embedding
5534 //fsmear->SetParameters(1,0,3.28);//for 2010 PbPb jets, R=0.4, ptmintrack = 0.15 GeV/c, cent 50-60%, delta-pt width estimated via single track embedding
5536 //fsmear->SetParameters(1,0,4.472208);// for 2010 PbPb jets, R=0.2, ptmintrack = 0.15 GeV/c, cent 00-10%
5538 /* //delta-pt width for anti-kt jet finder:
5541 if((cl == 1)&&(jetRadius == 0.4)&&(ptmintrack == 0.15)){
5542 fsmear->SetParameters(1,0,10.178069);//(max.,mean,sigma) of gaussian, needs to be adjusted for every combination of jet cone size, centrality and min. pt constituents cut
5544 if((cl == 2)&&(jetRadius == 0.4)&&(ptmintrack == 0.15)){
5545 fsmear->SetParameters(1,0,8.536195);
5547 if((cl == 3)&&(jetRadius == 0.4)&&(ptmintrack == 0.15)){
5548 fsmear->SetParameters(1,0,?);
5550 if((cl == 4)&&(jetRadius == 0.4)&&(ptmintrack == 0.15)){
5551 fsmear->SetParameters(1,0,5.229839);
5555 if((cl == 1)&&(jetRadius == 0.3)&&(ptmintrack == 0.15)){
5556 fsmear->SetParameters(1,0,7.145967);
5558 if((cl == 2)&&(jetRadius == 0.3)&&(ptmintrack == 0.15)){
5559 fsmear->SetParameters(1,0,5.844796);
5561 if((cl == 3)&&(jetRadius == 0.3)&&(ptmintrack == 0.15)){
5562 fsmear->SetParameters(1,0,?);
5564 if((cl == 4)&&(jetRadius == 0.3)&&(ptmintrack == 0.15)){
5565 fsmear->SetParameters(1,0,3.630751);
5569 if((cl == 1)&&(jetRadius == 0.2)&&(ptmintrack == 0.15)){
5570 fsmear->SetParameters(1,0,4.472208);
5572 if((cl == 2)&&(jetRadius == 0.2)&&(ptmintrack == 0.15)){
5573 fsmear->SetParameters(1,0,3.543938);
5575 if((cl == 3)&&(jetRadius == 0.2)&&(ptmintrack == 0.15)){
5576 fsmear->SetParameters(1,0,?);
5578 if((cl == 4)&&(jetRadius == 0.2)&&(ptmintrack == 0.15)){
5579 fsmear->SetParameters(1,0,1.037476);
5584 Double_t r = fsmear.GetRandom();
5585 jetPtSmear = jetPt + r;
5587 // std::cout<<"jetPt: "<<jetPt<<std::endl;
5588 // std::cout<<"jetPtSmear: "<<jetPtSmear<<std::endl;
5589 // std::cout<<"r: "<<r<<std::endl;
5596 //______________________________________________________________________________________________________________________
5597 //____________________________________________________________________________________________________________________
5599 Bool_t AliAnalysisTaskJetChem::IsParticleInCone(const AliVParticle* part1, const AliVParticle* part2, Double_t dRMax) const
5601 // decides whether a particle is inside a jet cone
5602 if (!part1 || !part2)
5605 TVector3 vecMom2(part2->Px(),part2->Py(),part2->Pz());
5606 TVector3 vecMom1(part1->Px(),part1->Py(),part1->Pz());
5607 Double_t dR = vecMom2.DeltaR(vecMom1); // = sqrt(dEta*dEta+dPhi*dPhi)
5608 if(dR<dRMax) // momentum vectors of part1 and part2 are closer than dRMax
5612 //__________________________________________________________________________________________________________________
5615 Bool_t AliAnalysisTaskJetChem::IsRCJCOverlap(TList* recjetlist, const AliVParticle* part, Double_t dDistance) const{
5617 if(!recjetlist) return kFALSE;
5618 if(!part) return kFALSE;
5619 if(!dDistance) return kFALSE;
5620 Int_t nRecJetsCuts = fJetsRecCuts->GetEntries();
5622 for(Int_t i=0; i<nRecJetsCuts; ++i){ //loop over all reconstructed jets in events
5623 AliAODJet* jet = (AliAODJet*) (recjetlist->At(i));
5624 if(!jet){if(fDebug>2)std::cout<<"AliAnalysisTaskJetChem::IsRCJCOverlap jet pointer invalid!"<<std::endl;continue;}
5625 if(IsParticleInCone(jet, part, dDistance) == kTRUE)return kTRUE;//RC and JC are overlapping
5627 }//end loop testing RC-JC overlap
5628 return kFALSE;//RC and JC are not overlapping -> good!
5631 //_______________________________________________________________________________________________________________________
5632 AliAODJet* AliAnalysisTaskJetChem::GetRandomCone(TList* jetlist, Double_t dEtaConeMax, Double_t dDistance) const
5634 TLorentzVector vecRdCone;
5635 AliAODJet* jetRC = 0;//random cone candidate
5636 Double_t dEta, dPhi; //random eta and phi value for RC
5637 Bool_t IsRCoutJC = kFALSE;//check whether RC is not overlapping with any selected jet cone in event
5638 Int_t iRCTrials = 10;//search at maximum 10 times for random cone that doesn't overlap with jet cone
5640 for(Int_t i=0; i<iRCTrials; iRCTrials++){
5642 dEta = dEtaConeMax*(2*fRandom->Rndm()-1.); //random eta value in range: [-dEtaConeMax,+dEtaConeMax]
5643 dPhi = TMath::TwoPi()*fRandom->Rndm(); //random phi value in range: [0,2*Pi]
5644 vecRdCone.SetPtEtaPhiM(1.,dEta,dPhi,0.);
5645 jetRC = new AliAODJet(vecRdCone);//new RC candidate
5647 if (!IsRCJCOverlap(jetlist,jetRC,dDistance))
5649 IsRCoutJC = kTRUE; //std::cout<<"RC and JC are not overlapping!!!"<<std::endl;
5653 delete jetRC; //RC is overlapping with JC, delete this RC candidate
5656 if(!IsRCoutJC) {jetRC = 0;}//in case no random cone was selected
5662 // _______________________________________________________________________________________________________________________
5663 AliAODJet* AliAnalysisTaskJetChem::GetMedianCluster()
5665 // fill tracks from bckgCluster branch,
5666 // using cluster with median density (odd number of clusters)
5667 // or picking randomly one of the two closest to median (even number)
5669 Int_t nBckgClusters = fBckgJetsRec->GetEntries(); // not 'recCuts': use all clusters in full eta range
5671 if(nBckgClusters<3) return 0; // need at least 3 clusters (skipping 2 highest)
5673 Double_t* bgrDensity = new Double_t[nBckgClusters];
5674 Int_t* indices = new Int_t[nBckgClusters];
5676 for(Int_t ij=0; ij<nBckgClusters; ++ij){
5678 AliAODJet* bgrCluster = (AliAODJet*)(fBckgJetsRec->At(ij));
5679 Double_t clusterPt = bgrCluster->Pt();
5680 Double_t area = bgrCluster->EffectiveAreaCharged();
5682 Double_t density = 0;
5683 if(area>0) density = clusterPt/area;
5685 bgrDensity[ij] = density;
5689 TMath::Sort(nBckgClusters, bgrDensity, indices);
5691 // get median cluster
5693 AliAODJet* medianCluster = 0;
5694 Double_t medianDensity = 0;
5696 if(TMath::Odd(nBckgClusters)){
5698 Int_t medianIndex = indices[(Int_t) (0.5*(nBckgClusters+1))];
5700 medianCluster = (AliAODJet*)(fBckgJetsRec->At(medianIndex));
5702 Double_t clusterPt = medianCluster->Pt();
5703 Double_t area = medianCluster->EffectiveAreaCharged();
5705 if(area>0) medianDensity = clusterPt/area;
5709 Int_t medianIndex1 = indices[(Int_t) (0.5*nBckgClusters)];
5710 Int_t medianIndex2 = indices[(Int_t) (0.5*nBckgClusters+1)];
5712 AliAODJet* medianCluster1 = (AliAODJet*)(fBckgJetsRec->At(medianIndex1));
5713 AliAODJet* medianCluster2 = (AliAODJet*)(fBckgJetsRec->At(medianIndex2));
5715 Double_t density1 = 0;
5716 Double_t clusterPt1 = medianCluster1->Pt();
5717 Double_t area1 = medianCluster1->EffectiveAreaCharged();
5718 if(area1>0) density1 = clusterPt1/area1;
5720 Double_t density2 = 0;
5721 Double_t clusterPt2 = medianCluster2->Pt();
5722 Double_t area2 = medianCluster2->EffectiveAreaCharged();
5723 if(area2>0) density2 = clusterPt2/area2;
5725 medianDensity = 0.5*(density1+density2);
5727 medianCluster = ( (gRandom->Rndm()>0.5) ? medianCluster1 : medianCluster2 ); // select one randomly to avoid adding areas
5730 delete[] bgrDensity;
5733 return medianCluster;
5736 //____________________________________________________________________________________________
5738 Double_t AliAnalysisTaskJetChem::AreaCircSegment(Double_t dRadius, Double_t dDistance) const
5740 // calculate area of a circular segment defined by the circle radius and the (oriented) distance between the secant line and the circle centre
5741 Double_t dEpsilon = 1e-2;
5742 Double_t dR = dRadius;
5743 Double_t dD = dDistance;
5744 if (TMath::Abs(dR)<dEpsilon)
5746 if(fDebug>0) printf("AliAnalysisTaskJetChem::AreaCircSegment: Error: Too small radius: %f < %f\n",dR,dEpsilon);
5752 return TMath::Pi()*dR*dR;
5753 return dR*dR*TMath::ACos(dD/dR)-dD*TMath::Sqrt(dR*dR-dD*dD);