From dff47c26e4e7916e6ed1a2b27ae34c6c669f2925 Mon Sep 17 00:00:00 2001 From: smaff92 <33285879+smaff92@users.noreply.github.com> Date: Wed, 7 Aug 2024 21:56:50 +0900 Subject: [PATCH] PWGJE: Add matching of phi's in gen-to-rec jets (#7193) * Added K*(892) functionality, as well as a host of optimizations and features * Please consider the following formatting changes * Megalinter fixes for pull 6940 * Please consider the following formatting changes * Fixed rand_r() bug * Please consider the following formatting changes * Bugfix to random evaluator * New update, properly creates the response matrix efficiency on the generator level * Bugfix on the MC-REC matching * Please consider the following formatting changes * Megalinter requests for if/else braces * Please consider the following formatting changes * Reorganized and optimized histogram booking * Please consider the following formatting changes * Please consider the following formatting changes * PWGJE: Added phi gen-to-rec matched finding * Please consider the following formatting changes * Bugfix to fix whitespace/tab issue * Bugfix to remove tab and replace it with whitespace, take 2 * Bugfix to remove tab and replace it with whitespace, take 3 --------- Co-authored-by: ALICE Action Bot --- PWGJE/Tasks/phiInJets.cxx | 215 +++++++++++++++++++++++--------------- 1 file changed, 128 insertions(+), 87 deletions(-) diff --git a/PWGJE/Tasks/phiInJets.cxx b/PWGJE/Tasks/phiInJets.cxx index cec6faae811..f04008495f3 100644 --- a/PWGJE/Tasks/phiInJets.cxx +++ b/PWGJE/Tasks/phiInJets.cxx @@ -197,6 +197,10 @@ struct phiInJets { JEhistos.add("2DGenToRec", "2DGenToRec", kTH2F, {PtAxis, axisPt}); JEhistos.add("2DGenToRec_constrained", "2DGenToRec_constrained", kTH2F, {PtAxis, axisPt}); + JEhistos.add("RespGen_Matrix_MATCHED", "RespGen_Matrix_MATCHED", HistType::kTHnSparseD, {PtAxis, axisPt, PtAxis, axisPt}); // REC(Phi,Jet), GEN(Phi,Jet) + JEhistos.add("RespGen_Matrix_MATCHED_rand0", "RespGen_Matrix_MATCHED_rand0", HistType::kTHnSparseD, {PtAxis, axisPt, PtAxis, axisPt}); // REC(Phi,Jet), GEN(Phi,Jet) + JEhistos.add("RespGen_Matrix_MATCHED_rand1", "RespGen_Matrix_MATCHED_rand1", HistType::kTHnSparseD, {PtAxis, axisPt, PtAxis, axisPt}); // REC(Phi,Jet), GEN(Phi,Jet) + JEhistos.add("hMCTrue_hUSS_INSIDE", "hMCTrue_hUSS_INSIDE", kTH3F, {dRAxis, PtAxis, MinvAxis}); JEhistos.add("hMCTrue_hUSS_INSIDE_1D", "hMCTrue_hUSS_INSIDE_1D", kTH1F, {MinvAxis}); JEhistos.add("hMCTrue_hUSS_INSIDE_1D_2_3", "hMCTrue_hUSS_INSIDE_1D_2_3", kTH1F, {MinvAxis}); @@ -515,30 +519,30 @@ struct phiInJets { ///////////////////////////////////////////////////////////////////////////// // Fill outside Jet - /*if (!jetFlag) { - if (trk1.sign() * trk2.sign() < 0) { - if (!IsMC) { - JEhistos.fill(HIST("hUSS_OUTSIDE_1D"), lResonance.M()); - if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) - JEhistos.fill(HIST("hUSS_OUTSIDE_1D_2_3"), lResonance.M()); - JEhistos.fill(HIST("hUSS_OUTSIDE"), jetpt, lResonance.Pt(), lResonance.M()); - } - - if (IsMC) { - JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_1D"), lResonance.M()); - if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) - JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_1D_2_3"), lResonance.M()); - JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE"), jetpt, lResonance.Pt(), lResonance.M()); - } - - } else if (trk1.sign() * trk2.sign() > 0) { - - JEhistos.fill(HIST("hLSS_OUTSIDE_1D"), lResonance.M()); - if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) - JEhistos.fill(HIST("hLSS_OUTSIDE_1D_2_3"), lResonance.M()); - JEhistos.fill(HIST("hLSS_OUTSIDE"), jetpt, lResonance.Pt(), lResonance.M()); - } - } //! jetflag*/ + // if (!jetFlag) { + // if (trk1.sign() * trk2.sign() < 0) { + // if (!IsMC) { + // JEhistos.fill(HIST("hUSS_OUTSIDE_1D"), lResonance.M()); + // if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) + // JEhistos.fill(HIST("hUSS_OUTSIDE_1D_2_3"), lResonance.M()); + // JEhistos.fill(HIST("hUSS_OUTSIDE"), jetpt, lResonance.Pt(), lResonance.M()); + // } + + // if (IsMC) { + // JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_1D"), lResonance.M()); + // if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) + // JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_1D_2_3"), lResonance.M()); + // JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE"), jetpt, lResonance.Pt(), lResonance.M()); + // } + + // } else if (trk1.sign() * trk2.sign() > 0) { + + // JEhistos.fill(HIST("hLSS_OUTSIDE_1D"), lResonance.M()); + // if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) + // JEhistos.fill(HIST("hLSS_OUTSIDE_1D_2_3"), lResonance.M()); + // JEhistos.fill(HIST("hLSS_OUTSIDE"), jetpt, lResonance.Pt(), lResonance.M()); + // } + // } //! jetflag ///////////////////////////////////////////////////////////////////////////// if (!cfgIsKstar) { if (lResonance.M() > 1.005 && lResonance.M() < 1.035) { @@ -835,24 +839,23 @@ struct phiInJets { if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_INSIDE_1D_2_3"), lResonance.M()); JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_INSIDE"), jetpt, lResonance.Pt(), lResonance.M()); + } + // else if (!jetFlag && mcd_pt.size() > 0) { + // JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE_TRIG_1D"), lResonance.M()); - } /* else if (!jetFlag && mcd_pt.size() > 0) { - JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE_TRIG_1D"), lResonance.M()); - - if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) - JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE_TRIG_1D_2_3"), lResonance.M()); + // if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) + // JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE_TRIG_1D_2_3"), lResonance.M()); - JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE_TRIG"), jetpt, lResonance.Pt(), lResonance.M()); + // JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE_TRIG"), jetpt, lResonance.Pt(), lResonance.M()); - } else if (!jetFlag) { - JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE_1D"), lResonance.M()); + // } else if (!jetFlag) { + // JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE_1D"), lResonance.M()); - if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) - JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE_1D_2_3"), lResonance.M()); + // if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) + // JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE_1D_2_3"), lResonance.M()); - JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE"), jetpt, lResonance.Pt(), lResonance.M()); - }*/ - //! jetflag + // JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE"), jetpt, lResonance.Pt(), lResonance.M()); + // } //! jetflag if (hasJets) { JEhistos.fill(HIST("ptJEHistogramPhi_JetTrigger"), lResonance.Pt()); @@ -1005,24 +1008,24 @@ struct phiInJets { if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_INSIDE_1D_2_3"), lResonance.M()); JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_INSIDE"), jetpt, lResonance.Pt(), lResonance.M()); + } + // else if (!jetFlag && mcp_pt.size() > 0) { + // JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_OUTSIDE_TRIG_1D"), lResonance.M()); - } /* else if (!jetFlag && mcp_pt.size() > 0) { - JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_OUTSIDE_TRIG_1D"), lResonance.M()); - - if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) - JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_OUTSIDE_TRIG_1D_2_3"), lResonance.M()); + // if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) + // JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_OUTSIDE_TRIG_1D_2_3"), lResonance.M()); - JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_OUTSIDE_TRIG"), jetpt, lResonance.Pt(), lResonance.M()); + // JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_OUTSIDE_TRIG"), jetpt, lResonance.Pt(), lResonance.M()); - } else if (!jetFlag) { - JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_OUTSIDE_1D"), lResonance.M()); + // } else if (!jetFlag) { + // JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_OUTSIDE_1D"), lResonance.M()); - if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) - JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_OUTSIDE_1D_2_3"), lResonance.M()); + // if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) + // JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_OUTSIDE_1D_2_3"), lResonance.M()); - JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_OUTSIDE"), jetpt, lResonance.Pt(), lResonance.M()); + // JEhistos.fill(HIST("hMCTrue_nonmatch_hUSS_OUTSIDE"), jetpt, lResonance.Pt(), lResonance.M()); - }*/ //! jetflag + // } //! jetflag ////////////////////////////Phi found if (hasJets) { @@ -1051,7 +1054,8 @@ struct phiInJets { JetMCPTable const& mcpjets, myCompleteJetTracks const& tracks, myCompleteTracks const&, - aod::JMcParticles const& mcParticles) + aod::JMcParticles const& mcParticles, + aod::McParticles const&) { if (cDebugLevel > 0) { @@ -1137,6 +1141,10 @@ struct phiInJets { double phi_dgth_px[2] = {0}; double phi_dgth_py[2] = {0}; double phi_dgth_pz[2] = {0}; + double TEMP_phi_dgth_px[2] = {0}; + double TEMP_phi_dgth_py[2] = {0}; + double TEMP_phi_dgth_pz[2] = {0}; + bool good_daughter[2] = {false}; int dgth_index = 0; // First we check for Forced BR @@ -1152,21 +1160,35 @@ struct phiInJets { auto trk = track.track_as(); if (!trackSelection(trk)) continue; + if (fabs(trk.eta()) > cfgtrkMaxEta) + continue; if (cfgSimPID) { if (!trackPID(trk, true)) continue; } - if (track.globalIndex() == dgth.globalIndex()) { - phi_dgth_px[dgth_index] = track.px(); - phi_dgth_py[dgth_index] = track.py(); - phi_dgth_pz[dgth_index] = track.pz(); + if (!trk.has_mcParticle()) + continue; + auto part = trk.mcParticle(); + if (part.globalIndex() == dgth.globalIndex()) { + TEMP_phi_dgth_px[dgth_index] = track.px(); + TEMP_phi_dgth_py[dgth_index] = track.py(); + TEMP_phi_dgth_pz[dgth_index] = track.pz(); good_daughter[dgth_index] = true; dgth_index++; - } - } - } - } - } else { + if (dgth_index == 2) { + phi_dgth_px[0] = TEMP_phi_dgth_px[0]; + phi_dgth_py[0] = TEMP_phi_dgth_py[0]; + phi_dgth_pz[0] = TEMP_phi_dgth_pz[0]; + phi_dgth_px[1] = TEMP_phi_dgth_px[1]; + phi_dgth_py[1] = TEMP_phi_dgth_py[1]; + phi_dgth_pz[1] = TEMP_phi_dgth_pz[1]; + break; + } + } // index check + } // track loop + } // mc daughter loop + } // check if particle has daughters + } else { // check for kstar if (mcParticle.has_daughters()) for (auto& dgth : mcParticle.daughters_as()) if (fabs(dgth.pdgCode()) != 321 || fabs(dgth.pdgCode()) != 211) @@ -1187,9 +1209,9 @@ struct phiInJets { lDecayDaughter1_REC.SetXYZM(phi_dgth_px[0], phi_dgth_py[0], phi_dgth_pz[0], massKa); lDecayDaughter2_REC.SetXYZM(phi_dgth_px[1], phi_dgth_py[1], phi_dgth_pz[1], massKa); lResonance_REC = lDecayDaughter1_REC + lDecayDaughter2_REC; - if (cDebugLevel > 0) - if (good_daughter[0] && good_daughter[1]) - std::cout << "Reconstructed level phi pT: " << lResonance_REC.Pt() << std::endl; + // if (cDebugLevel > 0) + // if (good_daughter[0] && good_daughter[1]) + // std::cout << "Reconstructed level phi pT: " << lResonance_REC.Pt() << std::endl; bool jetFlag = false; for (int i = 0; i < mcp_pt.size(); i++) { @@ -1211,6 +1233,25 @@ struct phiInJets { } if (jetFlag) { + if (cDebugLevel > 0) { + std::cout << "******************************************" << std::endl; + std::cout << "GEN TO REC LEVEL: " << std::endl; + std::cout << "Rec. Phi Pt: " << lResonance_REC.Pt() << std::endl; + std::cout << "Rec. Phi Eta: " << lResonance_REC.Eta() << std::endl; + std::cout << "Rec. Jet Pt: " << jetpt_mcd << std::endl; + std::cout << "Gen. Phi Pt: " << lResonance.Pt() << std::endl; + std::cout << "Gen. Phi Eta: " << lResonance.Eta() << std::endl; + std::cout << "Gen. Jet Pt: " << jetpt_mcp << std::endl; + std::cout << "******************************************" << std::endl; + } + + JEhistos.fill(HIST("RespGen_Matrix_MATCHED"), lResonance_REC.Pt(), jetpt_mcd, lResonance.Pt(), jetpt_mcp); + unsigned int seed = static_cast(std::chrono::system_clock::now().time_since_epoch().count()); + int dice = rand_r(&seed) % 2; + if (dice > 0) + JEhistos.fill(HIST("RespGen_Matrix_MATCHED_rand0"), lResonance_REC.Pt(), jetpt_mcd, lResonance.Pt(), jetpt_mcp); + else + JEhistos.fill(HIST("RespGen_Matrix_MATCHED_rand1"), lResonance_REC.Pt(), jetpt_mcd, lResonance.Pt(), jetpt_mcp); JEhistos.fill(HIST("2DGenToRec"), lResonance.Pt(), jetpt_mcp); // //check constrained eff @@ -1221,25 +1262,24 @@ struct phiInJets { if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) JEhistos.fill(HIST("hMCTrue_hUSS_INSIDE_1D_2_3"), lResonance.M()); JEhistos.fill(HIST("hMCTrue_hUSS_INSIDE"), jetpt_mcp, lResonance.Pt(), lResonance.M()); + } + // else if (!jetFlag && mcp_pt.size() > 0) { + // JEhistos.fill(HIST("hMCTrue_hUSS_OUTSIDE_TRIG_1D"), lResonance.M()); - } /* else if (!jetFlag && mcp_pt.size() > 0) { - JEhistos.fill(HIST("hMCTrue_hUSS_OUTSIDE_TRIG_1D"), lResonance.M()); - - if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) - JEhistos.fill(HIST("hMCTrue_hUSS_OUTSIDE_TRIG_1D_2_3"), lResonance.M()); + // if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) + // JEhistos.fill(HIST("hMCTrue_hUSS_OUTSIDE_TRIG_1D_2_3"), lResonance.M()); - JEhistos.fill(HIST("hMCTrue_hUSS_OUTSIDE_TRIG"), jetpt_mcp, lResonance.Pt(), lResonance.M()); + // JEhistos.fill(HIST("hMCTrue_hUSS_OUTSIDE_TRIG"), jetpt_mcp, lResonance.Pt(), lResonance.M()); - } else if (!jetFlag) { - JEhistos.fill(HIST("hMCTrue_hUSS_OUTSIDE_1D"), lResonance.M()); + // } else if (!jetFlag) { + // JEhistos.fill(HIST("hMCTrue_hUSS_OUTSIDE_1D"), lResonance.M()); - if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) - JEhistos.fill(HIST("hMCTrue_hUSS_OUTSIDE_1D_2_3"), lResonance.M()); + // if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) + // JEhistos.fill(HIST("hMCTrue_hUSS_OUTSIDE_1D_2_3"), lResonance.M()); - JEhistos.fill(HIST("hMCTrue_hUSS_OUTSIDE"), jetpt_mcp, lResonance.Pt(), lResonance.M()); + // JEhistos.fill(HIST("hMCTrue_hUSS_OUTSIDE"), jetpt_mcp, lResonance.Pt(), lResonance.M()); - }*/ - //! jetflag + // } //! jetflag } // chech for phi } // MC Particles } // main fcn @@ -1421,6 +1461,7 @@ struct phiInJets { if (jetFlag) { // Fill Resp. Matrix if (cDebugLevel > 0) { std::cout << "******************************************" << std::endl; + std::cout << "REC TO GEN LEVEL: " << std::endl; std::cout << "Rec. Phi Pt: " << lResonance.Pt() << std::endl; std::cout << "Rec. Jet Pt: " << jetpt_mcd << std::endl; std::cout << "Gen. Phi Pt: " << mothers1Pt[0] << std::endl; @@ -1447,24 +1488,24 @@ struct phiInJets { if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) JEhistos.fill(HIST("hMCRec_hUSS_INSIDE_1D_2_3"), lResonance.M()); JEhistos.fill(HIST("hMCRec_hUSS_INSIDE"), jetpt_mcd, lResonance.Pt(), lResonance.M()); + } + // else if (!jetFlag && mcd_pt.size() > 0) { + // JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_TRIG_1D"), lResonance.M()); - } /* else if (!jetFlag && mcd_pt.size() > 0) { - JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_TRIG_1D"), lResonance.M()); - - if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) - JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_TRIG_1D_2_3"), lResonance.M()); + // if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) + // JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_TRIG_1D_2_3"), lResonance.M()); - JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_TRIG"), jetpt_mcd, lResonance.Pt(), lResonance.M()); + // JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_TRIG"), jetpt_mcd, lResonance.Pt(), lResonance.M()); - } else if (!jetFlag) { - JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_1D"), lResonance.M()); + // } else if (!jetFlag) { + // JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_1D"), lResonance.M()); - if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) - JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_1D_2_3"), lResonance.M()); + // if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) + // JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE_1D_2_3"), lResonance.M()); - JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE"), jetpt_mcd, lResonance.Pt(), lResonance.M()); + // JEhistos.fill(HIST("hMCRec_hUSS_OUTSIDE"), jetpt_mcd, lResonance.Pt(), lResonance.M()); - }*/ //! jetflag + // } //! jetflag } // pass track cut } // has mc particle