Skip to content

Commit 670bac0

Browse files
omvazqueomvazque
authored
[PWGLF] Fix collisions loop in MC closure (#11124)
1 parent 1a5c088 commit 670bac0

File tree

1 file changed

+69
-97
lines changed

1 file changed

+69
-97
lines changed

‎PWGLF/Tasks/GlobalEventProperties/uccZdc.cxx‎

Lines changed: 69 additions & 97 deletions
Original file line numberDiff line numberDiff line change
@@ -738,50 +738,56 @@ struct UccZdc {
738738
void processMCclosure(aod::McCollisions::iterator const& mccollision, soa::SmallGroups<o2::aod::SimCollisions> const& collisions, o2::aod::BCsRun3 const& /*bcs*/, aod::McParticles const& mcParticles, TheFilteredSimTracks const& simTracks)
739739
{
740740

741-
float rndNum = randPointer->Uniform(0.0, 1.0);
742-
registry.fill(HIST("RandomNumber"), rndNum);
743-
744-
// Half of the statistics for MC closure
745-
if (rndNum >= zEro && rndNum < oneHalf) {
746-
registry.fill(HIST("EvtsDivided"), 0);
747-
//----- MC reconstructed -----//
748-
for (const auto& collision : collisions) {
749-
750-
// To use run-by-run efficiency
751-
const auto& foundBC = collision.foundBC_as<o2::aod::BCsRun3>();
752-
auto efficiency = ccdb->getForTimeStamp<TH1F>(paTH.value, foundBC.timestamp());
753-
// auto efficiency = ccdb->getForRun<TH1F>(paTH.value, foundBC.runNumber());
754-
if (!efficiency) {
755-
LOGF(fatal, "Efficiency object not found!");
756-
}
741+
for (const auto& collision : collisions) {
742+
743+
// To use run-by-run efficiency
744+
const auto& foundBC = collision.foundBC_as<o2::aod::BCsRun3>();
745+
auto efficiency = ccdb->getForTimeStamp<TH1F>(paTH.value, foundBC.timestamp());
746+
if (!efficiency) {
747+
LOGF(fatal, "Efficiency object not found!");
748+
continue;
749+
}
757750

758-
// Event selection
759-
if (!isEventSelected(collision)) {
760-
continue;
761-
}
762-
// MC collision?
763-
if (!collision.has_mcCollision()) {
764-
continue;
765-
}
751+
registry.fill(HIST("hEventCounterMC"), EvCutLabel::All);
752+
753+
// Event selection
754+
if (!isEventSelected(collision)) {
755+
continue;
756+
}
757+
// MC collision?
758+
if (!collision.has_mcCollision()) {
759+
continue;
760+
}
761+
762+
if (std::fabs(mccollision.posZ()) > posZcut) {
763+
continue;
764+
}
765+
registry.fill(HIST("hEventCounterMC"), EvCutLabel::VtxZ);
766+
registry.fill(HIST("zPosMC"), mccollision.posZ());
766767

767-
registry.fill(HIST("T0Ccent"), collision.centFT0C());
768-
registry.fill(HIST("zPos"), collision.posZ());
768+
const auto& cent{collision.centFT0C()};
769+
registry.fill(HIST("nRecColvsCent"), collisions.size(), cent);
770+
registry.fill(HIST("T0Ccent"), cent);
771+
registry.fill(HIST("zPos"), collision.posZ());
769772

770-
const auto& groupedTracks{simTracks.sliceBy(perCollision, collision.globalIndex())};
773+
float rndNum = randPointer->Uniform(0.0, 1.0);
774+
registry.fill(HIST("RandomNumber"), rndNum);
775+
const auto& groupedTracks{simTracks.sliceBy(perCollision, collision.globalIndex())};
771776

777+
if (rndNum >= zEro && rndNum < oneHalf) { // Half of the statistics for MC closure
778+
779+
registry.fill(HIST("EvtsDivided"), 0);
772780
std::vector<float> pTs;
773781
std::vector<float> wIs;
774782
// Calculates the event weight, W_k
775783
for (const auto& track : groupedTracks) {
776784
// Track Selection
777-
if (myTrackSel.IsSelected(track)) {
778-
if (passesDCAxyCut(track)) {
779-
float pt{track.pt()};
780-
double weight{efficiency->GetBinContent(efficiency->FindBin(pt))};
781-
if (weight > 0.) {
782-
pTs.emplace_back(pt);
783-
wIs.emplace_back(weight);
784-
}
785+
if (myTrackSel.IsSelected(track) && passesDCAxyCut(track)) {
786+
float pt{track.pt()};
787+
double weight{efficiency->GetBinContent(efficiency->FindBin(pt))};
788+
if (weight > 0.) {
789+
pTs.emplace_back(pt);
790+
wIs.emplace_back(weight);
785791
}
786792
}
787793
}
@@ -815,13 +821,6 @@ struct UccZdc {
815821
registry.fill(HIST("NchVsFourParCorr"), w1, fourParCorr, denFourParCorr);
816822

817823
//--------------------------- Generated MC ---------------------------
818-
registry.fill(HIST("hEventCounterMC"), EvCutLabel::All);
819-
if (std::fabs(mccollision.posZ()) > posZcut) {
820-
continue;
821-
}
822-
registry.fill(HIST("zPosMC"), mccollision.posZ());
823-
registry.fill(HIST("hEventCounterMC"), EvCutLabel::VtxZ);
824-
825824
std::vector<float> pTsMC;
826825
std::vector<float> wIsMC;
827826
// Calculates the event weight, W_k
@@ -867,69 +866,42 @@ struct UccZdc {
867866
registry.fill(HIST("NchvsTwoParCorrGen"), nchMC, twoParCorrMC, denTwoParCorrMC);
868867
registry.fill(HIST("NchvsThreeParCorrGen"), nchMC, threeParCorrMC, denThreeParCorrMC);
869868
registry.fill(HIST("NchvsFourParCorrGen"), nchMC, fourParCorrMC, denFourParCorrMC);
870-
}
871-
} else { // Correction with the remaining half of the sample
872-
registry.fill(HIST("EvtsDivided"), 1);
873-
//----- MC reconstructed -----//
874-
for (const auto& collision : collisions) {
875-
// Event selection
876-
if (!isEventSelected(collision)) {
877-
continue;
878-
}
879-
// MC collision?
880-
if (!collision.has_mcCollision()) {
881-
continue;
882-
}
883-
884-
registry.fill(HIST("zPos"), collision.posZ());
885-
registry.fill(HIST("nRecColvsCent"), collisions.size(), collision.centFT0C());
886-
887-
const auto& cent{collision.centFT0C()};
888-
registry.fill(HIST("T0Ccent"), cent);
889-
890-
const auto& groupedTracks{simTracks.sliceBy(perCollision, collision.globalIndex())};
869+
} else { // Correction with the remaining half of the sample
870+
registry.fill(HIST("EvtsDivided"), 1);
871+
//----- MC reconstructed -----//
891872
for (const auto& track : groupedTracks) {
892873
// Has MC particle?
893874
if (!track.has_mcParticle()) {
894875
continue;
895876
}
896877
// Track selection
897-
if (myTrackSel.IsSelected(track)) {
898-
if (passesDCAxyCut(track)) {
899-
registry.fill(HIST("Pt_all_ch"), cent, track.pt());
900-
registry.fill(HIST("EtaVsPhi"), track.eta(), track.phi());
901-
902-
const auto& particle{track.mcParticle()};
903-
if (!particle.isPhysicalPrimary()) {
904-
continue;
905-
}
906-
907-
registry.fill(HIST("Pt_ch"), cent, track.pt());
908-
if (particle.pdgCode() == PDG_t::kPiPlus || particle.pdgCode() == PDG_t::kPiMinus) {
909-
registry.fill(HIST("Pt_pi"), cent, track.pt());
910-
} else if (particle.pdgCode() == PDG_t::kKPlus || particle.pdgCode() == PDG_t::kKMinus) {
911-
registry.fill(HIST("Pt_ka"), cent, track.pt());
912-
} else if (particle.pdgCode() == PDG_t::kProton || particle.pdgCode() == PDG_t::kProtonBar) {
913-
registry.fill(HIST("Pt_pr"), cent, track.pt());
914-
} else if (particle.pdgCode() == PDG_t::kSigmaPlus || particle.pdgCode() == PDG_t::kSigmaBarMinus) {
915-
registry.fill(HIST("Pt_sigpos"), cent, track.pt());
916-
} else if (particle.pdgCode() == PDG_t::kSigmaMinus || particle.pdgCode() == PDG_t::kSigmaBarPlus) {
917-
registry.fill(HIST("Pt_signeg"), cent, track.pt());
918-
} else {
919-
registry.fill(HIST("Pt_re"), cent, track.pt());
920-
}
878+
if (myTrackSel.IsSelected(track) && passesDCAxyCut(track)) {
879+
registry.fill(HIST("Pt_all_ch"), cent, track.pt());
880+
registry.fill(HIST("EtaVsPhi"), track.eta(), track.phi());
881+
882+
const auto& particle{track.mcParticle()};
883+
if (!particle.isPhysicalPrimary()) {
884+
continue;
921885
}
922-
}
923-
}
924886

925-
// Generated MC
926-
registry.fill(HIST("hEventCounterMC"), EvCutLabel::All);
927-
if (std::fabs(mccollision.posZ()) > posZcut) {
928-
continue;
887+
registry.fill(HIST("Pt_ch"), cent, track.pt());
888+
if (particle.pdgCode() == PDG_t::kPiPlus || particle.pdgCode() == PDG_t::kPiMinus) {
889+
registry.fill(HIST("Pt_pi"), cent, track.pt());
890+
} else if (particle.pdgCode() == PDG_t::kKPlus || particle.pdgCode() == PDG_t::kKMinus) {
891+
registry.fill(HIST("Pt_ka"), cent, track.pt());
892+
} else if (particle.pdgCode() == PDG_t::kProton || particle.pdgCode() == PDG_t::kProtonBar) {
893+
registry.fill(HIST("Pt_pr"), cent, track.pt());
894+
} else if (particle.pdgCode() == PDG_t::kSigmaPlus || particle.pdgCode() == PDG_t::kSigmaBarMinus) {
895+
registry.fill(HIST("Pt_sigpos"), cent, track.pt());
896+
} else if (particle.pdgCode() == PDG_t::kSigmaMinus || particle.pdgCode() == PDG_t::kSigmaBarPlus) {
897+
registry.fill(HIST("Pt_signeg"), cent, track.pt());
898+
} else {
899+
registry.fill(HIST("Pt_re"), cent, track.pt());
900+
}
901+
}
929902
}
930-
registry.fill(HIST("zPosMC"), mccollision.posZ());
931-
registry.fill(HIST("hEventCounterMC"), EvCutLabel::VtxZ);
932903

904+
//----- MC generated -----//
933905
for (const auto& particle : mcParticles) {
934906
if (particle.eta() < minEta || particle.eta() > maxEta) {
935907
continue;
@@ -956,8 +928,8 @@ struct UccZdc {
956928
registry.fill(HIST("PtMC_re"), cent, particle.pt());
957929
}
958930
}
959-
}
960-
} // Half of statistics for corrections
931+
} // Half of statistics for corrections
932+
} // Collisions loop
961933
}
962934
PROCESS_SWITCH(UccZdc, processMCclosure, "Process MC closure", false);
963935

0 commit comments

Comments
 (0)