Skip to content

Commit 13bc202

Browse files
ananthapadmanabhan18ananthapadmanabhan18
committed
Modified MC and other histograms
1 parent 53e407a commit 13bc202

File tree

1 file changed

+253
-0
lines changed

1 file changed

+253
-0
lines changed

‎PWGUD/Tasks/exclusiveRhoTo4Pi.cxx‎

Lines changed: 253 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -645,7 +645,260 @@ struct ExclusiveRhoTo4Pi {
645645
using MCtracks = soa::Join<aod::UDTracks, aod::UDTracksPID, aod::UDTracksExtra, aod::UDTracksFlags, aod::UDTracksDCA, aod::UDMcTrackLabels>;
646646
using MCCollisions = soa::Join<aod::UDCollisions, aod::SGCollisions, aod::UDCollisionSelExtras, aod::UDCollisionsSels, aod::UDZdcsReduced, aod::UDMcCollsLabels>;
647647

648+
void processMCrec(soa::Filtered<MCCollisions>::iterator const& collision, soa::Filtered<MCtracks> const& tracks)
649+
{
650+
651+
if (collision.has_udMcCollision()) {
652+
return;
653+
}
654+
655+
// Check if the Event is reconstructed in UPC mode and RCT flag
656+
if ((collision.flags() != ifUPC) || (!sgSelector.isCBTHadronOk(collision))) {
657+
return;
658+
}
659+
660+
int runIndex = getRunNumberIndex(collision.runNumber());
661+
662+
histosQA.fill(HIST("Events/selected/UPCmode"), collision.flags());
663+
histosQA.fill(HIST("Events/selected/GapSide"), collision.gapSide());
664+
histosQA.fill(HIST("Events/selected/TrueGapSide"), sgSelector.trueGap(collision, fv0Cut, ft0aCut, ft0cCut, zdcCut));
665+
histosQA.fill(HIST("Events/selected/isCBTOk"), sgSelector.isCBTOk(collision));
666+
histosQA.fill(HIST("Events/selected/isCBTHadronOk"), sgSelector.isCBTHadronOk(collision));
667+
histosQA.fill(HIST("Events/selected/isCBTZdcOk"), sgSelector.isCBTZdcOk(collision));
668+
histosQA.fill(HIST("Events/selected/isCBTHadronZdcOk"), sgSelector.isCBTHadronZdcOk(collision));
669+
histosQA.fill(HIST("Events/selected/vertexX"), collision.posX());
670+
histosQA.fill(HIST("Events/selected/vertexY"), collision.posY());
671+
histosQA.fill(HIST("Events/selected/vertexZ"), collision.posZ());
672+
histosQA.fill(HIST("Events/selected/occupancy"), collision.occupancyInTime());
673+
histosQA.fill(HIST("Events/selected/FV0A"), collision.totalFV0AmplitudeA());
674+
histosQA.fill(HIST("Events/selected/FT0A"), collision.totalFT0AmplitudeA());
675+
histosQA.fill(HIST("Events/selected/FT0C"), collision.totalFT0AmplitudeC());
676+
histosQA.fill(HIST("Events/selected/ZDC_A"), collision.energyCommonZNA());
677+
histosQA.fill(HIST("Events/selected/ZDC_C"), collision.energyCommonZNC());
678+
histosQA.fill(HIST("Events/selected/FDDA"), collision.totalFDDAmplitudeA());
679+
histosQA.fill(HIST("Events/selected/FDDC"), collision.totalFDDAmplitudeC());
680+
681+
std::vector<decltype(tracks.begin())> selectedTracks;
682+
std::vector<decltype(tracks.begin())> selectedPionTracks;
683+
std::vector<decltype(tracks.begin())> selectedPionPlusTracks;
684+
std::vector<decltype(tracks.begin())> selectedPionMinusTracks;
685+
686+
for (const auto& t0 : tracks) {
687+
688+
PxPyPzMVector tVector(t0.px(), t0.py(), t0.pz(), o2::constants::physics::MassPionCharged);
689+
690+
// QA-Tracks before selection
691+
histosQA.fill(HIST("Tracks/all/dcaXY"), t0.tpcChi2NCl());
692+
histosQA.fill(HIST("Tracks/all/dcaZ"), t0.tpcChi2NCl());
693+
histosQA.fill(HIST("Tracks/all/itsChi2NCl"), t0.itsChi2NCl());
694+
histosQA.fill(HIST("Tracks/all/itsChi2"), t0.itsChi2NCl() * t0.itsNCls());
695+
histosQA.fill(HIST("Tracks/all/tpcChi2NCl"), t0.tpcChi2NCl());
696+
histosQA.fill(HIST("Tracks/all/tpcNClsFindable"), t0.tpcNClsFindable());
697+
698+
// PID before track selection
699+
histosPID.fill(HIST("all/tpcSignal"), tVector.P(), t0.tpcSignal());
700+
histosPID.fill(HIST("all/tpcNSigmaPi"), t0.tpcNSigmaPi(), tVector.Pt());
701+
histosPID.fill(HIST("all/tpcNSigmaKa"), t0.tpcNSigmaKa(), tVector.Pt());
702+
histosPID.fill(HIST("all/tpcNSigmaPr"), t0.tpcNSigmaPr(), tVector.Pt());
703+
histosPID.fill(HIST("all/tpcNSigmaEl"), t0.tpcNSigmaEl(), tVector.Pt());
704+
histosPID.fill(HIST("all/tpcNSigmaMu"), t0.tpcNSigmaMu(), tVector.Pt());
705+
histosPID.fill(HIST("all/tofBeta"), tVector.P(), t0.beta());
706+
histosPID.fill(HIST("all/tofNSigmaPi"), t0.tofNSigmaPi(), tVector.Pt());
707+
histosPID.fill(HIST("all/tofNSigmaKa"), t0.tofNSigmaKa(), tVector.Pt());
708+
histosPID.fill(HIST("all/tofNSigmaPr"), t0.tofNSigmaPr(), tVector.Pt());
709+
histosPID.fill(HIST("all/tofNSigmaEl"), t0.tofNSigmaEl(), tVector.Pt());
710+
histosPID.fill(HIST("all/tofNSigmaMu"), t0.tofNSigmaMu(), tVector.Pt());
711+
712+
// Kinematics for all particles before selection
713+
histosKin.fill(HIST("all"), tVector.Pt(), tVector.Eta(), tVector.Phi());
714+
715+
// Selecting good tracks
716+
if (!isSelectedTrack(t0, pTcut, etaCut, dcaXYcut, dcaZcut, useITStracksOnly, useTPCtracksOnly, itsChi2NClsCut, tpcChi2NClsCut, tpcNClsFindableCut)) {
717+
continue;
718+
}
719+
if (!t0.has_udMcParticle()) {
720+
continue;
721+
}
722+
723+
// QA-Tracks after selection
724+
histosQA.fill(HIST("Tracks/selected/dcaXY"), t0.dcaXY());
725+
histosQA.fill(HIST("Tracks/selected/dcaZ"), t0.dcaZ());
726+
histosQA.fill(HIST("Tracks/selected/itsChi2NCl"), t0.itsChi2NCl());
727+
histosQA.fill(HIST("Tracks/selected/itsChi2"), t0.itsChi2NCl() * t0.itsNCls());
728+
histosQA.fill(HIST("Tracks/selected/tpcChi2NCl"), t0.tpcChi2NCl());
729+
histosQA.fill(HIST("Tracks/selected/tpcNClsFindable"), t0.tpcNClsFindable());
730+
731+
// PID after track selection before selecting pions
732+
histosPID.fill(HIST("selected/tpcSignal"), tVector.P(), t0.tpcSignal());
733+
histosPID.fill(HIST("selected/tpcNSigmaPi"), t0.tpcNSigmaPi(), tVector.Pt());
734+
histosPID.fill(HIST("selected/tpcNSigmaKa"), t0.tpcNSigmaKa(), tVector.Pt());
735+
histosPID.fill(HIST("selected/tpcNSigmaPr"), t0.tpcNSigmaPr(), tVector.Pt());
736+
histosPID.fill(HIST("selected/tpcNSigmaEl"), t0.tpcNSigmaEl(), tVector.Pt());
737+
histosPID.fill(HIST("selected/tpcNSigmaMu"), t0.tpcNSigmaMu(), tVector.Pt());
738+
histosPID.fill(HIST("selected/tofBeta"), tVector.P(), t0.beta());
739+
histosPID.fill(HIST("selected/tofNSigmaPi"), t0.tofNSigmaPi(), tVector.Pt());
740+
histosPID.fill(HIST("selected/tofNSigmaKa"), t0.tofNSigmaKa(), tVector.Pt());
741+
histosPID.fill(HIST("selected/tofNSigmaPr"), t0.tofNSigmaPr(), tVector.Pt());
742+
histosPID.fill(HIST("selected/tofNSigmaEl"), t0.tofNSigmaEl(), tVector.Pt());
743+
histosPID.fill(HIST("selected/tofNSigmaMu"), t0.tofNSigmaMu(), tVector.Pt());
744+
745+
// Kinematics for all particles after track selection before selecting pions
746+
histosKin.fill(HIST("selected"), tVector.Pt(), tVector.Eta(), tVector.Phi());
747+
748+
selectedTracks.push_back(t0);
749+
if (ifPion(t0, useTOF, nSigmaTPCcut, nSigmaTOFcut)) {
750+
751+
selectedPionTracks.push_back(t0);
752+
753+
// QA-Tracks after selecting pions
754+
histosQA.fill(HIST("Tracks/pions/dcaXY"), t0.dcaXY());
755+
histosQA.fill(HIST("Tracks/pions/dcaZ"), t0.dcaZ());
756+
histosQA.fill(HIST("Tracks/pions/itsChi2NCl"), t0.itsChi2NCl());
757+
histosQA.fill(HIST("Tracks/pions/itsChi2"), t0.itsChi2NCl() * t0.itsNCls());
758+
histosQA.fill(HIST("Tracks/pions/tpcChi2NCl"), t0.tpcChi2NCl());
759+
histosQA.fill(HIST("Tracks/pions/tpcNClsFindable"), t0.tpcNClsFindable());
760+
761+
// PID after selecting pions
762+
histosPID.fill(HIST("pions/tpcSignal"), tVector.P(), t0.tpcSignal());
763+
histosPID.fill(HIST("pions/tpcNSigmaPi"), t0.tpcNSigmaPi(), tVector.Pt());
764+
histosPID.fill(HIST("pions/tpcNSigmaKa"), t0.tpcNSigmaKa(), tVector.Pt());
765+
histosPID.fill(HIST("pions/tpcNSigmaPr"), t0.tpcNSigmaPr(), tVector.Pt());
766+
histosPID.fill(HIST("pions/tpcNSigmaEl"), t0.tpcNSigmaEl(), tVector.Pt());
767+
histosPID.fill(HIST("pions/tpcNSigmaMu"), t0.tpcNSigmaMu(), tVector.Pt());
768+
histosPID.fill(HIST("pions/tofBeta"), tVector.P(), t0.beta());
769+
histosPID.fill(HIST("pions/tofNSigmaPi"), t0.tofNSigmaPi(), tVector.Pt());
770+
histosPID.fill(HIST("pions/tofNSigmaKa"), t0.tofNSigmaKa(), tVector.Pt());
771+
histosPID.fill(HIST("pions/tofNSigmaPr"), t0.tofNSigmaPr(), tVector.Pt());
772+
histosPID.fill(HIST("pions/tofNSigmaEl"), t0.tofNSigmaEl(), tVector.Pt());
773+
histosPID.fill(HIST("pions/tofNSigmaMu"), t0.tofNSigmaMu(), tVector.Pt());
774+
775+
// Kinematics for pions
776+
histosKin.fill(HIST("pions"), tVector.Pt(), tVector.Eta(), tVector.Phi());
777+
778+
if (t0.sign() == 1) {
779+
selectedPionPlusTracks.push_back(t0);
780+
}
781+
if (t0.sign() == -1) {
782+
selectedPionMinusTracks.push_back(t0);
783+
}
784+
} // End of Selection PID Pion
785+
} // End of loop over tracks
786+
787+
int numSelectedPionTracks = static_cast<int>(selectedPionTracks.size());
788+
int numPiPlusTracks = static_cast<int>(selectedPionPlusTracks.size());
789+
int numPionMinusTracks = static_cast<int>(selectedPionMinusTracks.size());
790+
791+
// event should have exactly 4 pions
792+
if (numSelectedPionTracks != numFourPionTracks) {
793+
return;
794+
}
795+
796+
// Selecting Events with net charge = 0
797+
if (numPionMinusTracks == numPiMinus && numPiPlusTracks == numPiPlus) {
798+
799+
// QA-Events-4pion
800+
histosQA.fill(HIST("Events/4pion/UPCmode"), collision.flags());
801+
histosQA.fill(HIST("Events/4pion/GapSide"), collision.gapSide());
802+
histosQA.fill(HIST("Events/4pion/TrueGapSide"), sgSelector.trueGap(collision, fv0Cut, ft0aCut, ft0cCut, zdcCut));
803+
histosQA.fill(HIST("Events/4pion/isCBTOk"), sgSelector.isCBTOk(collision));
804+
histosQA.fill(HIST("Events/4pion/isCBTHadronOk"), sgSelector.isCBTHadronOk(collision));
805+
histosQA.fill(HIST("Events/4pion/isCBTZdcOk"), sgSelector.isCBTZdcOk(collision));
806+
histosQA.fill(HIST("Events/4pion/isCBTHadronZdcOk"), sgSelector.isCBTHadronZdcOk(collision));
807+
histosQA.fill(HIST("Events/4pion/vertexX"), collision.posX());
808+
histosQA.fill(HIST("Events/4pion/vertexY"), collision.posY());
809+
histosQA.fill(HIST("Events/4pion/vertexZ"), collision.posZ());
810+
histosQA.fill(HIST("Events/4pion/occupancy"), collision.occupancyInTime());
811+
histosQA.fill(HIST("Events/4pion/FV0A"), collision.totalFV0AmplitudeA());
812+
histosQA.fill(HIST("Events/4pion/FT0A"), collision.totalFT0AmplitudeA());
813+
histosQA.fill(HIST("Events/4pion/FT0C"), collision.totalFT0AmplitudeC());
814+
histosQA.fill(HIST("Events/4pion/ZDC_A"), collision.energyCommonZNA());
815+
histosQA.fill(HIST("Events/4pion/ZDC_C"), collision.energyCommonZNC());
816+
histosQA.fill(HIST("Events/4pion/FDDA"), collision.totalFDDAmplitudeA());
817+
histosQA.fill(HIST("Events/4pion/FDDC"), collision.totalFDDAmplitudeC());
818+
819+
for (int i = 0; i < numFourPionTracks; i++) {
820+
PxPyPzMVector tVector(selectedPionTracks[i].px(), selectedPionTracks[i].py(), selectedPionTracks[i].pz(), o2::constants::physics::MassPionCharged);
821+
// Tracks QA for all four pions
822+
histosQA.fill(HIST("Tracks/pions-from-4pi/dcaXY"), selectedPionTracks[i].dcaXY());
823+
histosQA.fill(HIST("Tracks/pions-from-4pi/dcaZ"), selectedPionTracks[i].dcaZ());
824+
histosQA.fill(HIST("Tracks/pions-from-4pi/itsChi2NCl"), selectedPionTracks[i].itsChi2NCl());
825+
histosQA.fill(HIST("Tracks/pions-from-4pi/itsChi2"), selectedPionTracks[i].itsChi2NCl() * selectedPionTracks[i].itsNCls());
826+
histosQA.fill(HIST("Tracks/pions-from-4pi/tpcChi2NCl"), selectedPionTracks[i].tpcChi2NCl());
827+
histosQA.fill(HIST("Tracks/pions-from-4pi/tpcNClsFindable"), selectedPionTracks[i].tpcNClsFindable());
828+
// PID for all four pions
829+
histosPID.fill(HIST("pions-from-4pi/tpcSignal"), tVector.P(), selectedPionTracks[i].tpcSignal());
830+
histosPID.fill(HIST("pions-from-4pi/tpcNSigmaPi"), selectedPionTracks[i].tpcNSigmaPi(), tVector.Pt());
831+
histosPID.fill(HIST("pions-from-4pi/tpcNSigmaKa"), selectedPionTracks[i].tpcNSigmaKa(), tVector.Pt());
832+
histosPID.fill(HIST("pions-from-4pi/tpcNSigmaPr"), selectedPionTracks[i].tpcNSigmaPr(), tVector.Pt());
833+
histosPID.fill(HIST("pions-from-4pi/tpcNSigmaEl"), selectedPionTracks[i].tpcNSigmaEl(), tVector.Pt());
834+
histosPID.fill(HIST("pions-from-4pi/tpcNSigmaMu"), selectedPionTracks[i].tpcNSigmaMu(), tVector.Pt());
835+
histosPID.fill(HIST("pions-from-4pi/tofBeta"), tVector.P(), selectedPionTracks[i].beta());
836+
histosPID.fill(HIST("pions-from-4pi/tofNSigmaPi"), selectedPionTracks[i].tofNSigmaPi(), tVector.Pt());
837+
histosPID.fill(HIST("pions-from-4pi/tofNSigmaKa"), selectedPionTracks[i].tofNSigmaKa(), tVector.Pt());
838+
histosPID.fill(HIST("pions-from-4pi/tofNSigmaPr"), selectedPionTracks[i].tofNSigmaPr(), tVector.Pt());
839+
histosPID.fill(HIST("pions-from-4pi/tofNSigmaEl"), selectedPionTracks[i].tofNSigmaEl(), tVector.Pt());
840+
histosPID.fill(HIST("pions-from-4pi/tofNSigmaMu"), selectedPionTracks[i].tofNSigmaMu(), tVector.Pt());
841+
}
842+
843+
PxPyPzMVector p1(selectedPionPlusTracks[0].px(), selectedPionPlusTracks[0].py(), selectedPionPlusTracks[0].pz(), o2::constants::physics::MassPionCharged);
844+
PxPyPzMVector p2(selectedPionPlusTracks[1].px(), selectedPionPlusTracks[1].py(), selectedPionPlusTracks[1].pz(), o2::constants::physics::MassPionCharged);
845+
PxPyPzMVector p3(selectedPionMinusTracks[0].px(), selectedPionMinusTracks[0].py(), selectedPionMinusTracks[0].pz(), o2::constants::physics::MassPionCharged);
846+
PxPyPzMVector p4(selectedPionMinusTracks[1].px(), selectedPionMinusTracks[1].py(), selectedPionMinusTracks[1].pz(), o2::constants::physics::MassPionCharged);
847+
848+
// Kinematics for pions from 4 pion events
849+
histosKin.fill(HIST("pions-from-4pi"), p1.Pt(), p1.Eta(), p1.Phi(), p1.Rapidity());
850+
histosKin.fill(HIST("pions-from-4pi"), p2.Pt(), p2.Eta(), p2.Phi(), p2.Rapidity());
851+
histosKin.fill(HIST("pions-from-4pi"), p3.Pt(), p3.Eta(), p3.Phi(), p3.Rapidity());
852+
histosKin.fill(HIST("pions-from-4pi"), p4.Pt(), p4.Eta(), p4.Phi(), p4.Rapidity());
853+
854+
PxPyPzMVector p1234 = p1 + p2 + p3 + p4;
855+
PxPyPzMVector p13 = p1 + p3;
856+
PxPyPzMVector p14 = p1 + p4;
857+
PxPyPzMVector p23 = p2 + p3;
858+
PxPyPzMVector p24 = p2 + p4;
859+
860+
// Two Pion Mass combinations
861+
histos4piKin.fill(HIST("two-pion-mass"), p13.M(), p14.M(), p23.M(), p24.M());
862+
863+
double fourPiPhiPair1 = collinSoperPhi(p13, p1234);
864+
double fourPiPhiPair2 = collinSoperPhi(p14, p1234);
865+
double fourPiPhiPair3 = collinSoperPhi(p23, p1234);
866+
double fourPiPhiPair4 = collinSoperPhi(p24, p1234);
867+
868+
double fourPiCosThetaPair1 = collinSoperCosTheta(p13, p1234);
869+
double fourPiCosThetaPair2 = collinSoperCosTheta(p14, p1234);
870+
double fourPiCosThetaPair3 = collinSoperCosTheta(p23, p1234);
871+
double fourPiCosThetaPair4 = collinSoperCosTheta(p24, p1234);
872+
873+
double mDiff13 = std::abs((p13.M() - mRho0));
874+
double mDiff14 = std::abs((p14.M() - mRho0));
875+
double mDiff23 = std::abs((p23.M() - mRho0));
876+
double mDiff24 = std::abs((p24.M() - mRho0));
877+
if ((mDiff13 < mDiff14) && (mDiff13 < mDiff23) && (mDiff13 < mDiff24)) {
878+
histos4piKin.fill(HIST("zero-charge"), p1234.Pt(), p1234.Eta(), p1234.Phi(), p1234.Rapidity(), p1234.M(), fourPiCosThetaPair1, fourPiPhiPair1, runIndex);
879+
} else if ((mDiff14 < mDiff13) && (mDiff14 < mDiff23) && (mDiff14 < mDiff24)) {
880+
histos4piKin.fill(HIST("zero-charge"), p1234.Pt(), p1234.Eta(), p1234.Phi(), p1234.Rapidity(), p1234.M(), fourPiCosThetaPair2, fourPiPhiPair2, runIndex);
881+
} else if ((mDiff23 < mDiff13) && (mDiff23 < mDiff14) && (mDiff23 < mDiff24)) {
882+
histos4piKin.fill(HIST("zero-charge"), p1234.Pt(), p1234.Eta(), p1234.Phi(), p1234.Rapidity(), p1234.M(), fourPiCosThetaPair3, fourPiPhiPair3, runIndex);
883+
} else if ((mDiff24 < mDiff13) && (mDiff24 < mDiff14) && (mDiff24 < mDiff23)) {
884+
histos4piKin.fill(HIST("zero-charge"), p1234.Pt(), p1234.Eta(), p1234.Phi(), p1234.Rapidity(), p1234.M(), fourPiCosThetaPair4, fourPiPhiPair4, runIndex);
885+
}
886+
} // End of Analysis for 0 charge events
887+
888+
// Selecting Events with net charge != 0 for estimation of background
889+
if (numPionMinusTracks != numPiMinus && numPiPlusTracks != numPiPlus) {
890+
PxPyPzMVector p1(selectedPionTracks[0].px(), selectedPionTracks[0].py(), selectedPionTracks[0].pz(), o2::constants::physics::MassPionCharged);
891+
PxPyPzMVector p2(selectedPionTracks[1].px(), selectedPionTracks[1].py(), selectedPionTracks[1].pz(), o2::constants::physics::MassPionCharged);
892+
PxPyPzMVector p3(selectedPionTracks[2].px(), selectedPionTracks[2].py(), selectedPionTracks[2].pz(), o2::constants::physics::MassPionCharged);
893+
PxPyPzMVector p4(selectedPionTracks[3].px(), selectedPionTracks[3].py(), selectedPionTracks[3].pz(), o2::constants::physics::MassPionCharged);
894+
PxPyPzMVector p1234 = p1 + p2 + p3 + p4;
895+
// Kinematics for 4 pion system from non 0 charge events
896+
histos4piKin.fill(HIST("non-zero-charge"), p1234.Pt(), p1234.Eta(), p1234.Phi(), p1234.Rapidity(), p1234.M(), runIndex);
897+
} // End of Analysis for non 0 charge events
898+
} // End of 4 Pion Analysis Process function for Pass5 MC
899+
648900
PROCESS_SWITCH(ExclusiveRhoTo4Pi, processData, "Data Analysis Function", true);
901+
PROCESS_SWITCH(ExclusiveRhoTo4Pi, processMCrec, "MC reconstructed Analysis Function", false);
649902
PROCESS_SWITCH(ExclusiveRhoTo4Pi, processEventCounter, "Event Counter Function", true);
650903
PROCESS_SWITCH(ExclusiveRhoTo4Pi, processTrackCounter, "Track Counter Function", true);
651904

0 commit comments

Comments
 (0)