diff --git a/PWGEM/Dilepton/Core/DielectronCut.h b/PWGEM/Dilepton/Core/DielectronCut.h index 64bdd069f18..5c629c3f45d 100644 --- a/PWGEM/Dilepton/Core/DielectronCut.h +++ b/PWGEM/Dilepton/Core/DielectronCut.h @@ -170,7 +170,7 @@ class DielectronCut : public TNamed template bool IsSelectedTrack(TTrack const& track) const { - if (!track.hasITS()) { + if (!track.hasITS() || !track.hasTPC()) { return false; } @@ -225,34 +225,32 @@ class DielectronCut : public TNamed } } - if (!mIncludeITSsa && (!track.hasITS() || !track.hasTPC())) { // track has to be ITS-TPC matched track + // if (!mIncludeITSsa && (!track.hasITS() || !track.hasTPC())) { // track has to be ITS-TPC matched track + // return false; + // } + + // if ((track.hasITS() && !track.hasTPC() && !track.hasTRD() && !track.hasTOF()) && track.pt() > mMaxPtITSsa) { // ITSsa + // return false; + // } + + // TPC cuts + if (!IsSelectedTrack(track, DielectronCuts::kTPCNCls)) { return false; } - - if ((track.hasITS() && !track.hasTPC() && !track.hasTRD() && !track.hasTOF()) && track.pt() > mMaxPtITSsa) { // ITSsa + if (!IsSelectedTrack(track, DielectronCuts::kTPCCrossedRows)) { return false; } - - // TPC cuts - if (track.hasTPC()) { - if (!IsSelectedTrack(track, DielectronCuts::kTPCNCls)) { - return false; - } - if (!IsSelectedTrack(track, DielectronCuts::kTPCCrossedRows)) { - return false; - } - if (!IsSelectedTrack(track, DielectronCuts::kTPCCrossedRowsOverNCls)) { - return false; - } - if (!IsSelectedTrack(track, DielectronCuts::kTPCFracSharedClusters)) { - return false; - } - if (!IsSelectedTrack(track, DielectronCuts::kRelDiffPin)) { - return false; - } - if (!IsSelectedTrack(track, DielectronCuts::kTPCChi2NDF)) { - return false; - } + if (!IsSelectedTrack(track, DielectronCuts::kTPCCrossedRowsOverNCls)) { + return false; + } + if (!IsSelectedTrack(track, DielectronCuts::kTPCFracSharedClusters)) { + return false; + } + if (!IsSelectedTrack(track, DielectronCuts::kRelDiffPin)) { + return false; + } + if (!IsSelectedTrack(track, DielectronCuts::kTPCChi2NDF)) { + return false; } if (mApplyPF && !IsSelectedTrack(track, DielectronCuts::kPrefilter)) { diff --git a/PWGEM/Dilepton/Core/Dilepton.h b/PWGEM/Dilepton/Core/Dilepton.h index f0b0aa28156..9291314d00c 100644 --- a/PWGEM/Dilepton/Core/Dilepton.h +++ b/PWGEM/Dilepton/Core/Dilepton.h @@ -209,8 +209,8 @@ struct Dilepton { Configurable cfg_max_chi2tpc{"cfg_max_chi2tpc", 4.0, "max chi2/NclsTPC"}; Configurable cfg_max_chi2its{"cfg_max_chi2its", 5.0, "max chi2/NclsITS"}; Configurable cfg_max_chi2tof{"cfg_max_chi2tof", 1e+10, "max chi2 TOF"}; - Configurable cfg_max_dcaxy{"cfg_max_dcaxy", 0.2, "max dca XY for single track in cm"}; - Configurable cfg_max_dcaz{"cfg_max_dcaz", 0.2, "max dca Z for single track in cm"}; + Configurable cfg_max_dcaxy{"cfg_max_dcaxy", 1.0, "max dca XY for single track in cm"}; + Configurable cfg_max_dcaz{"cfg_max_dcaz", 1.0, "max dca Z for single track in cm"}; Configurable cfg_require_itsib_any{"cfg_require_itsib_any", false, "flag to require ITS ib any hits"}; Configurable cfg_require_itsib_1st{"cfg_require_itsib_1st", true, "flag to require ITS ib 1st hit"}; Configurable cfg_min_its_cluster_size{"cfg_min_its_cluster_size", 0.f, "min ITS cluster size"}; @@ -237,8 +237,6 @@ struct Dilepton { Configurable cfg_min_pin_pirejTPC{"cfg_min_pin_pirejTPC", 0.f, "min. pin for pion rejection in TPC"}; Configurable cfg_max_pin_pirejTPC{"cfg_max_pin_pirejTPC", 1e+10, "max. pin for pion rejection in TPC"}; Configurable enableTTCA{"enableTTCA", true, "Flag to enable or disable TTCA"}; - Configurable includeITSsa{"includeITSsa", false, "Flag to enable ITSsa tracks"}; - Configurable cfg_max_pt_track_ITSsa{"cfg_max_pt_track_ITSsa", 0.15, "max pt for ITSsa tracks"}; // configuration for PID ML Configurable> onnxFileNames{"onnxFileNames", std::vector{"filename"}, "ONNX file names for each bin (if not from CCDB full path)"}; @@ -730,7 +728,6 @@ struct Dilepton { fDielectronCut.RequireITSib1st(dielectroncuts.cfg_require_itsib_1st); fDielectronCut.SetChi2TOF(0, dielectroncuts.cfg_max_chi2tof); fDielectronCut.SetRelDiffPin(dielectroncuts.cfg_min_rel_diff_pin, dielectroncuts.cfg_max_rel_diff_pin); - fDielectronCut.IncludeITSsa(dielectroncuts.includeITSsa, dielectroncuts.cfg_max_pt_track_ITSsa); // for eID fDielectronCut.SetPIDScheme(dielectroncuts.cfg_pid_scheme); @@ -1143,7 +1140,8 @@ struct Dilepton { SliceCache cache; Preslice perCollision_electron = aod::emprimaryelectron::emeventId; - Filter trackFilter_electron = dielectroncuts.cfg_min_pt_track < o2::aod::track::pt && dielectroncuts.cfg_min_eta_track < o2::aod::track::eta && o2::aod::track::eta < dielectroncuts.cfg_max_eta_track && nabs(o2::aod::track::dcaXY) < dielectroncuts.cfg_max_dcaxy && nabs(o2::aod::track::dcaZ) < dielectroncuts.cfg_max_dcaz; + Filter trackFilter_electron = dielectroncuts.cfg_min_pt_track < o2::aod::track::pt && dielectroncuts.cfg_min_eta_track < o2::aod::track::eta && o2::aod::track::eta < dielectroncuts.cfg_max_eta_track && nabs(o2::aod::track::dcaXY) < dielectroncuts.cfg_max_dcaxy && nabs(o2::aod::track::dcaZ) < dielectroncuts.cfg_max_dcaz && o2::aod::track::itsChi2NCl < dielectroncuts.cfg_max_chi2its && o2::aod::track::tpcChi2NCl < dielectroncuts.cfg_max_chi2tpc; + Filter pidFilter_electron = dielectroncuts.cfg_min_TPCNsigmaEl < o2::aod::pidtpc::tpcNSigmaEl && o2::aod::pidtpc::tpcNSigmaEl < dielectroncuts.cfg_max_TPCNsigmaEl; Filter ttcaFilter_electron = ifnode(dielectroncuts.enableTTCA.node(), o2::aod::emprimaryelectron::isAssociatedToMPC == true || o2::aod::emprimaryelectron::isAssociatedToMPC == false, o2::aod::emprimaryelectron::isAssociatedToMPC == true); Filter prefilter_derived_electron = ifnode(dielectroncuts.cfg_apply_cuts_from_prefilter_derived.node() && dielectroncuts.cfg_prefilter_bits_derived.node() >= static_cast(1), ifnode((dielectroncuts.cfg_prefilter_bits_derived.node() & static_cast(1 << int(o2::aod::pwgem::dilepton::utils::pairutil::DileptonPrefilterBitDerived::kMee))) > static_cast(0), (o2::aod::emprimaryelectron::pfbderived & static_cast(1 << int(o2::aod::pwgem::dilepton::utils::pairutil::DileptonPrefilterBitDerived::kMee))) <= static_cast(0), true) && diff --git a/PWGEM/Dilepton/Core/DileptonMC.h b/PWGEM/Dilepton/Core/DileptonMC.h index e323fee60a9..b20a9be3820 100644 --- a/PWGEM/Dilepton/Core/DileptonMC.h +++ b/PWGEM/Dilepton/Core/DileptonMC.h @@ -211,8 +211,8 @@ struct DileptonMC { Configurable cfg_max_chi2tpc{"cfg_max_chi2tpc", 4.0, "max chi2/NclsTPC"}; Configurable cfg_max_chi2its{"cfg_max_chi2its", 5.0, "max chi2/NclsITS"}; Configurable cfg_max_chi2tof{"cfg_max_chi2tof", 1e+10, "max chi2 TOF"}; - Configurable cfg_max_dcaxy{"cfg_max_dcaxy", 0.2, "max dca XY for single track in cm"}; - Configurable cfg_max_dcaz{"cfg_max_dcaz", 0.2, "max dca Z for single track in cm"}; + Configurable cfg_max_dcaxy{"cfg_max_dcaxy", 1.0, "max dca XY for single track in cm"}; + Configurable cfg_max_dcaz{"cfg_max_dcaz", 1.0, "max dca Z for single track in cm"}; Configurable cfg_require_itsib_any{"cfg_require_itsib_any", false, "flag to require ITS ib any hits"}; Configurable cfg_require_itsib_1st{"cfg_require_itsib_1st", true, "flag to require ITS ib 1st hit"}; Configurable cfg_min_its_cluster_size{"cfg_min_its_cluster_size", 0.f, "min ITS cluster size"}; @@ -239,8 +239,6 @@ struct DileptonMC { Configurable cfg_min_pin_pirejTPC{"cfg_min_pin_pirejTPC", 0.f, "min. pin for pion rejection in TPC"}; Configurable cfg_max_pin_pirejTPC{"cfg_max_pin_pirejTPC", 1e+10, "max. pin for pion rejection in TPC"}; Configurable enableTTCA{"enableTTCA", true, "Flag to enable or disable TTCA"}; - Configurable includeITSsa{"includeITSsa", false, "Flag to enable ITSsa tracks"}; - Configurable cfg_max_pt_track_ITSsa{"cfg_max_pt_track_ITSsa", 0.15, "max pt for ITSsa tracks"}; // configuration for PID ML Configurable> onnxFileNames{"onnxFileNames", std::vector{"filename"}, "ONNX file names for each bin (if not from CCDB full path)"}; @@ -737,7 +735,6 @@ struct DileptonMC { fDielectronCut.RequireITSib1st(dielectroncuts.cfg_require_itsib_1st); fDielectronCut.SetChi2TOF(0.0, dielectroncuts.cfg_max_chi2tof); fDielectronCut.SetRelDiffPin(dielectroncuts.cfg_min_rel_diff_pin, dielectroncuts.cfg_max_rel_diff_pin); - fDielectronCut.IncludeITSsa(dielectroncuts.includeITSsa, dielectroncuts.cfg_max_pt_track_ITSsa); // for eID fDielectronCut.SetPIDScheme(dielectroncuts.cfg_pid_scheme); @@ -2484,7 +2481,8 @@ struct DileptonMC { SliceCache cache; Preslice perCollision_electron = aod::emprimaryelectron::emeventId; - Filter trackFilter_electron = nabs(o2::aod::track::dcaXY) < dielectroncuts.cfg_max_dcaxy && nabs(o2::aod::track::dcaZ) < dielectroncuts.cfg_max_dcaz; + Filter trackFilter_electron = nabs(o2::aod::track::dcaXY) < dielectroncuts.cfg_max_dcaxy && nabs(o2::aod::track::dcaZ) < dielectroncuts.cfg_max_dcaz && o2::aod::track::itsChi2NCl < dielectroncuts.cfg_max_chi2its && o2::aod::track::tpcChi2NCl < dielectroncuts.cfg_max_chi2tpc; + Filter pidFilter_electron = dielectroncuts.cfg_min_TPCNsigmaEl < o2::aod::pidtpc::tpcNSigmaEl && o2::aod::pidtpc::tpcNSigmaEl < dielectroncuts.cfg_max_TPCNsigmaEl; Filter ttcaFilter_electron = ifnode(dielectroncuts.enableTTCA.node(), o2::aod::emprimaryelectron::isAssociatedToMPC == true || o2::aod::emprimaryelectron::isAssociatedToMPC == false, o2::aod::emprimaryelectron::isAssociatedToMPC == true); Filter prefilter_derived_electron = ifnode(dielectroncuts.cfg_apply_cuts_from_prefilter_derived.node() && dielectroncuts.cfg_prefilter_bits_derived.node() >= static_cast(1), ifnode((dielectroncuts.cfg_prefilter_bits_derived.node() & static_cast(1 << int(o2::aod::pwgem::dilepton::utils::pairutil::DileptonPrefilterBitDerived::kMee))) > static_cast(0), (o2::aod::emprimaryelectron::pfbderived & static_cast(1 << int(o2::aod::pwgem::dilepton::utils::pairutil::DileptonPrefilterBitDerived::kMee))) <= static_cast(0), true) && diff --git a/PWGEM/Dilepton/TableProducer/createEMEventDilepton.cxx b/PWGEM/Dilepton/TableProducer/createEMEventDilepton.cxx index a4ab164a4b2..4c92285441c 100644 --- a/PWGEM/Dilepton/TableProducer/createEMEventDilepton.cxx +++ b/PWGEM/Dilepton/TableProducer/createEMEventDilepton.cxx @@ -14,7 +14,6 @@ /// \author Daiki Sekihata, daiki.sekihata@cern.ch #include "PWGEM/Dilepton/DataModel/dileptonTables.h" -#include "PWGEM/PhotonMeson/DataModel/gammaTables.h" #include "Common/CCDB/EventSelectionParams.h" #include "Common/DataModel/Centrality.h" @@ -294,16 +293,13 @@ struct CreateEMEventDilepton { PROCESS_SWITCH(CreateEMEventDilepton, processDummy, "processDummy", true); }; struct AssociateDileptonToEMEvent { - Produces v0kfeventid; Produces prmeleventid; Produces prmmueventid; Produces prmtrackeventid; - Preslice perCollision_pcm = aod::v0photonkf::collisionId; PresliceUnsorted perCollision_el = aod::emprimaryelectron::collisionId; PresliceUnsorted perCollision_mu = aod::emprimarymuon::collisionId; Preslice perCollision_track = aod::emprimarytrack::collisionId; - // Preslice perCollision_track = aod::track::collisionId; void init(o2::framework::InitContext&) {} @@ -323,11 +319,6 @@ struct AssociateDileptonToEMEvent { // This struct is for both data and MC. // Note that reconstructed collisions without mc collisions are already rejected in CreateEMEventDilepton in MC. - // void processPCM(aod::EMEvents const& collisions, aod::V0PhotonsKF const& photons) - // { - // fillEventId(collisions, photons, v0kfeventid, perCollision_pcm); - // } - void processElectron(aod::EMEvents const& collisions, aod::EMPrimaryElectrons const& tracks) { fillEventId(collisions, tracks, prmeleventid, perCollision_el); @@ -345,7 +336,6 @@ struct AssociateDileptonToEMEvent { void processDummy(aod::EMEvents const&) {} - // PROCESS_SWITCH(AssociateDileptonToEMEvent, processPCM, "process indexing for PCM", false); PROCESS_SWITCH(AssociateDileptonToEMEvent, processElectron, "process indexing for electrons", false); PROCESS_SWITCH(AssociateDileptonToEMEvent, processFwdMuon, "process indexing for forward muons", false); PROCESS_SWITCH(AssociateDileptonToEMEvent, processChargedTrack, "process indexing for charged tracks", false); diff --git a/PWGEM/Dilepton/Tasks/createResolutionMap.cxx b/PWGEM/Dilepton/Tasks/createResolutionMap.cxx index 5df5855d171..891667627a1 100644 --- a/PWGEM/Dilepton/Tasks/createResolutionMap.cxx +++ b/PWGEM/Dilepton/Tasks/createResolutionMap.cxx @@ -179,6 +179,7 @@ struct CreateResolutionMap { Configurable requireMFTHitMap{"requireMFTHitMap", false, "flag to require MFT hit map"}; Configurable> requiredMFTDisks{"requiredMFTDisks", std::vector{4}, "hit map on MFT disks [0,1,2,3,4]. logical-OR of each double-sided disk"}; Configurable matchingZ{"matchingZ", -77.5, "z position where matching is performed"}; + Configurable cfgApplyPreselectionInBestMatch{"cfgApplyPreselectionInBestMatch", false, "flag to apply preselection in find best match function"}; } muoncuts; HistogramRegistry registry{"registry", {}, OutputObjHandlingPolicy::AnalysisObject}; @@ -556,7 +557,7 @@ struct CreateResolutionMap { if (muon.chi2MatchMCHMID() < 0.f) { // this should never happen. only for protection. return; } - o2::dataformats::GlobalFwdTrack propmuonAtPV = propagateMuon(muon, muon, collision, propagationPoint::kToVertex, muoncuts.matchingZ, mBzMFT); + o2::dataformats::GlobalFwdTrack propmuonAtPV = propagateMuon(muon, muon, collision, propagationPoint::kToVertex, muoncuts.matchingZ, mBzMFT, 0.0); float pt = propmuonAtPV.getPt(); float eta = propmuonAtPV.getEta(); @@ -592,17 +593,17 @@ struct CreateResolutionMap { if (mcparticle.globalIndex() != mcparticle_MCHMID.globalIndex()) { // this should not happen. this is only for protection. return; } - if (cfg_reject_fake_match_mft_mch && mcparticle.globalIndex() != mcparticle_MFT.globalIndex()) { // evaluate mismatch + if (cfg_reject_fake_match_mft_mch && mcparticle_MCHMID.globalIndex() != mcparticle_MFT.globalIndex()) { // evaluate mismatch return; } - o2::dataformats::GlobalFwdTrack propmuonAtPV_Matched = propagateMuon(mchtrack, mchtrack, collision, propagationPoint::kToVertex, muoncuts.matchingZ, mBzMFT); + o2::dataformats::GlobalFwdTrack propmuonAtPV_Matched = propagateMuon(mchtrack, mchtrack, collision, propagationPoint::kToVertex, muoncuts.matchingZ, mBzMFT, 0.0); ptMatchedMCHMID = propmuonAtPV_Matched.getPt(); etaMatchedMCHMID = propmuonAtPV_Matched.getEta(); phiMatchedMCHMID = propmuonAtPV_Matched.getPhi(); o2::math_utils::bringTo02Pi(phiMatchedMCHMID); - // o2::dataformats::GlobalFwdTrack propmuonAtDCA_Matched = propagateMuon(mchtrack, mchtrack, collision, propagationPoint::kToDCA, muoncuts.matchingZ, mBzMFT); + // o2::dataformats::GlobalFwdTrack propmuonAtDCA_Matched = propagateMuon(mchtrack, mchtrack, collision, propagationPoint::kToDCA, muoncuts.matchingZ, mBzMFT, 0.0); // float dcaX_Matched = propmuonAtDCA_Matched.getX() - collision.posX(); // float dcaY_Matched = propmuonAtDCA_Matched.getY() - collision.posY(); // float dcaXY_Matched = std::sqrt(dcaX_Matched * dcaX_Matched + dcaY_Matched * dcaY_Matched); @@ -626,7 +627,7 @@ struct CreateResolutionMap { if constexpr (withMFTCov) { auto mfttrackcov = mftCovs.rawIteratorAt(map_mfttrackcovs[mfttrack.globalIndex()]); - auto muonAtMP = propagateMuon(mchtrack, mchtrack, collision, propagationPoint::kToMatchingPlane, muoncuts.matchingZ, mBzMFT); // propagated to matching plane + auto muonAtMP = propagateMuon(mchtrack, mchtrack, collision, propagationPoint::kToMatchingPlane, muoncuts.matchingZ, mBzMFT, 0.0); // propagated to matching plane o2::track::TrackParCovFwd mftsaAtMP = getTrackParCovFwd(mfttrack, mfttrackcov); // values at innermost update mftsaAtMP.propagateToZhelix(muoncuts.matchingZ, mBzMFT); // propagated to matching plane etaMatchedMFTatMP = mftsaAtMP.getEta(); @@ -680,12 +681,12 @@ struct CreateResolutionMap { } } else if (muon.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::MuonStandaloneTrack) { - o2::dataformats::GlobalFwdTrack propmuonAtRabs = propagateMuon(muon, muon, collision, propagationPoint::kToRabs, muoncuts.matchingZ, mBzMFT); // this is necessary only for MuonStandaloneTrack + o2::dataformats::GlobalFwdTrack propmuonAtRabs = propagateMuon(muon, muon, collision, propagationPoint::kToRabs, muoncuts.matchingZ, mBzMFT, 0.0); // this is necessary only for MuonStandaloneTrack float xAbs = propmuonAtRabs.getX(); float yAbs = propmuonAtRabs.getY(); rAtAbsorberEnd = std::sqrt(xAbs * xAbs + yAbs * yAbs); // Redo propagation only for muon tracks // propagation of MFT tracks alredy done in reconstruction - o2::dataformats::GlobalFwdTrack propmuonAtDCA = propagateMuon(muon, muon, collision, propagationPoint::kToDCA, muoncuts.matchingZ, mBzMFT); + o2::dataformats::GlobalFwdTrack propmuonAtDCA = propagateMuon(muon, muon, collision, propagationPoint::kToDCA, muoncuts.matchingZ, mBzMFT, 0.0); dcaX = propmuonAtDCA.getX() - collision.posX(); dcaY = propmuonAtDCA.getY() - collision.posY(); dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY); @@ -874,31 +875,80 @@ struct CreateResolutionMap { } std::vector> vec_min_chi2MatchMCHMFT; // std::pair -> chi2MatchMCHMFT; - template - void findBestMatchPerMCHMID(TMuons const& muons) + template + void findBestMatchPerMCHMID(TCollision const& collision, TFwdTrack const& fwdtrack, TFwdTracks const& fwdtracks, TMFTTracks const&) { - vec_min_chi2MatchMCHMFT.reserve(muons.size()); - for (const auto& muon : muons) { - if (muon.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::MuonStandaloneTrack) { - const auto& muons_per_MCHMID = muons.sliceBy(fwdtracksPerMCHTrack, muon.globalIndex()); - // LOGF(info, "stanadalone: muon.globalIndex() = %d, muon.chi2MatchMCHMFT() = %f", muon.globalIndex(), muon.chi2MatchMCHMFT()); - // LOGF(info, "muons_per_MCHMID.size() = %d", muons_per_MCHMID.size()); - - float min_chi2MatchMCHMFT = 1e+10; - std::tuple tupleIds_at_min; - for (const auto& muon_tmp : muons_per_MCHMID) { - if (muon_tmp.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack) { - // LOGF(info, "muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, muon_tmp.chi2MatchMCHMFT() = %f", muon_tmp.globalIndex(), muon_tmp.matchMCHTrackId(), muon_tmp.matchMFTTrackId(), muon_tmp.chi2MatchMCHMFT()); - if (0.f < muon_tmp.chi2MatchMCHMFT() && muon_tmp.chi2MatchMCHMFT() < min_chi2MatchMCHMFT) { - min_chi2MatchMCHMFT = muon_tmp.chi2MatchMCHMFT(); - tupleIds_at_min = std::make_tuple(muon_tmp.globalIndex(), muon_tmp.matchMCHTrackId(), muon_tmp.matchMFTTrackId()); - } + if (fwdtrack.trackType() != o2::aod::fwdtrack::ForwardTrackTypeEnum::MuonStandaloneTrack) { + return; + } + if (!fwdtrack.has_mcParticle()) { + return; + } + + o2::dataformats::GlobalFwdTrack propmuonAtPV_Matched = propagateMuon(fwdtrack, fwdtrack, collision, propagationPoint::kToVertex, muoncuts.matchingZ, mBzMFT, 0.0); + float etaMatchedMCHMID = propmuonAtPV_Matched.getEta(); + float phiMatchedMCHMID = propmuonAtPV_Matched.getPhi(); + o2::math_utils::bringTo02Pi(phiMatchedMCHMID); + + o2::dataformats::GlobalFwdTrack propmuonAtDCA_Matched = propagateMuon(fwdtrack, fwdtrack, collision, propagationPoint::kToDCA, muoncuts.matchingZ, mBzMFT, 0.0); + float dcaX_Matched = propmuonAtDCA_Matched.getX() - collision.posX(); + float dcaY_Matched = propmuonAtDCA_Matched.getY() - collision.posY(); + float dcaXY_Matched = std::sqrt(dcaX_Matched * dcaX_Matched + dcaY_Matched * dcaY_Matched); + float pDCA = fwdtrack.p() * dcaXY_Matched; + + auto muons_per_MCHMID = fwdtracks.sliceBy(fwdtracksPerMCHTrack, fwdtrack.globalIndex()); + + float min_chi2MatchMCHMFT = 1e+10; + std::tuple tupleIds_at_min; + for (const auto& muon_tmp : muons_per_MCHMID) { + if (muon_tmp.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack) { + // LOGF(info, "muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, muon_tmp.chi2MatchMCHMFT() = %f", muon_tmp.globalIndex(), muon_tmp.matchMCHTrackId(), muon_tmp.matchMFTTrackId(), muon_tmp.chi2MatchMCHMFT()); + auto mchtrack = muon_tmp.template matchMCHTrack_as(); // MCH-MID + auto mfttrack = muon_tmp.template matchMFTTrack_as(); // MFTsa + + if (muon_tmp.chi2() < 0.f || muon_tmp.chi2MatchMCHMFT() < 0.f || muon_tmp.chi2MatchMCHMID() < 0.f || mfttrack.chi2() < 0.f) { // reject negative chi2, i.e. wrong. + continue; + } + + o2::dataformats::GlobalFwdTrack propmuonAtPV = propagateMuon(muon_tmp, muon_tmp, collision, propagationPoint::kToVertex, muoncuts.matchingZ, mBzMFT, 0.0); + float pt = propmuonAtPV.getPt(); + float eta = propmuonAtPV.getEta(); + float phi = propmuonAtPV.getPhi(); + o2::math_utils::bringTo02Pi(phi); + + if (muoncuts.refitGlobalMuon) { + pt = propmuonAtPV_Matched.getP() * std::sin(2.f * std::atan(std::exp(-eta))); + } + + float deta = etaMatchedMCHMID - eta; + float dphi = phiMatchedMCHMID - phi; + o2::math_utils::bringToPMPi(dphi); + int ndf = 2 * (mchtrack.nClusters() + mfttrack.nClusters()) - 5; + + float dcaX = propmuonAtPV.getX() - collision.posX(); + float dcaY = propmuonAtPV.getY() - collision.posY(); + float dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY); + + if (muoncuts.cfgApplyPreselectionInBestMatch) { + if (!isSelectedMuon(pt, eta, muon_tmp.rAtAbsorberEnd(), pDCA, muon_tmp.chi2() / ndf, muon_tmp.trackType(), dcaXY)) { + continue; + } + if (std::sqrt(std::pow(deta / muoncuts.cfg_max_deta, 2) + std::pow(dphi / muoncuts.cfg_max_dphi, 2)) > 1.f) { + continue; + } + if (muon_tmp.chi2MatchMCHMFT() > muoncuts.cfg_max_matching_chi2_mftmch) { + continue; } } - vec_min_chi2MatchMCHMFT.emplace_back(tupleIds_at_min); - // LOGF(info, "min: muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, muon_tmp.chi2MatchMCHMFT() = %f", std::get<0>(tupleIds_at_min), std::get<1>(tupleIds_at_min), std::get<2>(tupleIds_at_min), min_chi2MatchMCHMFT); + + if (0.f < muon_tmp.chi2MatchMCHMFT() && muon_tmp.chi2MatchMCHMFT() < min_chi2MatchMCHMFT) { + min_chi2MatchMCHMFT = muon_tmp.chi2MatchMCHMFT(); + tupleIds_at_min = std::make_tuple(muon_tmp.globalIndex(), muon_tmp.matchMCHTrackId(), muon_tmp.matchMFTTrackId()); + } } - } // end of muon loop + } + vec_min_chi2MatchMCHMFT.emplace_back(tupleIds_at_min); + // LOGF(info, "min: muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, muon_tmp.chi2MatchMCHMFT() = %f", std::get<0>(tupleIds_at_min), std::get<1>(tupleIds_at_min), std::get<2>(tupleIds_at_min), min_chi2MatchMCHMFT); } void processElectronSA(MyCollisions const& collisions, aod::BCsWithTimestamps const&, MyTracks const& tracks, aod::McCollisions const&, aod::McParticles const&) @@ -976,9 +1026,19 @@ struct CreateResolutionMap { // Partition sa_muons = o2::aod::fwdtrack::trackType == uint8_t(o2::aod::fwdtrack::ForwardTrackTypeEnum::MuonStandaloneTrack); // MCH-MID // Partition global_muons = o2::aod::fwdtrack::trackType == uint8_t(o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack); // MFT-MCH-MID - void processMuonSA(MyCollisions const& collisions, aod::BCsWithTimestamps const&, MyFwdTracks const& fwdtracks, MyMFTTracks const&, aod::McCollisions const&, aod::McParticles const&) + void processMuonSA(MyCollisions const& collisions, aod::BCsWithTimestamps const&, MyFwdTracks const& fwdtracks, MyMFTTracks const& mfttracks, aod::McCollisions const&, aod::McParticles const&) { - findBestMatchPerMCHMID(fwdtracks); + vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); + + for (const auto& collision : collisions) { + auto bc = collision.template bc_as(); + initCCDB(bc); + const auto& fwdtracks_per_coll = fwdtracks.sliceBy(perCollision_fwd, collision.globalIndex()); + for (const auto& fwdtrack : fwdtracks_per_coll) { + findBestMatchPerMCHMID(collision, fwdtrack, fwdtracks, mfttracks); + } // end of fwdtrack loop + } // end of collision loop + for (const auto& collision : collisions) { auto bc = collision.template foundBC_as(); initCCDB(bc); @@ -1016,51 +1076,18 @@ struct CreateResolutionMap { PROCESS_SWITCH(CreateResolutionMap, processMuonSA, "create resolution map for muon at forward rapidity", true); Preslice fwdtrackIndicesPerCollision = aod::track_association::collisionId; - void processMuonTTCA(MyCollisions const& collisions, aod::BCsWithTimestamps const&, MyFwdTracks const& fwdtracks, MyMFTTracks const&, aod::FwdTrackAssoc const& fwdtrackIndices, aod::McCollisions const&, aod::McParticles const&) + void processMuonTTCA(MyCollisions const& collisions, aod::BCsWithTimestamps const&, MyFwdTracks const& fwdtracks, MyMFTTracks const& mfttracks, aod::FwdTrackAssoc const& fwdtrackIndices, aod::McCollisions const&, aod::McParticles const&) { - findBestMatchPerMCHMID(fwdtracks); + vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); for (const auto& collision : collisions) { - auto bc = collision.template foundBC_as(); + auto bc = collision.template bc_as(); initCCDB(bc); - - if (!isSelectedEvent(collision)) { - continue; - } - - if (!collision.has_mcCollision()) { - continue; - } - - float centrality = std::array{collision.centFT0M(), collision.centFT0A(), collision.centFT0C()}[cfgCentEstimator]; - // auto mccollision = collision.template mcCollision_as(); - // registry.fill(HIST("Event/Muon/hImpPar_Centrality"), mccollision.impactParameter(), centrality); - auto fwdtrackIdsThisCollision = fwdtrackIndices.sliceBy(fwdtrackIndicesPerCollision, collision.globalIndex()); for (const auto& fwdtrackId : fwdtrackIdsThisCollision) { - auto muon = fwdtrackId.template fwdtrack_as(); - if (!muon.has_mcParticle()) { - continue; - } - auto mctrack = muon.template mcParticle_as(); - auto mccollision_from_mctrack = mctrack.template mcCollision_as(); - if (cfgEventGeneratorType >= 0 && mccollision_from_mctrack.getSubGeneratorId() != cfgEventGeneratorType) { - continue; - } - fillMuon(collision, muon, nullptr, centrality); + auto fwdtrack = fwdtrackId.template fwdtrack_as(); + findBestMatchPerMCHMID(collision, fwdtrack, fwdtracks, mfttracks); } // end of fwdtrack loop } // end of collision loop - vec_min_chi2MatchMCHMFT.clear(); - vec_min_chi2MatchMCHMFT.shrink_to_fit(); - } - PROCESS_SWITCH(CreateResolutionMap, processMuonTTCA, "create resolution map for muon at forward rapidity", false); - - std::unordered_map map_mfttrackcovs; - void processMuonTTCA_withMFTCov(MyCollisions const& collisions, aod::BCsWithTimestamps const&, MyFwdTracks const& fwdtracks, MyMFTTracks const&, aod::MFTTracksCov const& mftCovs, aod::FwdTrackAssoc const& fwdtrackIndices, aod::McCollisions const&, aod::McParticles const&) - { - for (const auto& mfttrackConv : mftCovs) { - map_mfttrackcovs[mfttrackConv.matchMFTTrackId()] = mfttrackConv.globalIndex(); - } - findBestMatchPerMCHMID(fwdtracks); for (const auto& collision : collisions) { auto bc = collision.template foundBC_as(); @@ -1089,14 +1116,58 @@ struct CreateResolutionMap { if (cfgEventGeneratorType >= 0 && mccollision_from_mctrack.getSubGeneratorId() != cfgEventGeneratorType) { continue; } - fillMuon(collision, muon, mftCovs, centrality); + fillMuon(collision, muon, nullptr, centrality); } // end of fwdtrack loop } // end of collision loop - map_mfttrackcovs.clear(); vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); } - PROCESS_SWITCH(CreateResolutionMap, processMuonTTCA_withMFTCov, "create resolution map for muon at forward rapidity", false); + PROCESS_SWITCH(CreateResolutionMap, processMuonTTCA, "create resolution map for muon at forward rapidity", false); + + std::unordered_map map_mfttrackcovs; + + // void processMuonTTCA_withMFTCov(MyCollisions const& collisions, aod::BCsWithTimestamps const&, MyFwdTracks const& fwdtracks, MyMFTTracks const&, aod::MFTTracksCov const& mftCovs, aod::FwdTrackAssoc const& fwdtrackIndices, aod::McCollisions const&, aod::McParticles const&) + // { + // for (const auto& mfttrackConv : mftCovs) { + // map_mfttrackcovs[mfttrackConv.matchMFTTrackId()] = mfttrackConv.globalIndex(); + // } + // findBestMatchPerMCHMID(fwdtracks); + + // for (const auto& collision : collisions) { + // auto bc = collision.template foundBC_as(); + // initCCDB(bc); + + // if (!isSelectedEvent(collision)) { + // continue; + // } + + // if (!collision.has_mcCollision()) { + // continue; + // } + + // float centrality = std::array{collision.centFT0M(), collision.centFT0A(), collision.centFT0C()}[cfgCentEstimator]; + // // auto mccollision = collision.template mcCollision_as(); + // // registry.fill(HIST("Event/Muon/hImpPar_Centrality"), mccollision.impactParameter(), centrality); + + // auto fwdtrackIdsThisCollision = fwdtrackIndices.sliceBy(fwdtrackIndicesPerCollision, collision.globalIndex()); + // for (const auto& fwdtrackId : fwdtrackIdsThisCollision) { + // auto muon = fwdtrackId.template fwdtrack_as(); + // if (!muon.has_mcParticle()) { + // continue; + // } + // auto mctrack = muon.template mcParticle_as(); + // auto mccollision_from_mctrack = mctrack.template mcCollision_as(); + // if (cfgEventGeneratorType >= 0 && mccollision_from_mctrack.getSubGeneratorId() != cfgEventGeneratorType) { + // continue; + // } + // fillMuon(collision, muon, mftCovs, centrality); + // } // end of fwdtrack loop + // } // end of collision loop + // map_mfttrackcovs.clear(); + // vec_min_chi2MatchMCHMFT.clear(); + // vec_min_chi2MatchMCHMFT.shrink_to_fit(); + // } + // PROCESS_SWITCH(CreateResolutionMap, processMuonTTCA_withMFTCov, "create resolution map for muon at forward rapidity", false); void processGen(aod::McCollisions const& mcCollisions) { diff --git a/PWGEM/Dilepton/Tasks/matchingMFT.cxx b/PWGEM/Dilepton/Tasks/matchingMFT.cxx index 558d0f418d0..e38648e48d0 100644 --- a/PWGEM/Dilepton/Tasks/matchingMFT.cxx +++ b/PWGEM/Dilepton/Tasks/matchingMFT.cxx @@ -88,7 +88,6 @@ struct matchingMFT { Configurable cfgManualZShift{"cfgManualZShift", 0, "manual z-shift for propagation of global muon to PV"}; Configurable requireTrueAssociation{"requireTrueAssociation", false, "flag to require true mc collision association"}; - Configurable maxRelDPt{"maxRelDPt", 1e+10f, "max. relative dpt between MFT-MCH-MID and MCH-MID"}; Configurable maxDEta{"maxDEta", 1e+10f, "max. deta between MFT-MCH-MID and MCH-MID"}; Configurable maxDPhi{"maxDPhi", 1e+10f, "max. dphi between MFT-MCH-MID and MCH-MID"}; Configurable maxDXMP{"maxDXMP", 1e+10f, "max. dx between MFT and MCH-MID at matching plane Z"}; @@ -96,12 +95,7 @@ struct matchingMFT { Configurable requireMFTHitMap{"requireMFTHitMap", false, "flag to require MFT hit map"}; Configurable> requiredMFTDisks{"requiredMFTDisks", std::vector{0}, "hit map on MFT disks [0,1,2,3,4]. logical-OR of each double-sided disk"}; Configurable matchingZ{"matchingZ", -77.5, "z position where matching is performed"}; - Configurable cfgBestMatchFinder{"cfgBestMatchFinder", 0, "matching chi2:0, dr:1, 2d:2"}; Configurable cfgApplyPreselectionInBestMatch{"cfgApplyPreselectionInBestMatch", false, "flag to apply preselection in find best match function"}; - Configurable cfgSlope_dr_chi2MatchMFTMCH{"cfgSlope_dr_chi2MatchMFTMCH", -0.15 / 30, "slope of chiMatchMCHMFT vs. dR"}; - Configurable cfgIntercept_dr_chi2MatchMFTMCH{"cfgIntercept_dr_chi2MatchMFTMCH", 1e+10f, "intercept of chiMatchMCHMFT vs. dR"}; - Configurable cfgPeakPosition_chi2MatchMFTMCH{"cfgPeakPosition_chi2MatchMFTMCH", 2.f, "peak position of chiMatchMCHMFT distribution"}; - Configurable cfgPeakPosition_dr{"cfgPeakPosition_dr", 0.01, "peak position of dr distribution"}; Configurable cfgCentEstimator{"cfgCentEstimator", 2, "FT0M:0, FT0A:1, FT0C:2"}; Configurable cfgCentMin{"cfgCentMin", -1.f, "min. centrality"}; @@ -541,7 +535,6 @@ struct matchingMFT { float deta = etaMatchedMCHMID - eta; float dphi = phiMatchedMCHMID - phi; o2::math_utils::bringToPMPi(dphi); - float dr = std::sqrt(deta * deta + dphi * dphi); if (std::sqrt(std::pow(deta / maxDEta, 2) + std::pow(dphi / maxDPhi, 2)) > 1.f) { return; @@ -549,13 +542,6 @@ struct matchingMFT { if (std::sqrt(std::pow(dxMP / maxDXMP, 2) + std::pow(dyMP / maxDYMP, 2)) > 1.f) { return; } - if (std::fabs(dpt) > maxRelDPt) { - return; - } - - if (cfgSlope_dr_chi2MatchMFTMCH * fwdtrack.chi2MatchMCHMFT() + cfgIntercept_dr_chi2MatchMFTMCH < dr) { - return; - } if (!isSelected(pt, eta, rAtAbsorberEnd, pDCA, fwdtrack.chi2() / (2.f * (mchtrack.nClusters() + nClustersMFT) - 5.f), fwdtrack.trackType(), dcaXY)) { return; @@ -753,9 +739,7 @@ struct matchingMFT { } std::vector> vec_min_chi2MatchMCHMFT; // std::pair -> chi2MatchMCHMFT; - std::vector> vec_min_dr; // std::pair -> deta + dphi; - std::vector> vec_min_2d; // std::pair -> deta + dphi; - std::map, bool> mapCorrectMatch; + // std::map, bool> mapCorrectMatch; template void findBestMatchPerMCHMID(TCollision const& collision, TFwdTrack const& fwdtrack, TFwdTracks const& fwdtracks, TMFTTracks const&) @@ -768,11 +752,7 @@ struct matchingMFT { } std::tuple tupleIds_at_min_chi2mftmch; - std::tuple tupleIds_at_min_dr; - std::tuple tupleIds_at_min_distance_2d; float min_chi2MatchMCHMFT = 1e+10; - float min_dr = 1e+10; - float min_distance_2d = 1e+10; auto muons_per_MCHMID = fwdtracks.sliceBy(fwdtracksPerMCHTrack, fwdtrack.globalIndex()); // LOGF(info, "muons_per_MCHMID.size() = %d", muons_per_MCHMID.size()); @@ -781,12 +761,22 @@ struct matchingMFT { float phiMatchedMCHMID = propmuonAtPV_Matched.getPhi(); o2::math_utils::bringTo02Pi(phiMatchedMCHMID); + o2::dataformats::GlobalFwdTrack propmuonAtDCA_Matched = propagateMuon(fwdtrack, fwdtrack, collision, propagationPoint::kToDCA, matchingZ, mBz, mZShift); + float dcaX_Matched = propmuonAtDCA_Matched.getX() - collision.posX(); + float dcaY_Matched = propmuonAtDCA_Matched.getY() - collision.posY(); + float dcaXY_Matched = std::sqrt(dcaX_Matched * dcaX_Matched + dcaY_Matched * dcaY_Matched); + float pDCA = fwdtrack.p() * dcaXY_Matched; + for (const auto& muon_tmp : muons_per_MCHMID) { if (muon_tmp.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack) { auto tupleId = std::make_tuple(muon_tmp.globalIndex(), muon_tmp.matchMCHTrackId(), muon_tmp.matchMFTTrackId()); auto mchtrack = muon_tmp.template matchMCHTrack_as(); // MCH-MID auto mfttrack = muon_tmp.template matchMFTTrack_as(); + if (muon_tmp.chi2() < 0.f || muon_tmp.chi2MatchMCHMFT() < 0.f || muon_tmp.chi2MatchMCHMID() < 0.f || mfttrack.chi2() < 0.f) { // reject negative chi2, i.e. wrong. + continue; + } + // LOGF(info, "muon_tmp.has_mcParticle() = %d, mchtrack.has_mcParticle() = %d, mfttrack.has_mcParticle() = %d", muon_tmp.has_mcParticle(), mchtrack.has_mcParticle(), mfttrack.has_mcParticle()); if (!muon_tmp.has_mcParticle() || !mchtrack.has_mcParticle() || !mfttrack.has_mcParticle()) { @@ -801,23 +791,23 @@ struct matchingMFT { bool isMatched = (mcParticle_MFT.globalIndex() == mcParticle_MCHMID.globalIndex()) && (mcParticle_MFT.mcCollisionId() == mcParticle_MCHMID.mcCollisionId()); o2::dataformats::GlobalFwdTrack propmuonAtPV = propagateMuon(muon_tmp, muon_tmp, collision, propagationPoint::kToVertex, matchingZ, mBz, mZShift); - // float pt = propmuonAtPV.getPt(); + float pt = propmuonAtPV.getPt(); float eta = propmuonAtPV.getEta(); float phi = propmuonAtPV.getPhi(); o2::math_utils::bringTo02Pi(phi); - // if (refitGlobalMuon) { - // pt = propmuonAtPV_Matched.getP() * std::sin(2.f * std::atan(std::exp(-eta))); - // } + if (refitGlobalMuon) { + pt = propmuonAtPV_Matched.getP() * std::sin(2.f * std::atan(std::exp(-eta))); + } float deta = etaMatchedMCHMID - eta; float dphi = phiMatchedMCHMID - phi; o2::math_utils::bringToPMPi(dphi); float dr = std::sqrt(deta * deta + dphi * dphi); - // float chi2ndf = muon_tmp.chi2() / (2.f * (mchtrack.nClusters() + mfttrack.nClusters()) - 5.f); + float chi2ndf = muon_tmp.chi2() / (2.f * (mchtrack.nClusters() + mfttrack.nClusters()) - 5.f); - if (cfgApplyPreselectionInBestMatch && cfgSlope_dr_chi2MatchMFTMCH * muon_tmp.chi2MatchMCHMFT() + cfgIntercept_dr_chi2MatchMFTMCH < dr) { - continue; - } + float dcaX = propmuonAtPV.getX() - collision.posX(); + float dcaY = propmuonAtPV.getY() - collision.posY(); + float dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY); if (isPrimary) { if (isMatched) { @@ -833,11 +823,11 @@ struct matchingMFT { } } - if (mcParticle_MFT.globalIndex() == mcParticle_MCHMID.globalIndex()) { - mapCorrectMatch[tupleId] = true; - } else { - mapCorrectMatch[tupleId] = false; - } + // if (mcParticle_MFT.globalIndex() == mcParticle_MCHMID.globalIndex()) { + // mapCorrectMatch[tupleId] = true; + // } else { + // mapCorrectMatch[tupleId] = false; + // } // if (std::abs(mcParticle_MCHMID.pdgCode()) == 13 && mcParticle_MCHMID.isPhysicalPrimary()) { // if (mcParticle_MFT.globalIndex() == mcParticle_MCHMID.globalIndex()) { @@ -847,34 +837,31 @@ struct matchingMFT { // } // } - if (0.f < muon_tmp.chi2MatchMCHMFT() && std::sqrt(std::pow(muon_tmp.chi2MatchMCHMFT() - cfgPeakPosition_chi2MatchMFTMCH, 2)) < min_chi2MatchMCHMFT) { - min_chi2MatchMCHMFT = std::sqrt(std::pow(muon_tmp.chi2MatchMCHMFT() - cfgPeakPosition_chi2MatchMFTMCH, 2)); - tupleIds_at_min_chi2mftmch = tupleId; - } - - if (std::sqrt(std::pow(dr - cfgPeakPosition_dr, 2)) < min_dr) { - min_dr = std::sqrt(std::pow(dr - cfgPeakPosition_dr, 2)); - tupleIds_at_min_dr = tupleId; + if (cfgApplyPreselectionInBestMatch) { + if (!isSelected(pt, eta, muon_tmp.rAtAbsorberEnd(), pDCA, chi2ndf, fwdtrack.trackType(), dcaXY)) { + continue; + } + if (std::sqrt(std::pow(deta / maxDEta, 2) + std::pow(dphi / maxDPhi, 2)) > 1.f) { + continue; + } + if (muon_tmp.chi2MatchMCHMFT() > maxMatchingChi2MCHMFT) { + continue; + } } - float distance_2d = std::sqrt(std::pow(muon_tmp.chi2MatchMCHMFT() - cfgPeakPosition_chi2MatchMFTMCH, 2) + std::pow(dr - cfgPeakPosition_dr, 2)); - if (distance_2d < min_distance_2d) { - min_distance_2d = distance_2d; - tupleIds_at_min_distance_2d = tupleId; + if (0.f < muon_tmp.chi2MatchMCHMFT() && muon_tmp.chi2MatchMCHMFT() < min_chi2MatchMCHMFT) { + min_chi2MatchMCHMFT = muon_tmp.chi2MatchMCHMFT(); + tupleIds_at_min_chi2mftmch = tupleId; } } // end of if global muon } // end of candidates loop vec_min_chi2MatchMCHMFT.emplace_back(tupleIds_at_min_chi2mftmch); - vec_min_dr.emplace_back(tupleIds_at_min_dr); - vec_min_2d.emplace_back(tupleIds_at_min_distance_2d); // auto mcParticleTMP = fwdtrack.template mcParticle_as(); // this is identical to mcParticle_MFTMCHMID // if (std::abs(mcParticleTMP.pdgCode()) == 13 && mcParticleTMP.isPhysicalPrimary()) { // LOGF(info, "min chi2: muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, muon_tmp.chi2MatchMCHMFT() = %f, correct match = %d", std::get<0>(tupleIds_at_min_chi2mftmch), std::get<1>(tupleIds_at_min_chi2mftmch), std::get<2>(tupleIds_at_min_chi2mftmch), min_chi2MatchMCHMFT, mapCorrectMatch[tupleIds_at_min_chi2mftmch]); - // LOGF(info, "min dr: muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, correct match = %d", std::get<0>(tupleIds_at_min_dr), std::get<1>(tupleIds_at_min_dr), std::get<2>(tupleIds_at_min_dr), mapCorrectMatch[tupleIds_at_min_dr]); - // LOGF(info, "min 2d: muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, correct match = %d", std::get<0>(tupleIds_at_min_distance_2d), std::get<1>(tupleIds_at_min_distance_2d), std::get<2>(tupleIds_at_min_distance_2d), mapCorrectMatch[tupleIds_at_min_distance_2d]); // } } @@ -1022,8 +1009,6 @@ struct matchingMFT { void processWithoutFTTCA(FilteredMyCollisions const& collisions, MyFwdTracks const& fwdtracks, MyMFTTracks const& mfttracks, aod::BCsWithTimestamps const&, aod::McParticles const&) { vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); - vec_min_dr.reserve(fwdtracks.size()); - vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); initCCDB(bc); @@ -1060,18 +1045,8 @@ struct matchingMFT { continue; } - if (cfgBestMatchFinder == 0) { // chi2 - if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_chi2MatchMCHMFT.begin(), vec_min_chi2MatchMCHMFT.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_chi2MatchMCHMFT.end()) { - continue; - } - } else if (cfgBestMatchFinder == 1) { // dr - if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_dr.begin(), vec_min_dr.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_dr.end()) { - continue; - } - } else if (cfgBestMatchFinder == 2) { // 2d - if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_2d.begin(), vec_min_2d.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_2d.end()) { - continue; - } + if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_chi2MatchMCHMFT.begin(), vec_min_chi2MatchMCHMFT.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_chi2MatchMCHMFT.end()) { + continue; } fillHistograms(collision, fwdtrack, fwdtracks, mfttracks, nullptr); @@ -1080,18 +1055,12 @@ struct matchingMFT { vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); - vec_min_dr.clear(); - vec_min_dr.shrink_to_fit(); - vec_min_2d.clear(); - vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(matchingMFT, processWithoutFTTCA, "process without FTTCA", false); void processWithFTTCA(FilteredMyCollisions const& collisions, MyFwdTracks const& fwdtracks, MyMFTTracks const& mfttracks, aod::BCsWithTimestamps const&, aod::FwdTrackAssoc const& fwdtrackIndices, aod::McParticles const&) { vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); - vec_min_dr.reserve(fwdtracks.size()); - vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); initCCDB(bc); @@ -1129,18 +1098,8 @@ struct matchingMFT { continue; } - if (cfgBestMatchFinder == 0) { // chi2 - if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_chi2MatchMCHMFT.begin(), vec_min_chi2MatchMCHMFT.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_chi2MatchMCHMFT.end()) { - continue; - } - } else if (cfgBestMatchFinder == 1) { // dr - if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_dr.begin(), vec_min_dr.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_dr.end()) { - continue; - } - } else if (cfgBestMatchFinder == 2) { // 2d - if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_2d.begin(), vec_min_2d.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_2d.end()) { - continue; - } + if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_chi2MatchMCHMFT.begin(), vec_min_chi2MatchMCHMFT.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_chi2MatchMCHMFT.end()) { + continue; } fillHistograms(collision, fwdtrack, fwdtracks, mfttracks, nullptr); @@ -1149,10 +1108,6 @@ struct matchingMFT { vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); - vec_min_dr.clear(); - vec_min_dr.shrink_to_fit(); - vec_min_2d.clear(); - vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(matchingMFT, processWithFTTCA, "process with FTTCA", true); @@ -1165,8 +1120,6 @@ struct matchingMFT { } vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); - vec_min_dr.reserve(fwdtracks.size()); - vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); @@ -1205,18 +1158,8 @@ struct matchingMFT { continue; } - if (cfgBestMatchFinder == 0) { // chi2 - if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_chi2MatchMCHMFT.begin(), vec_min_chi2MatchMCHMFT.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_chi2MatchMCHMFT.end()) { - continue; - } - } else if (cfgBestMatchFinder == 1) { // dr - if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_dr.begin(), vec_min_dr.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_dr.end()) { - continue; - } - } else if (cfgBestMatchFinder == 2) { // 2d - if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_2d.begin(), vec_min_2d.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_2d.end()) { - continue; - } + if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_chi2MatchMCHMFT.begin(), vec_min_chi2MatchMCHMFT.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_chi2MatchMCHMFT.end()) { + continue; } fillHistograms(collision, fwdtrack, fwdtracks, mfttracks, mftCovs); } // end of fwdtrack loop @@ -1225,10 +1168,6 @@ struct matchingMFT { map_mfttrackcovs.clear(); vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); - vec_min_dr.clear(); - vec_min_dr.shrink_to_fit(); - vec_min_2d.clear(); - vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(matchingMFT, processWithFTTCA_withMFTCov, "process with FTTCA with MFTCov", false); }; diff --git a/PWGEM/Dilepton/Tasks/prefilterDielectron.cxx b/PWGEM/Dilepton/Tasks/prefilterDielectron.cxx index 2b8d6a12321..c384f47a749 100644 --- a/PWGEM/Dilepton/Tasks/prefilterDielectron.cxx +++ b/PWGEM/Dilepton/Tasks/prefilterDielectron.cxx @@ -50,13 +50,13 @@ using namespace o2::soa; using namespace o2::aod::pwgem::dilepton::utils::emtrackutil; using namespace o2::aod::pwgem::dilepton::utils::pairutil; -using MyCollisions = soa::Join; -using MyCollision = MyCollisions::iterator; +struct prefilterDielectron { + using MyCollisions = soa::Join; + using MyCollision = MyCollisions::iterator; -using MyTracks = soa::Join; -using MyTrack = MyTracks::iterator; + using MyTracks = soa::Join; + using MyTrack = MyTracks::iterator; -struct prefilterDielectron { Produces pfb_derived; // Configurables @@ -249,7 +249,7 @@ struct prefilterDielectron { fEMEventCut = EMEventCut("fEMEventCut", "fEMEventCut"); fEMEventCut.SetRequireSel8(eventcuts.cfgRequireSel8); fEMEventCut.SetRequireFT0AND(eventcuts.cfgRequireFT0AND); - fEMEventCut.SetZvtxRange(eventcuts.cfgZvtxMin, +eventcuts.cfgZvtxMax); + fEMEventCut.SetZvtxRange(eventcuts.cfgZvtxMin, eventcuts.cfgZvtxMax); fEMEventCut.SetRequireNoTFB(eventcuts.cfgRequireNoTFB); fEMEventCut.SetRequireNoITSROFB(eventcuts.cfgRequireNoITSROFB); fEMEventCut.SetRequireNoSameBunchPileup(eventcuts.cfgRequireNoSameBunchPileup); @@ -356,11 +356,11 @@ struct prefilterDielectron { int ndf = 0; void processPFB(FilteredMyCollisions const& collisions, MyTracks const& tracks) { - for (auto& track : tracks) { + for (const auto& track : tracks) { map_pfb[track.globalIndex()] = 0; } // end of track loop - for (auto& collision : collisions) { + for (const auto& collision : collisions) { initCCDB(collision); const float centralities[3] = {collision.centFT0M(), collision.centFT0A(), collision.centFT0C()}; bool is_cent_ok = true; @@ -372,10 +372,10 @@ struct prefilterDielectron { auto negTracks_per_coll = negTracks->sliceByCached(o2::aod::emprimaryelectron::emeventId, collision.globalIndex(), cache); if (!fEMEventCut.IsSelected(collision) || !is_cent_ok) { - for (auto& pos : posTracks_per_coll) { + for (const auto& pos : posTracks_per_coll) { map_pfb[pos.globalIndex()] = 0; } - for (auto& neg : negTracks_per_coll) { + for (const auto& neg : negTracks_per_coll) { map_pfb[neg.globalIndex()] = 0; } continue; @@ -383,7 +383,7 @@ struct prefilterDielectron { // LOGF(info, "centrality = %f , posTracks_per_coll.size() = %d, negTracks_per_coll.size() = %d", centralities[cfgCentEstimator], posTracks_per_coll.size(), negTracks_per_coll.size()); - for (auto& [pos, ele] : combinations(CombinationsFullIndexPolicy(posTracks_per_coll, negTracks_per_coll))) { // ULS + for (const auto& [pos, ele] : combinations(CombinationsFullIndexPolicy(posTracks_per_coll, negTracks_per_coll))) { // ULS if (!fDielectronCut.IsSelectedTrack(pos) || !fDielectronCut.IsSelectedTrack(ele)) { continue; } @@ -417,7 +417,7 @@ struct prefilterDielectron { } } // end of ULS pairing - for (auto& [pos1, pos2] : combinations(CombinationsStrictlyUpperIndexPolicy(posTracks_per_coll, posTracks_per_coll))) { // LS++ + for (const auto& [pos1, pos2] : combinations(CombinationsStrictlyUpperIndexPolicy(posTracks_per_coll, posTracks_per_coll))) { // LS++ if (!fDielectronCut.IsSelectedTrack(pos1) || !fDielectronCut.IsSelectedTrack(pos2)) { continue; } @@ -441,7 +441,7 @@ struct prefilterDielectron { } } // end of LS++ pairing - for (auto& [ele1, ele2] : combinations(CombinationsStrictlyUpperIndexPolicy(negTracks_per_coll, negTracks_per_coll))) { // LS-- + for (const auto& [ele1, ele2] : combinations(CombinationsStrictlyUpperIndexPolicy(negTracks_per_coll, negTracks_per_coll))) { // LS-- if (!fDielectronCut.IsSelectedTrack(ele1) || !fDielectronCut.IsSelectedTrack(ele2)) { continue; } @@ -467,13 +467,13 @@ struct prefilterDielectron { } // end of collision loop - for (auto& track : tracks) { + for (const auto& track : tracks) { // LOGF(info, "map_pfb[%d] = %d", track.globalIndex(), map_pfb[track.globalIndex()]); pfb_derived(map_pfb[track.globalIndex()]); } // end of track loop // check pfb. - for (auto& collision : collisions) { + for (const auto& collision : collisions) { const float centralities[3] = {collision.centFT0M(), collision.centFT0A(), collision.centFT0C()}; if (centralities[cfgCentEstimator] < cfgCentMin || cfgCentMax < centralities[cfgCentEstimator]) { continue; @@ -486,7 +486,7 @@ struct prefilterDielectron { continue; } - for (auto& [pos, ele] : combinations(CombinationsFullIndexPolicy(posTracks_per_coll, negTracks_per_coll))) { // ULS + for (const auto& [pos, ele] : combinations(CombinationsFullIndexPolicy(posTracks_per_coll, negTracks_per_coll))) { // ULS if (!fDielectronCut.IsSelectedTrack(pos) || !fDielectronCut.IsSelectedTrack(ele)) { continue; } @@ -507,7 +507,7 @@ struct prefilterDielectron { fRegistry.fill(HIST("Pair/after/uls/hDeltaEtaDeltaPhi"), dphi, deta); } - for (auto& [pos1, pos2] : combinations(CombinationsStrictlyUpperIndexPolicy(posTracks_per_coll, posTracks_per_coll))) { // LS++ + for (const auto& [pos1, pos2] : combinations(CombinationsStrictlyUpperIndexPolicy(posTracks_per_coll, posTracks_per_coll))) { // LS++ if (!fDielectronCut.IsSelectedTrack(pos1) || !fDielectronCut.IsSelectedTrack(pos2)) { continue; } @@ -528,7 +528,7 @@ struct prefilterDielectron { fRegistry.fill(HIST("Pair/after/lspp/hDeltaEtaDeltaPhi"), dphi, deta); } - for (auto& [ele1, ele2] : combinations(CombinationsStrictlyUpperIndexPolicy(negTracks_per_coll, negTracks_per_coll))) { // LS-- + for (const auto& [ele1, ele2] : combinations(CombinationsStrictlyUpperIndexPolicy(negTracks_per_coll, negTracks_per_coll))) { // LS-- if (!fDielectronCut.IsSelectedTrack(ele1) || !fDielectronCut.IsSelectedTrack(ele2)) { continue; }