@@ -498,6 +498,11 @@ struct Kshortlambda {
498498 return true ;
499499 }
500500
501+ double massK0s = o2::constants::physics::MassK0Short;
502+ double massLambda = o2::constants::physics::MassLambda;
503+ double massPr = o2::constants::physics::MassProton;
504+ double massPi = o2::constants::physics::MassPionCharged;
505+
501506 // Defining filters for events (event selection)
502507 // Filter eventFilter = (o2::aod::evsel::sel8 == true);
503508
@@ -516,14 +521,13 @@ struct Kshortlambda {
516521
517522 ROOT::Math::PxPyPzMVector daughter1, daughter2;
518523 ROOT::Math::PxPyPzMVector protonVec, pionVec, lambdaVec;
524+ ROOT::Math::PxPyPzMVector fourVecDau, fourVecMother, fourVecDauCM;
525+ ROOT::Math::XYZVector threeVecDauCM;
526+ ROOT::Math::XYZVector helicityVec, normalVec, randomVec, beamVec;
519527
520528 void processSE (EventCandidates::iterator const & collision, TrackCandidates const & /* tracks*/ , aod::V0Datas const & V0s)
521529 {
522530 hvzero.fill (HIST (" heventscheck" ), 0.5 );
523- const double massK0s = o2::constants::physics::MassK0Short;
524- const double massLambda = o2::constants::physics::MassLambda;
525- const double massPr = o2::constants::physics::MassProton;
526- const double massPi = o2::constants::physics::MassPionCharged;
527531
528532 float multiplicity = 0 .0f ;
529533 if (cfgMultFOTM) {
@@ -550,7 +554,7 @@ struct Kshortlambda {
550554 }
551555
552556 std::vector<int > v0indexes;
553-
557+ TLorentzVector lv1, lv3, lvLambda, lv4, lv5;
554558 for (const auto & [v1, v2] : combinations (CombinationsStrictlyUpperIndexPolicy (V0s, V0s))) {
555559 hvzero.fill (HIST (" heventscheck" ), 2.5 );
556560 if (v1.size () == 0 || v2.size () == 0 ) {
@@ -610,9 +614,9 @@ struct Kshortlambda {
610614 }
611615 lambdaVec = protonVec + pionVec;
612616
613- TLorentzVector lv1, lv2, lv3, lvLambda, lv4, lv5;
614617 lv1.SetPtEtaPhiM (v1.pt (), v1.eta (), v1.phi (), massK0s);
615618 lvLambda.SetPtEtaPhiM (lambdaVec.pt (), lambdaVec.eta (), lambdaVec.phi (), massLambda);
619+
616620 lv3 = lv1 + lvLambda;
617621
618622 if (qAv0daughters) {
@@ -628,61 +632,60 @@ struct Kshortlambda {
628632 daughter2 = ROOT::Math::PxPyPzMVector (v2.px (), v2.py (), v2.pz (), massLambda); // V02
629633
630634 // polarization calculations
631- ROOT::Math::PxPyPzMVector fourVecDau = ROOT::Math::PxPyPzMVector (daughter1.Px (), daughter1.Py (), daughter1.Pz (), massK0s); // Kshort
632- ROOT::Math::PxPyPzMVector fourVecMother = ROOT::Math::PxPyPzMVector (lv3.Px (), lv3.Py (), lv3.Pz (), lv3.M ()); // mass of KshortKshort pair
633- ROOT::Math::Boost boost{fourVecMother.BoostToCM ()}; // boost mother to center of mass frame
634- ROOT::Math::PxPyPzMVector fourVecDauCM = boost (fourVecDau); // boost the frame of daughter same as mother
635- ROOT::Math::XYZVector threeVecDauCM = fourVecDauCM.Vect (); // get the 3 vector of daughter in the frame of mother
635+ fourVecDau = ROOT::Math::PxPyPzMVector (daughter1.Px (), daughter1.Py (), daughter1.Pz (), massK0s); // Kshort
636+ fourVecMother = ROOT::Math::PxPyPzMVector (lv3.Px (), lv3.Py (), lv3.Pz (), lv3.M ()); // mass of KshortKshort pair
637+ ROOT::Math::Boost boost{fourVecMother.BoostToCM ()}; // boost mother to center of mass frame
638+ fourVecDauCM = boost (fourVecDau); // boost the frame of daughter same as mother
639+ threeVecDauCM = fourVecDauCM.Vect (); // get the 3 vector of daughter in the frame of mother
640+
636641 if (std::abs (lv3.Rapidity ()) < 0.5 ) {
637642 if (invMass1D) {
638643 hvzero.fill (HIST (" h1vzeropairInvMassRot" ), lv3.M ());
639644 }
640645
641646 if (activateTHnSparseCosThStarHelicity) {
642- ROOT::Math::XYZVector helicityVec = fourVecMother.Vect (); // 3 vector of mother in COM frame
647+ helicityVec = fourVecMother.Vect (); // 3 vector of mother in COM frame
643648 auto cosThetaStarHelicity = helicityVec.Dot (threeVecDauCM) / (std::sqrt (threeVecDauCM.Mag2 ()) * std::sqrt (helicityVec.Mag2 ()));
644649 hvzero.fill (HIST (" h3vzeropairInvMassDS" ), multiplicity, lv3.Pt (), lv3.M (), cosThetaStarHelicity, occupancyno);
645650
646651 for (int i = 0 ; i < cnofrotations; i++) {
647652 float theta2 = rn->Uniform (o2::constants::math::PI - o2::constants::math::PI / rotationalCut, o2::constants::math::PI + o2::constants::math::PI / rotationalCut);
648653 lv4.SetPtEtaPhiM (lambdaVec.pt (), lambdaVec.eta (), lambdaVec.phi () + theta2, massLambda); // for rotated background
649-
650654 lv5 = lv1 + lv4;
651-
652655 hvzero.fill (HIST (" h3vzeropairInvMassRot" ), multiplicity, lv5.Pt (), lv5.M (), cosThetaStarHelicity, occupancyno);
653656 }
654657
655658 } else if (activateTHnSparseCosThStarProduction) {
656- ROOT::Math::XYZVector normalVec = ROOT::Math::XYZVector (lv3.Py (), -lv3.Px (), 0 .f );
659+ normalVec = ROOT::Math::XYZVector (lv3.Py (), -lv3.Px (), 0 .f );
657660 auto cosThetaStarProduction = normalVec.Dot (threeVecDauCM) / (std::sqrt (threeVecDauCM.Mag2 ()) * std::sqrt (normalVec.Mag2 ()));
658661 hvzero.fill (HIST (" h3vzeropairInvMassDS" ), multiplicity, lv3.Pt (), lv3.M (), cosThetaStarProduction, occupancyno);
659662 for (int i = 0 ; i < cnofrotations; i++) {
660663
661- float theta2 = rn->Uniform (o2::constants::math::PI - o2::constants::math::PI / rotationalCut, o2::constants::math::PI + o2::constants::math::PI / rotationalCut);
664+ auto theta2 = rn->Uniform (o2::constants::math::PI - o2::constants::math::PI / rotationalCut, o2::constants::math::PI + o2::constants::math::PI / rotationalCut);
662665
663666 lv4.SetPtEtaPhiM (lambdaVec.pt (), lambdaVec.eta (), lambdaVec.phi () + theta2, massLambda); // for rotated background
664667 lv5 = lv1 + lv4;
665668 hvzero.fill (HIST (" h3vzeropairInvMassRot" ), multiplicity, lv5.Pt (), lv5.M (), cosThetaStarProduction, occupancyno);
666669 }
667670 } else if (activateTHnSparseCosThStarBeam) {
668- ROOT::Math::XYZVector beamVec = ROOT::Math::XYZVector (0 .f , 0 .f , 1 .f );
671+ beamVec = ROOT::Math::XYZVector (0 .f , 0 .f , 1 .f );
669672 auto cosThetaStarBeam = beamVec.Dot (threeVecDauCM) / std::sqrt (threeVecDauCM.Mag2 ());
670673 hvzero.fill (HIST (" h3vzeropairInvMassDS" ), multiplicity, lv3.Pt (), lv3.M (), cosThetaStarBeam, occupancyno);
671674 for (int i = 0 ; i < cnofrotations; i++) {
672- float theta2 = rn->Uniform (o2::constants::math::PI - o2::constants::math::PI / rotationalCut, o2::constants::math::PI + o2::constants::math::PI / rotationalCut);
675+ auto theta2 = rn->Uniform (o2::constants::math::PI - o2::constants::math::PI / rotationalCut, o2::constants::math::PI + o2::constants::math::PI / rotationalCut);
673676 lv4.SetPtEtaPhiM (lambdaVec.pt (), lambdaVec.eta (), lambdaVec.phi () + theta2, massLambda); // for rotated background
674677 lv5 = lv1 + lv4;
675678 hvzero.fill (HIST (" h3vzeropairInvMassRot" ), multiplicity, lv5.Pt (), lv5.M (), cosThetaStarBeam, occupancyno);
676679 }
677680 } else if (activateTHnSparseCosThStarRandom) {
678681 auto phiRandom = gRandom ->Uniform (0 .f , constants::math::TwoPI);
679682 auto thetaRandom = gRandom ->Uniform (0 .f , constants::math::PI);
680- ROOT::Math::XYZVector randomVec = ROOT::Math::XYZVector (std::sin (thetaRandom) * std::cos (phiRandom), std::sin (thetaRandom) * std::sin (phiRandom), std::cos (thetaRandom));
683+ randomVec = ROOT::Math::XYZVector (std::sin (thetaRandom) * std::cos (phiRandom), std::sin (thetaRandom) * std::sin (phiRandom), std::cos (thetaRandom));
681684 auto cosThetaStarRandom = randomVec.Dot (threeVecDauCM) / std::sqrt (threeVecDauCM.Mag2 ());
682685 hvzero.fill (HIST (" h3vzeropairInvMassDS" ), multiplicity, lv3.Pt (), lv3.M (), cosThetaStarRandom, occupancyno);
683686
684687 for (int i = 0 ; i < cnofrotations; i++) {
685- float theta2 = rn->Uniform (o2::constants::math::PI - o2::constants::math::PI / rotationalCut, o2::constants::math::PI + o2::constants::math::PI / rotationalCut);
688+ auto theta2 = rn->Uniform (o2::constants::math::PI - o2::constants::math::PI / rotationalCut, o2::constants::math::PI + o2::constants::math::PI / rotationalCut);
686689 lv4.SetPtEtaPhiM (lambdaVec.pt (), lambdaVec.eta (), lambdaVec.phi () + theta2, massLambda); // for rotated background
687690 lv5 = lv1 + lv4;
688691 hvzero.fill (HIST (" h3vzeropairInvMassRot" ), multiplicity, lv5.Pt (), lv5.M (), cosThetaStarRandom, occupancyno);
@@ -707,11 +710,6 @@ struct Kshortlambda {
707710
708711 void processME (EventCandidates const & collisions, TrackCandidates const & /* tracks*/ , V0TrackCandidate const & v0s)
709712 {
710- const double massK0s = o2::constants::physics::MassK0Short;
711- const double massLambda = o2::constants::physics::MassLambda;
712- const double massPr = o2::constants::physics::MassProton;
713- const double massPi = o2::constants::physics::MassPionCharged;
714-
715713 auto tracksTuple = std::make_tuple (v0s);
716714 BinningTypeVertexContributor binningOnPositions1{{mevz, memult}, true };
717715 BinningTypeCentralityM binningOnPositions2{{mevz, memult}, true };
@@ -720,6 +718,8 @@ struct Kshortlambda {
720718 SameKindPair<EventCandidates, V0TrackCandidate, BinningTypeCentralityM> pair2{binningOnPositions2, cfgNmixedEvents, -1 , collisions, tracksTuple, &cache}; // for pp
721719
722720 if (cfgMultFOTM) {
721+ TLorentzVector lv1, lv3, lvLambda, lv4, lv5;
722+
723723 for (const auto & [c1, tracks1, c2, tracks2] : pair2) // two different centrality c1 and c2 and tracks corresponding to them
724724 {
725725
@@ -733,7 +733,7 @@ struct Kshortlambda {
733733 auto occupancyno = c1.trackOccupancyInTimeRange ();
734734 auto occupancyno2 = c2.trackOccupancyInTimeRange ();
735735 if (applyOccupancyCut && (occupancyno < occupancyCut || occupancyno2 < occupancyCut)) {
736- return ;
736+ continue ;
737737 }
738738
739739 for (const auto & [t1, t2] : o2::soa::combinations (o2::soa::CombinationsFullIndexPolicy (tracks1, tracks2))) {
@@ -792,33 +792,32 @@ struct Kshortlambda {
792792 pionVec = ROOT::Math::PxPyPzMVector (postrack2.px (), postrack2.py (), postrack2.pz (), massPi);
793793 }
794794 lambdaVec = protonVec + pionVec;
795-
796- TLorentzVector lv1, lv2, lv3, lvLambda, lv4, lv5;
797795 lv1.SetPtEtaPhiM (t1.pt (), t1.eta (), t1.phi (), massK0s);
798796 lvLambda.SetPtEtaPhiM (lambdaVec.pt (), lambdaVec.eta (), lambdaVec.phi (), massLambda);
799797 lv3 = lv1 + lvLambda;
800- ROOT::Math::PxPyPzMVector fourVecDau = ROOT::Math::PxPyPzMVector (daughter1.Px (), daughter1.Py (), daughter1.Pz (), massK0s); // Kshort
801- ROOT::Math::PxPyPzMVector fourVecMother = ROOT::Math::PxPyPzMVector (lv3.Px (), lv3.Py (), lv3.Pz (), lv3.M ()); // mass of KshortKshort pair
798+ fourVecDau = ROOT::Math::PxPyPzMVector (daughter1.Px (), daughter1.Py (), daughter1.Pz (), massK0s); // Kshort
799+ fourVecMother = ROOT::Math::PxPyPzMVector (lv3.Px (), lv3.Py (), lv3.Pz (), lv3.M ()); // mass of KshortKshort pair
802800 ROOT::Math::Boost boost{fourVecMother.BoostToCM ()}; // boost mother to center of mass frame
803- ROOT::Math::PxPyPzMVector fourVecDauCM = boost (fourVecDau); // boost the frame of daughter same as mother
804- ROOT::Math::XYZVector threeVecDauCM = fourVecDauCM.Vect (); // get the 3 vector of daughter in the frame of mother
801+ fourVecDauCM = boost (fourVecDau); // boost the frame of daughter same as mother
802+ threeVecDauCM = fourVecDauCM.Vect (); // get the 3 vector of daughter in the frame of mother
803+
805804 if (std::abs (lv3.Rapidity ()) < 0.5 ) {
806805 if (activateTHnSparseCosThStarHelicity) {
807- ROOT::Math::XYZVector helicityVec = fourVecMother.Vect (); // 3 vector of mother in COM frame
806+ helicityVec = fourVecMother.Vect (); // 3 vector of mother in COM frame
808807 auto cosThetaStarHelicity = helicityVec.Dot (threeVecDauCM) / (std::sqrt (threeVecDauCM.Mag2 ()) * std::sqrt (helicityVec.Mag2 ()));
809808 hvzero.fill (HIST (" h3vzeropairInvMassME" ), multiplicity, lv3.Pt (), lv3.M (), cosThetaStarHelicity, occupancyno);
810809 } else if (activateTHnSparseCosThStarProduction) {
811- ROOT::Math::XYZVector normalVec = ROOT::Math::XYZVector (lv3.Py (), -lv3.Px (), 0 .f );
810+ normalVec = ROOT::Math::XYZVector (lv3.Py (), -lv3.Px (), 0 .f );
812811 auto cosThetaStarProduction = normalVec.Dot (threeVecDauCM) / (std::sqrt (threeVecDauCM.Mag2 ()) * std::sqrt (normalVec.Mag2 ()));
813812 hvzero.fill (HIST (" h3vzeropairInvMassME" ), multiplicity, lv3.Pt (), lv3.M (), cosThetaStarProduction, occupancyno);
814813 } else if (activateTHnSparseCosThStarBeam) {
815- ROOT::Math::XYZVector beamVec = ROOT::Math::XYZVector (0 .f , 0 .f , 1 .f );
814+ beamVec = ROOT::Math::XYZVector (0 .f , 0 .f , 1 .f );
816815 auto cosThetaStarBeam = beamVec.Dot (threeVecDauCM) / std::sqrt (threeVecDauCM.Mag2 ());
817816 hvzero.fill (HIST (" h3vzeropairInvMassME" ), multiplicity, lv3.Pt (), lv3.M (), cosThetaStarBeam, occupancyno);
818817 } else if (activateTHnSparseCosThStarRandom) {
819818 auto phiRandom = gRandom ->Uniform (0 .f , constants::math::TwoPI);
820819 auto thetaRandom = gRandom ->Uniform (0 .f , constants::math::PI);
821- ROOT::Math::XYZVector randomVec = ROOT::Math::XYZVector (std::sin (thetaRandom) * std::cos (phiRandom), std::sin (thetaRandom) * std::sin (phiRandom), std::cos (thetaRandom));
820+ randomVec = ROOT::Math::XYZVector (std::sin (thetaRandom) * std::cos (phiRandom), std::sin (thetaRandom) * std::sin (phiRandom), std::cos (thetaRandom));
822821 auto cosThetaStarRandom = randomVec.Dot (threeVecDauCM) / std::sqrt (threeVecDauCM.Mag2 ());
823822 hvzero.fill (HIST (" h3vzeropairInvMassME" ), multiplicity, lv3.Pt (), lv3.M (), cosThetaStarRandom, occupancyno);
824823 }
@@ -836,10 +835,11 @@ struct Kshortlambda {
836835 }
837836 auto occupancyno = c1.trackOccupancyInTimeRange ();
838837 auto occupancyno2 = c2.trackOccupancyInTimeRange ();
838+
839839 if (applyOccupancyCut && (occupancyno < occupancyCut || occupancyno2 < occupancyCut)) {
840- return ;
840+ continue ;
841841 }
842-
842+ TLorentzVector lv1, lv3, lvLambda, lv4, lv5;
843843 for (const auto & [t1, t2] : o2::soa::combinations (o2::soa::CombinationsFullIndexPolicy (tracks1, tracks2))) {
844844
845845 if (t1.size () == 0 || t2.size () == 0 ) {
@@ -897,34 +897,34 @@ struct Kshortlambda {
897897 }
898898 lambdaVec = protonVec + pionVec;
899899
900- TLorentzVector lv1, lv2, lv3, lvLambda, lv4, lv5;
901900 lv1.SetPtEtaPhiM (t1.pt (), t1.eta (), t1.phi (), massK0s);
902901 lvLambda.SetPtEtaPhiM (lambdaVec.pt (), lambdaVec.eta (), lambdaVec.phi (), massLambda);
903902 lv3 = lv1 + lvLambda;
904- ROOT::Math::PxPyPzMVector fourVecDau = ROOT::Math::PxPyPzMVector (daughter1.Px (), daughter1.Py (), daughter1.Pz (), massK0s); // Kshort
905903
906- ROOT::Math::PxPyPzMVector fourVecMother = ROOT::Math::PxPyPzMVector (lv3.Px (), lv3.Py (), lv3.Pz (), lv3.M ()); // mass of KshortKshort pair
904+ fourVecDau = ROOT::Math::PxPyPzMVector (daughter1.Px (), daughter1.Py (), daughter1.Pz (), massK0s); // Kshort
905+ fourVecMother = ROOT::Math::PxPyPzMVector (lv3.Px (), lv3.Py (), lv3.Pz (), lv3.M ()); // mass of KshortKshort pair
907906 ROOT::Math::Boost boost{fourVecMother.BoostToCM ()}; // boost mother to center of mass frame
908- ROOT::Math::PxPyPzMVector fourVecDauCM = boost (fourVecDau); // boost the frame of daughter same as mother
909- ROOT::Math::XYZVector threeVecDauCM = fourVecDauCM.Vect (); // get the 3 vector of daughter in the frame of mother
907+ fourVecDauCM = boost (fourVecDau); // boost the frame of daughter same as mother
908+ threeVecDauCM = fourVecDauCM.Vect (); // get the 3 vector of daughter in the frame of mother
909+
910910 if (std::abs (lv3.Rapidity ()) < 0.5 ) {
911911
912912 if (activateTHnSparseCosThStarHelicity) {
913- ROOT::Math::XYZVector helicityVec = fourVecMother.Vect (); // 3 vector of mother in COM frame
913+ helicityVec = fourVecMother.Vect (); // 3 vector of mother in COM frame
914914 auto cosThetaStarHelicity = helicityVec.Dot (threeVecDauCM) / (std::sqrt (threeVecDauCM.Mag2 ()) * std::sqrt (helicityVec.Mag2 ()));
915915 hvzero.fill (HIST (" h3vzeropairInvMassME" ), multiplicity, lv3.Pt (), lv3.M (), cosThetaStarHelicity, occupancyno);
916916 } else if (activateTHnSparseCosThStarProduction) {
917- ROOT::Math::XYZVector normalVec = ROOT::Math::XYZVector (lv3.Py (), -lv3.Px (), 0 .f );
917+ normalVec = ROOT::Math::XYZVector (lv3.Py (), -lv3.Px (), 0 .f );
918918 auto cosThetaStarProduction = normalVec.Dot (threeVecDauCM) / (std::sqrt (threeVecDauCM.Mag2 ()) * std::sqrt (normalVec.Mag2 ()));
919919 hvzero.fill (HIST (" h3vzeropairInvMassME" ), multiplicity, lv3.Pt (), lv3.M (), cosThetaStarProduction, occupancyno);
920920 } else if (activateTHnSparseCosThStarBeam) {
921- ROOT::Math::XYZVector beamVec = ROOT::Math::XYZVector (0 .f , 0 .f , 1 .f );
921+ beamVec = ROOT::Math::XYZVector (0 .f , 0 .f , 1 .f );
922922 auto cosThetaStarBeam = beamVec.Dot (threeVecDauCM) / std::sqrt (threeVecDauCM.Mag2 ());
923923 hvzero.fill (HIST (" h3vzeropairInvMassME" ), multiplicity, lv3.Pt (), lv3.M (), cosThetaStarBeam, occupancyno);
924924 } else if (activateTHnSparseCosThStarRandom) {
925925 auto phiRandom = gRandom ->Uniform (0 .f , constants::math::TwoPI);
926926 auto thetaRandom = gRandom ->Uniform (0 .f , constants::math::PI);
927- ROOT::Math::XYZVector randomVec = ROOT::Math::XYZVector (std::sin (thetaRandom) * std::cos (phiRandom), std::sin (thetaRandom) * std::sin (phiRandom), std::cos (thetaRandom));
927+ randomVec = ROOT::Math::XYZVector (std::sin (thetaRandom) * std::cos (phiRandom), std::sin (thetaRandom) * std::sin (phiRandom), std::cos (thetaRandom));
928928 auto cosThetaStarRandom = randomVec.Dot (threeVecDauCM) / std::sqrt (threeVecDauCM.Mag2 ());
929929 hvzero.fill (HIST (" h3vzeropairInvMassME" ), multiplicity, lv3.Pt (), lv3.M (), cosThetaStarRandom, occupancyno);
930930 }
0 commit comments