Skip to content

Commit 8f608c3

Browse files
dsekihatdsekihat
authored
[PWGEM/Dilepton] update on prefilter for photon conversion rejection (#8585)
1 parent a801731 commit 8f608c3

File tree

6 files changed

+618
-39
lines changed

6 files changed

+618
-39
lines changed

‎PWGEM/Dilepton/Core/Dilepton.h‎

Lines changed: 62 additions & 24 deletions
Original file line numberDiff line numberDiff line change
@@ -174,9 +174,9 @@ struct Dilepton {
174174
Configurable<bool> cfg_apply_detadphi{"cfg_apply_detadphi", false, "flag to apply deta-dphi elliptic cut"};
175175
Configurable<float> cfg_min_deta{"cfg_min_deta", 0.02, "min deta between 2 electrons (elliptic cut)"};
176176
Configurable<float> cfg_min_dphi{"cfg_min_dphi", 0.2, "min dphi between 2 electrons (elliptic cut)"};
177-
Configurable<bool> cfg_apply_detadphi_geom{"cfg_apply_detadphi_geom", false, "flag to apply generator deta-dphi elliptic cut"};
178-
Configurable<float> cfg_min_deta_geom{"cfg_min_deta_geom", 0.02, "geometrical min deta between 2 electrons (elliptic cut)"};
179-
Configurable<float> cfg_min_dphi_geom{"cfg_min_dphi_geom", 0.2, "geometrical min dphi between 2 electrons (elliptic cut)"};
177+
Configurable<bool> cfg_apply_dzrdphi_geom{"cfg_apply_dzrdphi_geom", false, "flag to apply generator dz-rdphi elliptic cut"};
178+
Configurable<float> cfg_min_dz_geom{"cfg_min_dz_geom", 1.2, "geometrical min dz between 2 electrons (elliptic cut) in cm"};
179+
Configurable<float> cfg_min_rdphi_geom{"cfg_min_rdphi_geom", 2.4, "geometrical min rdphi between 2 electrons (elliptic cut) in cm"};
180180
Configurable<float> cfg_min_opang{"cfg_min_opang", 0.0, "min opening angle"};
181181
Configurable<float> cfg_max_opang{"cfg_max_opang", 6.4, "max opening angle"};
182182
Configurable<bool> cfg_require_diff_sides{"cfg_require_diff_sides", false, "flag to require 2 tracks are from different sides."};
@@ -405,7 +405,7 @@ struct Dilepton {
405405
}
406406

407407
fRegistry.addClone("Event/before/hCollisionCounter", "Event/norm/hCollisionCounter");
408-
fRegistry.add("Pair/mix/hDiffBC", "diff. global BC in mixed event;|BC_{current} - BC_{mixed}|", kTH1D, {{1001, -0.5, 1000.5}}, true);
408+
fRegistry.add("Pair/mix/hDiffBC", "diff. global BC in mixed event;|BC_{current} - BC_{mixed}|", kTH1D, {{10001, -0.5, 10000.5}}, true);
409409
if (doprocessTriggerAnalysis) {
410410
fRegistry.add("Event/hNInspectedTVX", "N inspected TVX;run number;N_{TVX}", kTProfile, {{80000, 520000.5, 600000.5}}, true);
411411
}
@@ -525,9 +525,10 @@ struct Dilepton {
525525
if (cfgAnalysisType == static_cast<int>(o2::aod::pwgem::dilepton::utils::pairutil::DileptonAnalysisType::kQC)) {
526526
fRegistry.add("Pair/same/uls/hs", "dilepton", kTHnSparseD, {axis_mass, axis_pt, axis_dca}, true);
527527
fRegistry.add("Pair/same/uls/hsDeltaP", "difference of p between 2 tracks;|p_{T,1} - p_{T,2}|/|p_{T,1} + p_{T,2}|;#Delta#eta;#Delta#varphi (rad.);", kTHnSparseD, {{20, 0, 1}, {100, -0.5, +0.5}, {180, -M_PI, M_PI}}, true);
528-
fRegistry.add("Pair/same/uls/hGeomDeltaEtaDeltaPhi", Form("difference in #eta-#varphi plane between 2 tracks at X = %2.1f cm;#Delta#varphi (rad.);#Delta#eta;", dielectroncuts.cfg_x_to_go.value), kTH2D, {{180, -M_PI, M_PI}, {100, -0.5, +0.5}}, true);
528+
fRegistry.add("Pair/same/uls/hGeomDeltaZRDeltaPhi", Form("difference in z-r#varphi plane between 2 tracks at r = %2.1f cm;r#Delta#varphi (cm);#Deltaz (cm);", dielectroncuts.cfg_x_to_go.value), kTH2D, {{200, -50, 50}, {40, -10, 10}}, true);
529529
if constexpr (pairtype == o2::aod::pwgem::dilepton::utils::pairutil::DileptonPairType::kDielectron) {
530-
fRegistry.add("Pair/same/uls/hMvsPhiV", "m_{ee} vs. #varphi_{V};#varphi_{V} (rad.);m_{ee} (GeV/c^{2})", kTH2D, {{90, 0, M_PI}, {100, 0.0f, 1.0f}}, true); // phiv is only for dielectron
530+
fRegistry.add("Pair/same/uls/hMvsPhiV", "m_{ee} vs. #varphi_{V};#varphi_{V} (rad.);m_{ee} (GeV/c^{2})", kTH2D, {{90, 0, M_PI}, {100, 0.0f, 1.0f}}, true); // phiv is only for dielectron
531+
fRegistry.add("Pair/same/uls/hMvsPhiV_prop", "m_{ee} vs. #varphi_{V};#varphi_{V} (rad.);m_{ee} (GeV/c^{2})", kTH2D, {{90, 0, M_PI}, {100, 0.0f, 1.0f}}, true); // phiv is only for dielectron
531532
fRegistry.add("Pair/same/uls/hMvsOpAng", "m_{ee} vs. angle between 2 tracks;#omega (rad.);m_{ee} (GeV/c^{2})", kTH2D, {{100, 0, 2.0}, {20, 0.0f, 3.2}}, true);
532533
}
533534
fRegistry.addClone("Pair/same/uls/", "Pair/same/lspp/");
@@ -771,8 +772,11 @@ struct Dilepton {
771772
}
772773
}
773774

774-
std::map<std::tuple<int, int, int>, float> map_eta_geom; // map <dfId, collisionId, trackId> -> geometrical eta position at propagated point
775-
std::map<std::tuple<int, int, int>, float> map_phi_geom; // map <dfId, collisionId, trackId> -> geometrical phi position at propagated point
775+
std::map<std::tuple<int, int, int>, float> map_z_prop; // map <dfId, collisionId, trackId> -> geometrical z position at propagated point
776+
std::map<std::tuple<int, int, int>, float> map_phi_prop; // map <dfId, collisionId, trackId> -> geometrical phi position at propagated point
777+
std::map<std::tuple<int, int, int>, float> map_px_prop; // map <dfId, collisionId, trackId> -> px at propagated point
778+
std::map<std::tuple<int, int, int>, float> map_py_prop; // map <dfId, collisionId, trackId> -> py at propagated point
779+
std::map<std::tuple<int, int, int>, float> map_pz_prop; // map <dfId, collisionId, trackId> -> pz at propagated point
776780

777781
template <typename TTracks>
778782
void propagateElectron(TTracks const& tracks)
@@ -783,11 +787,22 @@ struct Dilepton {
783787
track_par_cov.setPID(o2::track::PID::Electron);
784788
o2::base::Propagator::Instance()->propagateToX(track_par_cov, dielectroncuts.cfg_x_to_go, d_bz, o2::base::PropagatorImpl<float>::MAX_SIN_PHI, o2::base::PropagatorImpl<float>::MAX_STEP, matCorr);
785789
auto xyz = track_par_cov.getXYZGlo();
786-
float eta = RecoDecay::eta(std::array{xyz.X(), xyz.Y(), xyz.Z()});
790+
// float eta = RecoDecay::eta(std::array{xyz.X(), xyz.Y(), xyz.Z()});
787791
float phi = RecoDecay::phi(std::array{xyz.X(), xyz.Y()});
788792
o2::math_utils::bringTo02Pi(phi);
789-
map_eta_geom[std::make_tuple(ndf, track.emeventId(), track.globalIndex())] = eta;
790-
map_phi_geom[std::make_tuple(ndf, track.emeventId(), track.globalIndex())] = phi;
793+
map_z_prop[std::make_tuple(ndf, track.emeventId(), track.globalIndex())] = xyz.Z();
794+
map_phi_prop[std::make_tuple(ndf, track.emeventId(), track.globalIndex())] = phi;
795+
796+
std::array<float, 3> pxpypz_prop = {0, 0, 0}; // px, py, pz
797+
getPxPyPz(track_par_cov, pxpypz_prop);
798+
map_px_prop[std::make_tuple(ndf, track.emeventId(), track.globalIndex())] = pxpypz_prop[0];
799+
map_py_prop[std::make_tuple(ndf, track.emeventId(), track.globalIndex())] = pxpypz_prop[1];
800+
map_pz_prop[std::make_tuple(ndf, track.emeventId(), track.globalIndex())] = pxpypz_prop[2];
801+
802+
// float r = sqrt(pow(xyz.X(),2) + pow(xyz.Y(), 2));
803+
// float theta = 2.f * std::atan(std::exp(-eta));
804+
// float z = r/std::tan(theta);
805+
// LOGF(info, "r = %f, z = %f , xyz.Z() = %f", r, z, xyz.Z());
791806
}
792807
}
793808

@@ -836,21 +851,39 @@ struct Dilepton {
836851
}
837852
}
838853

839-
float deta_geom = 999.f, dphi_geom = 999.f;
854+
float dphi_geom = 999.f, dz_geom = 999.f, rdphi_geom = 999.f;
855+
float phiv_prop = 999.f, mee_prop = 999.f;
840856
if constexpr (pairtype == o2::aod::pwgem::dilepton::utils::pairutil::DileptonPairType::kDielectron) {
841857
if (!cut.template IsSelectedPair(t1, t2, d_bz)) {
842858
return false;
843859
}
844860

845861
if constexpr (ev_id == 0) {
846-
deta_geom = t1.sign() * t1.pt() > t2.sign() * t2.pt() ? map_eta_geom[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())] - map_eta_geom[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] : map_eta_geom[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] - map_eta_geom[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())];
847-
dphi_geom = t1.sign() * t1.pt() > t2.sign() * t2.pt() ? map_phi_geom[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())] - map_phi_geom[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] : map_phi_geom[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] - map_phi_geom[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())];
862+
dz_geom = t1.sign() * t1.pt() > t2.sign() * t2.pt() ? map_z_prop[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())] - map_z_prop[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] : map_z_prop[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] - map_z_prop[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())];
863+
dphi_geom = t1.sign() * t1.pt() > t2.sign() * t2.pt() ? map_phi_prop[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())] - map_phi_prop[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] : map_phi_prop[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] - map_phi_prop[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())];
864+
o2::math_utils::bringToPMPi(dphi_geom);
865+
rdphi_geom = dielectroncuts.cfg_x_to_go * dphi_geom;
866+
867+
float px1 = map_px_prop[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())];
868+
float py1 = map_py_prop[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())];
869+
float pz1 = map_pz_prop[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())];
870+
float px2 = map_px_prop[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())];
871+
float py2 = map_py_prop[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())];
872+
float pz2 = map_pz_prop[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())];
873+
874+
ROOT::Math::PxPyPzMVector v1prop(px1, py1, pz1, leptonM1);
875+
ROOT::Math::PxPyPzMVector v2prop(px2, py2, pz2, leptonM2);
876+
ROOT::Math::PxPyPzMVector v12prop = v1prop + v2prop;
877+
mee_prop = v12prop.M();
878+
phiv_prop = o2::aod::pwgem::dilepton::utils::pairutil::getPhivPair(px1, py1, pz1, px2, py2, pz2, t1.sign(), t2.sign(), d_bz);
879+
848880
} else { // mixed event
849-
deta_geom = t1.sign() * t1.pt() > t2.sign() * t2.pt() ? map_eta_geom[std::make_tuple(t1.dfId(), t1.collisionId(), t1.globalIndex())] - map_eta_geom[std::make_tuple(t2.dfId(), t2.collisionId(), t2.globalIndex())] : map_eta_geom[std::make_tuple(t2.dfId(), t2.collisionId(), t2.globalIndex())] - map_eta_geom[std::make_tuple(t1.dfId(), t1.collisionId(), t1.globalIndex())];
850-
dphi_geom = t1.sign() * t1.pt() > t2.sign() * t2.pt() ? map_phi_geom[std::make_tuple(t1.dfId(), t1.collisionId(), t1.globalIndex())] - map_phi_geom[std::make_tuple(t2.dfId(), t2.collisionId(), t2.globalIndex())] : map_phi_geom[std::make_tuple(t2.dfId(), t2.collisionId(), t2.globalIndex())] - map_phi_geom[std::make_tuple(t1.dfId(), t1.collisionId(), t1.globalIndex())];
881+
dz_geom = t1.sign() * t1.pt() > t2.sign() * t2.pt() ? map_z_prop[std::make_tuple(t1.dfId(), t1.collisionId(), t1.globalIndex())] - map_z_prop[std::make_tuple(t2.dfId(), t2.collisionId(), t2.globalIndex())] : map_z_prop[std::make_tuple(t2.dfId(), t2.collisionId(), t2.globalIndex())] - map_z_prop[std::make_tuple(t1.dfId(), t1.collisionId(), t1.globalIndex())];
882+
dphi_geom = t1.sign() * t1.pt() > t2.sign() * t2.pt() ? map_phi_prop[std::make_tuple(t1.dfId(), t1.collisionId(), t1.globalIndex())] - map_phi_prop[std::make_tuple(t2.dfId(), t2.collisionId(), t2.globalIndex())] : map_phi_prop[std::make_tuple(t2.dfId(), t2.collisionId(), t2.globalIndex())] - map_phi_prop[std::make_tuple(t1.dfId(), t1.collisionId(), t1.globalIndex())];
883+
o2::math_utils::bringToPMPi(dphi_geom);
884+
rdphi_geom = dielectroncuts.cfg_x_to_go * dphi_geom;
851885
}
852-
o2::math_utils::bringToPMPi(dphi_geom);
853-
if (dielectroncuts.cfg_x_to_go > 0.f && dielectroncuts.cfg_apply_detadphi_geom && std::pow(deta_geom / dielectroncuts.cfg_min_deta_geom, 2) + std::pow(dphi_geom / dielectroncuts.cfg_min_dphi_geom, 2) < 1.f) {
886+
if (dielectroncuts.cfg_x_to_go > 0.f && dielectroncuts.cfg_apply_dzrdphi_geom && std::pow(dz_geom / dielectroncuts.cfg_min_dz_geom, 2) + std::pow(rdphi_geom / dielectroncuts.cfg_min_rdphi_geom, 2) < 1.f) {
854887
return false;
855888
}
856889
} else if constexpr (pairtype == o2::aod::pwgem::dilepton::utils::pairutil::DileptonPairType::kDimuon) {
@@ -901,25 +934,28 @@ struct Dilepton {
901934
if (t1.sign() * t2.sign() < 0) { // ULS
902935
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("uls/hs"), v12.M(), v12.Pt(), pair_dca, weight);
903936
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("uls/hsDeltaP"), dpt, deta, dphi, weight);
904-
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("uls/hGeomDeltaEtaDeltaPhi"), dphi_geom, deta_geom, weight);
937+
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("uls/hGeomDeltaZRDeltaPhi"), rdphi_geom, dz_geom, weight);
905938
if constexpr (pairtype == o2::aod::pwgem::dilepton::utils::pairutil::DileptonPairType::kDielectron) {
906939
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("uls/hMvsPhiV"), phiv, v12.M(), weight);
940+
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("uls/hMvsPhiV_prop"), phiv_prop, mee_prop, weight);
907941
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("uls/hMvsOpAng"), opAng, v12.M(), weight);
908942
}
909943
} else if (t1.sign() > 0 && t2.sign() > 0) { // LS++
910944
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lspp/hs"), v12.M(), v12.Pt(), pair_dca, weight);
911945
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lspp/hsDeltaP"), dpt, deta, dphi, weight);
912-
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lspp/hGeomDeltaEtaDeltaPhi"), dphi_geom, deta_geom, weight);
946+
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lspp/hGeomDeltaZRDeltaPhi"), rdphi_geom, dz_geom, weight);
913947
if constexpr (pairtype == o2::aod::pwgem::dilepton::utils::pairutil::DileptonPairType::kDielectron) {
914948
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lspp/hMvsPhiV"), phiv, v12.M(), weight);
949+
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lspp/hMvsPhiV_prop"), phiv_prop, mee_prop, weight);
915950
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lspp/hMvsOpAng"), opAng, v12.M(), weight);
916951
}
917952
} else if (t1.sign() < 0 && t2.sign() < 0) { // LS--
918953
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lsmm/hs"), v12.M(), v12.Pt(), pair_dca, weight);
919954
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lsmm/hsDeltaP"), dpt, deta, dphi, weight);
920-
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lsmm/hGeomDeltaEtaDeltaPhi"), dphi_geom, deta_geom, weight);
955+
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lsmm/hGeomDeltaZRDeltaPhi"), rdphi_geom, dz_geom, weight);
921956
if constexpr (pairtype == o2::aod::pwgem::dilepton::utils::pairutil::DileptonPairType::kDielectron) {
922957
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lsmm/hMvsPhiV"), phiv, v12.M(), weight);
958+
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lsmm/hMvsPhiV_prop"), phiv_prop, mee_prop, weight);
923959
fRegistry.fill(HIST("Pair/") + HIST(event_pair_types[ev_id]) + HIST("lsmm/hMvsOpAng"), opAng, v12.M(), weight);
924960
}
925961
}
@@ -1364,10 +1400,12 @@ struct Dilepton {
13641400
if (!cut.template IsSelectedPair(t1, t2, d_bz)) {
13651401
return false;
13661402
}
1367-
float deta_geom = t1.sign() * t1.pt() > t2.sign() * t2.pt() ? map_eta_geom[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())] - map_eta_geom[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] : map_eta_geom[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] - map_eta_geom[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())];
1368-
float dphi_geom = t1.sign() * t1.pt() > t2.sign() * t2.pt() ? map_phi_geom[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())] - map_phi_geom[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] : map_phi_geom[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] - map_phi_geom[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())];
1403+
float dz_geom = t1.sign() * t1.pt() > t2.sign() * t2.pt() ? map_z_prop[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())] - map_z_prop[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] : map_z_prop[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] - map_z_prop[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())];
1404+
float dphi_geom = t1.sign() * t1.pt() > t2.sign() * t2.pt() ? map_phi_prop[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())] - map_phi_prop[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] : map_phi_prop[std::make_tuple(ndf, t2.emeventId(), t2.globalIndex())] - map_phi_prop[std::make_tuple(ndf, t1.emeventId(), t1.globalIndex())];
13691405
o2::math_utils::bringToPMPi(dphi_geom);
1370-
if (dielectroncuts.cfg_x_to_go > 0.f && dielectroncuts.cfg_apply_detadphi_geom && std::pow(deta_geom / dielectroncuts.cfg_min_deta_geom, 2) + std::pow(dphi_geom / dielectroncuts.cfg_min_dphi_geom, 2) < 1.f) {
1406+
float rdphi_geom = dielectroncuts.cfg_x_to_go * dphi_geom;
1407+
1408+
if (dielectroncuts.cfg_x_to_go > 0.f && dielectroncuts.cfg_apply_dzrdphi_geom && std::pow(dz_geom / dielectroncuts.cfg_min_dz_geom, 2) + std::pow(rdphi_geom / dielectroncuts.cfg_min_rdphi_geom, 2) < 1.f) {
13711409
return false;
13721410
}
13731411
} else if constexpr (pairtype == o2::aod::pwgem::dilepton::utils::pairutil::DileptonPairType::kDimuon) {

0 commit comments

Comments
 (0)