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)
228 ,fhnInvMassEtaTrackPtK0s(0)
229 ,fhnInvMassEtaTrackPtLa(0)
230 ,fhnInvMassEtaTrackPtALa(0)
236 // ,fh2MCgenK0Cone(0)
237 // ,fh2MCgenLaCone(0)
238 // ,fh2MCgenALaCone(0)
239 // ,fh2MCEtagenK0Cone(0)
240 // ,fh2MCEtagenLaCone(0)
241 // ,fh2MCEtagenALaCone(0)
242 ,fh2CorrHijingLaProton(0)
243 ,fh2CorrInjectLaProton(0)
244 ,fh2CorrHijingALaAProton(0)
245 ,fh2CorrInjectALaAProton(0)
248 ,fh1IMALaConeSmear(0)
249 ,fh2MCEtaVsPtHijingLa(0)
250 ,fh2MCEtaVsPtInjectLa(0)
251 ,fh2MCEtaVsPtHijingALa(0)
252 ,fh2MCEtaVsPtInjectALa(0)
253 ,fhnrecMCHijingLaIncl(0)
254 ,fhnrecMCHijingLaCone(0)
255 ,fhnrecMCHijingALaIncl(0)
256 ,fhnrecMCHijingALaCone(0)
257 ,fhnrecMCInjectLaIncl(0)
258 ,fhnrecMCInjectLaCone(0)
259 ,fhnrecMCInjectALaIncl(0)
260 ,fhnrecMCInjectALaCone(0)
264 ,fhnMCrecK0ConeSmear(0)
265 ,fhnMCrecLaConeSmear(0)
266 ,fhnMCrecALaConeSmear(0)
267 ,fhnK0sSecContinCone(0)
268 ,fhnLaSecContinCone(0)
269 ,fhnALaSecContinCone(0)
294 ,fh1MCMultiplicityPrimary(0)
295 ,fh1MCMultiplicityTracks(0)
298 ,fhnFeedDownLaCone(0)
299 ,fhnFeedDownALaCone(0)
300 ,fh1MCProdRadiusK0s(0)
301 ,fh1MCProdRadiusLambda(0)
302 ,fh1MCProdRadiusAntiLambda(0)
306 ,fh1MCPtAntiLambda(0)
314 //,fh1MCRapAntiLambda(0)
318 ,fh1MCEtaAntiLambda(0)
321 // default constructor
324 //__________________________________________________________________________________________
325 AliAnalysisTaskJetChem::AliAnalysisTaskJetChem(const char *name)
326 : AliAnalysisTaskFragmentationFunction(name)
340 ,fCutV0cosPointAngle(0)
347 ,fCutDcaV0Daughters(0)
348 ,fCutDcaPosToPrimVertex(0)
349 ,fCutDcaNegToPrimVertex(0)
359 ,fFFHistosRecCutsK0Evt(0)
360 //,fFFHistosIMK0AllEvt(0)
361 //,fFFHistosIMK0Jet(0)
362 //,fFFHistosIMK0Cone(0)
366 //,fFFHistosIMLaAllEvt(0)
367 //,fFFHistosIMLaJet(0)
368 //,fFFHistosIMLaCone(0)
372 ,fListFeeddownLaCand(0)
373 ,fListFeeddownALaCand(0)
379 ,fListMCgenK0sCone(0)
381 ,fListMCgenALaCone(0)
382 ,IsArmenterosSelected(0)
383 //,fFFHistosIMALaAllEvt(0)
384 //,fFFHistosIMALaJet(0)
385 // ,fFFHistosIMALaCone(0)
401 ,fFFIMLaNBinsJetPt(0)
432 // ,fh1trackPosNCls(0)
433 // ,fh1trackNegNCls(0)
443 ,fh2ProperLifetimeK0sVsPtBeforeCut(0)
444 ,fh2ProperLifetimeK0sVsPtAfterCut(0)
446 ,fh1DcaV0Daughters(0)
447 ,fh1DcaPosToPrimVertex(0)
448 ,fh1DcaNegToPrimVertex(0)
449 ,fh2ArmenterosBeforeCuts(0)
450 ,fh2ArmenterosAfterCuts(0)
453 ,fh1PosDaughterCharge(0)
454 ,fh1NegDaughterCharge(0)
468 ,fhnInvMassEtaTrackPtK0s(0)
469 ,fhnInvMassEtaTrackPtLa(0)
470 ,fhnInvMassEtaTrackPtALa(0)
478 //,fh2MCgenALaCone(0)
479 //,fh2MCEtagenK0Cone(0)
480 //,fh2MCEtagenLaCone(0)
481 //,fh2MCEtagenALaCone(0)
482 ,fh2CorrHijingLaProton(0)
483 ,fh2CorrInjectLaProton(0)
484 ,fh2CorrHijingALaAProton(0)
485 ,fh2CorrInjectALaAProton(0)
488 ,fh1IMALaConeSmear(0)
489 ,fh2MCEtaVsPtHijingLa(0)
490 ,fh2MCEtaVsPtInjectLa(0)
491 ,fh2MCEtaVsPtHijingALa(0)
492 ,fh2MCEtaVsPtInjectALa(0)
493 ,fhnrecMCHijingLaIncl(0)
494 ,fhnrecMCHijingLaCone(0)
495 ,fhnrecMCHijingALaIncl(0)
496 ,fhnrecMCHijingALaCone(0)
497 ,fhnrecMCInjectLaIncl(0)
498 ,fhnrecMCInjectLaCone(0)
499 ,fhnrecMCInjectALaIncl(0)
500 ,fhnrecMCInjectALaCone(0)
504 ,fhnMCrecK0ConeSmear(0)
505 ,fhnMCrecLaConeSmear(0)
506 ,fhnMCrecALaConeSmear(0)
507 ,fhnK0sSecContinCone(0)
508 ,fhnLaSecContinCone(0)
509 ,fhnALaSecContinCone(0)
534 ,fh1MCMultiplicityPrimary(0)
535 ,fh1MCMultiplicityTracks(0)
538 ,fhnFeedDownLaCone(0)
539 ,fhnFeedDownALaCone(0)
540 ,fh1MCProdRadiusK0s(0)
541 ,fh1MCProdRadiusLambda(0)
542 ,fh1MCProdRadiusAntiLambda(0)
546 ,fh1MCPtAntiLambda(0)
554 //,fh1MCRapAntiLambda(0)
558 ,fh1MCEtaAntiLambda(0)
564 DefineOutput(1,TList::Class());
567 //__________________________________________________________________________________________________________________________
568 AliAnalysisTaskJetChem::AliAnalysisTaskJetChem(const AliAnalysisTaskJetChem ©)
569 : AliAnalysisTaskFragmentationFunction()
571 ,fRandom(copy.fRandom)
572 ,fAnalysisMC(copy.fAnalysisMC)
573 ,fDeltaVertexZ(copy.fDeltaVertexZ)
574 ,fCutjetEta(copy.fCutjetEta)
575 ,fCuttrackNegNcls(copy.fCuttrackNegNcls)
576 ,fCuttrackPosNcls(copy.fCuttrackPosNcls)
577 ,fCutPostrackRap(copy.fCutPostrackRap)
578 ,fCutNegtrackRap(copy.fCutNegtrackRap)
579 ,fCutRap(copy.fCutRap)
580 ,fCutPostrackEta(copy.fCutPostrackEta)
581 ,fCutNegtrackEta(copy.fCutNegtrackEta)
582 ,fCutEta(copy.fCutEta)
583 ,fCutV0cosPointAngle(copy.fCutV0cosPointAngle)
584 ,fKinkDaughters(copy.fKinkDaughters)
585 ,fRequireTPCRefit(copy.fRequireTPCRefit)
586 ,fCutArmenteros(copy.fCutArmenteros)
587 ,fCutV0DecayMin(copy.fCutV0DecayMin)
588 ,fCutV0DecayMax(copy.fCutV0DecayMax)
589 ,fCutV0totMom(copy.fCutV0totMom)
590 ,fCutDcaV0Daughters(copy.fCutDcaV0Daughters)
591 ,fCutDcaPosToPrimVertex(copy.fCutDcaPosToPrimVertex)
592 ,fCutDcaNegToPrimVertex(copy.fCutDcaNegToPrimVertex)
593 ,fCutV0RadiusMin(copy.fCutV0RadiusMin)
594 ,fCutV0RadiusMax(copy.fCutV0RadiusMax)
595 ,fCutBetheBloch(copy.fCutBetheBloch)
596 ,fCutRatio(copy.fCutRatio)
597 ,fK0Type(copy.fK0Type)
598 ,fFilterMaskK0(copy.fFilterMaskK0)
599 ,fListK0s(copy.fListK0s)
600 ,fPIDResponse(copy.fPIDResponse)
601 ,fV0QAK0(copy.fV0QAK0)
602 ,fFFHistosRecCutsK0Evt(copy.fFFHistosRecCutsK0Evt)
603 //,fFFHistosIMK0AllEvt(copy.fFFHistosIMK0AllEvt)
604 //,fFFHistosIMK0Jet(copy.fFFHistosIMK0Jet)
605 //,fFFHistosIMK0Cone(copy.fFFHistosIMK0Cone)
606 ,fLaType(copy.fLaType)
607 ,fFilterMaskLa(copy.fFilterMaskLa)
608 ,fListLa(copy.fListLa)
609 //,fFFHistosIMLaAllEvt(copy.fFFHistosIMLaAllEvt)
610 //,fFFHistosIMLaJet(copy.fFFHistosIMLaJet)
611 //,fFFHistosIMLaCone(copy.fFFHistosIMLaCone)
612 ,fALaType(copy.fALaType)
613 ,fFilterMaskALa(copy.fFilterMaskALa)
614 ,fListALa(copy.fListALa)
615 ,fListFeeddownLaCand(copy.fListFeeddownLaCand)
616 ,fListFeeddownALaCand(copy.fListFeeddownALaCand)
617 ,jetConeFDLalist(copy.jetConeFDLalist)
618 ,jetConeFDALalist(copy.jetConeFDALalist)
619 ,fListMCgenK0s(copy.fListMCgenK0s)
620 ,fListMCgenLa(copy.fListMCgenLa)
621 ,fListMCgenALa(copy.fListMCgenALa)
622 ,fListMCgenK0sCone(copy.fListMCgenK0sCone)
623 ,fListMCgenLaCone(copy.fListMCgenLaCone)
624 ,fListMCgenALaCone(copy.fListMCgenALaCone)
625 ,IsArmenterosSelected(copy.IsArmenterosSelected)
626 //,fFFHistosIMALaAllEvt(copy.fFFHistosIMALaAllEvt)
627 //,fFFHistosIMALaJet(copy.fFFHistosIMALaJet)
628 //,fFFHistosIMALaCone(copy.fFFHistosIMALaCone)
629 ,fFFIMNBinsJetPt(copy.fFFIMNBinsJetPt)
630 ,fFFIMJetPtMin(copy.fFFIMJetPtMin)
631 ,fFFIMJetPtMax(copy.fFFIMJetPtMax)
632 ,fFFIMNBinsInvM(copy.fFFIMNBinsInvM)
633 ,fFFIMInvMMin(copy.fFFIMInvMMin)
634 ,fFFIMInvMMax(copy.fFFIMInvMMax)
635 ,fFFIMNBinsPt(copy.fFFIMNBinsPt)
636 ,fFFIMPtMin(copy.fFFIMPtMin)
637 ,fFFIMPtMax(copy.fFFIMPtMax)
638 ,fFFIMNBinsXi(copy.fFFIMNBinsXi)
639 ,fFFIMXiMin(copy.fFFIMXiMin)
640 ,fFFIMXiMax(copy.fFFIMXiMax)
641 ,fFFIMNBinsZ(copy.fFFIMNBinsZ)
642 ,fFFIMZMin(copy.fFFIMZMin)
643 ,fFFIMZMax(copy.fFFIMZMax)
644 ,fFFIMLaNBinsJetPt(copy.fFFIMLaNBinsJetPt)
645 ,fFFIMLaJetPtMin(copy.fFFIMLaJetPtMin)
646 ,fFFIMLaJetPtMax(copy.fFFIMLaJetPtMax)
647 ,fFFIMLaNBinsInvM(copy.fFFIMLaNBinsInvM)
648 ,fFFIMLaInvMMin(copy.fFFIMLaInvMMin)
649 ,fFFIMLaInvMMax(copy.fFFIMLaInvMMax)
650 ,fFFIMLaNBinsPt(copy.fFFIMLaNBinsPt)
651 ,fFFIMLaPtMin(copy.fFFIMLaPtMin)
652 ,fFFIMLaPtMax(copy.fFFIMLaPtMax)
653 ,fFFIMLaNBinsXi(copy.fFFIMLaNBinsXi)
654 ,fFFIMLaXiMin(copy.fFFIMLaXiMin)
655 ,fFFIMLaXiMax(copy.fFFIMLaXiMax)
656 ,fFFIMLaNBinsZ(copy.fFFIMLaNBinsZ)
657 ,fFFIMLaZMin(copy.fFFIMLaZMin)
658 ,fFFIMLaZMax(copy.fFFIMLaZMax)
659 ,fh1EvtAllCent(copy.fh1EvtAllCent)
661 ,fh1K0Mult(copy.fh1K0Mult)
662 ,fh1dPhiJetK0(copy.fh1dPhiJetK0)
663 ,fh1LaMult(copy.fh1LaMult)
664 ,fh1dPhiJetLa(copy.fh1dPhiJetLa)
665 ,fh1ALaMult(copy.fh1ALaMult)
666 ,fh1dPhiJetALa(copy.fh1dPhiJetALa)
667 ,fh1JetEta(copy.fh1JetEta)
668 ,fh1JetPhi(copy.fh1JetPhi)
669 ,fh2JetEtaPhi(copy.fh2JetEtaPhi)
670 //,fh1V0JetPt(copy.fh1V0JetPt)
671 ,fh1IMK0Cone(copy.fh1IMK0Cone)
672 ,fh1IMLaCone(copy.fh1IMLaCone)
673 ,fh1IMALaCone(copy.fh1IMALaCone)
674 ,fh2FFJetTrackEta(copy.fh2FFJetTrackEta)
675 //,fh1trackPosNCls(copy.fh1trackPosNCls)
676 //,fh1trackNegNCls(copy.fh1trackNegNCls)
677 ,fh1trackPosRap(copy.fh1trackPosRap)
678 ,fh1trackNegRap(copy.fh1trackNegRap)
679 //,fh1V0Rap(copy.fh1V0Rap)
680 ,fh1trackPosEta(copy.fh1trackPosEta)
681 ,fh1trackNegEta(copy.fh1trackNegEta)
682 ,fh1V0Eta(copy.fh1V0Eta)
683 //,fh1V0totMom(copy.fh1V0totMom)
684 ,fh1CosPointAngle(copy.fh1CosPointAngle)
685 ,fh1DecayLengthV0(copy.fh1DecayLengthV0)
686 ,fh2ProperLifetimeK0sVsPtBeforeCut(copy.fh2ProperLifetimeK0sVsPtBeforeCut)
687 ,fh2ProperLifetimeK0sVsPtAfterCut(copy.fh2ProperLifetimeK0sVsPtAfterCut)
688 ,fh1V0Radius(copy.fh1V0Radius)
689 ,fh1DcaV0Daughters(copy.fh1DcaV0Daughters)
690 ,fh1DcaPosToPrimVertex(copy.fh1DcaPosToPrimVertex)
691 ,fh1DcaNegToPrimVertex(copy.fh1DcaNegToPrimVertex)
692 ,fh2ArmenterosBeforeCuts(copy.fh2ArmenterosBeforeCuts)
693 ,fh2ArmenterosAfterCuts(copy.fh2ArmenterosAfterCuts)
694 ,fh2BBLaPos(copy.fh2BBLaPos)
695 ,fh2BBLaNeg(copy.fh2BBLaPos)
696 ,fh1PosDaughterCharge(copy.fh1PosDaughterCharge)
697 ,fh1NegDaughterCharge(copy.fh1NegDaughterCharge)
698 ,fh1PtMCK0s(copy.fh1PtMCK0s)
699 ,fh1PtMCLa(copy.fh1PtMCLa)
700 ,fh1PtMCALa(copy.fh1PtMCALa)
701 ,fh1EtaK0s(copy.fh1EtaK0s)
702 ,fh1EtaLa(copy.fh1EtaLa)
703 ,fh1EtaALa(copy.fh1EtaALa)
705 ,fh1RCBiasK0(copy.fh1RCBiasK0)
706 ,fh1RCBiasLa(copy.fh1RCBiasLa)
707 ,fh1RCBiasALa(copy.fh1RCBiasALa)
711 ,fhnInvMassEtaTrackPtK0s(copy.fhnInvMassEtaTrackPtK0s)
712 ,fhnInvMassEtaTrackPtLa(copy.fhnInvMassEtaTrackPtLa)
713 ,fhnInvMassEtaTrackPtALa(copy.fhnInvMassEtaTrackPtALa)
714 ,fh1TrackMultCone(copy.fh1TrackMultCone)
715 ,fh2TrackMultCone(copy.fh2TrackMultCone)
716 ,fhnNJK0(copy.fhnNJK0)
717 ,fhnNJLa(copy.fhnNJLa)
718 ,fhnNJALa(copy.fhnNJALa)
719 //,fh2MCgenK0Cone(copy.fh2MCgenK0Cone)
720 //,fh2MCgenLaCone(copy.fh2MCgenLaCone)
721 //,fh2MCgenALaCone(copy.fh2MCgenALaCone)
722 //,fh2MCEtagenK0Cone(copy.fh2MCEtagenK0Cone)
723 //,fh2MCEtagenLaCone(copy.fh2MCEtagenLaCone)
724 //,fh2MCEtagenALaCone(copy.fh2MCEtagenALaCone)
725 ,fh2CorrHijingLaProton(copy.fh2CorrHijingLaProton)
726 ,fh2CorrInjectLaProton(copy.fh2CorrInjectLaProton)
727 ,fh2CorrHijingALaAProton(copy.fh2CorrHijingALaAProton)
728 ,fh2CorrInjectALaAProton(copy.fh2CorrInjectALaAProton)
729 ,fh1IMK0ConeSmear(copy.fh1IMK0ConeSmear)
730 ,fh1IMLaConeSmear(copy.fh1IMLaConeSmear)
731 ,fh1IMALaConeSmear(copy.fh1IMALaConeSmear)
732 ,fh2MCEtaVsPtHijingLa(copy.fh2MCEtaVsPtHijingLa)
733 ,fh2MCEtaVsPtInjectLa(copy.fh2MCEtaVsPtInjectLa)
734 ,fh2MCEtaVsPtHijingALa(copy.fh2MCEtaVsPtHijingALa)
735 ,fh2MCEtaVsPtInjectALa(copy.fh2MCEtaVsPtInjectALa)
736 ,fhnrecMCHijingLaIncl(copy.fhnrecMCHijingLaIncl)
737 ,fhnrecMCHijingLaCone(copy.fhnrecMCHijingLaCone)
738 ,fhnrecMCHijingALaIncl(copy.fhnrecMCHijingALaIncl)
739 ,fhnrecMCHijingALaCone(copy.fhnrecMCHijingALaCone)
740 ,fhnrecMCInjectLaIncl(copy.fhnrecMCInjectLaIncl)
741 ,fhnrecMCInjectLaCone(copy.fhnrecMCInjectLaCone)
742 ,fhnrecMCInjectALaIncl(copy.fhnrecMCInjectALaIncl)
743 ,fhnrecMCInjectALaCone(copy.fhnrecMCInjectALaCone)
744 ,fhnMCrecK0Cone(copy.fhnMCrecK0Cone)
745 ,fhnMCrecLaCone(copy.fhnMCrecLaCone)
746 ,fhnMCrecALaCone(copy.fhnMCrecALaCone)
747 ,fhnMCrecK0ConeSmear(copy.fhnMCrecK0ConeSmear)
748 ,fhnMCrecLaConeSmear(copy.fhnMCrecLaConeSmear)
749 ,fhnMCrecALaConeSmear(copy.fhnMCrecALaConeSmear)
750 ,fhnK0sSecContinCone(copy.fhnK0sSecContinCone)
751 ,fhnLaSecContinCone(copy.fhnLaSecContinCone)
752 ,fhnALaSecContinCone(copy.fhnALaSecContinCone)
753 ,fhnK0sIncl(copy.fhnK0sIncl)
754 ,fhnK0sCone(copy.fhnK0sCone)
755 ,fhnLaIncl(copy.fhnLaIncl)
756 ,fhnLaCone(copy.fhnLaCone)
757 ,fhnALaIncl(copy.fhnALaIncl)
758 ,fhnALaCone(copy.fhnALaCone)
759 ,fhnK0sPC(copy.fhnK0sPC)
760 ,fhnLaPC(copy.fhnLaPC)
761 ,fhnALaPC(copy.fhnALaPC)
762 ,fhnK0sMCC(copy.fhnK0sMCC)
763 ,fhnLaMCC(copy.fhnLaMCC)
764 ,fhnALaMCC(copy.fhnALaMCC)
765 ,fhnK0sRC(copy.fhnK0sRC)
766 ,fhnLaRC(copy.fhnLaRC)
767 ,fhnALaRC(copy.fhnALaRC)
768 ,fhnK0sRCBias(copy.fhnK0sRCBias)
769 ,fhnLaRCBias(copy.fhnLaRCBias)
770 ,fhnALaRCBias(copy.fhnALaRCBias)
771 ,fhnK0sOC(copy.fhnK0sOC)
772 ,fhnLaOC(copy.fhnLaOC)
773 ,fhnALaOC(copy.fhnALaOC)
774 ,fh1AreaExcluded(copy.fh1AreaExcluded)
775 ,fh1MedianEta(copy.fh1MedianEta)
776 ,fh1JetPtMedian(copy.fh1JetPtMedian)
777 ,fh1MCMultiplicityPrimary(copy.fh1MCMultiplicityPrimary)
778 ,fh1MCMultiplicityTracks(copy.fh1MCMultiplicityTracks)
779 ,fhnFeedDownLa(copy.fhnFeedDownLa)
780 ,fhnFeedDownALa(copy.fhnFeedDownALa)
781 ,fhnFeedDownLaCone(copy.fhnFeedDownLaCone)
782 ,fhnFeedDownALaCone(copy.fhnFeedDownALaCone)
783 ,fh1MCProdRadiusK0s(copy.fh1MCProdRadiusK0s)
784 ,fh1MCProdRadiusLambda(copy.fh1MCProdRadiusLambda)
785 ,fh1MCProdRadiusAntiLambda(copy.fh1MCProdRadiusAntiLambda)
786 ,fh1MCPtV0s(copy.fh1MCPtV0s)
787 ,fh1MCPtK0s(copy.fh1MCPtK0s)
788 ,fh1MCPtLambda(copy.fh1MCPtLambda)
789 ,fh1MCPtAntiLambda(copy.fh1MCPtAntiLambda)
790 ,fh1MCXiPt(copy.fh1MCXiPt)
791 ,fh1MCXibarPt(copy.fh1MCXibarPt)
792 ,fh2MCEtaVsPtK0s(copy.fh2MCEtaVsPtK0s)
793 ,fh2MCEtaVsPtLa(copy.fh2MCEtaVsPtLa)
794 ,fh2MCEtaVsPtALa(copy.fh2MCEtaVsPtALa)
795 //,fh1MCRapK0s(copy.fh1MCRapK0s)
796 //,fh1MCRapLambda(copy.fh1MCRapLambda)
797 //,fh1MCRapAntiLambda(copy.fh1MCRapAntiLambda)
798 ,fh1MCEtaAllK0s(copy.fh1MCEtaAllK0s)
799 ,fh1MCEtaK0s(copy.fh1MCEtaK0s)
800 ,fh1MCEtaLambda(copy.fh1MCEtaLambda)
801 ,fh1MCEtaAntiLambda(copy.fh1MCEtaAntiLambda)
808 // _________________________________________________________________________________________________________________________________
809 AliAnalysisTaskJetChem& AliAnalysisTaskJetChem::operator=(const AliAnalysisTaskJetChem& o)
814 AliAnalysisTaskFragmentationFunction::operator=(o);
817 fAnalysisMC = o.fAnalysisMC;
818 fDeltaVertexZ = o.fDeltaVertexZ;
819 fCutjetEta = o.fCutjetEta;
820 fCuttrackNegNcls = o.fCuttrackNegNcls;
821 fCuttrackPosNcls = o.fCuttrackPosNcls;
822 fCutPostrackRap = o.fCutPostrackRap;
823 fCutNegtrackRap = o.fCutNegtrackRap;
825 fCutPostrackEta = o.fCutPostrackEta;
826 fCutNegtrackEta = o.fCutNegtrackEta;
828 fCutV0cosPointAngle = o.fCutV0cosPointAngle;
829 fKinkDaughters = o.fKinkDaughters;
830 fRequireTPCRefit = o.fRequireTPCRefit;
831 fCutArmenteros = o.fCutArmenteros;
832 fCutV0DecayMin = o.fCutV0DecayMin;
833 fCutV0DecayMax = o.fCutV0DecayMax;
834 fCutV0totMom = o.fCutV0totMom;
835 fCutDcaV0Daughters = o.fCutDcaV0Daughters;
836 fCutDcaPosToPrimVertex = o.fCutDcaPosToPrimVertex;
837 fCutDcaNegToPrimVertex = o.fCutDcaNegToPrimVertex;
838 fCutV0RadiusMin = o.fCutV0RadiusMin;
839 fCutV0RadiusMax = o.fCutV0RadiusMax;
840 fCutBetheBloch = o.fCutBetheBloch;
841 fCutRatio = o.fCutRatio;
843 fFilterMaskK0 = o.fFilterMaskK0;
844 fListK0s = o.fListK0s;
845 fPIDResponse = o.fPIDResponse;
847 fFFHistosRecCutsK0Evt = o.fFFHistosRecCutsK0Evt;
848 //fFFHistosIMK0AllEvt = o.fFFHistosIMK0AllEvt;
849 //fFFHistosIMK0Jet = o.fFFHistosIMK0Jet;
850 //fFFHistosIMK0Cone = o.fFFHistosIMK0Cone;
852 fFilterMaskLa = o.fFilterMaskLa;
854 //fFFHistosIMLaAllEvt = o.fFFHistosIMLaAllEvt;
855 //fFFHistosIMLaJet = o.fFFHistosIMLaJet;
856 //fFFHistosIMLaCone = o.fFFHistosIMLaCone;
857 fALaType = o.fALaType;
858 fFilterMaskALa = o.fFilterMaskALa;
859 fListFeeddownLaCand = o.fListFeeddownLaCand;
860 fListFeeddownALaCand = o.fListFeeddownALaCand;
861 jetConeFDLalist = o.jetConeFDLalist;
862 jetConeFDALalist = o.jetConeFDALalist;
863 fListMCgenK0s = o.fListMCgenK0s;
864 fListMCgenLa = o.fListMCgenLa;
865 fListMCgenALa = o.fListMCgenALa;
866 fListMCgenK0sCone = o.fListMCgenK0sCone;
867 fListMCgenLaCone = o.fListMCgenLaCone;
868 fListMCgenALaCone = o.fListMCgenALaCone;
869 IsArmenterosSelected = o.IsArmenterosSelected;
870 // fFFHistosIMALaAllEvt = o.fFFHistosIMALaAllEvt;
871 // fFFHistosIMALaJet = o.fFFHistosIMALaJet;
872 // fFFHistosIMALaCone = o.fFFHistosIMALaCone;
873 fFFIMNBinsJetPt = o.fFFIMNBinsJetPt;
874 fFFIMJetPtMin = o.fFFIMJetPtMin;
875 fFFIMJetPtMax = o.fFFIMJetPtMax;
876 fFFIMNBinsPt = o.fFFIMNBinsPt;
877 fFFIMPtMin = o.fFFIMPtMin;
878 fFFIMPtMax = o.fFFIMPtMax;
879 fFFIMNBinsXi = o.fFFIMNBinsXi;
880 fFFIMXiMin = o.fFFIMXiMin;
881 fFFIMXiMax = o.fFFIMXiMax;
882 fFFIMNBinsZ = o.fFFIMNBinsZ;
883 fFFIMZMin = o.fFFIMZMin;
884 fFFIMZMax = o.fFFIMZMax;
885 fFFIMLaNBinsJetPt = o.fFFIMLaNBinsJetPt;
886 fFFIMLaJetPtMin = o.fFFIMLaJetPtMin;
887 fFFIMLaJetPtMax = o.fFFIMLaJetPtMax;
888 fFFIMLaNBinsPt = o.fFFIMLaNBinsPt;
889 fFFIMLaPtMin = o.fFFIMLaPtMin;
890 fFFIMLaPtMax = o.fFFIMLaPtMax;
891 fFFIMLaNBinsXi = o.fFFIMLaNBinsXi;
892 fFFIMLaXiMin = o.fFFIMLaXiMin;
893 fFFIMLaXiMax = o.fFFIMLaXiMax;
894 fFFIMLaNBinsZ = o.fFFIMLaNBinsZ;
895 fFFIMLaZMin = o.fFFIMLaZMin;
896 fFFIMLaZMax = o.fFFIMLaZMax;
897 fh1EvtAllCent = o.fh1EvtAllCent;
899 fh1K0Mult = o.fh1K0Mult;
900 fh1dPhiJetK0 = o.fh1dPhiJetK0;
901 fh1LaMult = o.fh1LaMult;
902 fh1dPhiJetLa = o.fh1dPhiJetLa;
903 fh1ALaMult = o.fh1ALaMult;
904 fh1dPhiJetALa = o.fh1dPhiJetALa;
905 fh1JetEta = o.fh1JetEta;
906 fh1JetPhi = o.fh1JetPhi;
907 fh2JetEtaPhi = o.fh2JetEtaPhi;
908 //fh1V0JetPt = o.fh1V0JetPt;
909 fh1IMK0Cone = o.fh1IMK0Cone;
910 fh1IMLaCone = o.fh1IMLaCone;
911 fh1IMALaCone = o.fh1IMALaCone;
912 fh2FFJetTrackEta = o.fh2FFJetTrackEta;
913 //fh1trackPosNCls = o.fh1trackPosNCls;
914 //fh1trackNegNCls = o.fh1trackNegNCls;
915 fh1trackPosRap = o.fh1trackPosRap;
916 fh1trackNegRap = o.fh1trackNegRap;
917 //fh1V0Rap = o.fh1V0Rap;
918 fh1trackPosEta = o.fh1trackPosEta;
919 fh1trackNegEta = o.fh1trackNegEta;
920 fh1V0Eta = o.fh1V0Eta;
921 // fh1V0totMom = o.fh1V0totMom;
922 fh1CosPointAngle = o.fh1CosPointAngle;
923 fh1DecayLengthV0 = o.fh1DecayLengthV0;
924 fh2ProperLifetimeK0sVsPtBeforeCut = o.fh2ProperLifetimeK0sVsPtBeforeCut;
925 fh2ProperLifetimeK0sVsPtAfterCut= o.fh2ProperLifetimeK0sVsPtAfterCut;
926 fh1V0Radius = o.fh1V0Radius;
927 fh1DcaV0Daughters = o.fh1DcaV0Daughters;
928 fh1DcaPosToPrimVertex = o.fh1DcaPosToPrimVertex;
929 fh1DcaNegToPrimVertex = o.fh1DcaNegToPrimVertex;
930 fh2ArmenterosBeforeCuts = o.fh2ArmenterosBeforeCuts;
931 fh2ArmenterosAfterCuts = o.fh2ArmenterosAfterCuts;
932 fh2BBLaPos = o.fh2BBLaPos;
933 fh2BBLaNeg = o.fh2BBLaPos;
934 fh1PosDaughterCharge = o.fh1PosDaughterCharge;
935 fh1NegDaughterCharge = o.fh1NegDaughterCharge;
936 fh1PtMCK0s = o.fh1PtMCK0s;
937 fh1PtMCLa = o.fh1PtMCLa;
938 fh1PtMCALa = o.fh1PtMCALa;
939 fh1EtaK0s = o.fh1EtaK0s;
940 fh1EtaLa = o.fh1EtaLa;
941 fh1EtaALa = o.fh1EtaALa;
943 fh1RCBiasK0 = o.fh1RCBiasK0;
944 fh1RCBiasLa = o.fh1RCBiasLa;
945 fh1RCBiasALa = o.fh1RCBiasALa;
949 fhnInvMassEtaTrackPtK0s = o.fhnInvMassEtaTrackPtK0s;
950 fhnInvMassEtaTrackPtLa = o.fhnInvMassEtaTrackPtLa;
951 fhnInvMassEtaTrackPtALa = o.fhnInvMassEtaTrackPtALa;
952 fh1TrackMultCone = o.fh1TrackMultCone;
953 fh2TrackMultCone = o.fh2TrackMultCone;
956 fhnNJALa = o.fhnNJALa;
957 //fh2MCgenK0Cone = o.fh2MCgenK0Cone;
958 //fh2MCgenLaCone = o.fh2MCgenLaCone;
959 //fh2MCgenALaCone = o.fh2MCgenALaCone;
960 //fh2MCEtagenK0Cone = o.fh2MCEtagenK0Cone;
961 //fh2MCEtagenLaCone = o.fh2MCEtagenLaCone;
962 //fh2MCEtagenALaCone = o.fh2MCEtagenALaCone;
963 fh1IMK0ConeSmear = o.fh1IMK0ConeSmear;
964 fh1IMLaConeSmear = o.fh1IMLaConeSmear;
965 fh1IMALaConeSmear = o.fh1IMALaConeSmear;
966 fh2MCEtaVsPtHijingLa = o.fh2MCEtaVsPtHijingLa;
967 fh2MCEtaVsPtInjectLa = o.fh2MCEtaVsPtInjectLa;
968 fh2MCEtaVsPtHijingALa = o.fh2MCEtaVsPtHijingALa;
969 fh2MCEtaVsPtInjectALa = o.fh2MCEtaVsPtInjectALa;
970 fhnrecMCHijingLaIncl = o.fhnrecMCHijingLaIncl;
971 fhnrecMCHijingLaCone = o.fhnrecMCHijingLaCone;
972 fhnrecMCHijingALaIncl = o.fhnrecMCHijingALaIncl;
973 fhnrecMCHijingALaCone = o.fhnrecMCHijingALaCone;
974 fhnrecMCInjectLaIncl = o.fhnrecMCInjectLaIncl;
975 fhnrecMCInjectLaCone = o.fhnrecMCInjectLaCone;
976 fhnrecMCInjectALaIncl = o.fhnrecMCInjectALaIncl;
977 fhnrecMCInjectALaCone = o.fhnrecMCInjectALaCone;
978 fhnMCrecK0Cone = o.fhnMCrecK0Cone;
979 fhnMCrecLaCone = o.fhnMCrecLaCone;
980 fhnMCrecALaCone = o.fhnMCrecALaCone;
981 fhnMCrecK0ConeSmear = o.fhnMCrecK0ConeSmear;
982 fhnMCrecLaConeSmear = o.fhnMCrecLaConeSmear;
983 fhnMCrecALaConeSmear = o.fhnMCrecALaConeSmear;
984 fhnK0sSecContinCone = o.fhnK0sSecContinCone;
985 fhnLaSecContinCone = o.fhnLaSecContinCone;
986 fhnALaSecContinCone = o.fhnALaSecContinCone;
987 fhnK0sIncl = o.fhnK0sIncl;
988 fhnK0sCone = o.fhnK0sCone;
989 fhnLaIncl = o.fhnLaIncl;
990 fhnLaCone = o.fhnLaCone;
991 fhnALaIncl = o.fhnALaIncl;
992 fhnALaCone = o.fhnALaCone;
993 fhnK0sPC = o.fhnK0sPC;
995 fhnALaPC = o.fhnALaPC;
996 fhnK0sRC = o.fhnK0sRC;
998 fhnALaRC = o.fhnALaRC;
999 fhnK0sRCBias = o.fhnK0sRCBias;
1000 fhnLaRCBias = o.fhnLaRCBias;
1001 fhnALaRCBias = o.fhnALaRCBias;
1002 fhnK0sOC = o.fhnK0sOC;
1003 fhnLaOC = o.fhnLaOC;
1004 fhnALaOC = o.fhnALaOC;
1005 fh1AreaExcluded = o.fh1AreaExcluded;
1006 fh1MedianEta = o.fh1MedianEta;
1007 fh1JetPtMedian = o.fh1JetPtMedian;
1008 fh1MCMultiplicityPrimary = o.fh1MCMultiplicityPrimary;
1009 fh1MCMultiplicityTracks = o.fh1MCMultiplicityTracks;
1010 fhnFeedDownLa = o.fhnFeedDownLa;
1011 fhnFeedDownALa = o.fhnFeedDownALa;
1012 fhnFeedDownLaCone = o.fhnFeedDownLaCone;
1013 fhnFeedDownALaCone = o.fhnFeedDownALaCone;
1014 fh1MCProdRadiusK0s = o.fh1MCProdRadiusK0s;
1015 fh1MCProdRadiusLambda = o.fh1MCProdRadiusLambda;
1016 fh1MCProdRadiusAntiLambda = o.fh1MCProdRadiusAntiLambda;
1017 fh1MCPtV0s = o.fh1MCPtV0s;
1018 fh1MCPtK0s = o.fh1MCPtK0s;
1019 fh1MCPtLambda = o.fh1MCPtLambda;
1020 fh1MCPtAntiLambda = o.fh1MCPtAntiLambda;
1021 fh1MCXiPt = o.fh1MCXiPt;
1022 fh1MCXibarPt = o.fh1MCXibarPt;
1023 fh2MCEtaVsPtK0s = o.fh2MCEtaVsPtK0s;
1024 fh2MCEtaVsPtLa = o.fh2MCEtaVsPtLa;
1025 fh2MCEtaVsPtALa = o.fh2MCEtaVsPtALa;
1026 //fh1MCRapK0s = o.fh1MCRapK0s;
1027 //fh1MCRapLambda = o.fh1MCRapLambda;
1028 //fh1MCRapAntiLambda = o.fh1MCRapAntiLambda;
1029 fh1MCEtaAllK0s = o.fh1MCEtaAllK0s;
1030 fh1MCEtaK0s = o.fh1MCEtaK0s;
1031 fh1MCEtaLambda = o.fh1MCEtaLambda;
1032 fh1MCEtaAntiLambda = o.fh1MCEtaAntiLambda;
1038 //_______________________________________________
1039 AliAnalysisTaskJetChem::~AliAnalysisTaskJetChem()
1044 if(fListK0s) delete fListK0s;
1045 if(fListLa) delete fListLa;
1046 if(fListALa) delete fListALa;
1047 if(fListFeeddownLaCand) delete fListFeeddownLaCand;
1048 if(fListFeeddownALaCand) delete fListFeeddownALaCand;
1049 if(jetConeFDLalist) delete jetConeFDLalist;
1050 if(jetConeFDALalist) delete jetConeFDALalist;
1051 if(fListMCgenK0s) delete fListMCgenK0s;
1052 if(fListMCgenLa) delete fListMCgenLa;
1053 if(fListMCgenALa) delete fListMCgenALa;
1054 if(fRandom) delete fRandom;
1057 //________________________________________________________________________________________________________________________________
1058 AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::AliFragFuncHistosInvMass(const char* name,
1059 Int_t nJetPt, Float_t jetPtMin, Float_t jetPtMax,
1060 Int_t nInvMass, Float_t invMassMin, Float_t invMassMax,
1061 Int_t nPt, Float_t ptMin, Float_t ptMax,
1062 Int_t nXi, Float_t xiMin, Float_t xiMax,
1063 Int_t nZ , Float_t zMin , Float_t zMax )
1065 ,fNBinsJetPt(nJetPt)
1066 ,fJetPtMin(jetPtMin)
1067 ,fJetPtMax(jetPtMax)
1068 ,fNBinsInvMass(nInvMass)
1069 ,fInvMassMin(invMassMin)
1070 ,fInvMassMax(invMassMax)
1086 // default constructor
1090 //______________________________________________________________________________________________________________
1091 AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::AliFragFuncHistosInvMass(const AliFragFuncHistosInvMass& copy)
1093 ,fNBinsJetPt(copy.fNBinsJetPt)
1094 ,fJetPtMin(copy.fJetPtMin)
1095 ,fJetPtMax(copy.fJetPtMax)
1096 ,fNBinsInvMass(copy.fNBinsInvMass)
1097 ,fInvMassMin(copy.fInvMassMin)
1098 ,fInvMassMax(copy.fInvMassMax)
1099 ,fNBinsPt(copy.fNBinsPt)
1100 ,fPtMin(copy.fPtMin)
1101 ,fPtMax(copy.fPtMax)
1102 ,fNBinsXi(copy.fNBinsXi)
1103 ,fXiMin(copy.fXiMin)
1104 ,fXiMax(copy.fXiMax)
1105 ,fNBinsZ(copy.fNBinsZ)
1108 ,fh3TrackPt(copy.fh3TrackPt)
1111 ,fh1JetPt(copy.fh1JetPt)
1112 ,fNameFF(copy.fNameFF)
1117 //______________________________________________________________________________________________________________________________________________________________________
1118 AliAnalysisTaskJetChem::AliFragFuncHistosInvMass& AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::operator=(const AliAnalysisTaskJetChem::AliFragFuncHistosInvMass& o)
1123 TObject::operator=(o);
1124 fNBinsJetPt = o.fNBinsJetPt;
1125 fJetPtMin = o.fJetPtMin;
1126 fJetPtMax = o.fJetPtMax;
1127 fNBinsInvMass = o.fNBinsInvMass;
1128 fInvMassMin = o.fInvMassMin;
1129 fInvMassMax = o.fInvMassMax;
1130 fNBinsPt = o.fNBinsPt;
1133 fNBinsXi = o.fNBinsXi;
1136 fNBinsZ = o.fNBinsZ;
1139 fh3TrackPt = o.fh3TrackPt;
1142 fh1JetPt = o.fh1JetPt;
1143 fNameFF = o.fNameFF;
1149 //___________________________________________________________________________
1150 AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::~AliFragFuncHistosInvMass()
1154 if(fh1JetPt) delete fh1JetPt;
1155 if(fh3TrackPt) delete fh3TrackPt;
1156 if(fh3Xi) delete fh3Xi;
1157 if(fh3Z) delete fh3Z;
1160 //_________________________________________________________________
1161 void AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::DefineHistos()
1165 fh1JetPt = new TH1F(Form("fh1FFJetPtIM%s", fNameFF.Data()),"",fNBinsJetPt,fJetPtMin,fJetPtMax);
1166 fh3TrackPt = new TH3F(Form("fh3FFTrackPtIM%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsInvMass, fInvMassMin, fInvMassMax, fNBinsPt, fPtMin, fPtMax);
1167 fh3Xi = new TH3F(Form("fh3FFXiIM%s", fNameFF.Data()),"", fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsInvMass, fInvMassMin, fInvMassMax, fNBinsXi, fXiMin, fXiMax);
1168 fh3Z = new TH3F(Form("fh3FFZIM%s", fNameFF.Data()),"", fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsInvMass, fInvMassMin, fInvMassMax, fNBinsZ, fZMin, fZMax);
1170 AliAnalysisTaskFragmentationFunction::SetProperties(fh1JetPt, "p_{t} (GeV/c)", "entries");
1171 AliAnalysisTaskJetChem::SetProperties(fh3TrackPt,"jet p_{t} (GeV/c)","inv Mass (GeV/c^2)","p_{t} (GeV/c)");
1172 AliAnalysisTaskJetChem::SetProperties(fh3Xi,"jet p_{t} (GeV/c)","inv Mass (GeV/c^2)","#xi");
1173 AliAnalysisTaskJetChem::SetProperties(fh3Z,"jet p_{t} (GeV/c)","inv Mass (GeV/c^2)","z");
1176 //________________________________________________________________________________________________________________________________
1177 void AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::FillFF(Float_t trackPt, Float_t invM, Float_t jetPt, Bool_t incrementJetPt)
1179 // fill FF, don't use TH3F anymore use THnSparse instead to save memory
1181 if(incrementJetPt) fh1JetPt->Fill(jetPt);
1182 //fh3TrackPt->Fill(jetPt,invM,trackPt);//Fill(x,y,z)
1185 if(jetPt>0) z = trackPt / jetPt;
1187 //if(z>0) xi = TMath::Log(1/z);
1189 //fh3Xi->Fill(jetPt,invM,xi);
1190 //fh3Z->Fill(jetPt,invM,z);
1193 //___________________________________________________________________________________
1194 void AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::AddToOutput(TList* list) const
1196 // add histos to list
1198 list->Add(fh1JetPt);
1199 //list->Add(fh3TrackPt);
1205 //____________________________________________________
1206 void AliAnalysisTaskJetChem::UserCreateOutputObjects()
1208 // create output objects
1210 fRandom = new TRandom3(0);
1211 fRandom->SetSeed(0);
1213 if(fDebug > 1) Printf("AliAnalysisTaskJetChem::UserCreateOutputObjects()");
1215 // create list of tracks and jets
1217 fTracksRecCuts = new TList();
1218 fTracksRecCuts->SetOwner(kFALSE); //objects in TList wont be deleted when TList is deleted
1219 fJetsRecCuts = new TList();
1220 fJetsRecCuts->SetOwner(kFALSE);
1221 fBckgJetsRec = new TList();
1222 fBckgJetsRec->SetOwner(kFALSE);
1223 fListK0s = new TList();
1224 fListK0s->SetOwner(kFALSE);
1225 fListLa = new TList();
1226 fListLa->SetOwner(kFALSE);
1227 fListALa = new TList();
1228 fListALa->SetOwner(kFALSE);
1229 fListFeeddownLaCand = new TList(); //feeddown Lambda candidates
1230 fListFeeddownLaCand->SetOwner(kFALSE);
1231 fListFeeddownALaCand = new TList(); //feeddown Antilambda candidates
1232 fListFeeddownALaCand->SetOwner(kFALSE);
1233 jetConeFDLalist = new TList();
1234 jetConeFDLalist->SetOwner(kFALSE); //feeddown Lambda candidates in jet cone
1235 jetConeFDALalist = new TList();
1236 jetConeFDALalist->SetOwner(kFALSE); //feeddown Antilambda candidates in jet cone
1237 fListMCgenK0s = new TList(); //MC generated K0s
1238 fListMCgenK0s->SetOwner(kFALSE);
1239 fListMCgenLa = new TList(); //MC generated Lambdas
1240 fListMCgenLa->SetOwner(kFALSE);
1241 fListMCgenALa = new TList(); //MC generated Antilambdas
1242 fListMCgenALa->SetOwner(kFALSE);
1245 // Create histograms / output container
1247 fCommonHistList = new TList();
1248 fCommonHistList->SetOwner();
1250 Bool_t oldStatus = TH1::AddDirectoryStatus();
1251 TH1::AddDirectory(kFALSE);//By default (fAddDirectory = kTRUE), histograms are automatically added to the list of objects in memory
1253 // histograms inherited from AliAnalysisTaskFragmentationFunction
1255 fh1EvtSelection = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
1256 fh1EvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
1257 fh1EvtSelection->GetXaxis()->SetBinLabel(2,"event trigger selection: rejected");
1258 fh1EvtSelection->GetXaxis()->SetBinLabel(3,"event class: rejected");
1259 fh1EvtSelection->GetXaxis()->SetBinLabel(4,"vertex Ncontr: rejected");
1260 fh1EvtSelection->GetXaxis()->SetBinLabel(5,"vertex z: rejected");
1261 fh1EvtSelection->GetXaxis()->SetBinLabel(6,"vertex type: rejected");
1264 fh1EvtCent = new TH1F("fh1EvtCent","centrality",100,0.,100.);
1265 fh1VertexNContributors = new TH1F("fh1VertexNContributors", "Vertex N contributors", 11,-.5, 10.5);
1266 fh1VertexZ = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
1267 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
1268 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
1269 fh1Trials = new TH1F("fh1Trials","trials from pyxsec.root",1,0,1);
1270 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
1271 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
1272 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
1273 fh1nRecJetsCuts = new TH1F("fh1nRecJetsCuts","reconstructed jets per event",100,-0.5,99.5);
1275 // histograms JetChem task
1277 fh1EvtAllCent = new TH1F("fh1EvtAllCent","before centrality selection",100,0.,100.);
1278 fh1Evt = new TH1F("fh1Evt", "All events runned over", 3, 0.,1.);
1279 fh1EvtMult = new TH1F("fh1EvtMult","multiplicity",240,0.,240.);
1280 fh1K0Mult = new TH1F("fh1K0Mult","K0 multiplicity",100,0.,100.);//500. all
1281 fh1dPhiJetK0 = new TH1F("fh1dPhiJetK0","",64,-1,5.4);
1282 fh1LaMult = new TH1F("fh1LaMult","La multiplicity",100,0.,100.);
1283 fh1dPhiJetLa = new TH1F("fh1dPhiJetLa","",64,-1,5.4);
1284 fh1ALaMult = new TH1F("fh1ALaMult","ALa multiplicity",100,0.,100.);
1285 fh1dPhiJetALa = new TH1F("fh1dPhiJetALa","",64,-1,5.4);
1286 fh1JetEta = new TH1F("fh1JetEta","#eta distribution of all jets",40,-2.,2.);
1287 fh1JetPhi = new TH1F("fh1JetPhi","#phi distribution of all jets",63,0.,6.3);
1288 fh2JetEtaPhi = new TH2F("fh2JetEtaPhi","#eta and #phi distribution of all jets",400,-2.,2.,63,0.,6.3);
1291 //fh1V0JetPt = new TH1F("fh1V0JetPt","p_{T} distribution of all jets containing v0s",200,0.,200.);
1292 fh1IMK0Cone = new TH1F("fh1IMK0Cone","p_{T} distribution of all jets containing K0s candidates",19,5.,100.);
1293 fh1IMLaCone = new TH1F("fh1IMLaCone","p_{T} distribution of all jets containing #Lambda candidates",19,5.,100.);
1294 fh1IMALaCone = new TH1F("fh1IMALaCone","p_{T} distribution of all jets containing #bar{#Lambda} candidates",19,5.,100.);
1296 fh2FFJetTrackEta = new TH2F("fh2FFJetTrackEta","charged track eta distr. in jet cone",200,-1.,1.,40,0.,200.);
1297 //fh1trackPosNCls = new TH1F("fh1trackPosNCls","NTPC clusters positive daughters",10,0.,100.);
1298 //fh1trackNegNCls = new TH1F("fh1trackNegNCls","NTPC clusters negative daughters",10,0.,100.);
1299 fh1trackPosEta = new TH1F("fh1trackPosEta","eta positive daughters",100,-2.,2.);
1300 fh1trackNegEta = new TH1F("fh1trackNegEta","eta negative daughters",100,-2.,2.);
1301 fh1V0Eta = new TH1F("fh1V0Eta","V0 eta",60,-1.5,1.5);
1302 //fh1V0totMom = new TH1F("fh1V0totMom","V0 tot mom",100,0.,20.);
1303 fh1CosPointAngle = new TH1F("fh1CosPointAngle", "Cosine of V0's pointing angle",50,0.99,1.0);
1304 fh1DecayLengthV0 = new TH1F("fh1DecayLengthV0", "V0s decay Length;decay length(cm)",1200,0.,120.);
1305 fh2ProperLifetimeK0sVsPtBeforeCut = new TH2F("fh2ProperLifetimeK0sVsPtBeforeCut"," K0s ProperLifetime vs Pt; p_{T} (GeV/#it{c})",150,0.,15.,250,0.,250.);
1306 fh2ProperLifetimeK0sVsPtAfterCut = new TH2F("fh2ProperLifetimeK0sVsPtAfterCut"," K0s ProperLifetime vs Pt; p_{T} (GeV/#it{c})",150,0.,15.,250,0.,250.);
1307 fh1V0Radius = new TH1F("fh1V0Radius", "V0s Radius;Radius(cm)",200,0.,40.);
1308 fh1DcaV0Daughters = new TH1F("fh1DcaV0Daughters", "DCA between daughters;dca(cm)",200,0.,2.);
1309 fh1DcaPosToPrimVertex = new TH1F("fh1DcaPosToPrimVertex", "Positive V0 daughter;dca(cm)",100,0.,10.);
1310 fh1DcaNegToPrimVertex = new TH1F("fh1DcaNegToPrimVertex", "Negative V0 daughter;dca(cm)",100,0.,10.);
1311 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);
1312 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);
1313 fh2BBLaPos = new TH2F("fh2BBLaPos","PID of the positive daughter of La candidates; P (GeV); -dE/dx (keV/cm ?)",100,0,10,200,0,200);
1314 fh2BBLaNeg = new TH2F("fh2BBLaNeg","PID of the negative daughter of La candidates; P (GeV); -dE/dx (keV/cm ?)",100,0,10,200,0,200);
1315 fh1PosDaughterCharge = new TH1F("fh1PosDaughterCharge","charge of V0 positive daughters; V0 daughters",3,-2.,2.);
1316 fh1NegDaughterCharge = new TH1F("fh1NegDaughterCharge","charge of V0 negative daughters; V0 daughters",3,-2.,2.);
1317 fh1PtMCK0s = new TH1F("fh1PtMCK0s","Pt of MC rec K0s; #it{p}_{T} (GeV/#it{c})",200,0.,20.);
1318 fh1PtMCLa = new TH1F("fh1PtMCLa","Pt of MC rec La; #it{p}_{T} (GeV/#it{c})",200,0.,20.);
1319 fh1PtMCALa = new TH1F("fh1PtMCALa","Pt of MC rec ALa; #it{p}_{T} (GeV/#it{c})",200,0.,20.);
1320 fh1EtaK0s = new TH1F("fh1EtaK0s","K^{0}_{s} entries ;#eta",200,-1.,1.);
1321 fh1EtaLa = new TH1F("fh1EtaLa","#Lambda entries ;#eta",200,-1.,1.);
1322 fh1EtaALa = new TH1F("fh1EtaALa","#bar{#Lambda} entries ;#eta",200,-1.,1.);
1324 //histos for normalisation of MCC, RC, OC and NJ
1326 fh1RC = new TH1F("fh1RC"," # random cones used",1,0.5,1.5);
1327 fh1RCBiasK0 = new TH1F("fh1RCBiasK0"," # random cones with K0s trigger particle",1,0.5,1.5);
1328 fh1RCBiasLa = new TH1F("fh1RCBiasLa"," # random cones with La trigger particle",1,0.5,1.5);
1329 fh1RCBiasALa = new TH1F("fh1RCBiasALa"," # random cones with ALa trigger particle",1,0.5,1.5);
1330 fh1MCC = new TH1F("fh1MCC","# median cluster cones used",1,0.5,1.5);
1331 fh1OC = new TH1F("fh1OC","# outside cones used, number of jet events",1,0.5,1.5);
1332 fh1NJ = new TH1F("fh1NJ","# non-jet events used",1,0.5,1.5);
1334 Int_t binsInvMassEtaTrackPtK0s[3] = {200, 200, 120};//eta,invM,trackPt
1335 Double_t xminInvMassEtaTrackPtK0s[3] = {-1.,0.3,0.};
1336 Double_t xmaxInvMassEtaTrackPtK0s[3] = {1.,0.7,12.};
1338 fhnInvMassEtaTrackPtK0s = new THnSparseF("fhnInvMassEtaTrackPtK0s","#eta; K0s invM (GeV/{#it{c}}^{2}); #it{p}_{T} (GeV/#it{c})",3,binsInvMassEtaTrackPtK0s,xminInvMassEtaTrackPtK0s,xmaxInvMassEtaTrackPtK0s);
1340 Int_t binsInvMassEtaTrackPtLa[3] = {200, 200, 120};//eta,invM,trackPt
1341 Double_t xminInvMassEtaTrackPtLa[3] = {-1.,1.05,0.};
1342 Double_t xmaxInvMassEtaTrackPtLa[3] = {1.,1.25,12.};
1344 fhnInvMassEtaTrackPtLa = new THnSparseF("fhnInvMassEtaTrackPtLa","#eta; #Lambda invM (GeV/{#it{c}}^{2}); #it{p}_{T} (GeV/#it{c})",3,binsInvMassEtaTrackPtLa,xminInvMassEtaTrackPtLa,xmaxInvMassEtaTrackPtLa);
1346 Int_t binsInvMassEtaTrackPtALa[3] = {200, 200, 120};//eta,invM,trackPt
1347 Double_t xminInvMassEtaTrackPtALa[3] = {-1.,1.05,0.};
1348 Double_t xmaxInvMassEtaTrackPtALa[3] = {1.,1.25,12.};
1350 fhnInvMassEtaTrackPtALa = new THnSparseF("fhnInvMassEtaTrackPtALa","#eta; #bar{#Lambda} invM (GeV/{#it{c}}^{2}); #it{p}_{T} (GeV/#it{c})",3,binsInvMassEtaTrackPtALa,xminInvMassEtaTrackPtALa,xmaxInvMassEtaTrackPtALa);
1352 Int_t binsK0sPC[4] = {19, 200, 200, 200};
1353 Double_t xminK0sPC[4] = {5.,0.3, 0., -1.};
1354 Double_t xmaxK0sPC[4] = {100.,0.7, 20., 1.};
1355 fhnK0sPC = new THnSparseF("fhnK0sPC","jet pT; K0s invM; particle pT; particle #eta",4,binsK0sPC,xminK0sPC,xmaxK0sPC);
1357 Int_t binsLaPC[4] = {19, 200, 200, 200};
1358 Double_t xminLaPC[4] = {5.,1.05, 0., -1.};
1359 Double_t xmaxLaPC[4] = {100.,1.25, 20., 1.};
1360 fhnLaPC = new THnSparseF("fhnLaPC","jet pT; #Lambda invM; particle pT; particle #eta",4,binsLaPC,xminLaPC,xmaxLaPC);
1362 Int_t binsALaPC[4] = {19, 200, 200, 200};
1363 Double_t xminALaPC[4] = {5.,1.05, 0., -1.};
1364 Double_t xmaxALaPC[4] = {100.,1.25, 20., 1.};
1365 fhnALaPC = new THnSparseF("fhnALaPC","jet pT; #bar#Lambda invM; particle pT; particle #eta",4,binsALaPC,xminALaPC,xmaxALaPC);
1367 Int_t binsK0sMCC[3] = {200, 200, 200};
1368 Double_t xminK0sMCC[3] = {0.3, 0., -1.};
1369 Double_t xmaxK0sMCC[3] = {0.7, 20., 1.};
1370 fhnK0sMCC = new THnSparseF("fhnK0sMCC","jet pT; K0s invM; particle pT; particle #eta",3,binsK0sMCC,xminK0sMCC,xmaxK0sMCC);
1372 Int_t binsLaMCC[3] = {200, 200, 200};
1373 Double_t xminLaMCC[3] = {1.05, 0., -1.};
1374 Double_t xmaxLaMCC[3] = {1.25, 20., 1.};
1375 fhnLaMCC = new THnSparseF("fhnLaMCC","jet pT; #Lambda invM; particle pT; particle #eta",3,binsLaMCC,xminLaMCC,xmaxLaMCC);
1377 Int_t binsALaMCC[3] = {200, 200, 200};
1378 Double_t xminALaMCC[3] = {1.05, 0., -1.};
1379 Double_t xmaxALaMCC[3] = {1.25, 20., 1.};
1380 fhnALaMCC = new THnSparseF("fhnALaMCC","jet pT; #bara#Lambda invM; particle pT; particle #eta",3,binsALaMCC,xminALaMCC,xmaxALaMCC);
1382 Int_t binsK0sRC[3] = {200, 200, 200};
1383 Double_t xminK0sRC[3] = {0.3, 0., -1.};
1384 Double_t xmaxK0sRC[3] = {0.7, 20., 1.};
1385 fhnK0sRC = new THnSparseF("fhnK0sRC","jet pT; K0s invM; particle pT; particle #eta",3,binsK0sRC,xminK0sRC,xmaxK0sRC);
1387 Int_t binsLaRC[3] = {200, 200, 200};
1388 Double_t xminLaRC[3] = {1.05, 0., -1.};
1389 Double_t xmaxLaRC[3] = {1.25, 20., 1.};
1390 fhnLaRC = new THnSparseF("fhnLaRC","jet pT; #Lambda invM; particle pT; particle #eta",3,binsLaRC,xminLaRC,xmaxLaRC);
1392 Int_t binsALaRC[3] = {200, 200, 200};
1393 Double_t xminALaRC[3] = {1.05, 0., -1.};
1394 Double_t xmaxALaRC[3] = {1.25, 20., 1.};
1395 fhnALaRC = new THnSparseF("fhnALaRC","jet pT; #bara#Lambda invM; particle pT; particle #eta",3,binsALaRC,xminALaRC,xmaxALaRC);
1397 Int_t binsK0sRCBias[3] = {200, 200, 200};
1398 Double_t xminK0sRCBias[3] = {0.3, 0., -1.};
1399 Double_t xmaxK0sRCBias[3] = {0.7, 20., 1.};
1400 fhnK0sRCBias = new THnSparseF("fhnK0sRCBias","jet pT; K0s invM; particle pT; particle #eta",3,binsK0sRCBias,xminK0sRCBias,xmaxK0sRCBias);
1402 Int_t binsLaRCBias[3] = {200, 200, 200};
1403 Double_t xminLaRCBias[3] = {1.05, 0., -1.};
1404 Double_t xmaxLaRCBias[3] = {1.25, 20., 1.};
1405 fhnLaRCBias = new THnSparseF("fhnLaRCBias","jet pT; #Lambda invM; particle pT; particle #eta",3,binsLaRCBias,xminLaRCBias,xmaxLaRCBias);
1407 Int_t binsALaRCBias[3] = {200, 200, 200};
1408 Double_t xminALaRCBias[3] = {1.05, 0., -1.};
1409 Double_t xmaxALaRCBias[3] = {1.25, 20., 1.};
1410 fhnALaRCBias = new THnSparseF("fhnALaRCBias","jet pT; #bara#Lambda invM; particle pT; particle #eta",3,binsALaRCBias,xminALaRCBias,xmaxALaRCBias);
1412 Int_t binsK0sOC[3] = {200, 200, 200};
1413 Double_t xminK0sOC[3] = {0.3, 0., -1.};
1414 Double_t xmaxK0sOC[3] = {0.7, 20., 1.};
1415 fhnK0sOC = new THnSparseF("fhnK0sOC","jet pT; K0s invM; particle pT; particle #eta",3,binsK0sOC,xminK0sOC,xmaxK0sOC);
1417 Int_t binsLaOC[3] = {200, 200, 200};
1418 Double_t xminLaOC[3] = {1.05, 0., -1.};
1419 Double_t xmaxLaOC[3] = {1.25, 20., 1.};
1420 fhnLaOC = new THnSparseF("fhnLaOC","jet pT; #Lambda invM; particle pT; particle #eta",3,binsLaOC,xminLaOC,xmaxLaOC);
1422 Int_t binsALaOC[3] = {200, 200, 200};
1423 Double_t xminALaOC[3] = {1.05, 0., -1.};
1424 Double_t xmaxALaOC[3] = {1.25, 20., 1.};
1426 fhnALaOC = new THnSparseF("fhnALaOC","jet pT; #bara#Lambda invM; particle pT; particle #eta",3,binsALaOC,xminALaOC,xmaxALaOC);
1428 fh1AreaExcluded = new TH1F("fh1AreaExcluded","area excluded for selected jets in event acceptance",50,0.,1.);
1430 fh1MedianEta = new TH1F("fh1MedianEta","Median cluster axis ;#eta",200,-1.,1.);
1431 fh1JetPtMedian = new TH1F("fh1JetPtMedian"," (selected) jet it{p}_{T} distribution for MCC method; #GeV/it{c}",19,5.,100.);
1433 fh1TrackMultCone = new TH1F("fh1TrackMultCone","track multiplicity in jet cone; number of tracks",20,0.,50.);
1435 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.);
1437 Int_t binsNJK0[3] = {200, 200, 200};
1438 Double_t xminNJK0[3] = {0.3, 0., -1.};
1439 Double_t xmaxNJK0[3] = {0.7, 20., 1.};
1440 fhnNJK0 = new THnSparseF("fhnNJK0","K0s candidates in events wo selected jets;",3,binsNJK0,xminNJK0,xmaxNJK0);
1442 Int_t binsNJLa[3] = {200, 200, 200};
1443 Double_t xminNJLa[3] = {1.05, 0., -1.};
1444 Double_t xmaxNJLa[3] = {1.25, 20., 1.};
1445 fhnNJLa = new THnSparseF("fhnNJLa","La candidates in events wo selected jets; ",3,binsNJLa,xminNJLa,xmaxNJLa);
1447 Int_t binsNJALa[3] = {200, 200, 200};
1448 Double_t xminNJALa[3] = {1.05, 0., -1.};
1449 Double_t xmaxNJALa[3] = {1.25, 20., 1.};
1450 fhnNJALa = new THnSparseF("fhnNJALa","ALa candidates in events wo selected jets; ",3,binsNJALa,xminNJALa,xmaxNJALa);
1452 fFFHistosRecCuts = new AliFragFuncHistos("RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1453 fFFNBinsPt, fFFPtMin, fFFPtMax,
1454 fFFNBinsXi, fFFXiMin, fFFXiMax,
1455 fFFNBinsZ , fFFZMin , fFFZMax);
1457 fV0QAK0 = new AliFragFuncQATrackHistos("V0QAK0",fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1458 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1459 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1460 fQATrackHighPtThreshold);
1462 fFFHistosRecCutsK0Evt = new AliFragFuncHistos("RecCutsK0Evt", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1463 fFFNBinsPt, fFFPtMin, fFFPtMax,
1464 fFFNBinsXi, fFFXiMin, fFFXiMax,
1465 fFFNBinsZ , fFFZMin , fFFZMax);
1468 fFFHistosIMK0AllEvt = new AliFragFuncHistosInvMass("K0AllEvt", fFFIMNBinsJetPt, fFFIMJetPtMin, fFFIMJetPtMax,
1469 fFFIMNBinsInvM,fFFIMInvMMin,fFFIMInvMMax,
1470 fFFIMNBinsPt, fFFIMPtMin, fFFIMPtMax,
1471 fFFIMNBinsXi, fFFIMXiMin, fFFIMXiMax,
1472 fFFIMNBinsZ , fFFIMZMin , fFFIMZMax);
1474 fFFHistosIMK0Jet = new AliFragFuncHistosInvMass("K0Jet", fFFIMNBinsJetPt, fFFIMJetPtMin, fFFIMJetPtMax,
1475 fFFIMNBinsInvM,fFFIMInvMMin,fFFIMInvMMax,
1476 fFFIMNBinsPt, fFFIMPtMin, fFFIMPtMax,
1477 fFFIMNBinsXi, fFFIMXiMin, fFFIMXiMax,
1478 fFFIMNBinsZ , fFFIMZMin , fFFIMZMax);
1480 fFFHistosIMK0Cone = new AliFragFuncHistosInvMass("K0Cone", fFFIMNBinsJetPt, fFFIMJetPtMin, fFFIMJetPtMax,
1481 fFFIMNBinsInvM,fFFIMInvMMin,fFFIMInvMMax,
1482 fFFIMNBinsPt, fFFIMPtMin, fFFIMPtMax,
1483 fFFIMNBinsXi, fFFIMXiMin, fFFIMXiMax,
1484 fFFIMNBinsZ , fFFIMZMin , fFFIMZMax);
1486 fFFHistosIMLaAllEvt = new AliFragFuncHistosInvMass("LaAllEvt", fFFIMLaNBinsJetPt, fFFIMLaJetPtMin, fFFIMLaJetPtMax,
1487 fFFIMLaNBinsInvM,fFFIMLaInvMMin,fFFIMLaInvMMax,
1488 fFFIMLaNBinsPt, fFFIMLaPtMin, fFFIMLaPtMax,
1489 fFFIMLaNBinsXi, fFFIMLaXiMin, fFFIMLaXiMax,
1490 fFFIMLaNBinsZ , fFFIMLaZMin , fFFIMLaZMax);
1492 fFFHistosIMLaJet = new AliFragFuncHistosInvMass("LaJet", fFFIMLaNBinsJetPt, fFFIMLaJetPtMin, fFFIMLaJetPtMax,
1493 fFFIMLaNBinsInvM,fFFIMLaInvMMin,fFFIMLaInvMMax,
1494 fFFIMLaNBinsPt, fFFIMLaPtMin, fFFIMLaPtMax,
1495 fFFIMLaNBinsXi, fFFIMLaXiMin, fFFIMLaXiMax,
1496 fFFIMLaNBinsZ , fFFIMLaZMin , fFFIMLaZMax);
1499 fFFHistosIMLaCone = new AliFragFuncHistosInvMass("LaCone", fFFIMLaNBinsJetPt, fFFIMLaJetPtMin, fFFIMLaJetPtMax,
1500 fFFIMLaNBinsInvM,fFFIMLaInvMMin,fFFIMLaInvMMax,
1501 fFFIMLaNBinsPt, fFFIMLaPtMin, fFFIMLaPtMax,
1502 fFFIMLaNBinsXi, fFFIMLaXiMin, fFFIMLaXiMax,
1503 fFFIMLaNBinsZ , fFFIMLaZMin , fFFIMLaZMax);
1506 fFFHistosIMALaAllEvt = new AliFragFuncHistosInvMass("ALaAllEvt", fFFIMLaNBinsJetPt, fFFIMLaJetPtMin, fFFIMLaJetPtMax,
1507 fFFIMLaNBinsInvM,fFFIMLaInvMMin,fFFIMLaInvMMax,
1508 fFFIMLaNBinsPt, fFFIMLaPtMin, fFFIMLaPtMax,
1509 fFFIMLaNBinsXi, fFFIMLaXiMin, fFFIMLaXiMax,
1510 fFFIMLaNBinsZ , fFFIMLaZMin , fFFIMLaZMax);
1512 fFFHistosIMALaJet = new AliFragFuncHistosInvMass("ALaJet", fFFIMLaNBinsJetPt, fFFIMLaJetPtMin, fFFIMLaJetPtMax,
1513 fFFIMLaNBinsInvM,fFFIMLaInvMMin,fFFIMLaInvMMax,
1514 fFFIMLaNBinsPt, fFFIMLaPtMin, fFFIMLaPtMax,
1515 fFFIMLaNBinsXi, fFFIMLaXiMin, fFFIMLaXiMax,
1516 fFFIMLaNBinsZ , fFFIMLaZMin , fFFIMLaZMax);
1518 fFFHistosIMALaCone = new AliFragFuncHistosInvMass("ALaCone", fFFIMLaNBinsJetPt, fFFIMLaJetPtMin, fFFIMLaJetPtMax,
1519 fFFIMLaNBinsInvM,fFFIMLaInvMMin,fFFIMLaInvMMax,
1520 fFFIMLaNBinsPt, fFFIMLaPtMin, fFFIMLaPtMax,
1521 fFFIMLaNBinsXi, fFFIMLaXiMin, fFFIMLaXiMax,
1522 fFFIMLaNBinsZ , fFFIMLaZMin , fFFIMLaZMax);
1529 //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.);
1530 //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.);
1531 //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.);
1533 //fh2MCgenK0Cone->GetYaxis()->SetTitle("MC gen K^{0}}^{s} #it{p}_{T}");
1534 //fh2MCgenLaCone->GetYaxis()->SetTitle("MC gen #Lambda #it{p}_{T}");
1535 //fh2MCgenALaCone->GetYaxis()->SetTitle("MC gen #Antilambda #it{p}_{T}");
1537 //fh2MCEtagenK0Cone = new TH2F("fh2MCEtagenK0Cone","MC gen {K^{0}}^{s} #it{p}_{T} #eta distribution in jet cone;#eta",19,5.,100.,200,-1.,1.);
1538 //fh2MCEtagenLaCone = new TH2F("fh2MCEtagenLaCone","MC gen #Lambda #it{p}_{T} #eta distribution in jet cone;#eta",19,5.,100.,200,-1.,1.);
1539 //fh2MCEtagenALaCone = new TH2F("fh2MCEtagenALaCone","MC gen #Antilambda #it{p}_{T} #eta distribution in jet cone;#eta",19,5.,100.,200,-1.,1.);
1540 fh1IMK0ConeSmear = new TH1F("fh1IMK0ConeSmear","Smeared jet pt study for K0s-in-cone-jets; smeared jet #it{p}_{T}", 19,5.,100.);
1541 fh1IMLaConeSmear = new TH1F("fh1IMLaConeSmear","Smeared jet pt study for La-in-cone-jets; smeared jet #it{p}_{T}", 19,5.,100.);
1542 fh1IMALaConeSmear = new TH1F("fh1IMALaConeSmear","Smeared jet pt study for ALa-in-cone-jets; smeared jet #it{p}_{T}", 19,5.,100.);
1544 fh2CorrHijingLaProton = new TH2F("fh2CorrHijingLaProton","#Lambda - proton pT correlation, Hijing;#it{p^{#Lambda}}_{T} (GeV/c);#it{p^{proton}}_{T} (GeV/c)",20,0.,20.,20,0.,20.);
1545 fh2CorrInjectLaProton = new TH2F("fh2CorrInjectLaProton","#Lambda - proton pT correlation, Injected;#it{p^{#Lambda}}_{T} (GeV/c);#it{p^{proton}}_{T} (GeV/c)",20,0.,20.,20,0.,20.);
1546 fh2CorrHijingALaAProton = new TH2F("fh2CorrHijingALaAProton","#bar{#Lambda} - proton pT correlation, Hijing;#it{p^{#Lambda}}_{T} (GeV/c);#it{p^{#bar{proton}}}_{T} (GeV/c)",20,0.,20.,20,0.,20.);
1547 fh2CorrInjectALaAProton = new TH2F("fh2CorrInjectALaAProton","#bar{#Lambda} - proton pT correlation, Injected;#it{p^{#Lambda}}_{T} (GeV/c);#it{p^{#bar{proton}}}_{T} (GeV/c)",20,0.,20.,20,0.,20.);
1548 //12 new histograms: Cone, Incl, Lambda, Antilambda, Hijing, Injected:
1550 fh2MCEtaVsPtHijingLa = new TH2F("fh2MCEtaVsPtHijingLa","MC Hijing gen. #Lambda #eta; #it{p}_{T}",200,0.,20.,200,-1.,1.);
1551 fh2MCEtaVsPtInjectLa = new TH2F("fh2MCEtaVsPtInjectLa","MC injected gen. #Lambda #eta; #it{p}_{T}",200,0.,20.,200,-1.,1.);
1552 fh2MCEtaVsPtHijingALa = new TH2F("fh2MCEtaVsPtHijingALa","MC gen. Hijing #bar{#Lambda} #eta; #it{p}_{T}",200,0.,20.,200,-1.,1.);
1553 fh2MCEtaVsPtInjectALa = new TH2F("fh2MCEtaVsPtInjectALa","MC gen. injected #bar{#Lambda} #eta; #it{p}_{T}",200,0.,20.,200,-1.,1.);
1555 Int_t binsrecMCHijingLaIncl[3] = {200, 200, 200};
1556 Double_t xminrecMCHijingLaIncl[3] = {1.05, 0., -1.};
1557 Double_t xmaxrecMCHijingLaIncl[3] = {1.25, 20., 1.};
1558 fhnrecMCHijingLaIncl = new THnSparseF("fhnrecMCHijingLaIncl","La inv. mass; particle pT; particle #eta",3,binsrecMCHijingLaIncl,xminrecMCHijingLaIncl,xmaxrecMCHijingLaIncl);
1560 Int_t binsrecMCHijingLaCone[4] = {19, 200, 200, 200};
1561 Double_t xminrecMCHijingLaCone[4] = {5., 1.05, 0., -1.};
1562 Double_t xmaxrecMCHijingLaCone[4] = {100., 1.25, 20., 1.};
1563 fhnrecMCHijingLaCone = new THnSparseF("fhnrecMCHijingLaCone","La inv. mass; particle pT; particle #eta",4,binsrecMCHijingLaCone,xminrecMCHijingLaCone,xmaxrecMCHijingLaCone);
1565 Int_t binsrecMCHijingALaIncl[3] = {200, 200, 200};
1566 Double_t xminrecMCHijingALaIncl[3] = {1.05, 0., -1.};
1567 Double_t xmaxrecMCHijingALaIncl[3] = {1.25, 20., 1.};
1568 fhnrecMCHijingALaIncl = new THnSparseF("fhnrecMCHijingALaIncl","ALa inv. mass; particle pT; particle #eta",3,binsrecMCHijingALaIncl,xminrecMCHijingALaIncl,xmaxrecMCHijingALaIncl);
1570 Int_t binsrecMCHijingALaCone[4] = {19, 200, 200, 200};
1571 Double_t xminrecMCHijingALaCone[4] = {5., 1.05, 0., -1.};
1572 Double_t xmaxrecMCHijingALaCone[4] = {100., 1.25, 20., 1.};
1573 fhnrecMCHijingALaCone = new THnSparseF("fhnrecMCHijingALaCone","ALa inv. mass; particle pT; particle #eta",4,binsrecMCHijingALaCone,xminrecMCHijingALaCone,xmaxrecMCHijingALaCone);
1575 Int_t binsrecMCInjectLaIncl[3] = {200, 200, 200};
1576 Double_t xminrecMCInjectLaIncl[3] = {1.05, 0., -1.};
1577 Double_t xmaxrecMCInjectLaIncl[3] = {1.25, 20., 1.};
1578 fhnrecMCInjectLaIncl = new THnSparseF("fhnrecMCInjectLaIncl","La inv. mass; particle pT; particle #eta",3,binsrecMCInjectLaIncl,xminrecMCInjectLaIncl,xmaxrecMCInjectLaIncl);
1580 Int_t binsrecMCInjectLaCone[4] = {19, 200, 200, 200};
1581 Double_t xminrecMCInjectLaCone[4] = {5., 1.05, 0., -1.};
1582 Double_t xmaxrecMCInjectLaCone[4] = {100., 1.25, 20., 1.};
1583 fhnrecMCInjectLaCone = new THnSparseF("fhnrecMCInjectLaCone","La jet pT;inv. mass; particle pT; particle #eta",4,binsrecMCInjectLaCone,xminrecMCInjectLaCone,xmaxrecMCInjectLaCone);
1585 Int_t binsrecMCInjectALaIncl[3] = {200, 200, 200};
1586 Double_t xminrecMCInjectALaIncl[3] = {1.05, 0., -1.};
1587 Double_t xmaxrecMCInjectALaIncl[3] = {1.25, 20., 1.};
1588 fhnrecMCInjectALaIncl = new THnSparseF("fhnrecMCInjectALaIncl","ALa inv. mass; particle pT; particle #eta",3,binsrecMCInjectALaIncl,xminrecMCInjectALaIncl,xmaxrecMCInjectALaIncl);
1590 Int_t binsrecMCInjectALaCone[4] = {19, 200, 200, 200};
1591 Double_t xminrecMCInjectALaCone[4] = {5., 1.05, 0., -1.};
1592 Double_t xmaxrecMCInjectALaCone[4] = {100., 1.25, 20., 1.};
1593 fhnrecMCInjectALaCone = new THnSparseF("fhnrecMCInjectALaCone","ALa inv. mass; particle pT; particle #eta",4,binsrecMCInjectALaCone,xminrecMCInjectALaCone,xmaxrecMCInjectALaCone);
1596 Int_t binsMCrecK0Cone[4] = {19, 200, 200, 200};
1597 Double_t xminMCrecK0Cone[4] = {5.,0.3, 0., -1.};
1598 Double_t xmaxMCrecK0Cone[4] = {100.,0.7, 20., 1.};
1599 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);
1601 Int_t binsMCrecLaCone[4] = {19, 200, 200, 200};
1602 Double_t xminMCrecLaCone[4] = {5.,0.3, 0., -1.};
1603 Double_t xmaxMCrecLaCone[4] = {100.,0.7, 20., 1.};
1604 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);
1606 Int_t binsMCrecALaCone[4] = {19, 200, 200, 200};
1607 Double_t xminMCrecALaCone[4] = {5.,0.3, 0., -1.};
1608 Double_t xmaxMCrecALaCone[4] = {100.,0.7, 20., 1.};
1609 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);
1611 Int_t binsMCrecK0ConeSmear[4] = {19, 200, 200, 200};
1612 Double_t xminMCrecK0ConeSmear[4] = {5.,0.3, 0., -1.};
1613 Double_t xmaxMCrecK0ConeSmear[4] = {100.,0.7, 20., 1.};
1614 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);
1616 Int_t binsMCrecLaConeSmear[4] = {19, 200, 200, 200};
1617 Double_t xminMCrecLaConeSmear[4] = {5.,1.05, 0., -1.};
1618 Double_t xmaxMCrecLaConeSmear[4] = {100.,1.25, 20., 1.};
1619 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);
1621 Int_t binsMCrecALaConeSmear[4] = {19, 200, 200, 200};
1622 Double_t xminMCrecALaConeSmear[4] = {5.,1.05, 0., -1.};
1623 Double_t xmaxMCrecALaConeSmear[4] = {100.,1.25, 20., 1.};
1624 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);
1626 Int_t binsK0sSecContinCone[3] = {19, 200, 200};
1627 Double_t xminK0sSecContinCone[3] = {5.,0., -1.};
1628 Double_t xmaxK0sSecContinCone[3] = {100.,20., 1.};
1629 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);
1631 Int_t binsLaSecContinCone[3] = {19, 200, 200};
1632 Double_t xminLaSecContinCone[3] = {5.,0., -1.};
1633 Double_t xmaxLaSecContinCone[3] = {100.,20., 1.};
1634 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);
1636 Int_t binsALaSecContinCone[3] = {19, 200, 200};
1637 Double_t xminALaSecContinCone[3] = {5.,0., -1.};
1638 Double_t xmaxALaSecContinCone[3] = {100.,20., 1.};
1639 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);
1641 Int_t binsK0sIncl[3] = {200, 200, 200};
1642 Double_t xminK0sIncl[3] = {0.3, 0., -1.};
1643 Double_t xmaxK0sIncl[3] = {0.7, 20., 1.};
1644 fhnK0sIncl = new THnSparseF("fhnK0sIncl","K0s inv. mass; particle pT; particle #eta",3,binsK0sIncl,xminK0sIncl,xmaxK0sIncl);
1646 Int_t binsK0sCone[4] = {19, 200, 200, 200};
1647 Double_t xminK0sCone[4] = {5.,0.3, 0., -1.};
1648 Double_t xmaxK0sCone[4] = {100.,0.7, 20., 1.};
1649 fhnK0sCone = new THnSparseF("fhnK0sCone","jet pT; K0s inv. mass; particle pT; particle #eta",4,binsK0sCone,xminK0sCone,xmaxK0sCone);
1651 Int_t binsLaIncl[3] = {200, 200, 200};
1652 Double_t xminLaIncl[3] = {1.05, 0., -1.};
1653 Double_t xmaxLaIncl[3] = {1.25, 20., 1.};
1654 fhnLaIncl = new THnSparseF("fhnLaIncl","La inv. mass; particle pT; particle #eta",3,binsLaIncl,xminLaIncl,xmaxLaIncl);
1656 Int_t binsLaCone[4] = {19, 200, 200, 200};
1657 Double_t xminLaCone[4] = {5.,1.05, 0., -1.};
1658 Double_t xmaxLaCone[4] = {100.,1.25, 20., 1.};
1659 fhnLaCone = new THnSparseF("fhnLaCone","jet pT; La inv. mass; particle pT; particle #eta",4,binsLaCone,xminLaCone,xmaxLaCone);
1661 Int_t binsALaIncl[3] = {200, 200, 200};
1662 Double_t xminALaIncl[3] = {1.05, 0., -1.};
1663 Double_t xmaxALaIncl[3] = {1.25, 20., 1.};
1664 fhnALaIncl = new THnSparseF("fhnALaIncl","ALa inv. mass; particle pT; particle #eta",3,binsALaIncl,xminALaIncl,xmaxALaIncl);
1666 Int_t binsALaCone[4] = {19, 200, 200, 200};
1667 Double_t xminALaCone[4] = {5.,1.05, 0., -1.};
1668 Double_t xmaxALaCone[4] = {100.,1.25, 20., 1.};
1669 fhnALaCone = new THnSparseF("fhnALaCone","jet pT; ALa inv. mass; particle pT; particle #eta",4,binsALaCone,xminALaCone,xmaxALaCone);
1671 fh1MCMultiplicityPrimary = new TH1F("fh1MCMultiplicityPrimary", "MC Primary Particles;NPrimary;Count", 201, -0.5, 200.5);
1672 fh1MCMultiplicityTracks = new TH1F("h1MCMultiplicityTracks", "MC Tracks;Ntracks;Count", 201, -0.5, 200.5);
1675 Int_t binsFeedDownLa[3] = {19, 200, 200};
1676 Double_t xminFeedDownLa[3] = {5.,1.05, 0.};
1677 Double_t xmaxFeedDownLa[3] = {100.,1.25, 20.};
1678 fhnFeedDownLa = new THnSparseF("fhnFeedDownLa","#Lambda stemming from feeddown from Xi(0/-)",3,binsFeedDownLa,xminFeedDownLa,xmaxFeedDownLa);
1680 Int_t binsFeedDownALa[3] = {19, 200, 200};
1681 Double_t xminFeedDownALa[3] = {5.,1.05, 0.};
1682 Double_t xmaxFeedDownALa[3] = {100.,1.25, 20.};
1683 fhnFeedDownALa = new THnSparseF("fhnFeedDownALa","#bar#Lambda stemming from feeddown from Xibar(0/+)",3,binsFeedDownALa,xminFeedDownALa,xmaxFeedDownALa);
1685 Int_t binsFeedDownLaCone[3] = {19, 200, 200};
1686 Double_t xminFeedDownLaCone[3] = {5.,1.05, 0.};
1687 Double_t xmaxFeedDownLaCone[3] = {100.,1.25, 20.};
1688 fhnFeedDownLaCone = new THnSparseF("fhnFeedDownLaCone","#Lambda stemming from feeddown from Xi(0/-) in jet cone",3,binsFeedDownLaCone,xminFeedDownLaCone,xmaxFeedDownLaCone);
1690 Int_t binsFeedDownALaCone[3] = {19, 200, 200};
1691 Double_t xminFeedDownALaCone[3] = {5.,1.05, 0.};
1692 Double_t xmaxFeedDownALaCone[3] = {100.,1.25, 20.};
1693 fhnFeedDownALaCone = new THnSparseF("fhnFeedDownALaCone","#bar#Lambda stemming from feeddown from Xibar(0/+) in jet cone",3,binsFeedDownALaCone,xminFeedDownALaCone,xmaxFeedDownALaCone);
1696 fh1MCProdRadiusK0s = new TH1F("fh1MCProdRadiusK0s","MC gen. MC K0s prod radius",200,0.,200.);
1697 fh1MCProdRadiusLambda = new TH1F("fh1MCProdRadiusLambda","MC gen. MC La prod radius",200,0.,200.);
1698 fh1MCProdRadiusAntiLambda = new TH1F("fh1MCProdRadiusAntiLambda","MC gen. MC ALa prod radius",200,0.,200.);
1700 // Pt and inv mass distributions
1702 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
1703 fh1MCPtK0s = new TH1F("fh1MCPtK0s", "MC gen. K^{0}_{s} in eta range;#it{p}_{T} (GeV/#it{c})",200,0.,20.);
1704 fh1MCPtLambda = new TH1F("fh1MCPtLambda", "MC gen. #Lambda in rap range;#it{p}_{T} (GeV/#it{c})",200,0.,20.);
1705 fh1MCPtAntiLambda = new TH1F("fh1MCPtAntiLambda", "MC gen. #AntiLambda in rap range;#it{p}_{T} (GeV/#it{c})",200,0.,20.);
1706 fh1MCXiPt = new TH1F("fh1MCXiPt", "MC gen. #Xi^{-/o};#it{p}_{T} (GeV/#it{c})",200,0.,20.);
1707 fh1MCXibarPt = new TH1F("fh1MCXibarPt", "MC gen. #bar{#Xi}^{+/o};#it{p}_{T} (GeV/#it{c})",200,0.,20.);
1708 fh2MCEtaVsPtK0s = new TH2F("fh2MCEtaVsPtK0s","MC gen. K^{0}_{s} #eta; #it{p}_{T}",200,0.,20.,200,-1.,1.);
1709 fh2MCEtaVsPtLa = new TH2F("fh2MCEtaVsPtLa","MC gen. #Lambda #eta; #it{p}_{T}",200,0.,20.,200,-1.,1.);
1710 fh2MCEtaVsPtALa = new TH2F("fh2MCEtaVsPtALa","MC gen. #bar{#Lambda} #eta; #it{p}_{T}",200,0.,20.,200,-1.,1.);
1713 //fh1MCRapK0s = new TH1F("fh1MCRapK0s", "MC gen. K0s;rap with cut",200,-10,10);
1714 //fh1MCRapLambda = new TH1F("fh1MCRapLambda", "MC gen. #Lambda;rap",200,-10,10);
1715 //fh1MCRapAntiLambda = new TH1F("fh1MCRapAntiLambda", "MC gen. #bar{#Lambda};rap",200,-10,10);
1716 fh1MCEtaAllK0s = new TH1F("fh1MCEtaAllK0s", "MC gen. K0s;#eta",200,-1.,1.);
1717 fh1MCEtaK0s = new TH1F("fh1MCEtaK0s", "MC gen. K0s;#eta with cut",200,-1.,1.);
1718 fh1MCEtaLambda = new TH1F("fh1MCEtaLambda", "MC gen. #Lambda;#eta",200,-1.,1.);
1719 fh1MCEtaAntiLambda = new TH1F("fh1MCEtaAntiLambda", "MC gen. #bar{#Lambda};#eta",200,-1.,1.);
1721 fV0QAK0->DefineHistos();
1722 fFFHistosRecCuts->DefineHistos();
1723 fFFHistosRecCutsK0Evt->DefineHistos();
1724 /* fFFHistosIMK0AllEvt->DefineHistos();
1725 fFFHistosIMK0Jet->DefineHistos();
1726 fFFHistosIMK0Cone->DefineHistos();
1727 fFFHistosIMLaAllEvt->DefineHistos();
1728 fFFHistosIMLaJet->DefineHistos();
1729 fFFHistosIMLaCone->DefineHistos();
1730 fFFHistosIMALaAllEvt->DefineHistos();
1731 fFFHistosIMALaJet->DefineHistos();
1732 fFFHistosIMALaCone->DefineHistos();
1735 const Int_t saveLevel = 5;
1738 fCommonHistList->Add(fh1EvtAllCent);
1739 fCommonHistList->Add(fh1Evt);
1740 fCommonHistList->Add(fh1EvtSelection);
1741 fCommonHistList->Add(fh1EvtCent);
1742 fCommonHistList->Add(fh1VertexNContributors);
1743 fCommonHistList->Add(fh1VertexZ);
1744 fCommonHistList->Add(fh1Xsec);
1745 fCommonHistList->Add(fh1Trials);
1746 fCommonHistList->Add(fh1PtHard);
1747 fCommonHistList->Add(fh1PtHardTrials);
1748 fCommonHistList->Add(fh1nRecJetsCuts);
1749 fCommonHistList->Add(fh1EvtMult);
1750 fCommonHistList->Add(fh1K0Mult);
1751 fCommonHistList->Add(fh1dPhiJetK0);
1752 fCommonHistList->Add(fh1LaMult);
1753 fCommonHistList->Add(fh1dPhiJetLa);
1754 fCommonHistList->Add(fh1ALaMult);
1755 fCommonHistList->Add(fh1dPhiJetALa);
1756 fCommonHistList->Add(fh1JetEta);
1757 fCommonHistList->Add(fh1JetPhi);
1758 fCommonHistList->Add(fh2JetEtaPhi);
1759 //fCommonHistList->Add(fh1V0JetPt);
1760 fCommonHistList->Add(fh1IMK0Cone);
1761 fCommonHistList->Add(fh1IMLaCone);
1762 fCommonHistList->Add(fh1IMALaCone);
1763 fCommonHistList->Add(fh2FFJetTrackEta);
1764 // fCommonHistList->Add(fh1trackPosNCls);
1765 //fCommonHistList->Add(fh1trackNegNCls);
1766 fCommonHistList->Add(fh1trackPosEta);
1767 fCommonHistList->Add(fh1trackNegEta);
1768 fCommonHistList->Add(fh1V0Eta);
1769 // fCommonHistList->Add(fh1V0totMom);
1770 fCommonHistList->Add(fh1CosPointAngle);
1771 fCommonHistList->Add(fh1DecayLengthV0);
1772 fCommonHistList->Add(fh2ProperLifetimeK0sVsPtBeforeCut);
1773 fCommonHistList->Add(fh2ProperLifetimeK0sVsPtAfterCut);
1774 fCommonHistList->Add(fh1V0Radius);
1775 fCommonHistList->Add(fh1DcaV0Daughters);
1776 fCommonHistList->Add(fh1DcaPosToPrimVertex);
1777 fCommonHistList->Add(fh1DcaNegToPrimVertex);
1778 fCommonHistList->Add(fh2ArmenterosBeforeCuts);
1779 fCommonHistList->Add(fh2ArmenterosAfterCuts);
1780 fCommonHistList->Add(fh2BBLaPos);
1781 fCommonHistList->Add(fh2BBLaNeg);
1782 fCommonHistList->Add(fh1PosDaughterCharge);
1783 fCommonHistList->Add(fh1NegDaughterCharge);
1784 fCommonHistList->Add(fh1PtMCK0s);
1785 fCommonHistList->Add(fh1PtMCLa);
1786 fCommonHistList->Add(fh1PtMCALa);
1787 fCommonHistList->Add(fh1EtaK0s);
1788 fCommonHistList->Add(fh1EtaLa);
1789 fCommonHistList->Add(fh1EtaALa);
1790 fCommonHistList->Add(fh1RC);
1791 fCommonHistList->Add(fh1RCBiasK0);
1792 fCommonHistList->Add(fh1RCBiasLa);
1793 fCommonHistList->Add(fh1RCBiasALa);
1794 fCommonHistList->Add(fh1MCC);
1795 fCommonHistList->Add(fh1OC);
1796 fCommonHistList->Add(fh1NJ);
1797 fCommonHistList->Add(fhnInvMassEtaTrackPtK0s);
1798 fCommonHistList->Add(fhnInvMassEtaTrackPtLa);
1799 fCommonHistList->Add(fhnInvMassEtaTrackPtALa);
1800 fCommonHistList->Add(fh1TrackMultCone);
1801 fCommonHistList->Add(fh2TrackMultCone);
1802 fCommonHistList->Add(fhnNJK0);
1803 fCommonHistList->Add(fhnNJLa);
1804 fCommonHistList->Add(fhnNJALa);
1805 //fCommonHistList->Add(fh2MCgenK0Cone);
1806 //fCommonHistList->Add(fh2MCgenLaCone);
1807 //fCommonHistList->Add(fh2MCgenALaCone);
1808 //fCommonHistList->Add(fh2MCEtagenK0Cone);
1809 //fCommonHistList->Add(fh2MCEtagenLaCone);
1810 //fCommonHistList->Add(fh2MCEtagenALaCone);
1811 fCommonHistList->Add(fh2CorrHijingLaProton);
1812 fCommonHistList->Add(fh2CorrInjectLaProton);
1813 fCommonHistList->Add(fh2CorrHijingALaAProton);
1814 fCommonHistList->Add(fh2CorrInjectALaAProton);
1815 fCommonHistList->Add(fh2MCEtaVsPtHijingLa);
1816 fCommonHistList->Add(fh2MCEtaVsPtInjectLa);
1817 fCommonHistList->Add(fh2MCEtaVsPtHijingALa);
1818 fCommonHistList->Add(fh2MCEtaVsPtInjectALa);
1819 fCommonHistList->Add(fh1IMK0ConeSmear);
1820 fCommonHistList->Add(fh1IMLaConeSmear);
1821 fCommonHistList->Add(fh1IMALaConeSmear);
1822 fCommonHistList->Add(fhnrecMCHijingLaIncl);
1823 fCommonHistList->Add(fhnrecMCHijingLaCone);
1824 fCommonHistList->Add(fhnrecMCHijingALaIncl);
1825 fCommonHistList->Add(fhnrecMCHijingALaCone);
1826 fCommonHistList->Add(fhnrecMCInjectLaIncl);
1827 fCommonHistList->Add(fhnrecMCInjectLaCone);
1828 fCommonHistList->Add(fhnrecMCInjectALaIncl);
1829 fCommonHistList->Add(fhnrecMCInjectALaCone);
1830 fCommonHistList->Add(fhnMCrecK0Cone);
1831 fCommonHistList->Add(fhnMCrecLaCone);
1832 fCommonHistList->Add(fhnMCrecALaCone);
1833 fCommonHistList->Add(fhnMCrecK0ConeSmear);
1834 fCommonHistList->Add(fhnMCrecLaConeSmear);
1835 fCommonHistList->Add(fhnMCrecALaConeSmear);
1836 fCommonHistList->Add(fhnK0sSecContinCone);
1837 fCommonHistList->Add(fhnLaSecContinCone);
1838 fCommonHistList->Add(fhnALaSecContinCone);
1839 fCommonHistList->Add(fhnK0sIncl);
1840 fCommonHistList->Add(fhnK0sCone);
1841 fCommonHistList->Add(fhnLaIncl);
1842 fCommonHistList->Add(fhnLaCone);
1843 fCommonHistList->Add(fhnALaIncl);
1844 fCommonHistList->Add(fhnALaCone);
1845 fCommonHistList->Add(fhnK0sPC);
1846 fCommonHistList->Add(fhnLaPC);
1847 fCommonHistList->Add(fhnALaPC);
1848 fCommonHistList->Add(fhnK0sMCC);
1849 fCommonHistList->Add(fhnLaMCC);
1850 fCommonHistList->Add(fhnALaMCC);
1851 fCommonHistList->Add(fhnK0sRC);
1852 fCommonHistList->Add(fhnLaRC);
1853 fCommonHistList->Add(fhnALaRC);
1854 fCommonHistList->Add(fhnK0sRCBias);
1855 fCommonHistList->Add(fhnLaRCBias);
1856 fCommonHistList->Add(fhnALaRCBias);
1857 fCommonHistList->Add(fhnK0sOC);
1858 fCommonHistList->Add(fhnLaOC);
1859 fCommonHistList->Add(fhnALaOC);
1860 fCommonHistList->Add(fh1AreaExcluded);
1861 fCommonHistList->Add(fh1MedianEta);
1862 fCommonHistList->Add(fh1JetPtMedian);
1863 fCommonHistList->Add(fh1MCMultiplicityPrimary);
1864 fCommonHistList->Add(fh1MCMultiplicityTracks);
1865 fCommonHistList->Add(fhnFeedDownLa);
1866 fCommonHistList->Add(fhnFeedDownALa);
1867 fCommonHistList->Add(fhnFeedDownLaCone);
1868 fCommonHistList->Add(fhnFeedDownALaCone);
1869 fCommonHistList->Add(fh1MCProdRadiusK0s);
1870 fCommonHistList->Add(fh1MCProdRadiusLambda);
1871 fCommonHistList->Add(fh1MCProdRadiusAntiLambda);
1872 fCommonHistList->Add(fh1MCPtV0s);
1873 fCommonHistList->Add(fh1MCPtK0s);
1874 fCommonHistList->Add(fh1MCPtLambda);
1875 fCommonHistList->Add(fh1MCPtAntiLambda);
1876 fCommonHistList->Add(fh1MCXiPt);
1877 fCommonHistList->Add(fh1MCXibarPt);
1878 fCommonHistList->Add(fh2MCEtaVsPtK0s);
1879 fCommonHistList->Add(fh2MCEtaVsPtLa);
1880 fCommonHistList->Add(fh2MCEtaVsPtALa);
1881 //fCommonHistList->Add(fh1MCRapK0s);
1882 //fCommonHistList->Add(fh1MCRapLambda);
1883 //fCommonHistList->Add(fh1MCRapAntiLambda);
1884 fCommonHistList->Add(fh1MCEtaAllK0s);
1885 fCommonHistList->Add(fh1MCEtaK0s);
1886 fCommonHistList->Add(fh1MCEtaLambda);
1887 fCommonHistList->Add(fh1MCEtaAntiLambda);
1891 fV0QAK0->AddToOutput(fCommonHistList);
1892 fFFHistosRecCuts->AddToOutput(fCommonHistList);
1893 fFFHistosRecCutsK0Evt->AddToOutput(fCommonHistList);
1894 // fFFHistosIMK0AllEvt->AddToOutput(fCommonHistList);
1895 // fFFHistosIMK0Jet->AddToOutput(fCommonHistList);
1896 // fFFHistosIMK0Cone->AddToOutput(fCommonHistList);
1897 // fFFHistosIMLaAllEvt->AddToOutput(fCommonHistList);
1898 // fFFHistosIMLaJet->AddToOutput(fCommonHistList);
1899 // fFFHistosIMLaCone->AddToOutput(fCommonHistList);
1900 // fFFHistosIMALaAllEvt->AddToOutput(fCommonHistList);
1901 // fFFHistosIMALaJet->AddToOutput(fCommonHistList);
1902 // fFFHistosIMALaCone->AddToOutput(fCommonHistList);
1907 // =========== Switch on Sumw2 for all histos ===========
1908 for (Int_t i=0; i<fCommonHistList->GetEntries(); ++i){
1910 TH1 *h1 = dynamic_cast<TH1*>(fCommonHistList->At(i));
1912 if (h1) h1->Sumw2();//The error per bin will be computed as sqrt(sum of squares of weight) for each bin
1914 THnSparse *hnSparse = dynamic_cast<THnSparse*>(fCommonHistList->At(i));
1915 if(hnSparse) hnSparse->Sumw2();
1919 TH1::AddDirectory(oldStatus);
1920 PostData(1, fCommonHistList);
1923 //_______________________________________________
1924 void AliAnalysisTaskJetChem::UserExec(Option_t *)
1927 // Called for each event
1929 if(fDebug > 1) Printf("AliAnalysisTaskJetChem::UserExec()");
1931 if(fDebug > 1) Printf("Analysis event #%5d", (Int_t) fEntry);
1933 // Trigger selection
1934 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
1935 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
1938 //for AliPIDResponse:
1939 //AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
1940 //AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
1941 fPIDResponse = inputHandler->GetPIDResponse();
1943 if (!fPIDResponse){if(fDebug > 1) Printf("AliAnalysisTaskJetChem::UserExec(): fPIDResponse does not exist!"); return;}
1945 //std::cout<<"inputHandler->IsEventSelected(): "<<inputHandler->IsEventSelected()<<std::endl;
1946 //std::cout<<"fEvtSelectionMask: "<<fEvtSelectionMask<<std::endl;
1948 if(!(inputHandler->IsEventSelected() & fEvtSelectionMask)){
1949 //std::cout<<"########event rejected!!############"<<std::endl;
1950 fh1EvtSelection->Fill(1.);
1951 if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
1952 PostData(1, fCommonHistList);
1956 fESD = dynamic_cast<AliESDEvent*>(InputEvent());//casting of pointers for inherited class, only for ESDs
1958 if(fDebug>3) Printf("%s:%d ESDEvent not found in the input", (char*)__FILE__,__LINE__);
1961 fMCEvent = MCEvent();
1963 if(fDebug>3) Printf("%s:%d MCEvent not found in the input", (char*)__FILE__,__LINE__);
1966 // get AOD event from input/output
1967 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
1968 if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
1969 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
1970 if(fUseAODInputJets) fAODJets = fAOD;
1971 if (fDebug > 1) Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
1974 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
1975 if( handler && handler->InheritsFrom("AliAODHandler") ) {
1976 fAOD = ((AliAODHandler*)handler)->GetAOD();
1978 if (fDebug > 1) Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
1982 if(!fAODJets && !fUseAODInputJets){ // case we have AOD in input & output and want jets from output
1983 TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
1984 if( outHandler && outHandler->InheritsFrom("AliAODHandler") ){
1985 fAODJets = ((AliAODHandler*)outHandler)->GetAOD();
1986 if (fDebug > 1) Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
1990 if(fNonStdFile.Length()!=0){
1991 // case we have an AOD extension - fetch the jets from the extended output
1993 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
1994 fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
1996 if(fDebug>1)Printf("AODExtension not found for %s",fNonStdFile.Data());
2001 Printf("%s:%d AODEvent not found", (char*)__FILE__,__LINE__);
2005 Printf("%s:%d AODEvent with jet branch not found", (char*)__FILE__,__LINE__);
2009 //primary vertex position:
2010 AliAODVertex *myPrimaryVertex = NULL;
2011 myPrimaryVertex = (AliAODVertex*)fAOD->GetPrimaryVertex();
2012 if (!myPrimaryVertex) return;
2013 fh1Evt->Fill(1.);//fill in every event that was accessed with InputHandler
2015 // event selection *****************************************
2017 // *** vertex cut ***
2018 AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
2019 Int_t nTracksPrim = primVtx->GetNContributors();
2020 fh1VertexNContributors->Fill(nTracksPrim);
2022 if (fDebug > 1) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
2024 if(nTracksPrim <= 2){
2025 if (fDebug > 1) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__);
2026 fh1EvtSelection->Fill(3.);
2027 PostData(1, fCommonHistList);
2031 fh1VertexZ->Fill(primVtx->GetZ());
2033 if(TMath::Abs(primVtx->GetZ())>fMaxVertexZ){
2034 if (fDebug > 1) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
2035 fh1EvtSelection->Fill(4.);
2036 PostData(1, fCommonHistList);
2040 // accepts only events that have same "primary" and SPD vertex, special issue of LHC11h PbPb data
2042 //fAOD: pointer to global primary vertex
2044 const AliAODVertex* spdVtx = fAOD->GetPrimaryVertexSPD();
2046 if (TMath::Abs(spdVtx->GetZ() - primVtx->GetZ())>fDeltaVertexZ) { if (fDebug > 1) Printf("deltaZVertex: event REJECTED..."); return;}
2049 //check for vertex radius to be smaller than 1 cm, (that was first applied by Vit Kucera in his analysis)
2051 Double_t vtxX = primVtx->GetX();
2052 Double_t vtxY = primVtx->GetY();
2054 if(TMath::Sqrt(vtxX*vtxX + vtxY*vtxY)>=1){
2055 if (fDebug > 1) Printf("%s:%d primary vertex r = %f: event REJECTED...",(char*)__FILE__,__LINE__,TMath::Sqrt(vtxX*vtxX + vtxY*vtxY));
2060 TString primVtxName(primVtx->GetName());
2062 if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
2063 if (fDebug > 1) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
2064 fh1EvtSelection->Fill(5.);
2065 PostData(1, fCommonHistList);
2069 Bool_t selectedHelper = AliAnalysisHelperJetTasks::Selected();
2070 if(!selectedHelper){
2071 fh1EvtSelection->Fill(6.);
2072 PostData(1, fCommonHistList);
2076 // event selection *****************************************
2078 Double_t centPercent = -1;
2082 if(handler && handler->InheritsFrom("AliAODInputHandler")){
2084 centPercent = fAOD->GetHeader()->GetCentrality();
2086 //std::cout<<"centPercent: "<<centPercent<<std::endl;
2088 fh1EvtAllCent->Fill(centPercent);
2090 if(centPercent>10) cl = 2; //standard PWG-JE binning
2091 if(centPercent>30) cl = 3;
2092 if(centPercent>50) cl = 4;
2096 if(centPercent < 0) cl = -1;
2097 if(centPercent >= 0) cl = 1;
2098 if(centPercent > 10) cl = 2; //standard PWG-JE binning
2099 if(centPercent > 30) cl = 3;
2100 if(centPercent > 50) cl = 4;
2101 if(centPercent > 80) cl = 5; //takes centralities higher than my upper edge of 80%, not to be used
2106 cl = AliAnalysisHelperJetTasks::EventClass();
2108 if(fESD) centPercent = fESD->GetCentrality()->GetCentralityPercentile("V0M"); //ESD JetServices Task has the centrality binning 0-10,10-30,30-50,50-80
2109 fh1EvtAllCent->Fill(centPercent);
2112 if(cl!=fEventClass){ // event not in selected event class, reject event#########################################
2114 if (fDebug > 1) Printf("%s:%d event not in selected event class: event REJECTED ...",(char*)__FILE__,__LINE__);
2115 fh1EvtSelection->Fill(2.);
2116 PostData(1, fCommonHistList);
2119 }//end if fEventClass > 0
2122 if (fDebug > 1) Printf("%s:%d event ACCEPTED ...",(char*)__FILE__,__LINE__);
2125 //Printf("Analysis event #%5d", (Int_t) fEntry);
2127 fh1EvtSelection->Fill(0.);
2128 fh1EvtCent->Fill(centPercent);
2130 //___ get MC information __________________________________________________________________
2133 Double_t ptHard = 0.; //parton energy bins -> energy of particle
2134 Double_t nTrials = 1; // trials for MC trigger weight for real data
2137 AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
2138 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);//check usage of Pythia (pp) or Hijing (PbPb)
2139 AliGenHijingEventHeader* hijingGenHeader = 0x0;
2141 if(pythiaGenHeader){
2142 if(fDebug>3) Printf("%s:%d pythiaGenHeader found", (char*)__FILE__,__LINE__);
2143 nTrials = pythiaGenHeader->Trials();
2144 ptHard = pythiaGenHeader->GetPtHard();
2146 fh1PtHard->Fill(ptHard);
2147 fh1PtHardTrials->Fill(ptHard,nTrials);
2150 } else { // no pythia, hijing?
2152 if(fDebug>3) Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
2154 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
2155 if(!hijingGenHeader){
2156 if(fDebug>3) Printf("%s:%d no pythiaGenHeader or hjingGenHeader found", (char*)__FILE__,__LINE__);
2158 if(fDebug>3) Printf("%s:%d hijingGenHeader found", (char*)__FILE__,__LINE__);
2162 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
2165 //____ fetch jets _______________________________________________________________
2167 Int_t nJCuts = GetListOfJets(fJetsRecCuts, kJetsRecAcceptance);//fetch list with jets that survived all jet cuts: fJetsRecCuts
2169 Int_t nRecJetsCuts = 0; //number of reconstructed jets after jet cuts
2170 if(nJCuts>=0) nRecJetsCuts = fJetsRecCuts->GetEntries();
2171 if(fDebug>2)Printf("%s:%d Selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2172 if(nRecJetsCuts != nJCuts) Printf("%s:%d Mismatch selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2173 fh1nRecJetsCuts->Fill(nRecJetsCuts);
2176 //____ fetch background clusters ___________________________________________________
2177 if(fBranchRecBckgClusters.Length() != 0){
2179 Int_t nBJ = GetListOfBckgJets(fBckgJetsRec, kJetsRec);
2180 Int_t nRecBckgJets = 0;
2181 if(nBJ>=0) nRecBckgJets = fBckgJetsRec->GetEntries();
2182 if(fDebug>2)Printf("%s:%d Selected Rec background jets: %d %d",(char*)__FILE__,__LINE__,nBJ,nRecBckgJets);
2183 if(nBJ != nRecBckgJets) Printf("%s:%d Mismatch Selected Rec background jets: %d %d",(char*)__FILE__,__LINE__,nBJ,nRecBckgJets);
2187 //____ fetch reconstructed particles __________________________________________________________
2189 Int_t nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODCuts);//all tracks of event
2190 if(fDebug>2)Printf("%s:%d selected reconstructed tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,fTracksRecCuts->GetEntries());
2191 if(fTracksRecCuts->GetEntries() != nTCuts)
2192 Printf("%s:%d Mismatch selected reconstructed tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,fTracksRecCuts->GetEntries());
2193 fh1EvtMult->Fill(fTracksRecCuts->GetEntries());
2195 Int_t nK0s = GetListOfV0s(fListK0s,fK0Type,kK0,myPrimaryVertex,fAOD);//all V0s in event with K0s assumption
2197 if(fDebug>5){std::cout<<"fK0Type: "<<fK0Type<<" kK0: "<<kK0<<" myPrimaryVertex: "<<myPrimaryVertex<<" fAOD: "<<fAOD<<std::endl;}
2199 //std::cout<< "nK0s: "<<nK0s<<std::endl;
2201 if(fDebug>2)Printf("%s:%d Selected V0s after cuts: %d %d",(char*)__FILE__,__LINE__,nK0s,fListK0s->GetEntries());
2202 if(nK0s != fListK0s->GetEntries()) Printf("%s:%d Mismatch selected K0s: %d %d",(char*)__FILE__,__LINE__,nK0s,fListK0s->GetEntries());
2203 fh1K0Mult->Fill(fListK0s->GetEntries());
2206 Int_t nLa = GetListOfV0s(fListLa,fLaType,kLambda,myPrimaryVertex,fAOD);//all V0s in event with Lambda particle assumption
2207 if(fDebug>2)Printf("%s:%d Selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nLa,fListLa->GetEntries());
2208 if(nLa != fListLa->GetEntries()) Printf("%s:%d Mismatch selected La: %d %d",(char*)__FILE__,__LINE__,nLa,fListLa->GetEntries());
2209 fh1LaMult->Fill(fListLa->GetEntries());
2211 Int_t nALa = GetListOfV0s(fListALa,fALaType,kAntiLambda,myPrimaryVertex,fAOD);//all V0s in event with Antilambda particle assumption
2212 if(fDebug>2)Printf("%s:%d Selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nALa,fListALa->GetEntries());
2213 if(nALa != fListALa->GetEntries()) Printf("%s:%d Mismatch selected ALa: %d %d",(char*)__FILE__,__LINE__,nALa,fListALa->GetEntries());
2214 fh1ALaMult->Fill(fListALa->GetEntries());
2218 //fetch MC gen particles_______________________________________________________
2220 if(fAnalysisMC){ // here
2222 //fill feeddown histo for associated particles
2224 // Access MC generated particles, fill TLists and histograms :
2226 Int_t nMCgenK0s = GetListOfMCParticles(fListMCgenK0s,kK0,fAOD); //fill TList with MC generated primary true K0s (list to fill, particletype, mc aod event)
2227 if(nMCgenK0s != fListMCgenK0s->GetEntries()) Printf("%s:%d Mismatch selected MCgenK0s: %d %d",(char*)__FILE__,__LINE__,nMCgenK0s,fListMCgenK0s->GetEntries());
2230 for(Int_t it=0; it<fListMCgenK0s->GetSize(); ++it){ // loop MC generated K0s, filling histograms
2232 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenK0s->At(it));
2237 //Double_t fRapCurrentPart = MyRapidity(mcp0->E(),mcp0->Pz());
2238 Double_t fEtaCurrentPart = mcp0->Eta();
2239 Double_t fPtCurrentPart = mcp0->Pt();
2241 fh1MCEtaK0s->Fill(fEtaCurrentPart);
2242 //fh1MCRapK0s->Fill(fRapCurrentPart);
2243 fh1MCPtK0s->Fill(fPtCurrentPart);
2245 fh2MCEtaVsPtK0s->Fill(fPtCurrentPart,fEtaCurrentPart); //eta cut, physical primary selection and decay mode considered
2249 Int_t nMCgenLa = GetListOfMCParticles(fListMCgenLa,kLambda,fAOD); //fill TList with MC generated primary true Lambdas (list to fill, particletype, mc aod event)
2250 if(nMCgenLa != fListMCgenLa->GetEntries()) Printf("%s:%d Mismatch selected MCgenLa: %d %d",(char*)__FILE__,__LINE__,nMCgenLa,fListMCgenLa->GetEntries());
2252 TList *mclist = fAOD->GetList();
2253 TClonesArray *stackMC = 0x0;
2254 stackMC = (TClonesArray*)mclist->FindObject(AliAODMCParticle::StdBranchName());
2256 Printf("ERROR: AliAnalysisTaskJetChem.cxx: loop over MC gen. particles: stackMC not available!");
2259 AliAODMCHeader *mcHdr=(AliAODMCHeader*)mclist->FindObject(AliAODMCHeader::StdBranchName());
2260 if(!mcHdr)Printf("ERROR: AliAnalysisTaskJetChem.cxx: loop over MC gen. particles: mcHdr not available!");
2262 for(Int_t it=0; it<fListMCgenLa->GetSize(); ++it){ // loop MC generated La, filling histograms
2264 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenLa->At(it));
2269 //Double_t fRapCurrentPart = MyRapidity(mcp0->E(),mcp0->Pz());
2270 Double_t fEtaCurrentPart = mcp0->Eta();
2271 Double_t fPtCurrentPart = mcp0->Pt();
2272 TString generatorName;
2274 fh1MCEtaLambda->Fill(fEtaCurrentPart);
2275 //fh1MCRapLambda->Fill(fRapCurrentPart);
2276 fh1MCPtLambda->Fill(fPtCurrentPart);
2277 fh2MCEtaVsPtLa->Fill(fPtCurrentPart,fEtaCurrentPart); //eta cut, physical primary selection and decay mode considered
2279 Int_t mcp0label = mcp0->GetLabel();
2280 Bool_t istrackInject = IsTrackInjected(mcp0label, mcHdr, stackMC, generatorName);
2282 //std::cout<<"generatorName: "<<generatorName<<std::endl;
2285 if(generatorName == "Hijing"){
2286 fh2MCEtaVsPtHijingLa->Fill(fPtCurrentPart,fEtaCurrentPart);
2289 if(istrackInject == kTRUE){
2290 fh2MCEtaVsPtHijingLa->Fill(fPtCurrentPart,fEtaCurrentPart);
2296 Int_t nMCgenALa = GetListOfMCParticles(fListMCgenALa,kAntiLambda,fAOD); //fill TList with MC generated primary true Antilambdas (list to fill, particletype, mc aod event)
2297 if(nMCgenALa != fListMCgenALa->GetEntries()) Printf("%s:%d Mismatch selected MCgenALa: %d %d",(char*)__FILE__,__LINE__,nMCgenALa,fListMCgenALa->GetEntries());
2300 for(Int_t it=0; it<fListMCgenALa->GetSize(); ++it){ // loop MC generated ALa, filling histograms
2302 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenALa->At(it));
2305 //MC gen Antilambdas
2307 //Double_t fRapCurrentPart = MyRapidity(mcp0->E(),mcp0->Pz());
2308 Double_t fEtaCurrentPart = mcp0->Eta();
2309 Double_t fPtCurrentPart = mcp0->Pt();
2311 fh1MCEtaAntiLambda->Fill(fEtaCurrentPart);
2312 //fh1MCRapAntiLambda->Fill(fRapCurrentPart);
2313 fh1MCPtAntiLambda->Fill(fPtCurrentPart);
2314 fh2MCEtaVsPtALa->Fill(fPtCurrentPart,fEtaCurrentPart); //eta cut, physical primary selection and decay mode considered
2320 //loop over MC feeddown candidates in TList
2325 } //end MCAnalysis part for gen particles
2328 // ___ V0 QA + K0s + La + ALa pt spectra all events _______________________________________________
2330 Double_t lPrimaryVtxPosition[3];
2331 Double_t lV0Position[3];
2332 lPrimaryVtxPosition[0] = primVtx->GetX();
2333 lPrimaryVtxPosition[1] = primVtx->GetY();
2334 lPrimaryVtxPosition[2] = primVtx->GetZ();
2335 Double_t dRadiusExcludeCone = 2*GetFFRadius(); //2 times jet radius
2336 //------------------------------------------
2337 for(Int_t it=0; it<fListK0s->GetSize(); ++it){
2339 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0s->At(it));
2342 // VO's main characteristics to check the reconstruction cuts
2344 //Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2347 Double_t fV0Radius = -999;
2348 Double_t fDcaV0Daughters = v0->DcaV0Daughters();
2349 Double_t fDcaPosToPrimVertex = v0->DcaPosToPrimVertex();//IP of positive charged daughter
2350 Double_t fDcaNegToPrimVertex = v0->DcaNegToPrimVertex();//IP of negative charged daughter
2351 Int_t negDaughterpdg = 0;
2352 Int_t posDaughterpdg = 0;
2353 Int_t motherType = 0;
2356 Bool_t fPhysicalPrimary = kFALSE;//don't use IsPhysicalPrimary() anymore for MC analysis, use instead 2D distance from primary to secondary vertex
2357 Int_t MCv0PdgCode = 0;
2359 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0));
2360 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1));
2362 Double_t PosEta = trackPos->AliAODTrack::Eta();//daughter track charge is sometimes wrong here, account for that!!!
2363 Double_t NegEta = trackNeg->AliAODTrack::Eta();
2365 //Double_t trackPosNcls = trackPos->GetTPCNcls();//Get number of clusters for positive charged tracks
2366 //Double_t trackNegNcls = trackNeg->GetTPCNcls();//Get number of clusters for negative charged tracks
2368 CalculateInvMass(v0, kK0, invMK0s, trackPt); //function to calculate invMass with TLorentzVector class
2371 //OUTSIDE CONES:########
2373 Double_t fEta = v0->PseudoRapV0();
2374 Bool_t bIsInCone = kFALSE;//init boolean, is not in any cone (OC)
2377 for(Int_t ij=0; ij<nRecJetsCuts; ++ij){ // loop over all jets in event
2379 AliAODJet* jet = (AliAODJet*) (fJetsRecCuts->At(ij));
2380 TList* jettracklist = new TList();
2381 Double_t sumPt = 0.;
2382 Bool_t isBadJet = kFALSE;
2384 if(GetFFRadius()<=0){
2385 GetJetTracksTrackrefs(jettracklist, jet, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);// list of jet tracks from trackrefs
2387 GetJetTracksPointing(fTracksRecCuts, jettracklist, jet, GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet); // fill list of tracks in cone around jet axis with cone Radius (= 0.4 standard)
2391 //leading track pt bias on jets inside this small jet loop
2392 if(isBadJet) continue;//all bad jets are rejected
2397 //if jet is selected, then check whether V0 is part of the jet cone:
2398 if(IsParticleInCone(jet, v0, dRadiusExcludeCone) == kTRUE) {bIsInCone = kTRUE;}
2400 delete jettracklist;
2404 if(bIsInCone==kFALSE){//K0s is not part of any selected jet in event
2405 Double_t vK0sOC[3] = {invMK0s,trackPt,fEta};
2406 fhnK0sOC->Fill(vK0sOC);
2409 //end of outside cone K0s
2411 Double_t fV0cosPointAngle = v0->CosPointingAngle(lPrimaryVtxPosition);
2413 lV0Position[0]= v0->DecayVertexV0X();
2414 lV0Position[1]= v0->DecayVertexV0Y();
2415 lV0Position[2]= v0->DecayVertexV0Z();
2417 Double_t fV0DecayLength = v0->DecayLengthV0(lPrimaryVtxPosition);
2418 fV0Radius = TMath::Sqrt(lV0Position[0]*lV0Position[0]+lV0Position[1]*lV0Position[1]);
2421 fV0QAK0->FillTrackQA(v0->Eta(), TVector2::Phi_0_2pi(v0->Phi()), v0->Pt());
2422 //fFFHistosIMK0AllEvt->FillFF(trackPt, invMK0s, jetPt, incrementJetPt);
2423 //fh1trackPosNCls->Fill(trackPosNcls);
2424 //fh1trackNegNCls->Fill(trackNegNcls);
2425 fh1EtaK0s->Fill(fEta);
2427 Double_t vK0sIncl[3] = {invMK0s,trackPt,fEta}; //fill all K0s in event into THnSparse of 3 dimensions
2428 fhnK0sIncl->Fill(vK0sIncl);
2432 TString generatorName;
2434 TList *listmc = fAOD->GetList();
2435 Bool_t mclabelcheck = MCLabelCheck(v0, kK0, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
2436 //if(fPhysicalPrimary == kFALSE)continue;
2437 //std::cout<<"mclabelcheck: "<<mclabelcheck<<std::endl;
2438 //std::cout<<"IsPhysicalPrimary: "<<fPhysicalPrimary<<std::endl;
2440 if(mclabelcheck == kFALSE)continue;
2442 Double_t vInvMassEtaTrackPtK0s[3] = {fEta,invMK0s,trackPt};
2443 fhnInvMassEtaTrackPtK0s->Fill(vInvMassEtaTrackPtK0s);//includes also feeddown particles, mainly phi particles whose decay products are considered here as primary
2446 fh1PtMCK0s->Fill(MCPt);
2450 fh1V0Eta->Fill(fEta);
2451 //fh1V0totMom->Fill(fV0TotalMomentum);
2452 fh1CosPointAngle->Fill(fV0cosPointAngle);
2453 fh1DecayLengthV0->Fill(fV0DecayLength);
2454 fh1V0Radius->Fill(fV0Radius);
2455 fh1DcaV0Daughters->Fill(fDcaV0Daughters);
2456 fh1DcaPosToPrimVertex->Fill(fDcaPosToPrimVertex);
2457 fh1DcaNegToPrimVertex->Fill(fDcaNegToPrimVertex);
2458 fh1trackPosEta->Fill(PosEta);
2459 fh1trackNegEta->Fill(NegEta);
2463 // __La pt spectra all events _______________________________________________
2466 for(Int_t it=0; it<fListLa->GetSize(); ++it){ // loop all La
2468 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListLa->At(it));
2471 // VO's main characteristics to check the reconstruction cuts
2472 // Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2475 Double_t fV0Radius = -999;
2476 Double_t fDcaV0Daughters = v0->DcaV0Daughters();
2477 Double_t fDcaPosToPrimVertex = v0->DcaPosToPrimVertex();//IP of positive charged daughter
2478 Double_t fDcaNegToPrimVertex = v0->DcaNegToPrimVertex();//IP of negative charged daughter
2479 Int_t negDaughterpdg = 0;
2480 Int_t posDaughterpdg = 0;
2481 Int_t motherType = 0;
2484 Bool_t fPhysicalPrimary = kFALSE;
2485 Int_t MCv0PdgCode = 0;
2486 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0));
2487 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1));
2489 //Double_t trackPosNcls = trackPos->GetTPCNcls();//Get number of clusters for positive charged tracks
2490 //Double_t trackNegNcls = trackNeg->GetTPCNcls();//Get number of clusters for negative charged tracks
2492 Double_t PosEta = trackPos->AliAODTrack::Eta();//daughter track charge is sometimes wrong here, account for that!!!
2493 Double_t NegEta = trackNeg->AliAODTrack::Eta();
2495 Double_t fEta = v0->PseudoRapV0();
2496 Bool_t bIsInCone = kFALSE;//init boolean, is not in any cone (OC)
2498 CalculateInvMass(v0, kLambda, invMLa, trackPt);//function to calculate invMass with TLorentzVector class
2501 for(Int_t ij=0; ij<nRecJetsCuts; ++ij){ // loop over all jets in event
2503 AliAODJet* jet = (AliAODJet*) (fJetsRecCuts->At(ij));
2504 TList* jettracklist = new TList();
2505 Double_t sumPt = 0.;
2506 Bool_t isBadJet = kFALSE;
2508 if(GetFFRadius()<=0){
2509 GetJetTracksTrackrefs(jettracklist, jet, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);// list of jet tracks from trackrefs
2511 GetJetTracksPointing(fTracksRecCuts, jettracklist, jet, GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet); // fill list of tracks in cone around jet axis with cone Radius (= 0.4 standard)
2515 //leading track pt bias on jets inside this small jet loop
2516 if(isBadJet) continue;
2518 if(IsParticleInCone(jet, v0, dRadiusExcludeCone) == kTRUE) {bIsInCone = kTRUE;}
2520 delete jettracklist;
2523 if(bIsInCone == kFALSE){//success! Lambda doesn't belong to any selected jet in event
2524 Double_t vLaOC[3] = {invMLa, trackPt,fEta};
2525 fhnLaOC->Fill(vLaOC);
2528 // Double_t jetPt = fFFIMJetPtMin; // assign pro forma jet energy
2529 // Double_t fRap = v0->Y(3122);
2531 Double_t fV0DecayLength = v0->DecayLengthV0(lPrimaryVtxPosition);
2532 Double_t fV0cosPointAngle = v0->CosPointingAngle(lPrimaryVtxPosition);
2533 lV0Position[0]= v0->DecayVertexV0X();
2534 lV0Position[1]= v0->DecayVertexV0Y();
2535 lV0Position[2]= v0->DecayVertexV0Z();
2537 fV0Radius = TMath::Sqrt(lV0Position[0]*lV0Position[0]+lV0Position[1]*lV0Position[1]);
2539 //fFFHistosIMLaAllEvt->FillFF(trackPt, invMLa, jetPt, incrementJetPt);
2540 //fh1trackPosNCls->Fill(trackPosNcls);
2541 //fh1trackNegNCls->Fill(trackNegNcls);
2542 fh1EtaLa->Fill(fEta);
2544 Double_t vLaIncl[3] = {invMLa,trackPt,fEta};
2545 fhnLaIncl->Fill(vLaIncl);
2549 TString generatorName;
2551 TList* listmc = fAOD->GetList();
2552 Bool_t mclabelcheck = MCLabelCheck(v0, kLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
2553 if(mclabelcheck == kFALSE)continue;
2554 //if(fPhysicalPrimary == kFALSE)continue;
2556 if(generatorName == "Hijing"){
2557 Double_t vrecMCHijingLaIncl[3] = {invMLa,trackPt,fEta};
2558 fhnrecMCHijingLaIncl->Fill(vrecMCHijingLaIncl);
2560 Double_t protonPt = trackPos->Pt();
2561 fh2CorrHijingLaProton->Fill(trackPt,protonPt);
2564 if(isinjected == kTRUE){
2565 Double_t vrecMCInjectLaIncl[3] = {invMLa,trackPt,fEta};
2566 fhnrecMCInjectLaIncl->Fill(vrecMCInjectLaIncl);
2568 Double_t protonPt = trackPos->Pt();
2569 fh2CorrInjectLaProton->Fill(trackPt,protonPt);
2572 Double_t vInvMassEtaTrackPtLa[3] = {fEta,invMLa,trackPt};
2573 fhnInvMassEtaTrackPtLa->Fill(vInvMassEtaTrackPtLa);//includes also feed-down particles
2574 fh1PtMCLa->Fill(MCPt);
2577 fh1PtMCLa->Fill(MCPt);
2581 fh1V0Eta->Fill(fEta);
2582 //fh1V0totMom->Fill(fV0TotalMomentum);
2583 fh1CosPointAngle->Fill(fV0cosPointAngle);
2584 fh1DecayLengthV0->Fill(fV0DecayLength);
2585 fh1V0Radius->Fill(fV0Radius);
2586 fh1DcaV0Daughters->Fill(fDcaV0Daughters);
2587 fh1DcaPosToPrimVertex->Fill(fDcaPosToPrimVertex);
2588 fh1DcaNegToPrimVertex->Fill(fDcaNegToPrimVertex);
2589 fh1trackPosEta->Fill(PosEta);
2590 fh1trackNegEta->Fill(NegEta);
2593 // __ALa pt spectra all events _______________________________________________
2595 for(Int_t it=0; it<fListALa->GetSize(); ++it){ // loop all ALa
2597 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListALa->At(it));
2601 //VO's main characteristics to check the reconstruction cuts
2602 Double_t invMALa =0;
2604 Double_t fV0Radius = -999;
2605 Double_t fDcaV0Daughters = v0->DcaV0Daughters();
2606 Double_t fDcaPosToPrimVertex = v0->DcaPosToPrimVertex();//IP of positive charged daughter
2607 Double_t fDcaNegToPrimVertex = v0->DcaNegToPrimVertex();//IP of negative charged daughter
2608 Int_t negDaughterpdg = 0;
2609 Int_t posDaughterpdg = 0;
2610 Int_t motherType = 0;
2613 Bool_t fPhysicalPrimary = kFALSE;
2614 Int_t MCv0PdgCode = 0;
2616 AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0));
2617 AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1));
2619 Double_t PosEta = trackPos->AliAODTrack::Eta();//daughter track charge is sometimes wrong here, account for that!!!
2620 Double_t NegEta = trackNeg->AliAODTrack::Eta();
2622 Double_t fEta = v0->PseudoRapV0();
2623 Bool_t bIsInCone = kFALSE;//init boolean for OC
2626 CalculateInvMass(v0, kAntiLambda, invMALa, trackPt); //function to calculate invMass with TLorentzVector class
2628 for(Int_t ij=0; ij<nRecJetsCuts; ++ij){ // loop over all jets in event
2630 AliAODJet* jet = (AliAODJet*) (fJetsRecCuts->At(ij));
2631 TList* jettracklist = new TList();
2632 Double_t sumPt = 0.;
2633 Bool_t isBadJet = kFALSE;
2635 if(GetFFRadius()<=0){
2636 GetJetTracksTrackrefs(jettracklist, jet, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);// list of jet tracks from trackrefs
2638 GetJetTracksPointing(fTracksRecCuts, jettracklist, jet, GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet); // fill list of tracks in cone around jet axis with cone Radius (= 0.4 standard)
2641 //leading track pt bias on jets inside this small jet loop
2642 if(isBadJet) continue;
2644 if(IsParticleInCone(jet, v0, dRadiusExcludeCone) == kTRUE){
2648 delete jettracklist;
2651 if(bIsInCone == kFALSE){//success!
2652 Double_t vALaOC[3] = {invMALa, trackPt,fEta};
2653 fhnALaOC->Fill(vALaOC);
2656 //Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2657 //Double_t jetPt = fFFIMJetPtMin; // assign pro forma jet energy
2658 // Double_t fRap = v0->Y(-3122);
2661 Double_t fV0cosPointAngle = v0->CosPointingAngle(lPrimaryVtxPosition);
2662 lV0Position[0]= v0->DecayVertexV0X();
2663 lV0Position[1]= v0->DecayVertexV0Y();
2664 lV0Position[2]= v0->DecayVertexV0Z();
2665 Double_t fV0DecayLength = v0->DecayLengthV0(lPrimaryVtxPosition);
2666 fV0Radius = TMath::Sqrt(lV0Position[0]*lV0Position[0]+lV0Position[1]*lV0Position[1]);
2668 //fFFHistosIMALaAllEvt->FillFF(trackPt, invMALa, jetPt, incrementJetPt);
2669 //fh1trackPosNCls->Fill(trackPosNcls);
2670 //fh1trackNegNCls->Fill(trackNegNcls);
2671 fh1EtaALa->Fill(fEta);
2673 Double_t vALaIncl[3] = {invMALa,trackPt,fEta};
2674 fhnALaIncl->Fill(vALaIncl);
2677 TString generatorName;
2679 TList* listmc = fAOD->GetList();
2680 Bool_t mclabelcheck = MCLabelCheck(v0, kAntiLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
2681 if(mclabelcheck == kFALSE)continue;
2682 //if(fPhysicalPrimary == kFALSE)continue;//take also feeddown particles into account
2684 if(generatorName == "Hijing"){
2685 Double_t vrecMCHijingALaIncl[3] = {invMALa,trackPt,fEta};
2686 fhnrecMCHijingALaIncl->Fill(vrecMCHijingALaIncl);
2688 Double_t aprotonPt = trackNeg->Pt();
2689 fh2CorrHijingALaAProton->Fill(trackPt,aprotonPt);
2693 if(isinjected == kTRUE){
2694 Double_t vrecMCInjectALaIncl[3] = {invMALa,trackPt,fEta};
2695 fhnrecMCInjectALaIncl->Fill(vrecMCInjectALaIncl);
2697 Double_t aprotonPt = trackNeg->Pt();
2698 fh2CorrInjectALaAProton->Fill(trackPt,aprotonPt);
2703 Double_t vInvMassEtaTrackPtALa[3] = {fEta,invMALa,trackPt};
2704 fhnInvMassEtaTrackPtALa->Fill(vInvMassEtaTrackPtALa);
2705 fh1PtMCALa->Fill(MCPt);
2708 fh1V0Eta->Fill(fEta);
2709 //fh1V0totMom->Fill(fV0TotalMomentum);
2710 fh1CosPointAngle->Fill(fV0cosPointAngle);
2711 fh1DecayLengthV0->Fill(fV0DecayLength);
2712 fh1V0Radius->Fill(fV0Radius);
2713 fh1DcaV0Daughters->Fill(fDcaV0Daughters);
2714 fh1DcaPosToPrimVertex->Fill(fDcaPosToPrimVertex);
2715 fh1DcaNegToPrimVertex->Fill(fDcaNegToPrimVertex);
2716 fh1trackPosEta->Fill(PosEta);
2717 fh1trackNegEta->Fill(NegEta);
2720 //_____no jets events______________________________________________________________________________________________________________________________________
2722 if(nRecJetsCuts == 0){//no jet events, before the remaining jet cuts are applied, the second part for the non-jet events comes inside the jet loop
2724 fh1NJ->Fill(1.);//for normalisation by number of NJ events for events in which no rec. jets are found right from the beginning and before even the leading track bias is applied
2726 if(fDebug>6) { std::cout<<"################## nRecJetsCuts == 0 ###################"<<std::endl;
2727 //std::cout<<"fListK0s->GetSize() in NJ event: "<<fListK0s->GetSize()<<std::endl;
2730 for(Int_t it=0; it<fListK0s->GetSize(); ++it){ // loop all K0s
2732 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0s->At(it));
2735 Double_t invMK0s =0;
2737 CalculateInvMass(v0, kK0, invMK0s, trackPt);
2738 Double_t fEta = v0->Eta();
2740 Double_t vNJK0[3] = {invMK0s,trackPt,fEta}; //fill all K0s in events wo selected jets
2741 fhnNJK0->Fill(vNJK0);
2745 for(Int_t it=0; it<fListLa->GetSize(); ++it){ // loop all La
2747 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListLa->At(it));
2752 CalculateInvMass(v0, kLambda, invMLa, trackPt);
2753 Double_t fEta = v0->Eta();
2755 Double_t vNJLa[3] = {invMLa,trackPt,fEta}; //fill all K0s in events wo selected jets
2756 fhnNJLa->Fill(vNJLa);
2760 for(Int_t it=0; it<fListALa->GetSize(); ++it){ // loop all ALa
2762 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListALa->At(it));
2765 Double_t invMALa =0;
2767 CalculateInvMass(v0, kAntiLambda, invMALa, trackPt);
2769 Double_t fEta = v0->Eta();
2771 Double_t vNJALa[3] = {invMALa,trackPt,fEta}; //fill all K0s in events wo selected jets
2772 fhnNJALa->Fill(vNJALa);
2779 //____ fill all jet related histos ________________________________________________________________________________________________________________________
2780 //##########################jet loop########################################################################################################################
2782 Int_t nSelJets = nRecJetsCuts; //init value
2783 Bool_t IsOCEvt = kFALSE; //init for this outside cones normalisation histo (total number of OC events)
2784 Bool_t IsRCEvt = kFALSE; //init for that the random cone is placed only once per event
2785 Bool_t IsMCCEvt = kFALSE; //init for that the median cluster cone is placed only once per event
2787 //fill jet histos in general
2788 for(Int_t ij=0; ij<nRecJetsCuts; ++ij){ // ij is an index running over the list of the reconstructed jets after most of the cuts, but not yet the leading track bias, all jets in event are looped
2790 AliAODJet* jet = (AliAODJet*) (fJetsRecCuts->At(ij));
2792 Double_t jetPt = jet->Pt();
2793 Double_t jetEta = jet->Eta();
2794 Double_t jetPhi = jet->Phi();
2796 //if(ij==0){ // loop over leading jets for ij = 0, for ij>= 0 look into all jets
2798 if(ij>=0){//all jets in event
2800 TList* jettracklist = new TList();
2801 Double_t sumPt = 0.;
2802 Bool_t isBadJet = kFALSE;
2803 Int_t njetTracks = 0;
2805 if(GetFFRadius()<=0){
2806 GetJetTracksTrackrefs(jettracklist, jet, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);// list of jet tracks from trackrefs
2808 GetJetTracksPointing(fTracksRecCuts, jettracklist, jet, GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet); // fill list of charged hybrid tracks in cone around jet axis with cone Radius (= 0.4 standard), application of leading track cut
2811 //not applied at the moment:
2812 if(GetFFMinNTracks()>0 && jettracklist->GetSize() <= GetFFMinNTracks()) isBadJet = kTRUE; // reject jets with less tracks than fFFMinNTracks
2814 //APPLICATION OF REMAINING JET CUTS (leading track pt bias etc..) + NJ events
2817 nSelJets = nSelJets-1;//remove one jet from nSelJets (was initialized with nRecJetsCuts)
2819 if(nSelJets == 0){//case that event doesn't contain no selected jets at all and there are no jets remaining to be looped over
2821 fh1NJ->Fill(1.);//for normalisation by number of NJ events
2823 for(Int_t it=0; it<fListK0s->GetSize(); ++it){ // loop all K0s
2825 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0s->At(it));
2828 Double_t invMK0s =0;
2830 CalculateInvMass(v0, kK0, invMK0s, trackPt);
2831 Double_t fEta = v0->Eta();
2833 Double_t vNJK0[3] = {invMK0s,trackPt,fEta}; //fill all K0s in events wo selected jets
2834 fhnNJK0->Fill(vNJK0);
2838 for(Int_t it=0; it<fListLa->GetSize(); ++it){ // loop all La
2840 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListLa->At(it));
2845 CalculateInvMass(v0, kLambda, invMLa, trackPt);
2846 Double_t fEta = v0->Eta();
2848 Double_t vNJLa[3] = {invMLa,trackPt,fEta}; //fill all K0s in events wo selected jets
2849 fhnNJLa->Fill(vNJLa);
2853 for(Int_t it=0; it<fListALa->GetSize(); ++it){ // loop all ALa
2855 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListALa->At(it));
2858 Double_t invMALa =0;
2860 CalculateInvMass(v0, kAntiLambda, invMALa, trackPt);
2862 Double_t fEta = v0->Eta();
2864 Double_t vNJALa[3] = {invMALa,trackPt,fEta}; //fill all K0s in events wo selected jets
2865 fhnNJALa->Fill(vNJALa);
2871 continue;//rejection of current jet
2872 } // rejects jets in which no track has a track pt higher than 5 GeV/c (see AddTask macro)
2874 if(IsOCEvt == kFALSE){IsOCEvt = kTRUE;fh1OC->Fill(1.);}//the first found jet triggers an OC event and is filled only once into normalisation histo
2876 //Float_t fJetAreaMin = 0.6*TMath::Pi()*GetFFRadius()*GetFFRadius(); // minimum jet area cut, already applied in JetListOfJets() in FF Task
2878 //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
2880 Double_t dAreaExcluded = TMath::Pi()*dRadiusExcludeCone*dRadiusExcludeCone; // area of the cone
2881 dAreaExcluded -= AreaCircSegment(dRadiusExcludeCone,fCutjetEta-jet->Eta()); // positive eta overhang
2882 dAreaExcluded -= AreaCircSegment(dRadiusExcludeCone,fCutjetEta+jet->Eta()); // negative eta overhang
2883 fh1AreaExcluded->Fill(dAreaExcluded);//histo contains all areas that are jet related and have to be excluded concerning OC UE pt spectrum normalisation by area
2885 fh1JetEta->Fill(jetEta);
2886 fh1JetPhi->Fill(jetPhi);
2887 fh2JetEtaPhi->Fill(jetEta,jetPhi);
2889 // printf("pT = %f, eta = %f, phi = %f, leadtr pt = %f\n, ",jetPt,jetEta,jetphi,leadtrack);
2891 for(Int_t it=0; it<jettracklist->GetSize(); ++it){//loop over all particles in jet
2893 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));//all tracks in jet cone
2894 if(!trackVP)continue;
2896 Float_t trackPt = trackVP->Pt();//transversal momentum of jet particle
2897 Float_t trackEta = trackVP->Eta();
2899 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
2901 fFFHistosRecCuts->FillFF(trackPt, jetPt, incrementJetPt);//histo with tracks/jets after cut selection, for all events
2902 if(nK0s>0) fFFHistosRecCutsK0Evt->FillFF(trackPt, jetPt, incrementJetPt);//only for K0s events
2903 fh2FFJetTrackEta->Fill(trackEta,jetPt);
2908 njetTracks = jettracklist->GetSize();
2910 //____________________________________________________________________________________________________________________
2911 //strangeness constribution to jet cone
2915 TList *list = fAOD->GetList();
2916 AliAODMCHeader *mcHeadr=(AliAODMCHeader*)list->FindObject(AliAODMCHeader::StdBranchName());
2917 if(!mcHeadr)continue;
2919 Double_t mcXv=0., mcYv=0., mcZv=0.;//MC primary vertex position
2921 mcXv=mcHeadr->GetVtxX(); mcYv=mcHeadr->GetVtxY(); mcZv=mcHeadr->GetVtxZ(); // position of the MC primary vertex
2923 for(Int_t it=0; it<jettracklist->GetSize(); ++it){//loop over all tracks in the jet
2925 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));//track in jet cone
2926 if(!trackVP)continue;
2927 AliAODTrack *tr = dynamic_cast<AliAODTrack*> (trackVP); //fetch one jet track from the TList
2930 //get MC label information
2931 TList *mclist = fAOD->GetList();
2933 //fetch the MC stack
2934 TClonesArray* stackMC = (TClonesArray*)mclist->FindObject(AliAODMCParticle::StdBranchName()); //get MCAOD branch in data
2935 if (!stackMC) {Printf("ERROR: stack not available");}
2939 Int_t trackLabel = TMath::Abs(tr->GetLabel()); //fetch jet track label in MC stack
2941 AliAODMCParticle *part = dynamic_cast<AliAODMCParticle*>(stackMC->At(trackLabel)); //fetch MC gen. particle for rec. jet track
2943 if(!part)continue; //skip non-existing objects
2946 //Bool_t IsPhysicalPrimary = part->IsPhysicalPrimary();//not recommended to check, better use distance between primary vertex and secondary vertex
2948 Float_t fDistPrimaryMax = 0.01;
2949 // Get the distance between production point of the MC mother particle and the primary vertex
2951 Double_t dx = mcXv-part->Xv();//mc primary vertex - mc gen. v0 vertex
2952 Double_t dy = mcYv-part->Yv();
2953 Double_t dz = mcZv-part->Zv();
2955 Float_t fDistPrimary = TMath::Sqrt(dx*dx + dy*dy + dz*dz);
2956 Bool_t fPhysicalPrimary = (fDistPrimary < fDistPrimaryMax);
2958 // std::cout<<"fDistPrimary"<<fDistPrimary<<std::endl;
2959 // std::cout<<"fPhysicalPrimary"<<fPhysicalPrimary<<std::endl;
2961 if(!fPhysicalPrimary)continue;//rejects Kstar and other strong decaying particles from Secondary Contamination
2963 Bool_t isFromStrange = kFALSE;// flag to check whether particle has strange mother
2965 Int_t iMother = part->GetMother(); //get mother MC gen. particle label
2968 AliAODMCParticle *partM = dynamic_cast<AliAODMCParticle*>(stackMC->At(iMother)); //fetch mother of MC gen. particle
2969 if(!partM) continue;
2971 Int_t codeM = TMath::Abs(partM->GetPdgCode()); //mothers pdg code
2973 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..)
2975 if (mfl == 3 && codeM != 3) isFromStrange = kTRUE;
2978 if(isFromStrange == kTRUE){
2980 Double_t trackPt = part->Pt();
2981 Double_t trackEta = part->Eta();
2982 //fh3StrContinCone->Fill(jetPt, trackPt, trackEta);//MC gen. particle parameters, but rec. jet pt
2984 }//isFromStrange is kTRUE
2986 }//end loop over jet tracks
2991 fh1TrackMultCone->Fill(njetTracks);
2992 fh2TrackMultCone->Fill(njetTracks,jetPt);
2996 // fQAJetHistosRecCutsLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
2998 for(Int_t it=0; it<fListK0s->GetSize(); ++it){ // loop all K0s
3000 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0s->At(it));
3001 if(!v0) continue;//rejection of events with no V0 vertex
3005 TVector3 v0MomVect(v0Mom);
3007 Double_t dPhiJetK0 = (jet->MomentumVector()->Vect()).DeltaPhi(v0MomVect);
3008 // Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3010 // if(incrementJetPt){fh1V0JetPt->Fill(jetPt);}
3012 Double_t invMK0s =0;
3014 CalculateInvMass(v0, kK0, invMK0s, trackPt); //function to calculate invMass with TLorentzVector class
3016 // fFFHistosIMK0Jet->FillFF(trackPt, invMK0s, jetPt, incrementJetPt);
3019 if(dPhiJetK0<fh1dPhiJetK0->GetXaxis()->GetXmin()) dPhiJetK0 += 2*TMath::Pi();
3020 fh1dPhiJetK0->Fill(dPhiJetK0);
3024 // if(fListK0s->GetSize() == 0){ // no K0: increment jet pt spectrum
3026 // Bool_t incrementJetPt = kTRUE;
3027 // fFFHistosIMK0Jet->FillFF(-1, -1, jetPt, incrementJetPt);
3030 //____fetch reconstructed K0s in cone around jet axis:_______________________________________________________________________________
3032 TList* jetConeK0list = new TList();
3034 Double_t sumPtK0 = 0.;
3036 Bool_t isBadJetK0 = kFALSE; // dummy, do not use
3038 GetTracksInCone(fListK0s, jetConeK0list, jet, GetFFRadius(), sumPtK0, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetK0); //reconstructed K0s in cone around jet axis
3040 if(fDebug>2)Printf("%s:%d nK0s total: %d, in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nK0s,jetConeK0list->GetEntries(),GetFFRadius());
3043 for(Int_t it=0; it<jetConeK0list->GetSize(); ++it){ // loop for K0s in jet cone
3045 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetConeK0list->At(it));
3048 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3049 Double_t invMK0s =0;
3054 CalculateInvMass(v0, kK0, invMK0s, trackPt); //function to calculate invMass with TLorentzVector class
3058 Double_t jetPtSmear = -1;
3059 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
3060 if(incrementJetPt == kTRUE){fh1IMK0ConeSmear->Fill(jetPtSmear);} //fill TH1F for normalization purposes
3063 if(incrementJetPt==kTRUE){
3064 fh1IMK0Cone->Fill(jetPt);}//normalisation by number of selected jets
3066 //fFFHistosIMK0Cone->FillFF(trackPt, invMK0s, jetPt, incrementJetPt);
3068 Double_t vK0sCone[4] = {jetPt, invMK0s,trackPt,fEta};
3069 fhnK0sCone->Fill(vK0sCone);
3073 if(jetConeK0list->GetSize() == 0){ // no K0: increment jet pt spectrum
3076 Bool_t incrementJetPt = kTRUE;//jets without K0s will be only filled in TH1F only once, so no increment needed
3077 //fFFHistosIMK0Cone->FillFF(-1, -1, jetPt, incrementJetPt);
3078 Double_t vK0sCone[4] = {jetPt, -1, -1, -1};
3079 fhnK0sCone->Fill(vK0sCone);
3081 if(incrementJetPt==kTRUE){
3082 fh1IMK0Cone->Fill(jetPt);}//normalisation by number of selected jets
3085 Double_t jetPtSmear = -1;
3086 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
3087 if(incrementJetPt == kTRUE){fh1IMK0ConeSmear->Fill(jetPtSmear);} //fill TH1F for normalization purposes
3091 //Random cones________________________________________________________________________
3094 if(IsRCEvt == kFALSE){//fetch random cone V0s only once per event
3097 IsRCEvt = kTRUE;//set boolean to kTRUE once a random cone is placed per event
3099 AliAODJet* jetRC = 0;
3100 jetRC = GetRandomCone(fJetsRecCuts, fCutjetEta, 2*GetFFRadius());//fetch one random cone for each event
3102 TList* fListK0sRC = new TList();//list for K0s in random cone (RC), one RC per event
3103 TList* fListLaRC = new TList();
3104 TList* fListALaRC = new TList();
3106 Double_t sumPtK0sRC = 0;
3107 Double_t sumPtLaRC = 0;
3108 Double_t sumPtALaRC = 0;
3109 Bool_t isBadJetK0sRC = kFALSE;
3110 Bool_t isBadJetLaRC = kFALSE;
3111 Bool_t isBadJetALaRC = kFALSE;
3114 if(jetRC != 0) {//if random cone was selected properly and fullfilling all the requirements
3117 fh1RC->Fill(1.);//for normalisation purposes
3119 GetTracksInCone(fListK0s, fListK0sRC, jetRC, GetFFRadius(), sumPtK0sRC, 0, 0, isBadJetK0sRC);
3121 //________________fill RC with all V0s__________________
3122 for(Int_t it=0; it<fListK0sRC->GetSize(); ++it){ // loop for K0s in random cone
3124 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0sRC->At(it));
3127 Double_t invMK0s =0;
3132 CalculateInvMass(v0, kK0, invMK0s, trackPt); //function to calculate invMass with TLorentzVector class
3134 Double_t vK0sRC[3] = {invMK0s,trackPt,fEta};
3135 fhnK0sRC->Fill(vK0sRC);
3140 GetTracksInCone(fListLa, fListLaRC, jetRC, GetFFRadius(), sumPtLaRC, 0, 0, isBadJetLaRC);
3142 for(Int_t it=0; it<fListLaRC->GetSize(); ++it){ // loop for Lambdas in random cone
3144 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListLaRC->At(it));
3152 CalculateInvMass(v0, kLambda, invMLa, trackPt); //function to calculate invMass with TLorentzVector class
3154 Double_t vLaRC[3] = {invMLa,trackPt,fEta};
3155 fhnLaRC->Fill(vLaRC);
3160 GetTracksInCone(fListALa, fListALaRC, jetRC, GetFFRadius(), sumPtALaRC, 0, 0, isBadJetALaRC);
3162 for(Int_t it=0; it<fListALaRC->GetSize(); ++it){ // loop for Lambdas in random cone
3164 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListALaRC->At(it));
3167 Double_t invMALa =0;
3172 CalculateInvMass(v0, kAntiLambda, invMALa, trackPt); //function to calculate invMass with TLorentzVector class
3174 Double_t vALaRC[3] = {invMALa,trackPt,fEta};
3175 fhnALaRC->Fill(vALaRC);
3179 if(isBadJetK0sRC == kFALSE){ //in case RC contains at least one K0s with minimum pT
3180 fh1RCBiasK0->Fill(1.);//for normalisation purposes
3182 //________________fill RC (with trigger particle bias)_____________
3183 for(Int_t it=0; it<fListK0sRC->GetSize(); ++it){ // loop for K0s in random cone
3185 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0sRC->At(it));
3188 Double_t invMK0s =0;
3193 CalculateInvMass(v0, kK0, invMK0s, trackPt); //function to calculate invMass with TLorentzVector class
3195 //Double_t vK0sRC[3] = {invMK0s,trackPt,fEta};
3196 //fhnK0sRCBias->Fill(vK0sRC);
3201 if(isBadJetLaRC == kFALSE){ //in case RC contains at least one Lambda with minimum pT
3202 fh1RCBiasLa->Fill(1.);//for normalisation purposes
3203 for(Int_t it=0; it<fListLaRC->GetSize(); ++it){ // loop for Lambdas in random cone
3205 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListLaRC->At(it));
3213 CalculateInvMass(v0, kLambda, invMLa, trackPt); //function to calculate invMass with TLorentzVector class
3215 //Double_t vLaRC[3] = {invMLa,trackPt,fEta};
3216 //fhnLaRCBias->Fill(vLaRC);
3222 if(isBadJetALaRC == kFALSE){ //in case RC contains at least one Antilambda with minimum pT
3223 fh1RCBiasALa->Fill(1.);//for normalisation purposes
3224 for(Int_t it=0; it<fListALaRC->GetSize(); ++it){ // loop for Lambdas in random cone
3226 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListALaRC->At(it));
3229 Double_t invMALa =0;
3234 CalculateInvMass(v0, kAntiLambda, invMALa, trackPt); //function to calculate invMass with TLorentzVector class
3236 //Double_t vALaRC[3] = {invMALa,trackPt,fEta};
3237 //fhnALaRCBias->Fill(vALaRC);
3250 //fetch particles in perpendicular cone to estimate UE event contribution to particle spectrum
3251 //these perpendicular cone particle spectra serve to subtract the particles in jet cones, that are stemming from the Underlying event, on a statistical basis
3252 //for normalization the common jet pT spectrum is used: fh1IMK0Cone, fh1IMLaCone and fh1IMALaCone
3254 //____fetch reconstructed K0s in cone perpendicular to jet axis:_______________________________________________________________________________
3256 TList* jetPerpConeK0list = new TList();
3258 Double_t sumPerpPtK0 = 0.;
3260 GetTracksInPerpCone(fListK0s, jetPerpConeK0list, jet, GetFFRadius(), sumPerpPtK0); //reconstructed K0s in cone around jet axis
3262 if(fDebug>2)Printf("%s:%d nK0s total: %d, in perp jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nK0s,jetPerpConeK0list->GetEntries(),GetFFRadius());
3264 for(Int_t it=0; it<jetPerpConeK0list->GetSize(); ++it){ // loop for K0s in perpendicular cone
3266 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetPerpConeK0list->At(it));
3269 Double_t invMPerpK0s =0;
3274 CalculateInvMass(v0, kK0, invMPerpK0s, trackPt); //function to calculate invMass with TLorentzVector class
3275 Double_t vK0sPC[4] = {jetPt, invMPerpK0s,trackPt,fEta};
3277 fhnK0sPC->Fill(vK0sPC); //(x,y,z) //pay attention, this histogram contains the V0 content of both (+/- 90 degrees) perp. cones!!
3282 if(jetPerpConeK0list->GetSize() == 0){ // no K0s in jet cone
3284 Double_t vK0sPC[4] = {jetPt, -1, -1 , -999};//default values for case: no K0s is found in PC
3285 fhnK0sPC->Fill(vK0sPC);
3290 if(IsMCCEvt == kFALSE){//median cluster only once for event
3296 AliAODJet* medianCluster = GetMedianCluster();
3299 // ____ rec K0s in median cluster___________________________________________________________________________________________________________
3301 TList* jetMedianConeK0list = new TList();
3302 TList* jetMedianConeLalist = new TList();
3303 TList* jetMedianConeALalist = new TList();
3306 Double_t medianEta = medianCluster->Eta();
3308 if(TMath::Abs(medianEta)<=fCutjetEta){
3310 fh1MedianEta->Fill(medianEta);
3311 fh1JetPtMedian->Fill(jetPt);
3312 fh1MCC->Fill(1.);//for normalisation by total number of median cluster jets
3313 Double_t sumMedianPtK0 = 0.;
3315 Bool_t isBadJetK0Median = kFALSE; // dummy, do not use
3317 GetTracksInCone(fListK0s, jetMedianConeK0list, medianCluster, GetFFRadius(), sumMedianPtK0, 0., 0., isBadJetK0Median); //reconstructed K0s in median cone around jet axis
3318 //GetTracksInCone(fListK0s, jetConeK0list, jet, GetFFRadius(), sumPtK0, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetK0); //original use of function
3320 //cut parameters from Fragmentation Function task:
3321 //Float_t fFFMinLTrackPt; // reject jets with leading track with pt smaller than this value, use GetFFMinLTrackPt()
3322 //Float_t fFFMaxTrackPt; // reject jetscontaining any track with pt larger than this value, use GetFFMaxTrackPt()
3324 for(Int_t it=0; it<jetMedianConeK0list->GetSize(); ++it){ // loop for K0s in median cone
3326 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetMedianConeK0list->At(it));
3329 Double_t invMMedianK0s =0;
3334 CalculateInvMass(v0, kK0, invMMedianK0s, trackPt); //function to calculate invMass with TLorentzVector class
3335 Double_t vK0sMCC[3] = {invMMedianK0s,trackPt,fEta};
3336 fhnK0sMCC->Fill(vK0sMCC);
3340 if(jetMedianConeK0list->GetSize() == 0){ // no K0s in median cluster cone
3342 Double_t vK0sMCC[3] = {-1, -1, -999};
3343 fhnK0sMCC->Fill(vK0sMCC);
3347 //__________________________________________________________________________________________________________________________________________
3348 // ____ rec Lambdas in median cluster___________________________________________________________________________________________________________
3350 Double_t sumMedianPtLa = 0.;
3351 Bool_t isBadJetLaMedian = kFALSE; // dummy, do not use
3353 GetTracksInCone(fListLa, jetMedianConeLalist, medianCluster, GetFFRadius(), sumMedianPtLa, 0, 0, isBadJetLaMedian); //reconstructed Lambdas in median cone around jet axis
3355 //cut parameters from Fragmentation Function task:
3356 //Float_t fFFMinLTrackPt; // reject jets with leading track with pt smaller than this value, use GetFFMinLTrackPt()
3357 //Float_t fFFMaxTrackPt; // reject jets containing any track with pt larger than this value, use GetFFMaxTrackPt()
3359 for(Int_t it=0; it<jetMedianConeLalist->GetSize(); ++it){ // loop for Lambdas in perpendicular cone
3361 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetMedianConeLalist->At(it));
3364 Double_t invMMedianLa =0;
3369 CalculateInvMass(v0, kLambda, invMMedianLa, trackPt); //function to calculate invMass with TLorentzVector class
3371 Double_t vLaMCC[3] = {invMMedianLa,trackPt,fEta};
3372 fhnLaMCC->Fill(vLaMCC);
3375 if(jetMedianConeLalist->GetSize() == 0){ // no Lambdas in median cluster cone
3377 Double_t vLaMCC[4] = {jetPt, -1, -1, -999};
3378 fhnLaMCC->Fill(vLaMCC);
3383 // ____ rec Antilambdas in median cluster___________________________________________________________________________________________________________
3386 Double_t sumMedianPtALa = 0.;
3388 Bool_t isBadJetALaMedian = kFALSE; // dummy, do not use
3390 GetTracksInCone(fListALa, jetMedianConeALalist, medianCluster, GetFFRadius(), sumMedianPtALa, 0, 0, isBadJetALaMedian); //reconstructed Antilambdas in median cone around jet axis
3393 //cut parameters from Fragmentation Function task:
3394 //Float_t fFFMinLTrackPt; // reject jets with leading track with pt smaller than this value, use GetFFMinLTrackPt()
3395 //Float_t fFFMaxTrackPt; // reject jets containing any track with pt larger than this value, use GetFFMaxTrackPt()
3397 for(Int_t it=0; it<jetMedianConeALalist->GetSize(); ++it){ // loop for Antilambdas in median cluster cone
3399 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetMedianConeALalist->At(it));
3402 Double_t invMMedianALa =0;
3408 CalculateInvMass(v0, kAntiLambda, invMMedianALa, trackPt); //function to calculate invMass with TLorentzVector class
3409 Double_t vALaMCC[3] = {invMMedianALa,trackPt,fEta};
3410 fhnALaMCC->Fill(vALaMCC);
3414 if(jetMedianConeALalist->GetSize() == 0){ // no Antilambdas in median cluster cone
3416 Double_t vALaMCC[4] = {jetPt, -1, -1, -999};
3417 fhnALaMCC->Fill(vALaMCC);
3420 }//median cluster eta cut
3422 delete jetMedianConeK0list;
3423 delete jetMedianConeLalist;
3424 delete jetMedianConeALalist;
3426 }//if mediancluster is existing
3427 }//end (IsMCCEvt == kFALSE)
3428 //_________________________________________________________________________________________________________________________________________
3430 //____fetch reconstructed Lambdas in cone perpendicular to jet axis:__________________________________________________________________________
3432 TList* jetPerpConeLalist = new TList();
3434 Double_t sumPerpPtLa = 0.;
3436 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!!
3438 if(fDebug>2)Printf("%s:%d nLa total: %d, in perp jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nLa,jetPerpConeLalist->GetEntries(),GetFFRadius());
3440 for(Int_t it=0; it<jetPerpConeLalist->GetSize(); ++it){ // loop for Lambdas in perpendicular cone
3442 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetPerpConeLalist->At(it));
3445 Double_t invMPerpLa =0;
3450 CalculateInvMass(v0, kLambda, invMPerpLa, trackPt); //function to calculate invMass with TLorentzVector class
3451 Double_t vLaPC[4] = {jetPt, invMPerpLa,trackPt,fEta};
3452 fhnLaPC->Fill(vLaPC); //(x,y,z) //pay attention, this histogram contains the V0 content of both (+/- 90 degrees) perp. cones!!
3457 if(jetPerpConeLalist->GetSize() == 0){ // no Lambdas in jet
3459 Double_t vLaPC[4] = {jetPt, -1, -1 , -999};//default values for case: no K0s is found in PC
3460 fhnLaPC->Fill(vLaPC);
3466 //____fetch reconstructed Antilambdas in cone perpendicular to jet axis:___________________________________________________________________
3468 TList* jetPerpConeALalist = new TList();
3470 Double_t sumPerpPtALa = 0.;
3472 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!!
3474 if(fDebug>2)Printf("%s:%d nALa total: %d, in perp jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nALa,jetPerpConeALalist->GetEntries(),GetFFRadius());
3476 for(Int_t it=0; it<jetPerpConeALalist->GetSize(); ++it){ // loop for ALa in perpendicular cone
3478 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetPerpConeALalist->At(it));
3481 Double_t invMPerpALa =0;
3486 CalculateInvMass(v0, kAntiLambda, invMPerpALa, trackPt); //function to calculate invMass with TLorentzVector class
3487 Double_t vALaPC[4] = {jetPt, invMPerpALa,trackPt,fEta};
3488 fhnALaPC->Fill(vALaPC);
3493 if(jetPerpConeALalist->GetSize() == 0){ // no Antilambda
3495 Double_t vALaPC[4] = {jetPt, -1, -1, -999};
3496 fhnALaPC->Fill(vALaPC);
3502 //###########################################################################################################
3504 //__________________________________________________________________________________________________________________________________________
3508 //fill feeddown candidates from TList
3509 //std::cout<<"fListFeeddownLaCand entries: "<<fListFeeddownLaCand->GetSize()<<std::endl;
3511 Double_t sumPtFDLa = 0.;
3512 Bool_t isBadJetFDLa = kFALSE; // dummy, do not use
3514 GetTracksInCone(fListFeeddownLaCand, jetConeFDLalist, jet, GetFFRadius(), sumPtFDLa, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetFDLa);
3516 Double_t sumPtFDALa = 0.;
3517 Bool_t isBadJetFDALa = kFALSE; // dummy, do not use
3519 GetTracksInCone(fListFeeddownALaCand, jetConeFDALalist, jet, GetFFRadius(), sumPtFDALa, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetFDALa);
3521 //_________________________________________________________________
3522 for(Int_t it=0; it<fListFeeddownLaCand->GetSize(); ++it){
3524 AliAODv0* mcfd = dynamic_cast<AliAODv0*>(fListFeeddownLaCand->At(it));
3527 Double_t invMLaFDcand = 0;
3528 Double_t trackPt = 0;//pt of ass. particle, not used for the histos
3530 CalculateInvMass(mcfd, kLambda, invMLaFDcand, trackPt);
3532 //Get MC gen. Lambda transverse momentum
3533 TClonesArray *st = 0x0;
3536 TList *lt = fAOD->GetList();
3539 st = (TClonesArray*)lt->FindObject(AliAODMCParticle::StdBranchName()); //get MCAOD branch in data
3542 AliAODTrack *daughtertrack = (AliAODTrack *) (mcfd->GetSecondaryVtx()->GetDaughter(0));//fetch the first of the two daughter tracks
3543 Int_t AssLabel = TMath::Abs(daughtertrack->GetLabel());
3545 AliAODMCParticle *mcDaughterPart =(AliAODMCParticle*)st->UncheckedAt(AssLabel);
3547 Int_t v0lab = mcDaughterPart->GetMother();
3549 // Int_t v0lab= TMath::Abs(mcfd->GetLabel());//GetLabel doesn't work for AliAODv0 class!!! Only for AliAODtrack
3551 if((!v0lab) || (v0lab<0) || (v0lab > st->GetEntriesFast()))continue;//validity checks
3553 AliAODMCParticle *mcp=(AliAODMCParticle*)st->UncheckedAt(v0lab);
3555 Double_t genLaPt = mcp->Pt();
3557 //std::cout<<"Incl FD, genLaPt:"<<genLaPt<<std::endl;
3559 Double_t vFeedDownLa[3] = {5., invMLaFDcand, genLaPt};
3560 fhnFeedDownLa->Fill(vFeedDownLa);
3563 }//end loop over feeddown candidates for Lambda particles in jet cone
3564 //fetch MC truth in jet cones, denominator of rec. efficiency in jet cones
3565 //_________________________________________________________________
3566 for(Int_t it=0; it<jetConeFDLalist->GetSize(); ++it){
3568 AliAODv0* mcfd = dynamic_cast<AliAODv0*>(jetConeFDLalist->At(it));
3571 //std::cout<<"Cone, recLaPt:"<<mcfd->Pt()<<std::endl;
3573 Double_t invMLaFDcand = 0;
3574 Double_t trackPt = mcfd->Pt();//pt of ass. particle, not used for the histos
3576 CalculateInvMass(mcfd, kLambda, invMLaFDcand, trackPt);
3578 //Get MC gen. Lambda transverse momentum
3579 TClonesArray *st = 0x0;
3581 TList *lt = fAOD->GetList();
3584 st = (TClonesArray*)lt->FindObject(AliAODMCParticle::StdBranchName());
3586 AliAODTrack *daughtertrack = (AliAODTrack *) (mcfd->GetSecondaryVtx()->GetDaughter(0));//fetch the first of the two daughter tracks
3587 Int_t AssLabel = TMath::Abs(daughtertrack->GetLabel());
3589 AliAODMCParticle *mcDaughterPart =(AliAODMCParticle*)st->UncheckedAt(AssLabel);
3591 Int_t v0lab = mcDaughterPart->GetMother();
3593 //std::cout<<"v0lab: "<<v0lab<<std::endl;
3595 if((!v0lab) || (v0lab<0) || (v0lab > st->GetEntriesFast()))continue;//validity checks
3597 AliAODMCParticle *mcp=(AliAODMCParticle*)st->UncheckedAt(v0lab);
3599 Double_t genLaPt = mcp->Pt();
3602 //std::cout<<"Cone FD, genLaPt:"<<genLaPt<<std::endl;
3604 Double_t vFeedDownLaCone[3] = {jetPt, invMLaFDcand, genLaPt};
3605 fhnFeedDownLaCone->Fill(vFeedDownLaCone);
3608 }//end loop over feeddown candidates for Lambda particles in jet cone
3610 //_________________________________________________________________
3611 for(Int_t it=0; it<fListFeeddownALaCand->GetSize(); ++it){
3613 AliAODv0* mcfd = dynamic_cast<AliAODv0*>(fListFeeddownALaCand->At(it));
3616 Double_t invMALaFDcand = 0;
3617 Double_t trackPt = 0;//pt of ass. particle, not used for the histos
3619 CalculateInvMass(mcfd, kAntiLambda, invMALaFDcand, trackPt);
3621 //Get MC gen. Antilambda transverse momentum
3622 TClonesArray *st = 0x0;
3624 TList *lt = fAOD->GetList();
3627 st = (TClonesArray*)lt->FindObject(AliAODMCParticle::StdBranchName());
3629 AliAODTrack *daughtertrack = (AliAODTrack *) (mcfd->GetSecondaryVtx()->GetDaughter(0));//fetch the first of the two daughter tracks
3630 Int_t AssLabel = TMath::Abs(daughtertrack->GetLabel());
3632 AliAODMCParticle *mcDaughterPart =(AliAODMCParticle*)st->UncheckedAt(AssLabel);
3634 Int_t v0lab = mcDaughterPart->GetMother();
3637 if((!v0lab) || (v0lab<0) || (v0lab > st->GetEntriesFast()))continue;//validity checks
3639 AliAODMCParticle *mcp=(AliAODMCParticle*)st->UncheckedAt(v0lab);
3641 Double_t genALaPt = mcp->Pt();
3643 Double_t vFeedDownALa[3] = {5., invMALaFDcand, genALaPt};
3644 fhnFeedDownALa->Fill(vFeedDownALa);
3647 }//end loop over feeddown candidates for Antilambda particles
3649 //_________________________________________________________________
3650 //feeddown for Antilambdas from Xi(bar)+ and Xi(bar)0 in jet cone:
3652 for(Int_t it=0; it<jetConeFDALalist->GetSize(); ++it){
3654 AliAODv0* mcfd = dynamic_cast<AliAODv0*>(jetConeFDALalist->At(it));
3657 Double_t invMALaFDcand = 0;
3658 Double_t trackPt = 0;//pt of ass. particle, not used for the histos
3660 CalculateInvMass(mcfd, kAntiLambda, invMALaFDcand, trackPt);
3662 //Get MC gen. Antilambda transverse momentum
3663 TClonesArray *st = 0x0;
3665 TList *lt = fAOD->GetList();
3668 st = (TClonesArray*)lt->FindObject(AliAODMCParticle::StdBranchName());
3670 AliAODTrack *daughtertrack = (AliAODTrack *) (mcfd->GetSecondaryVtx()->GetDaughter(0));//fetch the first of the two daughter tracks
3671 Int_t AssLabel = TMath::Abs(daughtertrack->GetLabel());
3673 AliAODMCParticle *mcDaughterPart =(AliAODMCParticle*)st->UncheckedAt(AssLabel);
3675 Int_t v0lab = mcDaughterPart->GetMother();
3677 if((!v0lab) || (v0lab<0) || (v0lab > st->GetEntriesFast()))continue;//validity checks
3679 AliAODMCParticle *mcp=(AliAODMCParticle*)st->UncheckedAt(v0lab);
3681 Double_t genALaPt = mcp->Pt();
3683 Double_t vFeedDownALaCone[3] = {jetPt, invMALaFDcand, genALaPt};
3684 fhnFeedDownALaCone->Fill(vFeedDownALaCone);
3687 }//end loop over feeddown candidates for Antilambda particles in jet cone
3691 //____fetch MC generated K0s in cone around jet axis__(note: particles can stem from fragmentation but also from underlying event)________
3693 Double_t sumPtMCgenK0s = 0.;
3694 Bool_t isBadJetMCgenK0s = kFALSE; // dummy, do not use
3697 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!!)
3699 //first: sampling MC gen K0s
3701 GetTracksInCone(fListMCgenK0s, fListMCgenK0sCone, jet, GetFFRadius(), sumPtMCgenK0s, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetMCgenK0s); //MC generated K0s in cone around jet axis
3703 if(fDebug>2)Printf("%s:%d nMCgenK0s in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,fListMCgenK0sCone->GetEntries(),GetFFRadius());
3706 /* for(Int_t it=0; it<fListMCgenK0sCone->GetSize(); ++it){ // loop MC generated K0s in cone around jet axis
3708 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenK0sCone->At(it));
3711 //Double_t fRapMCgenK0s = MyRapidity(mcp0->E(),mcp0->Pz());//get rec. particle in cone information
3712 Double_t fEtaMCgenK0s = mcp0->Eta();
3713 Double_t fPtMCgenK0s = mcp0->Pt();
3715 //fh2MCgenK0Cone->Fill(jetPt,fPtMCgenK0s);
3716 // fh2MCEtagenK0Cone->Fill(jetPt,fEtaMCgenK0s);
3720 //check whether the reconstructed K0s in jet cone are stemming from MC gen K0s (on MCgenK0s list):__________________________________________________
3722 for(Int_t ic=0; ic<jetConeK0list->GetSize(); ++ic){ //loop over all reconstructed K0s in jet cone
3724 //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
3726 Int_t negDaughterpdg;
3727 Int_t posDaughterpdg;
3730 Double_t fPtMCrecK0Match;
3731 Double_t invMK0Match;
3735 Bool_t fPhysicalPrimary = -1;
3736 Int_t MCv0PDGCode =0;
3737 Double_t jetPtSmear = -1;
3739 AliAODv0* v0c = dynamic_cast<AliAODv0*>(jetConeK0list->At(ic));//pointer to reconstructed K0s inside jet cone (cone is placed around reconstructed jet axis)
3741 //AliAODv0* v0c = dynamic_cast<AliAODv0*>(fListK0s->At(ic));//pointer to reconstructed K0s
3744 Bool_t daughtercheck = DaughterTrackCheck(v0c, nnum, pnum);//check daughter tracks have proper sign
3745 if(daughtercheck == kFALSE)continue;
3747 const AliAODTrack *trackMCNeg=(AliAODTrack *)(v0c->GetDaughter(nnum));
3748 const AliAODTrack *trackMCPos=(AliAODTrack *)(v0c->GetDaughter(pnum));
3750 TString generatorName;
3751 TList *listmc = fAOD->GetList();
3753 Bool_t mclabelcheck = MCLabelCheck(v0c, kK0, trackMCNeg, trackMCPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PDGCode, generatorName, isinjected);
3755 if(mclabelcheck == kFALSE)continue;
3756 if(fPhysicalPrimary == kFALSE)continue; //requirements for rec. V0 associated to MC true primary particle
3758 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
3760 //for(Int_t it=0; it<fListMCgenK0sCone->GetSize(); ++it){//belongs to previous definition of rec. eff. of V0s within jet cone
3762 //Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3763 //AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenK0sCone->At(it));
3764 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenK0s->At(it));
3767 Bool_t particleMatching = IsParticleMatching(mcp0, v0Label);
3769 if(particleMatching == kFALSE)continue; //if reconstructed V0 particle doesn't match to the associated MC particle go to next stack entry
3770 CalculateInvMass(v0c, kK0, invMK0Match, fPtMCrecK0Match);
3771 Double_t fEta = v0c->Eta();
3772 Double_t fPtMCgenK0s = mcp0->Pt();//pt has to be always MC truth value!
3774 Double_t vMCrecK0Cone[4] = {jetPt, invMK0Match,fPtMCgenK0s,fEta};
3775 fhnMCrecK0Cone->Fill(vMCrecK0Cone); //fill matching rec. K0s in 3D histogram
3777 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear); //jetPt, cent, jetRadius, ptmintrack, &jetPtSmear
3779 Double_t vMCrecK0ConeSmear[4] = {jetPtSmear, invMK0Match,fPtMCgenK0s,fEta};
3780 fhnMCrecK0ConeSmear->Fill(vMCrecK0ConeSmear);
3782 //fill matching rec. K0s in 3D histogram, jet pT smeared according to deltaptjet distribution width
3785 } // end MCgenK0s / MCgenK0sCone loop
3788 //check the K0s daughters contamination of the jet tracks:
3790 TClonesArray *stackMC = 0x0;
3792 for(Int_t it=0; it<jettracklist->GetSize(); ++it){//loop over all particles in the jet
3794 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));//all tracks in jet cone
3795 if(!trackVP)continue;
3796 AliAODTrack *tr = dynamic_cast<AliAODTrack*> (trackVP); //fetch one jet track from the TList
3799 //get MC label information
3800 TList *mclist = fAOD->GetList(); //fetch the MC stack
3801 if(!mclist)continue;
3803 stackMC = (TClonesArray*)mclist->FindObject(AliAODMCParticle::StdBranchName()); //get MCAOD branch in data
3804 if (!stackMC) {Printf("ERROR: stack not available");}
3807 Int_t particleLabel = TMath::Abs(tr->GetLabel()); //fetch jet track label in MC stack
3809 //v0c is pointer to K0s candidate, is fetched already above, here it is just checked again whether daughters are properly ordered by their charge
3811 Bool_t daughterchecks = DaughterTrackCheck(v0c, nnum, pnum);
3813 if(daughterchecks == kFALSE)continue; //make sure that daughters are properly ordered
3815 const AliAODTrack *trackNeg=(AliAODTrack *)(v0c->GetDaughter(nnum)); //fetch v0 daughters of reconstructed K0s
3816 const AliAODTrack *trackPos=(AliAODTrack *)(v0c->GetDaughter(pnum));
3818 if(!trackNeg)continue;
3819 if(!trackPos)continue;
3821 Int_t negAssLabel = TMath::Abs(trackNeg->GetLabel()); //negative (reconstructed) charged track label in MC stack
3822 Int_t posAssLabel = TMath::Abs(trackPos->GetLabel()); //positive (reconstructed) charged track label in MC stack
3825 if(particleLabel == posAssLabel){ //check whether jet track and each of the rec. K0s daughters have same MC stack label -> are identical
3826 AliAODMCParticle* mctrackPos = dynamic_cast<AliAODMCParticle*>(stackMC->At(posAssLabel));
3827 if(!mctrackPos) continue;
3828 Double_t trackPosPt = mctrackPos->Pt();
3829 Double_t trackPosEta = mctrackPos->Eta();
3831 Double_t vK0sSecContinCone[3] = {jetPt, trackPosPt, trackPosEta};
3832 fhnK0sSecContinCone->Fill(vK0sSecContinCone);} //if it's the case, fill jet pt, daughter track pt and track eta in histo
3834 if(particleLabel == negAssLabel){
3835 AliAODMCParticle* mctrackNeg = dynamic_cast<AliAODMCParticle*>(stackMC->At(negAssLabel));
3836 if(!mctrackNeg) continue;
3837 Double_t trackNegPt = mctrackNeg->Pt();
3838 Double_t trackNegEta = mctrackNeg->Eta();
3840 Double_t vK0sSecContinCone[3] = {jetPt, trackNegPt, trackNegEta};
3841 fhnK0sSecContinCone->Fill(vK0sSecContinCone);} //if it's the case, fill jet pt, daughter track pt and track eta in histo
3849 } //end rec-K0-in-cone loop
3851 //________________________________________________________________________________________________________________________________________________________
3853 delete fListMCgenK0sCone;
3858 delete jetConeK0list;
3859 delete jetPerpConeK0list;
3860 delete jetPerpConeLalist;
3861 delete jetPerpConeALalist;
3864 //---------------La--------------------------------------------------------------------------------------------------------------------------------------------
3866 // fQAJetHistosRecCutsLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
3868 for(Int_t it=0; it<fListLa->GetSize(); ++it){ // loop all La
3870 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListLa->At(it));
3875 TVector3 v0MomVect(v0Mom);
3877 Double_t dPhiJetLa = (jet->MomentumVector()->Vect()).DeltaPhi(v0MomVect);
3882 CalculateInvMass(v0, kLambda, invMLa, trackPt); //function to calculate invMass with TLorentzVector class
3883 // Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3885 //if(incrementJetPt){fh1V0JetPt->Fill(jetPt);}
3887 //fFFHistosIMLaJet->FillFF(trackPt, invMLa, jetPt, incrementJetPt);
3889 if(dPhiJetLa<fh1dPhiJetLa->GetXaxis()->GetXmin()) dPhiJetLa += 2*TMath::Pi();
3890 fh1dPhiJetLa->Fill(dPhiJetLa);
3893 /* if(fListLa->GetSize() == 0){ // no La: increment jet pt spectrum
3895 Bool_t incrementJetPt = kTRUE;
3896 fFFHistosIMLaJet->FillFF(-1, -1, jetPt, incrementJetPt);
3900 // ____fetch rec. Lambdas in cone around jet axis_______________________________________________________________________________________
3902 TList* jetConeLalist = new TList();
3903 Double_t sumPtLa = 0.;
3904 Bool_t isBadJetLa = kFALSE; // dummy, do not use
3906 GetTracksInCone(fListLa, jetConeLalist, jet, GetFFRadius(), sumPtLa, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetLa);//method inherited from FF
3908 if(fDebug>2)Printf("%s:%d nLa total: %d, in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nLa,jetConeLalist->GetEntries(),GetFFRadius());
3910 for(Int_t it=0; it<jetConeLalist->GetSize(); ++it){ // loop La in jet cone
3912 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetConeLalist->At(it));
3918 Bool_t daughtercheck = DaughterTrackCheck(v0, nnum, pnum);
3919 if(daughtercheck == kFALSE)continue;
3926 CalculateInvMass(v0, kLambda, invMLa, trackPt); //function to calculate invMass with TLorentzVector class
3928 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;//needed for all histos, which serve for normalisation
3932 Int_t negDaughterpdg;
3933 Int_t posDaughterpdg;
3936 Double_t jetPtSmear = -1;
3938 Bool_t fPhysicalPrimary = -1;
3939 Int_t MCv0PDGCode =0;
3940 TString generatorName;
3942 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
3943 if(incrementJetPt == kTRUE){fh1IMLaConeSmear->Fill(jetPtSmear);
3945 const AliAODTrack *trackNeg=(AliAODTrack *)(v0->GetDaughter(nnum));
3946 const AliAODTrack *trackPos=(AliAODTrack *)(v0->GetDaughter(pnum));
3948 TList *listmc = fAOD->GetList();
3950 Bool_t mclabelcheck = MCLabelCheck(v0, kLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PDGCode, generatorName, isinjected);
3951 if(mclabelcheck == kFALSE)continue;
3953 //std::cout<<"generatorName: "<<generatorName<<std::endl;
3955 if(generatorName == "Hijing"){
3956 Double_t vrecMCHijingLaCone[4] = {jetPt, invMLa,trackPt,fEta};
3957 fhnrecMCHijingLaCone->Fill(vrecMCHijingLaCone);
3960 if(isinjected == kTRUE){
3961 Double_t vrecMCInjectLaCone[4] = {jetPt, invMLa,trackPt,fEta};
3962 fhnrecMCInjectLaCone->Fill(vrecMCInjectLaCone);
3965 }//fill TH1F for normalization purposes
3966 }//end MC analysis part
3968 if(incrementJetPt==kTRUE){
3969 fh1IMLaCone->Fill(jetPt);}//normalisation by number of selected jets
3971 //fFFHistosIMLaCone->FillFF(trackPt, invMLa, jetPt, incrementJetPt);
3972 Double_t vLaCone[4] = {jetPt, invMLa,trackPt,fEta};
3973 fhnLaCone->Fill(vLaCone);
3976 if(jetConeLalist->GetSize() == 0){ // no La: increment jet pt spectrum
3978 Bool_t incrementJetPt = kTRUE;
3979 // fFFHistosIMLaCone->FillFF(-1, -1, jetPt, incrementJetPt);
3980 Double_t vLaCone[4] = {jetPt, -1, -1, -1};
3981 fhnLaCone->Fill(vLaCone);
3983 if(incrementJetPt==kTRUE){
3984 fh1IMLaCone->Fill(jetPt);}//normalisation by number of selected jets
3987 Double_t jetPtSmear;
3988 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
3989 if(incrementJetPt == kTRUE){
3990 fh1IMLaConeSmear->Fill(jetPtSmear);
3999 //____fetch MC generated Lambdas in cone around jet axis__(particles can stem from fragmentation but also from underlying event)_____________
4001 Double_t sumPtMCgenLa = 0.;
4002 Bool_t isBadJetMCgenLa = kFALSE; // dummy, do not use
4004 //sampling MC gen. Lambdas in cone around reconstructed jet axis
4005 fListMCgenLaCone = new TList();
4007 GetTracksInCone(fListMCgenLa, fListMCgenLaCone, jet, GetFFRadius(), sumPtMCgenLa, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetMCgenLa);//fetch MC generated Lambdas in cone of resolution parameter R around jet axis
4009 if(fDebug>2)Printf("%s:%d nMCgenLa in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,fListMCgenLaCone->GetEntries(),GetFFRadius());
4011 /* for(Int_t it=0; it<fListMCgenLaCone->GetSize(); ++it){ // loop MC generated La in cone around jet axis
4013 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenLaCone->At(it));
4016 //Double_t fRapMCgenLa = MyRapidity(mcp0->E(),mcp0->Pz());
4017 Double_t fEtaMCgenLa = mcp0->Eta();
4018 Double_t fPtMCgenLa = mcp0->Pt();
4020 // fh2MCgenLaCone->Fill(jetPt,fPtMCgenLa);
4021 //fh2MCEtagenLaCone->Fill(jetPt,fEtaMCgenLa);
4025 //check whether the reconstructed La are stemming from MC gen La on fListMCgenLa List:__________________________________________________
4027 for(Int_t ic=0; ic<jetConeLalist->GetSize(); ++ic){//loop over all reconstructed La within jet cone, new definition
4029 Int_t negDaughterpdg;
4030 Int_t posDaughterpdg;
4033 Double_t fPtMCrecLaMatch;
4034 Double_t invMLaMatch;
4038 Bool_t fPhysicalPrimary = -1;
4039 Int_t MCv0PDGCode =0;
4040 Double_t jetPtSmear = -1;
4041 TString generatorName;
4043 AliAODv0* v0c = dynamic_cast<AliAODv0*>(jetConeLalist->At(ic));//new definition
4048 Bool_t daughtercheck = DaughterTrackCheck(v0c, nnum, pnum);
4049 if(daughtercheck == kFALSE)continue;
4051 const AliAODTrack *trackMCNeg=(AliAODTrack *)(v0c->GetDaughter(nnum));
4052 const AliAODTrack *trackMCPos=(AliAODTrack *)(v0c->GetDaughter(pnum));
4054 TList *listmc = fAOD->GetList();
4056 Bool_t mclabelcheck = MCLabelCheck(v0c, kLambda, trackMCNeg, trackMCPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PDGCode, generatorName, isinjected);
4058 if(mclabelcheck == kFALSE)continue;
4059 if(fPhysicalPrimary == kFALSE)continue;
4061 for(Int_t it=0; it<fListMCgenLa->GetSize(); ++it){//new definition // loop over MC generated K0s in cone around jet axis
4064 //Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4066 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenLa->At(it));//new definition
4067 //AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenLaCone->At(it));//old definition
4071 Bool_t particleMatching = IsParticleMatching(mcp0, v0Label);
4074 if(particleMatching == kFALSE)continue; //particle doesn't match on any associated MC gen particle in cone around rec jet axis
4076 CalculateInvMass(v0c, kLambda, invMLaMatch, fPtMCrecLaMatch);
4078 Double_t fPtMCgenLa = mcp0->Pt();
4079 Double_t fEta = v0c->Eta();//rec. MC particle
4080 Double_t vMCrecLaCone[4] = {jetPt, invMLaMatch,fPtMCgenLa,fEta};
4081 fhnMCrecLaCone->Fill(vMCrecLaCone);
4083 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
4085 Double_t vMCrecLaConeSmear[4] = {jetPtSmear, invMLaMatch,fPtMCgenLa,fEta};
4086 fhnMCrecLaConeSmear->Fill(vMCrecLaConeSmear); //fill matching rec. Lambdas in 3D histogram, jet pT smeared according to deltaptjet distribution width
4089 } // end MCgenLa loop
4091 //check the Lambda daughters contamination of the jet tracks://///////////////////////////////////////////////////////////////////////////////////////////
4093 TClonesArray *stackMC = 0x0;
4095 for(Int_t it=0; it<jettracklist->GetSize(); ++it){//loop over all particles in the jet
4097 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));//all tracks in jet cone
4098 if(!trackVP)continue;
4099 AliAODTrack *tr = dynamic_cast<AliAODTrack*> (trackVP); //fetch one jet track from the TList
4102 //get MC label information
4103 TList *mclist = fAOD->GetList(); //fetch the MC stack
4105 stackMC = (TClonesArray*)mclist->FindObject(AliAODMCParticle::StdBranchName()); //get MCAOD branch in data
4106 if (!stackMC) {Printf("ERROR: stack not available");}
4109 Int_t particleLabel = TMath::Abs(tr->GetLabel()); //fetch jet track label in MC stack
4111 Bool_t daughterchecks = DaughterTrackCheck(v0c, nnum, pnum);
4113 if(daughterchecks == kFALSE)continue; //make sure that daughters are properly ordered
4115 const AliAODTrack *trackNeg=(AliAODTrack *)(v0c->GetDaughter(nnum)); //fetch v0 daughters of reconstructed K0s
4116 const AliAODTrack *trackPos=(AliAODTrack *)(v0c->GetDaughter(pnum));
4118 Int_t negAssLabel = TMath::Abs(trackNeg->GetLabel()); //negative (reconstructed) charged track label in MC stack
4119 Int_t posAssLabel = TMath::Abs(trackPos->GetLabel()); //positive (reconstructed) charged track label in MC stack
4122 if(particleLabel == posAssLabel){ //check whether jet track and each of the rec. K0s daughters have same MC stack label -> are identical
4124 AliAODMCParticle* mctrackPos = dynamic_cast<AliAODMCParticle*>(stackMC->At(posAssLabel));
4125 if(!mctrackPos) continue;
4127 Double_t trackPosPt = trackPos->Pt();
4128 Double_t trackPosEta = trackPos->Eta();
4129 Double_t vLaSecContinCone[3] = {jetPt, trackPosPt, trackPosEta};
4130 fhnLaSecContinCone->Fill(vLaSecContinCone);
4132 } //if it's the case, fill jet pt, daughter track pt and track eta in histo
4135 if(particleLabel == negAssLabel){
4137 AliAODMCParticle* mctrackNeg = dynamic_cast<AliAODMCParticle*>(stackMC->At(negAssLabel));
4138 if(!mctrackNeg) continue;
4140 Double_t trackNegPt = trackNeg->Pt();
4141 Double_t trackNegEta = trackNeg->Eta();
4143 Double_t vLaSecContinCone[3] = {jetPt, trackNegPt, trackNegEta};
4144 fhnLaSecContinCone->Fill(vLaSecContinCone);
4147 } //if it's the case, fill jet pt, daughter track pt and track eta in histo
4152 } //end rec-La-in-cone loop
4153 //________________________________________________________________________________________________________________________________________________________
4155 delete fListMCgenLaCone;
4159 delete jetConeLalist;
4163 //---------------ALa-----------
4166 // fQAJetHistosRecCutsLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
4168 for(Int_t it=0; it<fListALa->GetSize(); ++it){ // loop all ALa
4170 AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListALa->At(it));
4175 TVector3 v0MomVect(v0Mom);
4177 Double_t dPhiJetALa = (jet->MomentumVector()->Vect()).DeltaPhi(v0MomVect);
4179 Double_t invMALa =0;
4182 CalculateInvMass(v0, kAntiLambda, invMALa, trackPt); //function to calculate invMass with TLorentzVector class
4183 //Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4185 //if(incrementJetPt){fh1V0JetPt->Fill(jetPt);}
4187 //fFFHistosIMALaJet->FillFF(trackPt, invMALa, jetPt, incrementJetPt);
4189 if(dPhiJetALa<fh1dPhiJetALa->GetXaxis()->GetXmin()) dPhiJetALa += 2*TMath::Pi();
4190 fh1dPhiJetALa->Fill(dPhiJetALa);
4193 // if(fListALa->GetSize() == 0){ // no ALa: increment jet pt spectrum
4195 // Bool_t incrementJetPt = kTRUE;
4196 //fFFHistosIMALaJet->FillFF(-1, -1, jetPt, incrementJetPt);
4200 // ____fetch rec. Antilambdas in cone around jet axis_______________________________________________________________________________________
4202 TList* jetConeALalist = new TList();
4203 Double_t sumPtALa = 0.;
4204 Bool_t isBadJetALa = kFALSE; // dummy, do not use
4206 GetTracksInCone(fListALa, jetConeALalist, jet, GetFFRadius(), sumPtALa, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetALa);//method inherited from FF
4208 if(fDebug>2)Printf("%s:%d nALa total: %d, in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nALa,jetConeALalist->GetEntries(),GetFFRadius());
4210 for(Int_t it=0; it<jetConeALalist->GetSize(); ++it){ // loop ALa in jet cone
4212 AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetConeALalist->At(it));
4219 Bool_t daughtercheck = DaughterTrackCheck(v0, nnum, pnum);
4220 if(daughtercheck == kFALSE)continue;
4223 Double_t invMALa =0;
4229 CalculateInvMass(v0, kAntiLambda, invMALa, trackPt); //function to calculate invMass with TLorentzVector class
4231 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4233 if(fAnalysisMC){ //jet pt smearing study for Antilambdas
4235 Int_t negDaughterpdg;
4236 Int_t posDaughterpdg;
4239 Double_t jetPtSmear = -1;
4241 Bool_t fPhysicalPrimary = -1;
4242 Int_t MCv0PDGCode =0;
4243 TString generatorName;
4245 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
4246 const AliAODTrack *trackNeg=(AliAODTrack *)(v0->GetDaughter(nnum));
4247 const AliAODTrack *trackPos=(AliAODTrack *)(v0->GetDaughter(pnum));
4249 TList *listmc = fAOD->GetList();
4251 Bool_t mclabelcheck = MCLabelCheck(v0, kAntiLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PDGCode, generatorName, isinjected);
4252 if(mclabelcheck == kFALSE)continue;
4254 //std::cout<<"generatorName: "<<generatorName<<std::endl;
4256 if(generatorName == "Hijing"){
4257 Double_t vrecMCHijingALaCone[4] = {jetPt, invMALa,trackPt,fEta};
4258 fhnrecMCHijingALaCone->Fill(vrecMCHijingALaCone);
4261 if(isinjected == kTRUE){
4262 Double_t vrecMCInjectALaCone[4] = {jetPt, invMALa,trackPt,fEta};
4263 fhnrecMCInjectALaCone->Fill(vrecMCInjectALaCone);
4266 if(incrementJetPt == kTRUE){fh1IMALaConeSmear->Fill(jetPtSmear);} //fill TH1F for normalization purposes
4269 if(incrementJetPt==kTRUE){
4270 fh1IMALaCone->Fill(jetPt);}//normalisation by number of selected jets
4272 //fFFHistosIMALaCone->FillFF(trackPt, invMALa, jetPt, incrementJetPt);
4273 Double_t vALaCone[4] = {jetPt, invMALa,trackPt,fEta};
4274 fhnALaCone->Fill(vALaCone);
4277 if(jetConeALalist->GetSize() == 0){ // no ALa: increment jet pt spectrum
4279 Bool_t incrementJetPt = kTRUE;
4281 if(incrementJetPt==kTRUE){
4282 fh1IMALaCone->Fill(jetPt);}//normalisation by number of selected jets
4284 //fFFHistosIMALaCone->FillFF(-1, -1, jetPt, incrementJetPt);
4285 Double_t vALaCone[4] = {jetPt, -1, -1, -1};
4286 fhnALaCone->Fill(vALaCone);
4289 Double_t jetPtSmear;
4290 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
4291 if(incrementJetPt == kTRUE)fh1IMALaConeSmear->Fill(jetPtSmear);}
4297 //____fetch MC generated Antilambdas in cone around jet axis__(particles can stem from fragmentation but also from underlying event)_____________
4299 Double_t sumPtMCgenALa = 0.;
4300 Bool_t isBadJetMCgenALa = kFALSE; // dummy, do not use
4302 //sampling MC gen Antilambdas in cone around reconstructed jet axis
4303 fListMCgenALaCone = new TList();
4305 GetTracksInCone(fListMCgenALa, fListMCgenALaCone, jet, GetFFRadius(), sumPtMCgenALa, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetMCgenALa);//MC generated K0s in cone around jet axis
4307 if(fDebug>2)Printf("%s:%d nMCgenALa in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,fListMCgenALaCone->GetEntries(),GetFFRadius());
4309 /* for(Int_t it=0; it<fListMCgenALaCone->GetSize(); ++it){ // loop MC generated La in cone around jet axis
4311 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenALaCone->At(it));
4314 //Double_t fRapMCgenALa = MyRapidity(mcp0->E(),mcp0->Pz());
4315 Double_t fEtaMCgenALa = mcp0->Eta();
4316 Double_t fPtMCgenALa = mcp0->Pt();
4318 //fh2MCgenALaCone->Fill(jetPt,fPtMCgenALa);
4319 //fh2MCEtagenALaCone->Fill(jetPt,fEtaMCgenALa);
4323 //check whether the reconstructed ALa are stemming from MC gen ALa on MCgenALa List:__________________________________________________
4325 for(Int_t ic=0; ic<jetConeALalist->GetSize(); ++ic){//loop over all reconstructed ALa
4327 Int_t negDaughterpdg;
4328 Int_t posDaughterpdg;
4331 Double_t fPtMCrecALaMatch;
4332 Double_t invMALaMatch;
4336 Bool_t fPhysicalPrimary = -1;
4337 Int_t MCv0PDGCode =0;
4338 Double_t jetPtSmear = -1;
4339 TString generatorName;
4341 AliAODv0* v0c = dynamic_cast<AliAODv0*>(jetConeALalist->At(ic));
4344 Bool_t daughtercheck = DaughterTrackCheck(v0c, nnum, pnum);
4345 if(daughtercheck == kFALSE)continue;
4347 const AliAODTrack *trackMCNeg=(AliAODTrack *)(v0c->GetDaughter(nnum));
4348 const AliAODTrack *trackMCPos=(AliAODTrack *)(v0c->GetDaughter(pnum));
4350 TList *listmc = fAOD->GetList();
4351 if(!listmc)continue;
4353 Bool_t mclabelcheck = MCLabelCheck(v0c, kAntiLambda, trackMCNeg, trackMCPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PDGCode, generatorName, isinjected);
4355 if(mclabelcheck == kFALSE)continue;
4356 if(fPhysicalPrimary == kFALSE)continue;
4358 for(Int_t it=0; it<fListMCgenALa->GetSize(); ++it){ // loop over MC generated Antilambdas in cone around jet axis
4360 //Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
4362 AliAODMCParticle* mcp0 = dynamic_cast<AliAODMCParticle*>(fListMCgenALa->At(it));
4365 Bool_t particleMatching = IsParticleMatching(mcp0, v0Label);
4367 if(particleMatching == kFALSE)continue; //particle doesn't match on any associated MC gen particle in cone around rec jet axis
4369 CalculateInvMass(v0c, kAntiLambda, invMALaMatch, fPtMCrecALaMatch);
4371 Double_t fPtMCgenALa = mcp0->Pt();
4372 Double_t fEta = v0c->Eta();
4373 Double_t vMCrecALaCone[4] = {jetPt, invMALaMatch,fPtMCgenALa,fEta};
4374 fhnMCrecALaCone->Fill(vMCrecALaCone); //fill matching rec. Antilambda in 3D histogram
4376 SmearJetPt(jetPt,cl,GetFFRadius(),GetFFMinLTrackPt(),jetPtSmear);
4378 Double_t vMCrecALaConeSmear[4] = {jetPtSmear, invMALaMatch,fPtMCgenALa,fEta};
4379 fhnMCrecALaConeSmear->Fill(vMCrecALaConeSmear); //fill matching rec. Antilambda in 3D histogram
4381 } // end MCgenALa loop
4385 //check the Antilambda daughters contamination of the jet tracks:
4387 TClonesArray *stackMC = 0x0;
4389 for(Int_t it=0; it<jettracklist->GetSize(); ++it){//loop over all particles in the jet
4391 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));//all tracks in jet cone
4392 if(!trackVP)continue;
4393 AliAODTrack *tr = dynamic_cast<AliAODTrack*> (trackVP); //fetch one jet track from the TList
4396 //get MC label information
4397 TList *mclist = fAOD->GetList(); //fetch the MC stack
4398 if(!mclist)continue;
4400 stackMC = (TClonesArray*)mclist->FindObject(AliAODMCParticle::StdBranchName()); //get MCAOD branch in data
4401 if (!stackMC) {Printf("ERROR: stack not available");}
4404 Int_t particleLabel = TMath::Abs(tr->GetLabel()); //fetch jet track label in MC stack
4406 Bool_t daughterchecks = DaughterTrackCheck(v0c, nnum, pnum);
4408 if(daughterchecks == kFALSE)continue; //make sure that daughters are properly ordered
4410 const AliAODTrack *trackNeg=(AliAODTrack *)(v0c->GetDaughter(nnum)); //fetch v0 daughters of reconstructed K0s
4411 const AliAODTrack *trackPos=(AliAODTrack *)(v0c->GetDaughter(pnum));
4412 if(!trackPos)continue;
4413 if(!trackNeg)continue;
4415 Int_t negAssLabel = TMath::Abs(trackNeg->GetLabel()); //negative (reconstructed) charged track label in MC stack
4416 Int_t posAssLabel = TMath::Abs(trackPos->GetLabel()); //positive (reconstructed) charged track label in MC stack
4418 if(!negAssLabel)continue;
4419 if(!posAssLabel)continue;
4421 if(particleLabel == posAssLabel){ //check whether jet track and each of the rec. K0s daughters have same MC stack label -> are identical
4422 AliAODMCParticle* mctrackPos = dynamic_cast<AliAODMCParticle*>(stackMC->At(posAssLabel));
4423 if(!mctrackPos) continue;
4425 Double_t trackPosPt = trackPos->Pt();
4426 Double_t trackPosEta = trackPos->Eta();
4427 if(!trackPosPt)continue;
4428 if(!trackPosEta)continue;
4430 Double_t vLaSecContinCone[3] = {jetPt, trackPosPt, trackPosEta};
4431 fhnLaSecContinCone->Fill(vLaSecContinCone);
4435 //fh3SecContinCone->Fill(jetPt, trackPosPt, trackPosEta);
4436 } //if it's the case, fill jet pt, daughter track pt and track eta in histo
4438 if(particleLabel == negAssLabel){
4440 AliAODMCParticle* mctrackNeg = dynamic_cast<AliAODMCParticle*>(stackMC->At(negAssLabel));
4441 if(!mctrackNeg) continue;
4443 Double_t trackNegPt = trackNeg->Pt();
4444 Double_t trackNegEta = trackNeg->Eta();
4446 if(!trackNegPt)continue;
4447 if(!trackNegEta)continue;
4449 Double_t vLaSecContinCone[3] = {jetPt, trackNegPt, trackNegEta};
4450 fhnLaSecContinCone->Fill(vLaSecContinCone);
4452 //fh3SecContinCone->Fill(jetPt, trackNegPt, trackNegEta);
4453 } //if it's the case, fill jet pt, daughter track pt and track eta in histo
4457 } //end rec-ALa-in-cone loop
4458 //________________________________________________________________________________________________________________________________________________________
4460 delete fListMCgenALaCone;
4464 delete jetConeALalist;
4465 delete jettracklist; //had been initialised at jet loop beginning
4467 }//end of if 'leading' or 'all jet' requirement
4473 fTracksRecCuts->Clear();
4474 fJetsRecCuts->Clear();
4475 fBckgJetsRec->Clear();
4479 fListFeeddownLaCand->Clear();
4480 fListFeeddownALaCand->Clear();
4481 jetConeFDLalist->Clear();
4482 jetConeFDALalist->Clear();
4483 fListMCgenK0s->Clear();
4484 fListMCgenLa->Clear();
4485 fListMCgenALa->Clear();
4489 PostData(1, fCommonHistList);
4494 // ____________________________________________________________________________________________
4495 void AliAnalysisTaskJetChem::SetProperties(TH3F* h,const char* x, const char* y, const char* z)
4497 //Set properties of histos (x,y and z title)
4502 h->GetXaxis()->SetTitleColor(1);
4503 h->GetYaxis()->SetTitleColor(1);
4504 h->GetZaxis()->SetTitleColor(1);
4508 //________________________________________________________________________________________________________________________________________
4509 Bool_t AliAnalysisTaskJetChem::AcceptBetheBloch(AliAODv0 *v0, AliPIDResponse *PIDResponse, const Int_t particletype) //dont use for MC Analysis
4515 const AliAODTrack *ntracktest=(AliAODTrack *)v0->GetDaughter(nnum);
4516 if(ntracktest->Charge() > 0){nnum = 0; pnum = 1;}
4518 const AliAODTrack *trackNeg=(AliAODTrack *)(v0->GetDaughter(nnum));
4519 const AliAODTrack *trackPos=(AliAODTrack *)(v0->GetDaughter(pnum));
4521 //Check if both tracks are available
4522 if (!trackPos || !trackNeg) {
4523 Printf("strange analysis::UserExec:: Error:Could not retrieve one of the daughter tracks\n");
4527 //remove like sign V0s
4528 if ( trackPos->Charge() == trackNeg->Charge() ){
4529 //if(fDebug>1) Printf("%s:%d found like-sign V0", (char*)__FILE__,__LINE__);
4534 Double_t nsig_p = 0; //number of sigmas that positive daughter track has got in TPC pid information
4535 Double_t nsig_n = 0;
4537 const AliAODPid *pid_p=trackPos->GetDetPid(); // returns fDetPID, more detailed or detector specific pid information
4538 const AliAODPid *pid_n=trackNeg->GetDetPid();
4540 if(!pid_p)return kFALSE;
4541 if(!pid_n)return kFALSE;
4545 if(particletype == 1) //PID cut on positive charged Lambda daughters (only those with pt < 1 GeV/c)
4548 nsig_p=PIDResponse->NumberOfSigmasTPC(trackPos,AliPID::kProton);
4549 Double_t protonPt = trackPos->Pt();
4550 if ((TMath::Abs(nsig_p) >= fCutBetheBloch) && (fCutBetheBloch >0) && (protonPt < 1)) return kFALSE;
4559 if(particletype == 2)
4561 nsig_n=PIDResponse->NumberOfSigmasTPC(trackNeg,AliPID::kProton);
4562 Double_t antiprotonPt = trackNeg->Pt();
4563 if ((TMath::Abs(nsig_n) >= fCutBetheBloch) && (fCutBetheBloch >0) && (antiprotonPt < 1)) return kFALSE;
4571 //___________________________________________________________________
4572 Bool_t AliAnalysisTaskJetChem::IsK0InvMass(const Double_t mass) const
4574 // K0 mass ? Use FF histo limits
4576 if(fFFIMInvMMin <= mass && mass < fFFIMInvMMax) return kTRUE;
4580 //___________________________________________________________________
4581 Bool_t AliAnalysisTaskJetChem::IsLaInvMass(const Double_t mass) const
4583 // La mass ? Use FF histo limits
4586 if(fFFIMLaInvMMin <= mass && mass < fFFIMLaInvMMax) return kTRUE;
4591 //_____________________________________________________________________________________
4592 Int_t AliAnalysisTaskJetChem::GetListOfV0s(TList *list, const Int_t type, const Int_t particletype, AliAODVertex* primVertex, AliAODEvent* aod)
4594 // fill list of V0s selected according to type
4597 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
4602 if(fDebug>5){std::cout<<"AliAnalysisTaskJetChem::GetListOfV0s(): type: "<<type<<" particletype: "<<particletype<<"aod: "<<aod<<std::endl;
4603 if(type==kTrackUndef){std::cout<<"AliAnalysisTaskJetChem::GetListOfV0s(): kTrackUndef!! "<<std::endl;}
4607 if(type==kTrackUndef) return 0;
4609 if(!primVertex) return 0;
4611 Double_t lPrimaryVtxPosition[3];
4612 Double_t lV0Position[3];
4613 lPrimaryVtxPosition[0] = primVertex->GetX();
4614 lPrimaryVtxPosition[1] = primVertex->GetY();
4615 lPrimaryVtxPosition[2] = primVertex->GetZ();
4617 if(fDebug>5){ std::cout<<"AliAnalysisTaskJetChem::GetListOfV0s(): aod->GetNumberOfV0s: "<<aod->GetNumberOfV0s()<<std::endl; }
4620 for(int i=0; i<aod->GetNumberOfV0s(); i++){ // loop over V0s
4623 AliAODv0* v0 = aod->GetV0(i);
4627 std::cout << std::endl
4628 << "Warning in AliAnalysisTaskJetChem::GetListOfV0s:" << std::endl
4629 << "v0 = " << v0 << std::endl;
4633 Bool_t isOnFly = v0->GetOnFlyStatus();
4635 if(!isOnFly && (type == kOnFly || type == kOnFlyPID || type == kOnFlydEdx || type == kOnFlyPrim)) continue;
4636 if( isOnFly && (type == kOffl || type == kOfflPID || type == kOffldEdx || type == kOfflPrim)) continue;
4638 Int_t motherType = -1;
4639 //Double_t v0CalcMass = 0; //mass of MC v0
4640 Double_t MCPt = 0; //pt of MC v0
4642 Double_t pp[3]={0,0,0}; //3-momentum positive charged track
4643 Double_t pm[3]={0,0,0}; //3-momentum negative charged track
4644 Double_t v0mom[3]={0,0,0};
4655 Bool_t daughtercheck = DaughterTrackCheck(v0, nnum, pnum);
4657 if(daughtercheck == kFALSE)continue;
4659 const AliAODTrack *trackNeg=(AliAODTrack *)(v0->GetDaughter(nnum));
4660 const AliAODTrack *trackPos=(AliAODTrack *)(v0->GetDaughter(pnum));
4663 ///////////////////////////////////////////////////////////////////////////////////
4665 //calculate InvMass for every V0 particle assumption (Kaon=1,Lambda=2,Antilambda=3)
4666 switch(particletype){
4668 CalculateInvMass(v0, kK0, invM, trackPt); //function to calculate invMass with TLorentzVector class
4672 CalculateInvMass(v0, kLambda, invM, trackPt);
4676 CalculateInvMass(v0, kAntiLambda, invM, trackPt);
4680 std::cout<<"***NO VALID PARTICLETYPE***"<<std::endl;
4685 /////////////////////////////////////////////////////////////
4686 //V0 and track Cuts:
4689 if(fDebug>7){if(!(IsK0InvMass(invMK0s)) && !(IsLaInvMass(invMLa)) && !(IsLaInvMass(invMALa))){std::cout<<"AliAnalysisTaskJetChem::GetListOfV0s: invM not in selected mass window "<<std::endl;}}
4691 if(!(IsK0InvMass(invMK0s)) && !(IsLaInvMass(invMLa)) && !(IsLaInvMass(invMALa)))continue;
4693 // Double_t PosEta = trackPos->AliAODTrack::Eta();//daughter track charge is sometimes wrong here, account for that!!!
4694 // Double_t NegEta = trackNeg->AliAODTrack::Eta();
4696 Double_t PosEta = trackPos->Eta();//daughter track charge is sometimes wrong here, account for that!!!
4697 Double_t NegEta = trackNeg->Eta();
4699 Double_t PosCharge = trackPos->Charge();
4700 Double_t NegCharge = trackNeg->Charge();
4702 if((trackPos->Charge() == 1) && (trackNeg->Charge() == -1)) //Fill daughters charge into histo to check if they are symmetric distributed
4703 { fh1PosDaughterCharge->Fill(PosCharge);
4704 fh1NegDaughterCharge->Fill(NegCharge);
4707 //DistOverTotMom_in_2D___________
4709 Float_t fMassK0s = TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass();
4710 Float_t fMassLambda = TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass();
4713 AliAODVertex* primVtx = fAOD->GetPrimaryVertex(); // get the primary vertex
4714 Double_t dPrimVtxPos[3]; // primary vertex position {x,y,z}
4715 primVtx->GetXYZ(dPrimVtxPos);
4717 Float_t fPtV0 = TMath::Sqrt(v0->Pt2V0()); // transverse momentum of V0
4718 Double_t dSecVtxPos[3]; // V0 vertex position {x,y,z}
4719 v0->GetSecondaryVtx(dSecVtxPos);
4720 Double_t dDecayPath[3];
4721 for (Int_t iPos = 0; iPos < 3; iPos++)
4722 dDecayPath[iPos] = dSecVtxPos[iPos]-dPrimVtxPos[iPos]; // vector of the V0 path
4723 Float_t fDecLen2D = TMath::Sqrt(dDecayPath[0]*dDecayPath[0]+dDecayPath[1]*dDecayPath[1]); //transverse path length R
4724 Float_t fROverPt = fDecLen2D/fPtV0; // R/pT
4726 Float_t fMROverPtK0s = fMassK0s*fROverPt; // m*R/pT
4727 Float_t fMROverPtLambda = fMassLambda*fROverPt; // m*R/pT
4729 //___________________
4730 //Double_t fRap = -999;//init values
4731 Double_t fEta = -999;
4732 Double_t fV0cosPointAngle = -999;
4733 Double_t fV0DecayLength = v0->DecayLengthV0(lPrimaryVtxPosition);
4737 fV0mom[0]=v0->MomV0X();
4738 fV0mom[1]=v0->MomV0Y();
4739 fV0mom[2]=v0->MomV0Z();
4741 Double_t fV0TotalMomentum = TMath::Sqrt(fV0mom[0]*fV0mom[0]+fV0mom[1]*fV0mom[1]+fV0mom[2]*fV0mom[2]);
4742 // const Double_t K0sPDGmass = 0.497614;
4743 // const Double_t LambdaPDGmass = 1.115683;
4745 const Double_t K0sPDGmass = TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass();
4746 const Double_t LambdaPDGmass = TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass();
4748 Double_t fDistOverTotMomK0s = 0;
4749 Double_t fDistOverTotMomLa = 0;
4751 //calculate proper lifetime of particles in 3D (not recommended anymore)
4753 if(particletype == kK0){
4755 fDistOverTotMomK0s = fV0DecayLength * K0sPDGmass;
4756 fDistOverTotMomK0s /= (fV0TotalMomentum+1e-10);
4759 if((particletype == kLambda)||(particletype == kAntiLambda)){
4761 fDistOverTotMomLa = fV0DecayLength * LambdaPDGmass;
4762 fDistOverTotMomLa /= (fV0TotalMomentum+1e-10);
4765 //TPC cluster (not used anymore) and TPCRefit cuts
4767 //Double_t trackPosNcls = trackPos->GetTPCNcls();//Get number of clusters for positive charged tracks
4768 //Double_t trackNegNcls = trackNeg->GetTPCNcls();//Get number of clusters for negative charged tracks
4770 if(fRequireTPCRefit==kTRUE){//if kTRUE: accept only if daughter track is refitted in TPC!!
4771 Bool_t isPosTPCRefit = (trackPos->AliAODTrack::IsOn(AliESDtrack::kTPCrefit));
4772 Bool_t isNegTPCRefit = (trackNeg->AliAODTrack::IsOn(AliESDtrack::kTPCrefit));
4773 if (!isPosTPCRefit)continue;
4774 if (!isNegTPCRefit)continue;
4777 if(fKinkDaughters==kFALSE){//if kFALSE: no acceptance of kink daughters
4778 AliAODVertex* ProdVtxDaughtersPos = (AliAODVertex*) (trackPos->AliAODTrack::GetProdVertex());
4779 Char_t isAcceptKinkDaughtersPos = ProdVtxDaughtersPos->GetType();
4780 if(isAcceptKinkDaughtersPos==AliAODVertex::kKink)continue;
4782 AliAODVertex* ProdVtxDaughtersNeg = (AliAODVertex*) (trackNeg->AliAODTrack::GetProdVertex());
4783 Char_t isAcceptKinkDaughtersNeg = ProdVtxDaughtersNeg->GetType();
4784 if(isAcceptKinkDaughtersNeg==AliAODVertex::kKink)continue;
4788 Double_t fV0Radius = -999;
4789 Double_t fDcaV0Daughters = v0->DcaV0Daughters();
4790 Double_t fDcaPosToPrimVertex = v0->DcaPosToPrimVertex();//IP of positive charged daughter
4791 Double_t fDcaNegToPrimVertex = v0->DcaNegToPrimVertex();//IP of negative charged daughter
4792 Double_t avDecayLengthK0s = 2.6844;
4793 Double_t avDecayLengthLa = 7.89;
4795 //Float_t fCTauK0s = 2.6844; // [cm] c tau of K0S
4796 //Float_t fCTauLambda = 7.89; // [cm] c tau of Lambda and Antilambda
4798 fV0cosPointAngle = v0->CosPointingAngle(lPrimaryVtxPosition);
4799 lV0Position[0]= v0->DecayVertexV0X();
4800 lV0Position[1]= v0->DecayVertexV0Y();
4801 lV0Position[2]= v0->DecayVertexV0Z();
4803 fV0Radius = TMath::Sqrt(lV0Position[0]*lV0Position[0]+lV0Position[1]*lV0Position[1]);
4805 if(particletype == kK0) {//fRap = v0->RapK0Short();
4806 fEta = v0->PseudoRapV0();}
4807 if(particletype == kLambda) {//fRap = v0->RapLambda();
4808 fEta = v0->PseudoRapV0();}
4809 if(particletype == kAntiLambda) {//fRap = v0->Y(-3122);
4810 fEta = v0->PseudoRapV0();}
4813 //cut on 3D DistOverTotMom: (not used anymore)
4814 //if((particletype == kLambda)||(particletype == kAntiLambda)){if(fDistOverTotMomLa >= (fCutV0DecayMax * avDecayLengthLa)) continue;}
4816 //cut on K0s applied below after all other cuts for histo fill purposes..
4818 //cut on 2D DistOverTransMom: (recommended from Iouri)
4819 if((particletype == kLambda)||(particletype == kAntiLambda)){if(fMROverPtLambda > (fCutV0DecayMax * avDecayLengthLa))continue;}//fCutV0DecayMax set to 5 in AddTask macro
4821 //Armenteros Podolanski Plot for K0s:////////////////////////////
4823 Double_t ArmenterosAlpha=-999;
4824 Double_t ArmenterosPt=-999;
4830 if(particletype == kK0){
4832 pp[0]=v0->MomPosX();
4833 pp[1]=v0->MomPosY();
4834 pp[2]=v0->MomPosZ();
4836 pm[0]=v0->MomNegX();
4837 pm[1]=v0->MomNegY();
4838 pm[2]=v0->MomNegZ();
4841 v0mom[0]=v0->MomV0X();
4842 v0mom[1]=v0->MomV0Y();
4843 v0mom[2]=v0->MomV0Z();
4845 TVector3 v0Pos(pp[0],pp[1],pp[2]);
4846 TVector3 v0Neg(pm[0],pm[1],pm[2]);
4847 TVector3 v0totMom(v0mom[0], v0mom[1], v0mom[2]); //vector for tot v0 momentum
4849 //PosPt = v0Pos.Perp(v0totMom); //longitudinal momentum of positive charged daughter track
4850 PosPl = v0Pos.Dot(v0totMom)/v0totMom.Mag(); //transversal momentum of positive charged daughter track
4852 //NegPt = v0Neg.Perp(v0totMom); //longitudinal momentum of negative charged daughter track
4853 NegPl = v0Neg.Dot(v0totMom)/v0totMom.Mag(); //transversal momentum of nergative charged daughter track
4855 ArmenterosAlpha = 1.-2./(1+(PosPl/NegPl));
4856 ArmenterosPt= v0->PtArmV0();
4860 if(particletype == kK0){//only cut on K0s histos
4861 if(IsArmenterosSelected == 1){// Armenteros Cut to reject Lambdas contamination in K0s inv. massspectrum
4862 fh2ArmenterosBeforeCuts->Fill(ArmenterosAlpha,ArmenterosPt);
4866 //some more cuts on v0s and daughter tracks:
4869 if((TMath::Abs(PosEta)>fCutPostrackEta) || (TMath::Abs(NegEta)>fCutNegtrackEta))continue; //Daughters pseudorapidity cut
4870 if (fV0cosPointAngle < fCutV0cosPointAngle) continue; //cospointangle cut
4872 //if(TMath::Abs(fRap) > fCutRap)continue; //V0 Rapidity Cut
4873 if(TMath::Abs(fEta) > fCutEta) continue; //V0 Eta Cut
4874 if (fDcaV0Daughters > fCutDcaV0Daughters)continue;
4875 if ((fDcaPosToPrimVertex < fCutDcaPosToPrimVertex) || (fDcaNegToPrimVertex < fCutDcaNegToPrimVertex))continue;
4876 if ((fV0Radius < fCutV0RadiusMin) || (fV0Radius > fCutV0RadiusMax))continue;
4878 const AliAODPid *pid_p1=trackPos->GetDetPid();
4879 const AliAODPid *pid_n1=trackNeg->GetDetPid();
4882 if(particletype == kLambda){
4883 // if(AcceptBetheBloch(v0, fPIDResponse, 1) == kFALSE){std::cout<<"******PID cut rejects Lambda!!!************"<<std::endl;}
4884 if(AcceptBetheBloch(v0, fPIDResponse, 1) == kFALSE)continue;
4885 fh2BBLaPos->Fill(pid_p1->GetTPCmomentum(),pid_p1->GetTPCsignal());//positive lambda daughter
4886 fh2BBLaNeg->Fill(pid_n1->GetTPCmomentum(),pid_n1->GetTPCsignal());//negative lambda daughter
4888 //Double_t phi = v0->Phi();
4889 //Double_t massLa = v0->MassLambda();
4891 //printf("La: i = %d, m = %f, pT = %f, eta = %f, phi = %f\n, ",i,massLa,trackPt,fEta,phi);
4895 if(particletype == kAntiLambda){
4897 if(AcceptBetheBloch(v0, fPIDResponse, 2) == kFALSE)continue;
4898 fh2BBLaPos->Fill(pid_p1->GetTPCmomentum(),pid_p1->GetTPCsignal());//positive antilambda daughter
4899 fh2BBLaNeg->Fill(pid_n1->GetTPCmomentum(),pid_n1->GetTPCsignal());//negative antilambda daughter
4904 //Armenteros cut on K0s:
4905 if(particletype == kK0){
4906 if(IsArmenterosSelected == 1){// Armenteros Cut to reject Lambdas contamination in K0s inv. massspectrum
4908 if((ArmenterosPt<=(TMath::Abs(fCutArmenteros*ArmenterosAlpha))) && (fCutArmenteros!=-999))continue; //Cuts out Lambda contamination in K0s histos
4909 fh2ArmenterosAfterCuts->Fill(ArmenterosAlpha,ArmenterosPt);
4913 //not used anymore in 3D, z component of total momentum has bad resolution, cut instead in 2D and use pT
4914 //Proper Lifetime Cut: DecayLength3D * PDGmass / |p_tot| < 3*2.68cm (ctau(betagamma=1)) ; |p|/mass = beta*gamma
4915 //////////////////////////////////////////////
4918 //cut on 2D DistOverTransMom
4919 if(particletype == kK0){//the cut on Lambdas you can find above
4921 fh2ProperLifetimeK0sVsPtBeforeCut->Fill(trackPt,fMROverPtK0s); //fill these histos after all other cuts
4922 if(fMROverPtK0s > (fCutV0DecayMax * avDecayLengthK0s))continue;
4923 fh2ProperLifetimeK0sVsPtAfterCut->Fill(trackPt,fMROverPtK0s);
4925 //Double_t phi = v0->Phi();
4926 // Double_t massK0s = v0->MassK0Short();
4927 //printf("K0S: i = %d, m = %f, pT = %f, eta = %f, phi = %f\n",i,invMK0s,trackPt,fEta,phi);
4929 //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;
4932 //MC Associated V0 particles: (reconstructed particles associated with MC truth (MC truth: true primary MC generated particle))
4935 if(fAnalysisMC){// begin MC part
4937 Int_t negDaughterpdg = 0;
4938 Int_t posDaughterpdg = 0;
4940 Bool_t fPhysicalPrimary = -1; //v0 physical primary check
4941 Int_t MCv0PdgCode = 0;
4942 Bool_t mclabelcheck = kFALSE;
4944 TList *listmc = aod->GetList(); //AliAODEvent* is inherited from AliVEvent*, listmc is pointer to reconstructed event in MC list, member of AliAODEvent
4946 if(!listmc)continue;
4948 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
4950 //feeddown-correction for Lambda/Antilambda particles
4951 //feedddown comes mainly from charged and neutral Xi particles
4952 //feeddown from Sigma decays so quickly that it's not possible to distinguish from primary Lambdas with detector
4953 //feeddown for K0s from phi decays is neglectible
4954 //TH2F* fh2FeedDownMatrix = 0x0; //histo for feeddown already decleared above
4957 //first for all Lambda and Antilambda candidates____________________________________________________________________
4958 TString generatorName;
4960 if(particletype == kLambda){
4962 mclabelcheck = MCLabelCheck(v0, kLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
4965 if((motherType == 3312)||(motherType == 3322)){//mother of v0 is neutral or negative Xi
4967 fListFeeddownLaCand->Add(v0); //fill TList with ass. particles, stemming from feeddown from Xi(bar) decays
4971 if(particletype == kAntiLambda){
4973 mclabelcheck = MCLabelCheck(v0, kAntiLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
4975 if((motherType == -3312)||(motherType == -3322)){
4976 fListFeeddownALaCand->Add(v0); //fill TList with ass. particles, stemming from feeddown from Xi(bar) decays
4982 //_only true primary particles survive the following checks_______________________________________________________________________________________________
4983 TString generatorName;
4985 if(particletype == kK0){
4987 mclabelcheck = MCLabelCheck(v0, kK0, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
4988 if(mclabelcheck == kFALSE)continue;
4990 if(particletype == kLambda){
4992 mclabelcheck = MCLabelCheck(v0, kLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
4993 if(mclabelcheck == kFALSE)continue;
4995 if(particletype == kAntiLambda){
4997 mclabelcheck = MCLabelCheck(v0, kAntiLambda, trackNeg, trackPos, listmc, negDaughterpdg, posDaughterpdg, motherType, v0Label, MCPt, fPhysicalPrimary, MCv0PdgCode, generatorName, isinjected);
4998 if(mclabelcheck == kFALSE)continue;
5001 if(fPhysicalPrimary != 1)continue; //V0 candidate (K0s, Lambda or Antilambda) must be physical primary, this means there is no mother particle existing
5010 Int_t nPart=list->GetSize();
5013 } // end GetListOfV0s()
5015 // -------------------------------------------------------------------------------------------------------
5017 void AliAnalysisTaskJetChem::CalculateInvMass(AliAODv0* v0vtx, const Int_t particletype, Double_t& invM, Double_t& trackPt){
5027 Double_t pp[3]={0,0,0}; //3-momentum positive charged track
5028 Double_t pm[3]={0,0,0}; //3-momentum negative charged track
5030 const Double_t massPi = 0.13957018; //better use PDG code at this point
5031 const Double_t massP = 0.93827203;
5036 TLorentzVector vector; //lorentzvector V0 particle
5037 TLorentzVector fourmom1;//lorentzvector positive daughter
5038 TLorentzVector fourmom2;//lorentzvector negative daughter
5040 //--------------------------------------------------------------
5042 AliAODTrack *trackPos = (AliAODTrack *) (v0vtx->GetSecondaryVtx()->GetDaughter(0));//index 0 defined as positive charged track in AliESDFilter
5044 if( trackPos->Charge() == 1 ){
5046 pp[0]=v0vtx->MomPosX();
5047 pp[1]=v0vtx->MomPosY();
5048 pp[2]=v0vtx->MomPosZ();
5050 pm[0]=v0vtx->MomNegX();
5051 pm[1]=v0vtx->MomNegY();
5052 pm[2]=v0vtx->MomNegZ();
5055 if( trackPos->Charge() == -1 ){
5057 pm[0]=v0vtx->MomPosX();
5058 pm[1]=v0vtx->MomPosY();
5059 pm[2]=v0vtx->MomPosZ();
5061 pp[0]=v0vtx->MomNegX();
5062 pp[1]=v0vtx->MomNegY();
5063 pp[2]=v0vtx->MomNegZ();
5066 if (particletype == kK0){ // case K0s
5067 mass1 = massPi;//positive particle
5068 mass2 = massPi;//negative particle
5069 } else if (particletype == kLambda){ // case Lambda
5070 mass1 = massP;//positive particle
5071 mass2 = massPi;//negative particle
5072 } else if (particletype == kAntiLambda){ //case AntiLambda
5073 mass1 = massPi;//positive particle
5074 mass2 = massP; //negative particle
5077 fourmom1.SetXYZM(pp[0],pp[1],pp[2],mass1);//positive track
5078 fourmom2.SetXYZM(pm[0],pm[1],pm[2],mass2);//negative track
5079 vector=fourmom1 + fourmom2;
5082 trackPt = vector.Pt();
5084 /*// 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
5086 if(particletype == kK0){
5087 std::cout << "invMK0s: " << invM <<std::endl;
5088 std::cout << "v0vtx->MassK0Short(): " << v0vtx->MassK0Short() << std::endl;
5089 std::cout << " " <<std::endl;
5090 //invM = v0vtx->MassK0Short();
5093 if(particletype == kLambda){
5094 std::cout << "invMLambda: " << invM <<std::endl;
5095 std::cout << "v0vtx->MassMassLambda(): " << v0vtx->MassLambda() << std::endl;
5096 std::cout << " " <<std::endl;
5097 //invM = v0vtx->MassLambda();
5100 if(particletype == kAntiLambda){
5101 std::cout << "invMAntiLambda: " << invM <<std::endl;
5102 std::cout << "v0vtx->MassAntiLambda(): " << v0vtx->MassAntiLambda() << std::endl;
5103 std::cout << " " <<std::endl;
5104 //invM = v0vtx->MassAntiLambda();
5112 //_____________________________________________________________________________________
5113 Int_t AliAnalysisTaskJetChem::GetListOfMCParticles(TList *outputlist, const Int_t particletype, AliAODEvent *mcaodevent) //(list to fill here e.g. fListMCgenK0s, particle species to search for)
5116 outputlist->Clear();
5118 TClonesArray *stack = 0x0;
5119 Double_t mcXv=0., mcYv=0., mcZv=0.;//MC vertex position
5122 // get MC generated particles
5124 Int_t fPdgcodeCurrentPart = 0; //pdg code current particle
5125 //Double_t fRapCurrentPart = 0; //get rapidity
5126 //Double_t fPtCurrentPart = 0; //get transverse momentum
5127 Double_t fEtaCurrentPart = 0; //get pseudorapidity
5129 //variable for check: physical primary particle
5130 //Bool_t IsPhysicalPrimary = -1;
5131 //Int_t index = 0; //check number of injected particles
5132 //****************************
5133 // Start loop over MC particles
5135 TList *lst = mcaodevent->GetList();
5138 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
5142 stack = (TClonesArray*)lst->FindObject(AliAODMCParticle::StdBranchName());
5144 Printf("ERROR: stack not available");
5148 AliAODMCHeader *mcHdr=(AliAODMCHeader*)lst->FindObject(AliAODMCHeader::StdBranchName());
5149 if(!mcHdr)return -1;
5151 mcXv=mcHdr->GetVtxX(); mcYv=mcHdr->GetVtxY(); mcZv=mcHdr->GetVtxZ(); // position of the MC primary vertex
5154 ntrk=stack->GetEntriesFast();
5156 //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...
5159 for (Int_t iMc = 0; iMc < ntrk; iMc++) { //loop over mc generated particles
5162 AliAODMCParticle *p0=(AliAODMCParticle*)stack->UncheckedAt(iMc);
5164 //Printf("ERROR: particle with label %d not found in stack (mc loop)", iMc);
5167 fPdgcodeCurrentPart = p0->GetPdgCode();
5169 // Keep only K0s, Lambda and AntiLambda, Xi and Phi:
5170 //if ( (fPdgcodeCurrentPart != 310 ) && (fPdgcodeCurrentPart != 3122 ) && (fPdgcodeCurrentPart != -3122 ) && (fPdgcodeCurrentPart != 3312 ) && (fPdgcodeCurrentPart != -3312) && (fPdgcodeCurrentPart != -333) ) continue;
5174 //Rejection of Pythia injected particles with David Chinellatos method - not the latest method, better Method with TString from MC generator in IsInjected() function below!
5176 /* if( (p0->GetStatus()==21) ||
5177 ((p0->GetPdgCode() == 443) &&
5178 (p0->GetMother() == -1) &&
5179 (p0->GetDaughter(0) == (iMc))) ){ index++; }
5181 if(p0->GetStatus()==21){std::cout<< "hello !!!!" <<std::endl;}
5183 std::cout<< "MC particle status: " << p0->GetStatus() <<std::endl;
5187 //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
5190 //injected particles could be from GenBox (single high pt tracks) or jet related tracks, both generated from PYTHIA MC generator
5192 //Check: MC particle mother
5194 //for feed-down checks
5195 /* //MC gen particles
5196 Int_t iMother = p0->GetMother(); //Motherparticle of V0 candidate (e.g. phi particle,..)
5198 AliAODMCParticle *partM = (AliAODMCParticle*)stack->UncheckedAt(iMother);
5200 if(partM) codeM = TMath::Abs(partM->GetPdgCode());
5203 3312 Xi- -3312 Xibar+
5204 3322 Xi0 -3322 Xibar0
5207 if((codeM == 3312)||(codeM == 3322))// feeddown for Lambda coming from Xi- and Xi0
5213 /* //Check: MC gen. particle decays via 2-pion decay? -> only to be done for the rec. particles !! (-> branching ratio ~ 70 % for K0s -> pi+ pi-)
5215 Int_t daughter0Label = p0->GetDaughter(0);
5216 AliAODMCParticle *mcDaughter0 = (AliAODMCParticle *)stack->UncheckedAt(daughter0Label);
5217 if(daughter0Label >= 0)
5218 {daughter0Type = mcDaughter0->GetPdgCode();}
5220 Int_t daughter1Label = p0->GetDaughter(1);
5221 AliAODMCParticle *mcDaughter1 = (AliAODMCParticle *)stack->UncheckedAt(daughter1Label);
5223 if(daughter1Label >= 1)
5224 {daughter1Type = mcDaughter1->GetPdgCode();} //requirement that daughters are pions is only done for the reconstructed V0s in GetListofV0s() below
5229 // Keep only K0s, Lambda and AntiLambda:
5230 if ( (fPdgcodeCurrentPart != 310 ) && (fPdgcodeCurrentPart != 3122 ) && (fPdgcodeCurrentPart != -3122 )) continue;
5231 // Check: Is physical primary
5233 //do not use anymore: //IsPhysicalPrimary = p0->IsPhysicalPrimary();
5234 //if(!IsPhysicalPrimary)continue;
5236 Float_t fDistPrimaryMax = 0.01; // [cm] max distance of production point to the primary vertex (criterion for choice of MC particles considered as primary)
5238 // Get the distance between production point of the MC mother particle and the primary vertex
5240 Double_t dx = mcXv-p0->Xv();//mc primary vertex - mc gen. v0 vertex
5241 Double_t dy = mcYv-p0->Yv();
5242 Double_t dz = mcZv-p0->Zv();
5244 Double_t fDistPrimary = TMath::Sqrt(dx*dx + dy*dy + dz*dz);
5245 Bool_t fPhysicalPrimary = (fDistPrimary < fDistPrimaryMax);
5247 if(!fPhysicalPrimary)continue;
5249 //if(fPhysicalPrimary){std::cout<<"hello**********************"<<std::endl;}
5251 /* std::cout<<"dx: "<<dx<<std::endl;
5252 std::cout<<"dy: "<<dy<<std::endl;
5253 std::cout<<"dz: "<<dz<<std::endl;
5255 std::cout<<"start: "<<std::endl;
5256 std::cout<<"mcXv: "<<mcXv<<std::endl;
5257 std::cout<<"mcYv: "<<mcYv<<std::endl;
5258 std::cout<<"mcZv: "<<mcZv<<std::endl;
5260 std::cout<<"p0->Xv(): "<<p0->Xv()<<std::endl;
5261 std::cout<<"p0->Yv(): "<<p0->Yv()<<std::endl;
5262 std::cout<<"p0->Zv(): "<<p0->Zv()<<std::endl;
5263 std::cout<<" "<<std::endl;
5264 std::cout<<"fDistPrimary"<<fDistPrimary<<std::endl;
5265 std::cout<<"fPhysicalPrimary"<<fPhysicalPrimary<<std::endl;
5268 //Is close enough to primary vertex to be considered as primary-like?
5270 //fRapCurrentPart = MyRapidity(p0->E(),p0->Pz());
5271 fEtaCurrentPart = p0->Eta();
5272 //fPtCurrentPart = p0->Pt();
5274 if (TMath::Abs(fEtaCurrentPart) < fCutEta){
5275 // if (TMath::Abs(fRapCurrentPart) > fCutRap)continue; //rap cut for crosschecks
5277 if(particletype == kK0){ //MC gen. K0s
5278 if (fPdgcodeCurrentPart==310){
5279 outputlist->Add(p0);
5283 if(particletype == kLambda){ //MC gen. Lambdas
5284 if (fPdgcodeCurrentPart==3122) {
5285 outputlist->Add(p0);
5289 if(particletype == kAntiLambda){
5290 if (fPdgcodeCurrentPart==-3122) { //MC gen. Antilambdas
5291 outputlist->Add(p0);
5296 }//end loop over MC generated particle
5298 Int_t nMCPart=outputlist->GetSize();
5305 //---------------------------------------------------------------------------
5307 Bool_t AliAnalysisTaskJetChem::FillFeeddownMatrix(TList* fListFeeddownCand, Int_t particletype)
5310 // Define Feeddown matrix
5311 Double_t lFeedDownMatrix [100][100];
5312 // FeedDownMatrix [Lambda Bin][Xi Bin];
5314 //Initialize entries of matrix:
5315 for(Int_t ilb = 0; ilb<100; ilb++){
5316 for(Int_t ixb = 0; ixb<100; ixb++){
5317 lFeedDownMatrix[ilb][ixb]=0; //first lambda bins, xi bins
5322 //----------------------------------------------------------------------------
5324 Double_t AliAnalysisTaskJetChem::MyRapidity(Double_t rE, Double_t rPz) const
5326 // Local calculation for rapidity
5327 return 0.5*TMath::Log((rE+rPz)/(rE-rPz+1.e-13));
5329 //----------------------------------------------------------------------------
5332 void AliAnalysisTaskJetChem::GetTracksInCone(TList* inputlist, TList* outputlist, const AliAODJet* jet,
5333 const Double_t radius, Double_t& sumPt, const Double_t minPt, const Double_t maxPt, Bool_t& isBadPt)
5335 // fill list of V0 tracks in cone around jet axis
5338 Bool_t isBadMaxPt = kFALSE;
5339 Bool_t isBadMinPt = kTRUE;
5343 jet->PxPyPz(jetMom);
5344 TVector3 jet3mom(jetMom);
5346 //if(jetets < jetetscutr)continue;
5348 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){//loop over all K0s found in event
5350 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
5352 Double_t trackMom[3];
5353 track->PxPyPz(trackMom);
5354 TVector3 track3mom(trackMom);
5356 Double_t dR = jet3mom.DeltaR(track3mom);
5358 if(dR<radius){//fill all the V0s inside cone into outputlist, radius is reutrn value of GetFFRadius()
5360 outputlist->Add(track);
5362 sumPt += track->Pt();
5364 if(maxPt>0 && track->Pt()>maxPt) isBadMaxPt = kTRUE; // reject jets containing any track with pt larger than this value, use GetFFMaxTrackPt()
5365 if(minPt>0 && track->Pt()>minPt) isBadMinPt = kFALSE; // reject jets with leading track with pt smaller than this value, use GetFFMinLTrackPt()
5371 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)
5372 if(maxPt>0 && isBadMaxPt) isBadPt = kTRUE; //..or because of leading track with too high pt (could be fake track)
5378 //____________________________________________________________________________________________________________________
5381 void AliAnalysisTaskJetChem::GetTracksInPerpCone(TList* inputlist, TList* outputlist, const AliAODJet* jet,
5382 const Double_t radius, Double_t& sumPerpPt)
5384 // fill list of tracks in two cones around jet axis rotated in phi +/- 90 degrees
5386 Double_t jetMom[3]; //array for entries in TVector3
5387 Double_t perpjetplusMom[3]; //array for entries in TVector3
5388 Double_t perpjetnegMom[3];
5392 jet->PxPyPz(jetMom); //get 3D jet momentum
5393 Double_t jetPerpPt = jet->Pt(); //original jet pt, invariant under rotations
5394 Double_t jetPhi = jet->Phi(); //original jet phi
5396 Double_t jetPerpposPhi = jetPhi + ((TMath::Pi())*0.5);//get new perp. jet axis phi clockwise
5397 Double_t jetPerpnegPhi = jetPhi - ((TMath::Pi())*0.5);//get new perp. jet axis phi counterclockwise
5399 TVector3 jet3mom(jetMom); //3-Vector for original rec. jet axis
5401 //Double_t phitest = jet3mom.Phi();
5403 perpjetplusMom[0]=(TMath::Cos(jetPerpposPhi)*jetPerpPt); //x coordinate (sidewards - when looking in beam direction)
5404 perpjetplusMom[1]=(TMath::Sin(jetPerpposPhi)*jetPerpPt); //y coordinate (upwards - when looking in beam direction)
5405 perpjetplusMom[2]=jetMom[2]; //z coordinate (along beam axis), invariant under azimuthal rotation
5407 perpjetnegMom[0]=(TMath::Cos(jetPerpnegPhi)*jetPerpPt); //x coordinate (sidewards - when looking in beam direction)
5408 perpjetnegMom[1]=(TMath::Sin(jetPerpnegPhi)*jetPerpPt); //y coordinate (upwards - when looking in beam direction)
5409 perpjetnegMom[2]=jetMom[2]; //z coordinate (along beam axis), invariant under azimuthal rotation
5412 TVector3 perpjetplus3mom(perpjetplusMom); //3-Vector for new perp. jet axis, clockwise rotated
5413 TVector3 perpjetneg3mom(perpjetnegMom); //3-Vector for new perp. jet axis, counterclockwise rotated
5415 //for crosscheck TVector3 rotation method
5417 //Double_t jetMomplusTest[3];
5418 //Double_t jetMomminusTest[3];
5420 //jet3mom.RotateZ(TMath::Pi()*0.5);//rotate original jet axis around +90 degrees in phi
5422 //perpjetminus3momTest = jet3mom.RotateZ((-1)*TMath::Pi()*0.5);
5424 // jet3mom.RotateZ(TMath::Pi()*0.5);
5425 // jet3mom.RotateZ((-1)*TMath::Pi()*0.5);
5427 //jetMomplusTest[0] = jet3mom.X(); //fetching perp. axis coordinates
5428 //jetMomplusTest[1] = jet3mom.Y();
5429 //jetMomplusTest[2] = jet3mom.Z();
5431 //TVector3 perpjetplus3momTest(jetMomplusTest); //new TVector3 for +90deg rotated jet axis with rotation method from ROOT
5432 //TVector3 perpjetminus3momTest(jetMomminusTest); //new TVector3 for -90deg rotated jet axis with rotation method from ROOT
5435 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){ //collect V0 content in perp cone, rotated clockwise
5437 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack)); //inputlist is fListK0s, all reconstructed K0s in event
5438 if(!track){std::cout<<"K0s track not found!!!"<<std::endl; continue;}
5440 Double_t trackMom[3];//3-mom of V0 particle
5441 track->PxPyPz(trackMom);
5442 TVector3 track3mom(trackMom);
5444 Double_t dR = perpjetplus3mom.DeltaR(track3mom);
5448 outputlist->Add(track); // output list is jetPerpConeK0list
5450 sumPerpPt += track->Pt();
5457 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){//collect V0 content in perp cone, rotated counterclockwise
5459 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack)); //inputlist is fListK0s, all reconstructed K0s in event
5460 if(!track){std::cout<<"K0s track not found!!!"<<std::endl; continue;}
5462 Double_t trackMom[3];//3-mom of V0 particle
5463 track->PxPyPz(trackMom);
5464 TVector3 track3mom(trackMom);
5466 Double_t dR = perpjetneg3mom.DeltaR(track3mom);
5470 outputlist->Add(track); // output list is jetPerpConeK0list
5472 sumPerpPt += track->Pt();
5478 // pay attention: this list contains the double amount of V0s, found in both cones
5479 // before using it, devide spectra by 2!!!
5480 sumPerpPt = sumPerpPt*0.5; //correct to do this?
5488 // _______________________________________________________________________________________________________________________________________________________
5490 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){
5492 if(!v0)return kFALSE;
5494 TClonesArray *stackmc = 0x0;
5495 stackmc = (TClonesArray*)listmc->FindObject(AliAODMCParticle::StdBranchName()); //get MCAOD branch in data
5498 Printf("ERROR: stack not available");
5503 Int_t negAssLabel = TMath::Abs(trackNeg->GetLabel()); //negative (reconstructed) charged track label in MC stack
5504 Int_t posAssLabel = TMath::Abs(trackPos->GetLabel()); //positive (reconstructed) charged track label in MC stack
5506 //injected particle checks
5512 AliAODMCHeader *header=(AliAODMCHeader*)listmc->FindObject(AliAODMCHeader::StdBranchName());
5513 if(!header)return kFALSE;
5515 mcXv=header->GetVtxX(); mcYv=header->GetVtxY(); mcZv=header->GetVtxZ();
5521 if(negAssLabel>=0 && negAssLabel < stackmc->GetEntriesFast() && posAssLabel>=0 && posAssLabel < stackmc->GetEntriesFast()){//safety check if label has valid value of stack
5523 AliAODMCParticle *mcNegPart =(AliAODMCParticle*)stackmc->UncheckedAt(negAssLabel);//fetch the, with one MC truth track associated (reconstructed), negative charged track
5524 v0Label = mcNegPart->GetMother();// mother of negative charged particle is v0, get v0 label here
5525 negDaughterpdg = mcNegPart->GetPdgCode();
5526 AliAODMCParticle *mcPosPart =(AliAODMCParticle*)stackmc->UncheckedAt(posAssLabel);//fetch the, with one MC truth track associated (reconstructed), positive charged track
5527 Int_t v0PosLabel = mcPosPart->GetMother(); //get mother label of positive charged track label
5528 posDaughterpdg = mcPosPart->GetPdgCode();
5530 if(v0Label >= 0 && v0Label < stackmc->GetEntriesFast() && v0Label == v0PosLabel){//first v0 mc label check, then: check if both daughters are stemming from same particle
5532 AliAODMCParticle *mcv0 = (AliAODMCParticle *)stackmc->UncheckedAt(v0Label); //fetch MC ass. particle to v0 (mother of the both charged daughter tracks)
5534 Float_t fDistPrimaryMax = 0.01; // [cm] max distance of production point to the primary vertex (criterion for choice of MC particles considered as primary)
5536 // Get the distance between production point of the MC mother particle and the primary vertex
5538 Double_t dx = mcXv-mcv0->Xv();//mc primary vertex - mc particle production vertex
5539 Double_t dy = mcYv-mcv0->Yv();
5540 Double_t dz = mcZv-mcv0->Zv();
5542 Float_t fDistPrimary = TMath::Sqrt(dx*dx + dy*dy + dz*dz);
5543 fPhysicalPrimary = kFALSE;//init
5545 fPhysicalPrimary = (fDistPrimary < fDistPrimaryMax);
5546 MCv0PDGCode = mcv0->GetPdgCode();
5548 //if(fPhysicalPrimary == kTRUE){//look only at physical primary particles
5550 isinjected = IsTrackInjected(v0Label, header, stackmc, generatorName);
5552 //trackinjected is kFALSE if it is either Hijing or has no generator name
5553 // std::cout<<" "<<std::endl;
5554 // std::cout<<"#### next particle: ####"<<std::endl;
5555 //std::cout<<"Is track injected: "<< trackinjected <<std::endl;
5556 // std::cout<<"pdg code: "<<MCv0PDGCode<<std::endl;
5557 // std::cout<<"v0Label: "<<v0Label<<std::endl;
5559 MCPt = mcv0->Pt();//for MC data, always use MC gen. pt for any pt distributions, also for the spectra, used for normalisation
5560 //for feed-down checks later
5562 Int_t motherLabel = mcv0->GetMother(); //get mother particle label of v0 particle
5563 // std::cout<<"motherLabel: "<<motherLabel<<std::endl;
5565 if(motherLabel >= 0 && v0Label < stackmc->GetEntriesFast()) //do safety check for mother label
5567 AliAODMCParticle *mcMother = (AliAODMCParticle *)stackmc->UncheckedAt(motherLabel); //get mother particle
5568 motherType = mcMother->GetPdgCode(); //get PDG code of mother
5571 Double_t XibarPt = 0.;
5573 if(particletype == kLambda){
5574 if((motherType == 3312)||(motherType == 3322)){ //if v0 mother is Xi0 or Xi- fill MC gen. pt in FD La histogram
5575 XiPt = mcMother->Pt();
5576 fh1MCXiPt->Fill(XiPt);
5579 if(particletype == kAntiLambda){
5580 if((motherType == -3312)||(motherType == -3322)){ //if v0 mother is Xibar0 or Xibar+ fill MC gen. pt in FD ALa histogram
5581 XibarPt = mcMother->Pt();
5582 fh1MCXibarPt->Fill(XibarPt);
5588 //pdg code checks etc..
5590 if(particletype == kK0){
5592 if(TMath::Abs(posDaughterpdg) != 211){return kFALSE;}//one or both of the daughters are not a pion
5593 if(TMath::Abs(negDaughterpdg) != 211){return kFALSE;}
5595 if(MCv0PDGCode != 310) {return kFALSE;}
5598 if(particletype == kLambda){
5599 if(MCv0PDGCode != 3122)return kFALSE;//if particle is not Antilambda, v0 is rejected
5600 if(posDaughterpdg != 2212)return kFALSE;
5601 if(negDaughterpdg != -211)return kFALSE; //pdg code check for Lambda daughters
5603 //{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 //}
5606 if(particletype == kAntiLambda){
5607 if(MCv0PDGCode != -3122)return kFALSE;
5608 if(posDaughterpdg != 211)return kFALSE;
5609 if(negDaughterpdg !=-2212)return kFALSE; //pdg code check for Antilambda daughters
5612 //{if((motherType == -3312)||(motherType == -3322)){continue;}//if bar{Xi0} and Xi+ is motherparticle of Antilambda, particle is rejected
5616 return kTRUE; //check was successful
5617 }//end mc v0 label check
5618 }// end of stack label check
5623 return kFALSE; //check wasn't successful
5625 //________________________________________________________________________________________________________________________________________________________
5628 Bool_t AliAnalysisTaskJetChem::IsParticleMatching(const AliAODMCParticle* mcp0, const Int_t v0Label){
5630 const Int_t mcp0label = mcp0->GetLabel();
5632 if(v0Label == mcp0label)return kTRUE;
5637 //_______________________________________________________________________________________________________________________________________________________
5639 Bool_t AliAnalysisTaskJetChem::DaughterTrackCheck(AliAODv0* v0, Int_t& nnum, Int_t& pnum){
5642 if(v0->GetNDaughters() != 2) return kFALSE;//case v0 has more or less than 2 daughters, avoids seg. break at some AOD files //reason?
5645 // safety check of input parameters
5648 if(fDebug > 1){std::cout << std::endl
5649 << "Warning in AliAnalysisTaskJetChem::DaughterTrackCheck:" << std::endl
5650 << "v0 = " << v0 << std::endl;}
5656 //Daughters track check: its Luke Hanrattys method to check daughters charge
5662 AliAODTrack *ntracktest =(AliAODTrack*)(v0->GetDaughter(nnum));
5664 if(ntracktest == NULL)
5666 if(fDebug > 1){std::cout << std::endl
5667 << "Warning in AliAnalysisTaskJetChem::DaughterTrackCheck:" << std::endl
5668 << "ntracktest = " << ntracktest << std::endl;}
5673 if(ntracktest->Charge() > 0)
5679 const AliAODTrack *trackNeg=(AliAODTrack *)(v0->GetDaughter(nnum));
5680 const AliAODTrack *trackPos=(AliAODTrack *)(v0->GetDaughter(pnum));
5682 //Check if both tracks are available
5683 if (!trackPos || !trackNeg) {
5684 if(fDebug > 1) Printf("strange analysis::UserExec:: Error:Could not retrieve one of the daughter tracks\n");
5689 //remove like sign V0s
5690 if ( trackPos->Charge() == trackNeg->Charge() ){
5691 //if(fDebug>1) Printf("%s:%d found like-sign V0", (char*)__FILE__,__LINE__);
5699 //______________________________________________________________________
5700 TString AliAnalysisTaskJetChem::GetGenerator(Int_t label, AliAODMCHeader* header){
5701 Int_t nsumpart=0;//number of particles
5702 TList *lh=header->GetCocktailHeaders();//TList with all generator headers
5703 Int_t nh=lh->GetEntries();//number of entries in TList with all headers
5705 for(Int_t i=0;i<nh;i++){
5706 AliGenEventHeader* gh=(AliGenEventHeader*)lh->At(i);
5707 TString genname=gh->GetName();//name of particle generator
5708 Int_t npart=gh->NProduced();//number of stable or undecayed particles in MC stack block (?)
5709 if(label>=nsumpart && label<(nsumpart+npart)) return genname;
5716 //_____________________________________________________________________
5717 void AliAnalysisTaskJetChem::GetTrackPrimaryGenerator(Int_t lab, AliAODMCHeader *header,TClonesArray *arrayMC,TString &nameGen){
5719 // method to check if a particle is stemming from a given generator
5721 nameGen=GetGenerator(lab,header);
5723 // Int_t countControl=0;
5725 while(nameGen.IsWhitespace()){
5726 AliAODMCParticle *mcpart= (AliAODMCParticle*)arrayMC->At(lab);//get MC generated particle for particle MC label
5728 printf("AliAnalysisTaskJetChem::IsTrackInjected - BREAK: No valid AliAODMCParticle at label %i\n",lab);
5731 Int_t mother = mcpart->GetMother();
5734 printf("AliAnalysisTaskJetChem::IsTrackInjected - BREAK: Reached primary particle without valid mother\n");
5738 nameGen=GetGenerator(mother,header);
5741 // if(countControl>=10){ // 10 = arbitrary number; protection from infinite loops
5742 // printf("AliVertexingHFUtils::IsTrackInjected - BREAK: Protection from infinite loop active\n");
5752 //---------------------------------------------------------------------------------------------------------------------
5753 Bool_t AliAnalysisTaskJetChem::IsTrackInjected(Int_t lab, AliAODMCHeader *header,TClonesArray *arrayMC, TString& nameGen){
5754 // method to check if a v0 particle comes from the signal event or from the underlying Hijing event
5757 GetTrackPrimaryGenerator(lab, header, arrayMC, nameGen);
5759 if(nameGen.IsWhitespace() || nameGen.Contains("ijing")) return kFALSE;//particle has either no info about generator or is Hijing particle, so it is not injected
5761 //std::cout<<"generator name: "<<nameGen<<std::endl;
5766 //_________________________________________________________________________________________________________________________________________
5767 Double_t AliAnalysisTaskJetChem::SmearJetPt(Double_t jetPt, Int_t /*cent*/, Double_t /*jetRadius*/, Double_t /*ptmintrack*/, Double_t& jetPtSmear){
5769 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
5773 /* if(cent>10) cl = 2;
5778 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
5779 //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
5781 //fsmear->SetParameters(1,0,4.472208);// for 2010 PbPb jets, R=0.2, ptmintrack = 0.15 GeV/c, cent 00-10%
5783 /* //delta-pt width for anti-kt jet finder:
5786 if((cl == 1)&&(jetRadius == 0.4)&&(ptmintrack == 0.15)){
5787 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
5789 if((cl == 2)&&(jetRadius == 0.4)&&(ptmintrack == 0.15)){
5790 fsmear->SetParameters(1,0,8.536195);
5792 if((cl == 3)&&(jetRadius == 0.4)&&(ptmintrack == 0.15)){
5793 fsmear->SetParameters(1,0,?);
5795 if((cl == 4)&&(jetRadius == 0.4)&&(ptmintrack == 0.15)){
5796 fsmear->SetParameters(1,0,5.229839);
5800 if((cl == 1)&&(jetRadius == 0.3)&&(ptmintrack == 0.15)){
5801 fsmear->SetParameters(1,0,7.145967);
5803 if((cl == 2)&&(jetRadius == 0.3)&&(ptmintrack == 0.15)){
5804 fsmear->SetParameters(1,0,5.844796);
5806 if((cl == 3)&&(jetRadius == 0.3)&&(ptmintrack == 0.15)){
5807 fsmear->SetParameters(1,0,?);
5809 if((cl == 4)&&(jetRadius == 0.3)&&(ptmintrack == 0.15)){
5810 fsmear->SetParameters(1,0,3.630751);
5814 if((cl == 1)&&(jetRadius == 0.2)&&(ptmintrack == 0.15)){
5815 fsmear->SetParameters(1,0,4.472208);
5817 if((cl == 2)&&(jetRadius == 0.2)&&(ptmintrack == 0.15)){
5818 fsmear->SetParameters(1,0,3.543938);
5820 if((cl == 3)&&(jetRadius == 0.2)&&(ptmintrack == 0.15)){
5821 fsmear->SetParameters(1,0,?);
5823 if((cl == 4)&&(jetRadius == 0.2)&&(ptmintrack == 0.15)){
5824 fsmear->SetParameters(1,0,1.037476);
5829 Double_t r = fsmear.GetRandom();
5830 jetPtSmear = jetPt + r;
5832 // std::cout<<"jetPt: "<<jetPt<<std::endl;
5833 // std::cout<<"jetPtSmear: "<<jetPtSmear<<std::endl;
5834 // std::cout<<"r: "<<r<<std::endl;
5841 //______________________________________________________________________________________________________________________
5842 //____________________________________________________________________________________________________________________
5844 Bool_t AliAnalysisTaskJetChem::IsParticleInCone(const AliVParticle* part1, const AliVParticle* part2, Double_t dRMax) const
5846 // decides whether a particle is inside a jet cone
5847 if (!part1 || !part2)
5850 TVector3 vecMom2(part2->Px(),part2->Py(),part2->Pz());
5851 TVector3 vecMom1(part1->Px(),part1->Py(),part1->Pz());
5852 Double_t dR = vecMom2.DeltaR(vecMom1); // = sqrt(dEta*dEta+dPhi*dPhi)
5853 if(dR<dRMax) // momentum vectors of part1 and part2 are closer than dRMax
5857 //__________________________________________________________________________________________________________________
5860 Bool_t AliAnalysisTaskJetChem::IsRCJCOverlap(TList* recjetlist, const AliVParticle* part, Double_t dDistance) const{
5862 if(!recjetlist) return kFALSE;
5863 if(!part) return kFALSE;
5864 if(!dDistance) return kFALSE;
5865 Int_t nRecJetsCuts = fJetsRecCuts->GetEntries();
5867 for(Int_t i=0; i<nRecJetsCuts; ++i){ //loop over all reconstructed jets in events
5868 AliAODJet* jet = (AliAODJet*) (recjetlist->At(i));
5869 if(!jet){if(fDebug>2)std::cout<<"AliAnalysisTaskJetChem::IsRCJCOverlap jet pointer invalid!"<<std::endl;continue;}
5870 if(IsParticleInCone(jet, part, dDistance) == kTRUE)return kTRUE;//RC and JC are overlapping
5872 }//end loop testing RC-JC overlap
5873 return kFALSE;//RC and JC are not overlapping -> good!
5876 //_______________________________________________________________________________________________________________________
5877 AliAODJet* AliAnalysisTaskJetChem::GetRandomCone(TList* jetlist, Double_t dEtaConeMax, Double_t dDistance) const
5879 TLorentzVector vecRdCone;
5880 AliAODJet* jetRC = 0;//random cone candidate
5881 Double_t dEta, dPhi; //random eta and phi value for RC
5882 Bool_t IsRCoutJC = kFALSE;//check whether RC is not overlapping with any selected jet cone in event
5883 Int_t iRCTrials = 10;//search at maximum 10 times for random cone that doesn't overlap with jet cone
5885 for(Int_t i=0; i<iRCTrials; iRCTrials++){
5887 dEta = dEtaConeMax*(2*fRandom->Rndm()-1.); //random eta value in range: [-dEtaConeMax,+dEtaConeMax]
5888 dPhi = TMath::TwoPi()*fRandom->Rndm(); //random phi value in range: [0,2*Pi]
5889 vecRdCone.SetPtEtaPhiM(1.,dEta,dPhi,0.);
5890 jetRC = new AliAODJet(vecRdCone);//new RC candidate
5892 if (!IsRCJCOverlap(jetlist,jetRC,dDistance))
5894 IsRCoutJC = kTRUE; //std::cout<<"RC and JC are not overlapping!!!"<<std::endl;
5898 delete jetRC; //RC is overlapping with JC, delete this RC candidate
5901 if(!IsRCoutJC) {jetRC = 0;}//in case no random cone was selected
5907 // _______________________________________________________________________________________________________________________
5908 AliAODJet* AliAnalysisTaskJetChem::GetMedianCluster()
5910 // fill tracks from bckgCluster branch,
5911 // using cluster with median density (odd number of clusters)
5912 // or picking randomly one of the two closest to median (even number)
5914 Int_t nBckgClusters = fBckgJetsRec->GetEntries(); // not 'recCuts': use all clusters in full eta range
5916 if(nBckgClusters<3) return 0; // need at least 3 clusters (skipping 2 highest)
5918 Double_t* bgrDensity = new Double_t[nBckgClusters];
5919 Int_t* indices = new Int_t[nBckgClusters];
5921 for(Int_t ij=0; ij<nBckgClusters; ++ij){
5923 AliAODJet* bgrCluster = (AliAODJet*)(fBckgJetsRec->At(ij));
5924 Double_t clusterPt = bgrCluster->Pt();
5925 Double_t area = bgrCluster->EffectiveAreaCharged();
5927 Double_t density = 0;
5928 if(area>0) density = clusterPt/area;
5930 bgrDensity[ij] = density;
5935 TMath::Sort(nBckgClusters, bgrDensity, indices);
5937 // get median cluster
5939 AliAODJet* medianCluster = 0;
5941 if(TMath::Odd(nBckgClusters)){
5943 Int_t medianIndex = indices[(Int_t) (0.5*(nBckgClusters+1))];
5945 medianCluster = (AliAODJet*)(fBckgJetsRec->At(medianIndex));
5947 //Double_t clusterPt = medianCluster->Pt();
5948 //Double_t area = medianCluster->EffectiveAreaCharged();
5952 Int_t medianIndex1 = indices[(Int_t) (0.5*nBckgClusters)];
5953 Int_t medianIndex2 = indices[(Int_t) (0.5*nBckgClusters+1)];
5955 AliAODJet* medianCluster1 = (AliAODJet*)(fBckgJetsRec->At(medianIndex1));
5956 AliAODJet* medianCluster2 = (AliAODJet*)(fBckgJetsRec->At(medianIndex2));
5958 // Double_t density1 = 0;
5959 //Double_t clusterPt1 = medianCluster1->Pt();
5960 //Double_t area1 = medianCluster1->EffectiveAreaCharged();
5961 //if(area1>0) Double_t density1 = clusterPt1/area1;
5963 // Double_t density2 = 0;
5964 //Double_t clusterPt2 = medianCluster2->Pt();
5965 //Double_t area2 = medianCluster2->EffectiveAreaCharged();
5966 // if(area2>0) Double_t density2 = clusterPt2/area2;
5968 medianCluster = ( (gRandom->Rndm()>0.5) ? medianCluster1 : medianCluster2 ); // select one randomly to avoid adding areas
5971 delete[] bgrDensity;
5974 return medianCluster;
5977 //____________________________________________________________________________________________
5979 Double_t AliAnalysisTaskJetChem::AreaCircSegment(Double_t dRadius, Double_t dDistance) const
5981 // calculate area of a circular segment defined by the circle radius and the (oriented) distance between the secant line and the circle centre
5982 Double_t dEpsilon = 1e-2;
5983 Double_t dR = dRadius;
5984 Double_t dD = dDistance;
5985 if (TMath::Abs(dR)<dEpsilon)
5987 if(fDebug>0) printf("AliAnalysisTaskJetChem::AreaCircSegment: Error: Too small radius: %f < %f\n",dR,dEpsilon);
5993 return TMath::Pi()*dR*dR;
5994 return dR*dR*TMath::ACos(dD/dR)-dD*TMath::Sqrt(dR*dR-dD*dD);