From 16aed78472168d573e35bac8287a0d70f8a6c6b2 Mon Sep 17 00:00:00 2001 From: Shunsuke-Kurita Date: Mon, 13 Apr 2026 17:37:04 +0900 Subject: [PATCH 1/4] Add electron-muon mixing --- PWGDQ/Tasks/tableReader_withAssoc.cxx | 333 +++++++++++++------------- 1 file changed, 166 insertions(+), 167 deletions(-) diff --git a/PWGDQ/Tasks/tableReader_withAssoc.cxx b/PWGDQ/Tasks/tableReader_withAssoc.cxx index 8fa36edcc28..e899a5844d1 100644 --- a/PWGDQ/Tasks/tableReader_withAssoc.cxx +++ b/PWGDQ/Tasks/tableReader_withAssoc.cxx @@ -26,48 +26,42 @@ #include "Common/CCDB/EventSelectionParams.h" #include "Common/Core/TableHelper.h" -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - +#include "CCDB/BasicCCDBManager.h" +#include "DataFormatsParameters/GRPMagField.h" +#include "DataFormatsParameters/GRPObject.h" +#include "DetectorsBase/GeometryManager.h" +#include "DetectorsBase/Propagator.h" +#include "Field/MagneticField.h" +#include "Framework/ASoAHelpers.h" +#include "Framework/AnalysisDataModel.h" +#include "Framework/AnalysisHelpers.h" +#include "Framework/AnalysisTask.h" +#include "Framework/Configurable.h" +#include "Framework/OutputObjHeader.h" +#include "Framework/runDataProcessing.h" +#include "ITSMFTBase/DPLAlpideParam.h" + +#include "TGeoGlobalMagField.h" #include +#include +#include #include #include -#include -#include #include +#include #include -#include - #include -#include -#include #include #include -#include #include #include -#include #include #include +#include +#include #include #include -#include #include using std::cout; @@ -78,7 +72,6 @@ using namespace o2; using namespace o2::framework; using namespace o2::framework::expressions; using namespace o2::aod; -using namespace o2::common::core; // Some definitions namespace o2::aod @@ -202,7 +195,6 @@ DECLARE_SOA_TABLE(JPsieeCandidates, "AOD", "DQPSEUDOPROPER", dqanalysisflags::Ma // Declarations of various short names using MyEvents = soa::Join; -using MyEventsBasic = soa::Join; using MyEventsMultExtra = soa::Join; using MyEventsMultExtraQVector = soa::Join; using MyEventsZdc = soa::Join; @@ -284,7 +276,6 @@ struct AnalysisEventSelection { // TODO: Provide the mixing variables and binning directly via configurables (e.g. vectors of float) Configurable fConfigMixingVariables{"cfgMixingVars", "", "Mixing configs separated by a comma, default no mixing"}; - Configurable fConfigMixingVariablesJson{"cfgMixingVarsJSON", "", "Mixing configs in JSON format"}; Configurable fConfigEventCuts{"cfgEventCuts", "eventStandard", "Event selection"}; Configurable fConfigEventCutsJSON{"cfgEventCutsJSON", "", "Additional event cuts specified in JSON format"}; Configurable fConfigAddEventHistogram{"cfgAddEventHistogram", "", "Comma separated list of histograms"}; @@ -316,7 +307,7 @@ struct AnalysisEventSelection { void init(o2::framework::InitContext& context) { - bool isAnyProcessEnabled = context.mOptions.get("processSkimmed") || context.mOptions.get("processSkimmedBasic") || context.mOptions.get("processSkimmedWithZdc") || context.mOptions.get("processSkimmedWithMultExtra") || context.mOptions.get("processSkimmedWithMultExtraZdc") || context.mOptions.get("processSkimmedWithMultExtraZdcFit") || context.mOptions.get("processSkimmedWithQvectorCentr"); + bool isAnyProcessEnabled = context.mOptions.get("processSkimmed") || context.mOptions.get("processSkimmedWithZdc") || context.mOptions.get("processSkimmedWithMultExtra") || context.mOptions.get("processSkimmedWithMultExtraZdc") || context.mOptions.get("processSkimmedWithMultExtraZdcFit") || context.mOptions.get("processSkimmedWithQvectorCentr"); bool isDummyEnabled = context.mOptions.get("processDummy"); if (isDummyEnabled) { @@ -365,18 +356,12 @@ struct AnalysisEventSelection { } TString mixVarsString = fConfigMixingVariables.value; - TString mixVarsJsonString = fConfigMixingVariablesJson.value; std::unique_ptr objArray(mixVarsString.Tokenize(",")); - if (objArray->GetEntries() > 0 || mixVarsJsonString != "") { + if (objArray->GetEntries() > 0) { fMixHandler = new MixingHandler("mixingHandler", "mixing handler"); fMixHandler->Init(); - if (objArray->GetEntries() > 0) { - for (int iVar = 0; iVar < objArray->GetEntries(); ++iVar) { - dqmixing::SetUpMixing(fMixHandler, objArray->At(iVar)->GetName()); - } - } - if (mixVarsJsonString != "") { - dqmixing::SetUpMixingFromJSON(fMixHandler, mixVarsJsonString.Data()); + for (int iVar = 0; iVar < objArray->GetEntries(); ++iVar) { + dqmixing::SetUpMixing(fMixHandler, objArray->At(iVar)->GetName()); } } @@ -543,11 +528,6 @@ struct AnalysisEventSelection { runEventSelection(events); publishSelections(events); } - void processSkimmedBasic(MyEventsBasic const& events) - { - runEventSelection(events); - publishSelections(events); - } void processSkimmedWithZdc(MyEventsZdc const& events) { runEventSelection(events); @@ -573,13 +553,12 @@ struct AnalysisEventSelection { runEventSelection(events); publishSelections(events); } - void processDummy(MyEventsBasic&) + void processDummy(MyEvents&) { // do nothing } PROCESS_SWITCH(AnalysisEventSelection, processSkimmed, "Run event selection on DQ skimmed events", false); - PROCESS_SWITCH(AnalysisEventSelection, processSkimmedBasic, "Run event selection on DQ skimmed events with basic tables", false); PROCESS_SWITCH(AnalysisEventSelection, processSkimmedWithZdc, "Run event selection on DQ skimmed events, with ZDC", false); PROCESS_SWITCH(AnalysisEventSelection, processSkimmedWithMultExtra, "Run event selection on DQ skimmed events, with mult extra", false); PROCESS_SWITCH(AnalysisEventSelection, processSkimmedWithMultExtraZdc, "Run event selection on DQ skimmed events, with mult extra and ZDC", false); @@ -835,7 +814,7 @@ struct AnalysisTrackSelection { { runTrackSelection(assocs, events, tracks); } - void processDummy(MyEventsBasic&) + void processDummy(MyEvents&) { // do nothing } @@ -1046,7 +1025,7 @@ struct AnalysisMuonSelection { { runMuonSelection(assocs, events, muons); } - void processDummy(MyEventsBasic&) + void processDummy(MyEvents&) { // do nothing } @@ -1195,7 +1174,7 @@ struct AnalysisPrefilterSelection { } // end loop over combinations } - void processBarrelSkimmed(MyEventsBasic const& events, soa::Join const& assocs, MyBarrelTracks const& tracks) + void processBarrelSkimmed(MyEvents const& events, soa::Join const& assocs, MyBarrelTracks const& tracks) { fPrefilterMap.clear(); @@ -1228,7 +1207,7 @@ struct AnalysisPrefilterSelection { } } - void processDummy(MyEventsBasic&) + void processDummy(MyEvents&) { // do nothing } @@ -1333,14 +1312,13 @@ struct AnalysisSameEventPairing { int fNCutsBarrel; int fNCutsMuon; int fNPairCuts; - int fNPairPerEvent; bool fEnableBarrelMixingHistos; bool fEnableBarrelHistos; bool fEnableMuonHistos; bool fEnableMuonMixingHistos; bool fEnableBarrelMuonHistos; - // bool fEnableBarrelMuonMixingHistos; + bool fEnableBarrelMuonMixingHistos; NoBinningPolicy hashBin; @@ -1355,9 +1333,10 @@ struct AnalysisSameEventPairing { fEnableMuonHistos = context.mOptions.get("processAllSkimmed") || context.mOptions.get("processMuonOnlySkimmed") || context.mOptions.get("processMuonOnlySkimmedMultExtra") || context.mOptions.get("processMuonOnlySkimmedFlow") || context.mOptions.get("processMixingMuonSkimmed"); fEnableMuonMixingHistos = context.mOptions.get("processMixingAllSkimmed") || context.mOptions.get("processMixingMuonSkimmed"); fEnableBarrelMuonHistos = context.mOptions.get("processElectronMuonSkimmed"); + fEnableBarrelMuonMixingHistos = context.mOptions.get("processMixingElectronMuonSkimmed"); if (context.mOptions.get("processDummy")) { - if (fEnableBarrelHistos || fEnableBarrelMixingHistos || fEnableMuonHistos || fEnableMuonMixingHistos || fEnableBarrelMuonHistos) { + if (fEnableBarrelHistos || fEnableBarrelMixingHistos || fEnableMuonHistos || fEnableMuonMixingHistos || fEnableBarrelMuonHistos || fEnableBarrelMuonMixingHistos) { LOG(fatal) << "No other processing tasks should be enabled if the processDummy is enabled!!"; } return; @@ -1599,7 +1578,7 @@ struct AnalysisSameEventPairing { VarManager::SetupMatLUTFwdDCAFitter(fLUT); } - if (fEnableBarrelMuonHistos) { + if (fEnableBarrelMuonHistos || fEnableBarrelMuonMixingHistos) { for (int iTrack = 0; iTrack < fNCutsBarrel; ++iTrack) { TString trackCutName = fTrackCuts[iTrack]; if (objArrayTrackCuts->FindObject(trackCutName.Data()) == nullptr) @@ -1616,15 +1595,14 @@ struct AnalysisSameEventPairing { Form("PairsEleMuSEMM_%s_%s", trackCutName.Data(), muonCutName.Data())}; histNames += Form("%s;%s;%s;", names[0].Data(), names[1].Data(), names[2].Data()); int index = iTrack * fNCutsMuon + iMuon; + if (fEnableBarrelMuonMixingHistos) { + names.push_back(Form("PairsEleMuMEPM_%s_%s", trackCutName.Data(), muonCutName.Data())); + names.push_back(Form("PairsEleMuMEPP_%s_%s", trackCutName.Data(), muonCutName.Data())); + names.push_back(Form("PairsEleMuMEMM_%s_%s", trackCutName.Data(), muonCutName.Data())); + histNames += Form("%s;%s;%s;", names[3].Data(), names[4].Data(), names[5].Data()); + } fTrackMuonHistNames[index] = names; - // if (fEnableBarrelMuonMixingHistos) { - // names.push_back(Form("PairsBarrelMuonMEPM_%s_%s", trackCutName.Data(), muonCutName.Data())); - // names.push_back(Form("PairsBarrelMuonMEPP_%s_%s", trackCutName.Data(), muonCutName.Data())); - // names.push_back(Form("PairsBarrelMuonMEMM_%s_%s", trackCutName.Data(), muonCutName.Data())); - // histNames += Form("%s;%s;%s;", names[3].Data(), names[4].Data(), names[5].Data()); - // } - TString cutNamesStr = fConfigCuts.pair.value; if (!cutNamesStr.IsNull()) { std::unique_ptr objArrayPair(cutNamesStr.Tokenize(",")); @@ -1635,6 +1613,12 @@ struct AnalysisSameEventPairing { Form("PairsEleMuSEPP_%s_%s_%s", trackCutName.Data(), muonCutName.Data(), objArrayPair->At(iPairCut)->GetName()), Form("PairsEleMuSEMM_%s_%s_%s", trackCutName.Data(), muonCutName.Data(), objArrayPair->At(iPairCut)->GetName())}; histNames += Form("%s;%s;%s;", names[0].Data(), names[1].Data(), names[2].Data()); + if (fEnableBarrelMuonMixingHistos) { + names.push_back(Form("PairsEleMuMEPM_%s_%s_%s", trackCutName.Data(), muonCutName.Data(), objArrayPair->At(iPairCut)->GetName())); + names.push_back(Form("PairsEleMuMEPP_%s_%s_%s", trackCutName.Data(), muonCutName.Data(), objArrayPair->At(iPairCut)->GetName())); + names.push_back(Form("PairsEleMuMEMM_%s_%s_%s", trackCutName.Data(), muonCutName.Data(), objArrayPair->At(iPairCut)->GetName())); + histNames += Form("%s;%s;%s;", names[3].Data(), names[4].Data(), names[5].Data()); + } index = iTrack * (fNCutsMuon * nPairCuts) + iMuon * nPairCuts + iPairCut; fTrackMuonHistNames[index] = names; } @@ -1649,14 +1633,8 @@ struct AnalysisSameEventPairing { fHistMan->SetDefaultVarNames(VarManager::fgVariableNames, VarManager::fgVariableUnits); VarManager::SetCollisionSystem((TString)fConfigOptions.collisionSystem, fConfigOptions.centerMassEnergy); // set collision system and center of mass energy DefineHistograms(fHistMan, histNames.Data(), fConfigAddSEPHistogram.value.data()); // define all histograms - if (fEnableBarrelHistos) { - DefineHistograms(fHistMan, "PairingSEQA", "sameevent-pairing"); // histograms for QA of the pairing - }; - if (fEnableBarrelMixingHistos) { - DefineHistograms(fHistMan, "PairingMEQA", "mixedevent-pairing"); // histograms for QA of the pairing - }; - dqhistograms::AddHistogramsFromJSON(fHistMan, fConfigAddJSONHistograms.value.c_str()); // ad-hoc histograms via JSON - VarManager::SetUseVars(fHistMan->GetUsedVars()); // provide the list of required variables so that VarManager knows what to fill + dqhistograms::AddHistogramsFromJSON(fHistMan, fConfigAddJSONHistograms.value.c_str()); // ad-hoc histograms via JSON + VarManager::SetUseVars(fHistMan->GetUsedVars()); // provide the list of required variables so that VarManager knows what to fill fOutputList.setObject(fHistMan->GetMainHistogramList()); } } @@ -1763,7 +1741,6 @@ struct AnalysisSameEventPairing { constexpr bool eventHasQvectorCentr = ((TEventFillMap & VarManager::ObjTypes::CollisionQvect) > 0); constexpr bool trackHasCov = ((TTrackFillMap & VarManager::ObjTypes::TrackCov) > 0 || (TTrackFillMap & VarManager::ObjTypes::ReducedTrackBarrelCov) > 0); bool isSelectedBDT = false; - fNPairPerEvent = 0; for (auto& event : events) { if (!event.isEventSelected_bit(0)) { @@ -1811,7 +1788,6 @@ struct AnalysisSameEventPairing { twoTrackFilter |= (static_cast(1) << 31); } - fNPairPerEvent++; VarManager::FillPair(t1, t2); // compute quantities which depend on the associated collision, such as DCA if (fConfigOptions.propTrack) { @@ -2117,10 +2093,6 @@ struct AnalysisSameEventPairing { } } // end loop (cuts) } // end loop over pairs of track associations - VarManager::fgValues[VarManager::kNPairsPerEvent] = fNPairPerEvent; - if (fEnableBarrelHistos && fConfigQA) { - fHistMan->FillHistClass("PairingSEQA", VarManager::fgValues); - } } // end loop over events } @@ -2140,7 +2112,6 @@ struct AnalysisSameEventPairing { } auto t1 = a1.template reducedtrack_as(); auto t2 = a2.template reducedtrack_as(); - fNPairPerEvent++; VarManager::FillPairME(t1, t2); if constexpr ((TEventFillMap & VarManager::ObjTypes::ReducedEventQvector) > 0) { VarManager::FillPairVn(t1, t2); @@ -2306,12 +2277,102 @@ struct AnalysisSameEventPairing { auto assocs2 = assocs.sliceBy(preSlice, event2.globalIndex()); assocs2.bindExternalIndices(&events); - fNPairPerEvent = 0; runMixedPairing(assocs1, assocs2, tracks, tracks); - VarManager::fgValues[VarManager::kNPairsPerEvent] = fNPairPerEvent; - if (fEnableBarrelMixingHistos && fConfigQA) { - fHistMan->FillHistClass("PairingMEQA", VarManager::fgValues); + } // end event loop + } + + template + void runEmuMixedPairing(TAssoc1 const& assocs1, TAssoc2 const& assocs2, TTracks1 const& /*tracks1*/, TTracks2 const& /*tracks2*/) + { + const auto& histNames = fTrackMuonHistNames; + int sign1 = 0; + int sign2 = 0; + int nPairCuts = (fPairCuts.size() > 0) ? fPairCuts.size() : 1; + constexpr bool eventHasQvector = ((TEventFillMap & VarManager::ObjTypes::ReducedEventQvector) > 0); + constexpr bool eventHasQvectorCentr = ((TEventFillMap & VarManager::ObjTypes::CollisionQvect) > 0); + + for (auto& a1 : assocs1) { + if (!(a1.isBarrelSelected_raw() & fTrackFilterMask)) { + continue; + } + for (auto& a2 : assocs2) { + if (!(a2.isMuonSelected_raw() & fMuonFilterMask)) { + continue; + } + + auto t1 = a1.template reducedtrack_as(); + auto t2 = a2.template reducedmuon_as(); + sign1 = t1.sign(); + sign2 = t2.sign(); + + VarManager::FillPairME(t1, t2); + if constexpr (eventHasQvector) { + VarManager::FillPairVn(t1, t2); + } + if constexpr (eventHasQvectorCentr) { + VarManager::FillPairVn(t1, t2); + } + + for (int iTrack = 0; iTrack < fNCutsBarrel; ++iTrack) { + if (!(a1.isBarrelSelected_raw() & (1u << iTrack))) { + continue; + } + for (int iMuon = 0; iMuon < fNCutsMuon; ++iMuon) { + if (!(a2.isMuonSelected_raw() & (1u << iMuon))) { + continue; + } + for (int iPairCut = 0; iPairCut < nPairCuts; ++iPairCut) { + if (!fPairCuts.empty()) { + AnalysisCompositeCut cut = fPairCuts.at(iPairCut); + if (!cut.IsSelected(VarManager::fgValues)) { + continue; + } + } + int index = iTrack * (fNCutsMuon * nPairCuts) + iMuon * nPairCuts + iPairCut; + auto itHist = histNames.find(index); + if (itHist == histNames.end() || itHist->second.size() < 6) { + continue; + } + if (sign1 * sign2 < 0) { + fHistMan->FillHistClass(itHist->second[3].Data(), VarManager::fgValues); + } else { + if (sign1 > 0) { + fHistMan->FillHistClass(itHist->second[4].Data(), VarManager::fgValues); + } else { + fHistMan->FillHistClass(itHist->second[5].Data(), VarManager::fgValues); + } + } + } // end pair cut loop + } // end muon cut loop + } // end barrel cut loop } + } + } + + template + void runEmuSameSideMixing(TEvents& events, Preslice& preslice1, TTrackAssocs const& assocs1, TTracks const& tracks1, + Preslice& preslice2, TMuonAssocs const& assocs2, TMuons const& tracks2) + { + events.bindExternalIndices(&assocs1); + events.bindExternalIndices(&assocs2); + int mixingDepth = fConfigMixingDepth.value; + for (auto& [event1, event2] : selfCombinations(hashBin, mixingDepth, -1, events, events)) { + VarManager::ResetValues(0, VarManager::kNVars); + VarManager::FillEvent(event1, VarManager::fgValues); + + auto groupedAssocs1 = assocs1.sliceBy(preslice1, event1.globalIndex()); + groupedAssocs1.bindExternalIndices(&events); + if (groupedAssocs1.size() == 0) { + continue; + } + + auto groupedAssocs2 = assocs2.sliceBy(preslice2, event2.globalIndex()); + groupedAssocs2.bindExternalIndices(&events); + if (groupedAssocs2.size() == 0) { + continue; + } + + runEmuMixedPairing(groupedAssocs1, groupedAssocs2, tracks1, tracks2); } // end event loop } @@ -2540,6 +2601,7 @@ struct AnalysisSameEventPairing { runSameSideMixing(events, muonAssocs, muons, muonAssocsPerCollision); } +<<<<<<< HEAD void processMixingMuonSkimmedFlow(soa::Filtered& events, soa::Join const& muonAssocs, MyMuonTracksWithCovWithAmbiguities const& muons) { @@ -2547,6 +2609,16 @@ struct AnalysisSameEventPairing { } void processDummy(MyEventsBasic&) +======= + void processMixingElectronMuonSkimmed(soa::Filtered& events, + soa::Join const& barrelAssocs, aod::ReducedTracks const& barrelTracks, + soa::Join const& muonAssocs, MyMuonTracksWithCovWithAmbiguities const& muons) + { + runEmuSameSideMixing(events, trackEmuAssocsPerCollision, barrelAssocs, barrelTracks, muonAssocsPerCollision, muonAssocs, muons); + } + + void processDummy(MyEvents&) +>>>>>>> d32d3fbdc (Add electron-muon mixing) { // do nothing } @@ -2567,7 +2639,11 @@ struct AnalysisSameEventPairing { PROCESS_SWITCH(AnalysisSameEventPairing, processMixingBarrelSkimmedFlow, "Run barrel type mixing pairing, with flow, with skimmed tracks", false); PROCESS_SWITCH(AnalysisSameEventPairing, processMixingBarrelWithQvectorCentrSkimmedNoCov, "Run barrel type mixing pairing, with skimmed tracks and with Qvector from central framework", false); PROCESS_SWITCH(AnalysisSameEventPairing, processMixingMuonSkimmed, "Run muon type mixing pairing, with skimmed muons", false); +<<<<<<< HEAD PROCESS_SWITCH(AnalysisSameEventPairing, processMixingMuonSkimmedFlow, "Run muon type mixing pairing, with skimmed muons and flow", false); +======= + PROCESS_SWITCH(AnalysisSameEventPairing, processMixingElectronMuonSkimmed, "Run electron-muon mixing pairing, with skimmed tracks/muons", false); +>>>>>>> d32d3fbdc (Add electron-muon mixing) PROCESS_SWITCH(AnalysisSameEventPairing, processDummy, "Dummy function, enabled only if none of the others are enabled", true); }; @@ -3293,7 +3369,7 @@ struct AnalysisAsymmetricPairing { runThreeProng(events, trackAssocsPerCollision, barrelAssocs, barrelTracks, VarManager::kTripleCandidateToPKPi); } - void processDummy(MyEventsBasic&) + void processDummy(MyEvents&) { // do nothing } @@ -3340,10 +3416,6 @@ struct AnalysisDileptonTrack { Configurable> fConfigFitmassEC{"cfgTFitmassEC", std::vector{-0.541438, 2.8, 3.2}, "parameter from the fit fuction and fit range"}; Configurable> fConfigTransRange{"cfgTransRange", std::vector{0.333333, 0.666667}, "Transverse region for the energy correlstor analysis"}; - Configurable fConfigApplyEfficiency{"cfgApplyEfficiency", false, "If true, apply efficiency correction for the energy correlator study"}; - Configurable fConfigApplyEfficiencyME{"cfgApplyEfficiencyME", false, "If true, apply efficiency correction for the energy correlator study"}; - Configurable fConfigAccCCDBPath{"AccCCDBPath", "Users/y/yalin/pptest/test2", "Path of the efficiency corrections"}; - int fCurrentRun; // needed to detect if the run changed and trigger update of calibrations etc. int fNCuts; // number of dilepton leg cuts int fNLegCuts; @@ -3377,12 +3449,6 @@ struct AnalysisDileptonTrack { TF1* fMassBkg = nullptr; - TH2F* hAcceptance_rec; - TH2F* hAcceptance_gen; - TH1F* hEfficiency_dilepton; - TH1F* hEfficiency_hadron; - TH1F* hMasswindow; - void init(o2::framework::InitContext& context) { bool isBarrel = context.mOptions.get("processBarrelSkimmed"); @@ -3640,22 +3706,6 @@ struct AnalysisDileptonTrack { } } - void initAccFromCCDB(uint64_t timestamp) - { - TList* listAccs = fCCDB->getForTimeStamp(fConfigAccCCDBPath, timestamp); - if (!listAccs) { - LOG(fatal) << "Problem getting TList object with efficiencies!"; - } - hEfficiency_dilepton = static_cast(listAccs->FindObject("hEfficiency_dilepton")); - hEfficiency_hadron = static_cast(listAccs->FindObject("hEfficiency_hadron")); - hAcceptance_rec = static_cast(listAccs->FindObject("hAcceptance_rec")); - hAcceptance_gen = static_cast(listAccs->FindObject("hAcceptance_gen")); - hMasswindow = static_cast(listAccs->FindObject("hMasswindow")); - if (!hAcceptance_rec || !hAcceptance_gen || !hEfficiency_dilepton || !hEfficiency_hadron || !hMasswindow) { - LOG(fatal) << "Problem getting histograms from the TList object with efficiencies!"; - } - } - // Template function to run pair - hadron combinations template void runDileptonHadron(TEvent const& event, TTrackAssocs const& assocs, TTracks const& tracks, TDileptons const& dileptons) @@ -3675,7 +3725,7 @@ struct AnalysisDileptonTrack { } // dilepton rap cut float rap = dilepton.rap(); - if (fConfigUseRapcut && std::abs(rap) > fConfigDileptonRapCutAbs) + if (fConfigUseRapcut && abs(rap) > fConfigDileptonRapCutAbs) continue; VarManager::FillTrack(dilepton, fValuesDilepton); @@ -3732,21 +3782,8 @@ struct AnalysisDileptonTrack { VarManager::FillDileptonTrackVertexing(event, lepton1, lepton2, track, fValuesHadron); // for the energy correlator analysis - float Effweight_rec = 1.0f; - if (fConfigApplyEfficiency) { - float dilepton_eta = dilepton.eta(); - float dilepton_phi = dilepton.phi(); - float hadron_eta = track.eta(); - float hadron_phi = track.phi(); - float deltaphi = RecoDecay::constrainAngle(dilepton_phi - hadron_phi, -0.5 * o2::constants::math::PI); - Effweight_rec = hAcceptance_rec->Interpolate(dilepton_eta - hadron_eta, deltaphi); - float Effdilepton = hEfficiency_dilepton->Interpolate(dilepton.pt()); - float Masswindow = hMasswindow->Interpolate(dilepton.pt()); - float Effhadron = hEfficiency_hadron->Interpolate(track.pt()); - Effweight_rec = Effweight_rec * Effdilepton * Effhadron * Masswindow; - } std::vector fTransRange = fConfigTransRange; - VarManager::FillEnergyCorrelatorTriple(lepton1, lepton2, track, fValuesHadron, fTransRange[0], fTransRange[1], fConfigApplyMassEC, fMassBkg->GetRandom(), 1. / Effweight_rec); + VarManager::FillEnergyCorrelator(dilepton, track, fValuesHadron, fTransRange[0], fTransRange[1], fConfigApplyMassEC, fMassBkg->GetRandom()); // table to be written out for ML analysis BmesonsTable(event.runNumber(), event.globalIndex(), event.timestamp(), fValuesHadron[VarManager::kPairMass], dilepton.mass(), fValuesHadron[VarManager::kDeltaMass], fValuesHadron[VarManager::kPairPt], fValuesHadron[VarManager::kPairEta], fValuesHadron[VarManager::kPairPhi], fValuesHadron[VarManager::kPairRap], @@ -3851,9 +3888,6 @@ struct AnalysisDileptonTrack { } if (fCurrentRun != events.begin().runNumber()) { // start: runNumber initParamsFromCCDB(events.begin().timestamp()); - if (fConfigApplyEfficiency) { - initAccFromCCDB(events.begin().timestamp()); - } fCurrentRun = events.begin().runNumber(); } // end: runNumber for (auto& event : events) { @@ -3906,19 +3940,11 @@ struct AnalysisDileptonTrack { void processBarrelMixedEvent(soa::Filtered& events, soa::Filtered> const& assocs, - MyBarrelTracksWithCov const& tracks, soa::Filtered const& dileptons) + MyBarrelTracksWithCov const&, soa::Filtered const& dileptons) { if (events.size() == 0) { return; } - - if (fCurrentRun != events.begin().runNumber()) { // start: runNumber - if (fConfigApplyEfficiency) { - initAccFromCCDB(events.begin().timestamp()); - } - fCurrentRun = events.begin().runNumber(); - } // end: runNumber - events.bindExternalIndices(&dileptons); events.bindExternalIndices(&assocs); @@ -3950,41 +3976,18 @@ struct AnalysisDileptonTrack { // loop over dileptons for (auto dilepton : evDileptons) { - // get full track info of tracks based on the index - auto lepton1 = tracks.rawIteratorAt(dilepton.index0Id()); - auto lepton2 = tracks.rawIteratorAt(dilepton.index1Id()); - // Check that the dilepton has zero charge - if (dilepton.sign() != 0) { - continue; - } + // dilepton rap cut float rap = dilepton.rap(); - if (fConfigUseRapcut && std::abs(rap) > fConfigDileptonRapCutAbs) + if (fConfigUseRapcut && abs(rap) > fConfigDileptonRapCutAbs) continue; // compute dilepton - track quantities VarManager::FillDileptonHadron(dilepton, track, VarManager::fgValues); // for the energy correlator analysis - float Effweight_rec = 1.0f; - if (fConfigApplyEfficiency) { - float dilepton_eta = dilepton.eta(); - float dilepton_phi = dilepton.phi(); - float hadron_eta = track.eta(); - float hadron_phi = track.phi(); - float deltaphi = RecoDecay::constrainAngle(dilepton_phi - hadron_phi, -0.5 * o2::constants::math::PI); - Effweight_rec = hAcceptance_rec->Interpolate(dilepton_eta - hadron_eta, deltaphi); - float Effdilepton = hEfficiency_dilepton->Interpolate(dilepton.pt()); - float Masswindow = hMasswindow->Interpolate(dilepton.pt()); - float Effhadron = hEfficiency_hadron->Interpolate(track.pt()); - if (fConfigApplyEfficiencyME) { - Effweight_rec = Effdilepton * Effhadron * Masswindow; // for the moment, apply the efficiency correction also for the mixed event pairs, but this can be changed in case we want to apply it only for the same event pairs - } else { - Effweight_rec = Effweight_rec * Effdilepton * Effhadron * Masswindow; // apply acceptance and efficiency correction for the real pairs - } - } std::vector fTransRange = fConfigTransRange; - VarManager::FillEnergyCorrelatorTriple(lepton1, lepton2, track, fValuesHadron, fTransRange[0], fTransRange[1], fConfigApplyMassEC, fMassBkg->GetRandom(), 1. / Effweight_rec); + VarManager::FillEnergyCorrelator(dilepton, track, VarManager::fgValues, fTransRange[0], fTransRange[1], fConfigApplyMassEC, fMassBkg->GetRandom()); // loop over dilepton leg cuts and track cuts and fill histograms separately for each combination for (int icut = 0; icut < fNCuts; icut++) { @@ -4050,7 +4053,7 @@ struct AnalysisDileptonTrack { } // end event loop } - void processDummy(MyEventsBasic&) + void processDummy(MyEvents&) { // do nothing } @@ -4305,7 +4308,7 @@ struct AnalysisDileptonTrackTrack { } } - void processDummy(MyEventsBasic&) + void processDummy(MyEvents&) { // do nothing } @@ -4374,10 +4377,6 @@ void DefineHistograms(HistogramManager* histMan, TString histClasses, const char dqhistograms::DefineHistograms(histMan, objArray->At(iclass)->GetName(), "pair", histName); } - if (classStr.Contains("Pairing")) { - dqhistograms::DefineHistograms(histMan, objArray->At(iclass)->GetName(), "event", histName); - } - if (classStr.Contains("Triplets")) { dqhistograms::DefineHistograms(histMan, objArray->At(iclass)->GetName(), "pair", histName); } From 01878861dc256cb9399f4db750f2c16ea8e91576 Mon Sep 17 00:00:00 2001 From: Shunsuke-Kurita Date: Mon, 13 Apr 2026 20:49:22 +0900 Subject: [PATCH 2/4] Add e-mu mixing & solve conflicts --- PWGDQ/Tasks/tableReader_withAssoc.cxx | 394 +++++++++++++++++--------- 1 file changed, 258 insertions(+), 136 deletions(-) diff --git a/PWGDQ/Tasks/tableReader_withAssoc.cxx b/PWGDQ/Tasks/tableReader_withAssoc.cxx index e899a5844d1..90c612bd44b 100644 --- a/PWGDQ/Tasks/tableReader_withAssoc.cxx +++ b/PWGDQ/Tasks/tableReader_withAssoc.cxx @@ -26,42 +26,48 @@ #include "Common/CCDB/EventSelectionParams.h" #include "Common/Core/TableHelper.h" -#include "CCDB/BasicCCDBManager.h" -#include "DataFormatsParameters/GRPMagField.h" -#include "DataFormatsParameters/GRPObject.h" -#include "DetectorsBase/GeometryManager.h" -#include "DetectorsBase/Propagator.h" -#include "Field/MagneticField.h" -#include "Framework/ASoAHelpers.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/AnalysisHelpers.h" -#include "Framework/AnalysisTask.h" -#include "Framework/Configurable.h" -#include "Framework/OutputObjHeader.h" -#include "Framework/runDataProcessing.h" -#include "ITSMFTBase/DPLAlpideParam.h" - -#include "TGeoGlobalMagField.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + #include -#include -#include #include #include +#include +#include #include -#include #include +#include + #include +#include +#include #include #include +#include #include #include +#include #include #include -#include -#include #include #include +#include #include using std::cout; @@ -72,6 +78,7 @@ using namespace o2; using namespace o2::framework; using namespace o2::framework::expressions; using namespace o2::aod; +using namespace o2::common::core; // Some definitions namespace o2::aod @@ -195,6 +202,7 @@ DECLARE_SOA_TABLE(JPsieeCandidates, "AOD", "DQPSEUDOPROPER", dqanalysisflags::Ma // Declarations of various short names using MyEvents = soa::Join; +using MyEventsBasic = soa::Join; using MyEventsMultExtra = soa::Join; using MyEventsMultExtraQVector = soa::Join; using MyEventsZdc = soa::Join; @@ -276,6 +284,7 @@ struct AnalysisEventSelection { // TODO: Provide the mixing variables and binning directly via configurables (e.g. vectors of float) Configurable fConfigMixingVariables{"cfgMixingVars", "", "Mixing configs separated by a comma, default no mixing"}; + Configurable fConfigMixingVariablesJson{"cfgMixingVarsJSON", "", "Mixing configs in JSON format"}; Configurable fConfigEventCuts{"cfgEventCuts", "eventStandard", "Event selection"}; Configurable fConfigEventCutsJSON{"cfgEventCutsJSON", "", "Additional event cuts specified in JSON format"}; Configurable fConfigAddEventHistogram{"cfgAddEventHistogram", "", "Comma separated list of histograms"}; @@ -307,7 +316,7 @@ struct AnalysisEventSelection { void init(o2::framework::InitContext& context) { - bool isAnyProcessEnabled = context.mOptions.get("processSkimmed") || context.mOptions.get("processSkimmedWithZdc") || context.mOptions.get("processSkimmedWithMultExtra") || context.mOptions.get("processSkimmedWithMultExtraZdc") || context.mOptions.get("processSkimmedWithMultExtraZdcFit") || context.mOptions.get("processSkimmedWithQvectorCentr"); + bool isAnyProcessEnabled = context.mOptions.get("processSkimmed") || context.mOptions.get("processSkimmedBasic") || context.mOptions.get("processSkimmedWithZdc") || context.mOptions.get("processSkimmedWithMultExtra") || context.mOptions.get("processSkimmedWithMultExtraZdc") || context.mOptions.get("processSkimmedWithMultExtraZdcFit") || context.mOptions.get("processSkimmedWithQvectorCentr"); bool isDummyEnabled = context.mOptions.get("processDummy"); if (isDummyEnabled) { @@ -356,12 +365,18 @@ struct AnalysisEventSelection { } TString mixVarsString = fConfigMixingVariables.value; + TString mixVarsJsonString = fConfigMixingVariablesJson.value; std::unique_ptr objArray(mixVarsString.Tokenize(",")); - if (objArray->GetEntries() > 0) { + if (objArray->GetEntries() > 0 || mixVarsJsonString != "") { fMixHandler = new MixingHandler("mixingHandler", "mixing handler"); fMixHandler->Init(); - for (int iVar = 0; iVar < objArray->GetEntries(); ++iVar) { - dqmixing::SetUpMixing(fMixHandler, objArray->At(iVar)->GetName()); + if (objArray->GetEntries() > 0) { + for (int iVar = 0; iVar < objArray->GetEntries(); ++iVar) { + dqmixing::SetUpMixing(fMixHandler, objArray->At(iVar)->GetName()); + } + } + if (mixVarsJsonString != "") { + dqmixing::SetUpMixingFromJSON(fMixHandler, mixVarsJsonString.Data()); } } @@ -528,6 +543,11 @@ struct AnalysisEventSelection { runEventSelection(events); publishSelections(events); } + void processSkimmedBasic(MyEventsBasic const& events) + { + runEventSelection(events); + publishSelections(events); + } void processSkimmedWithZdc(MyEventsZdc const& events) { runEventSelection(events); @@ -553,12 +573,13 @@ struct AnalysisEventSelection { runEventSelection(events); publishSelections(events); } - void processDummy(MyEvents&) + void processDummy(MyEventsBasic&) { // do nothing } PROCESS_SWITCH(AnalysisEventSelection, processSkimmed, "Run event selection on DQ skimmed events", false); + PROCESS_SWITCH(AnalysisEventSelection, processSkimmedBasic, "Run event selection on DQ skimmed events with basic tables", false); PROCESS_SWITCH(AnalysisEventSelection, processSkimmedWithZdc, "Run event selection on DQ skimmed events, with ZDC", false); PROCESS_SWITCH(AnalysisEventSelection, processSkimmedWithMultExtra, "Run event selection on DQ skimmed events, with mult extra", false); PROCESS_SWITCH(AnalysisEventSelection, processSkimmedWithMultExtraZdc, "Run event selection on DQ skimmed events, with mult extra and ZDC", false); @@ -814,7 +835,7 @@ struct AnalysisTrackSelection { { runTrackSelection(assocs, events, tracks); } - void processDummy(MyEvents&) + void processDummy(MyEventsBasic&) { // do nothing } @@ -1025,7 +1046,7 @@ struct AnalysisMuonSelection { { runMuonSelection(assocs, events, muons); } - void processDummy(MyEvents&) + void processDummy(MyEventsBasic&) { // do nothing } @@ -1174,7 +1195,7 @@ struct AnalysisPrefilterSelection { } // end loop over combinations } - void processBarrelSkimmed(MyEvents const& events, soa::Join const& assocs, MyBarrelTracks const& tracks) + void processBarrelSkimmed(MyEventsBasic const& events, soa::Join const& assocs, MyBarrelTracks const& tracks) { fPrefilterMap.clear(); @@ -1207,7 +1228,7 @@ struct AnalysisPrefilterSelection { } } - void processDummy(MyEvents&) + void processDummy(MyEventsBasic&) { // do nothing } @@ -1312,6 +1333,7 @@ struct AnalysisSameEventPairing { int fNCutsBarrel; int fNCutsMuon; int fNPairCuts; + int fNPairPerEvent; bool fEnableBarrelMixingHistos; bool fEnableBarrelHistos; @@ -1595,13 +1617,14 @@ struct AnalysisSameEventPairing { Form("PairsEleMuSEMM_%s_%s", trackCutName.Data(), muonCutName.Data())}; histNames += Form("%s;%s;%s;", names[0].Data(), names[1].Data(), names[2].Data()); int index = iTrack * fNCutsMuon + iMuon; + fTrackMuonHistNames[index] = names; + if (fEnableBarrelMuonMixingHistos) { names.push_back(Form("PairsEleMuMEPM_%s_%s", trackCutName.Data(), muonCutName.Data())); names.push_back(Form("PairsEleMuMEPP_%s_%s", trackCutName.Data(), muonCutName.Data())); names.push_back(Form("PairsEleMuMEMM_%s_%s", trackCutName.Data(), muonCutName.Data())); histNames += Form("%s;%s;%s;", names[3].Data(), names[4].Data(), names[5].Data()); } - fTrackMuonHistNames[index] = names; TString cutNamesStr = fConfigCuts.pair.value; if (!cutNamesStr.IsNull()) { @@ -1633,8 +1656,14 @@ struct AnalysisSameEventPairing { fHistMan->SetDefaultVarNames(VarManager::fgVariableNames, VarManager::fgVariableUnits); VarManager::SetCollisionSystem((TString)fConfigOptions.collisionSystem, fConfigOptions.centerMassEnergy); // set collision system and center of mass energy DefineHistograms(fHistMan, histNames.Data(), fConfigAddSEPHistogram.value.data()); // define all histograms - dqhistograms::AddHistogramsFromJSON(fHistMan, fConfigAddJSONHistograms.value.c_str()); // ad-hoc histograms via JSON - VarManager::SetUseVars(fHistMan->GetUsedVars()); // provide the list of required variables so that VarManager knows what to fill + if (fEnableBarrelHistos) { + DefineHistograms(fHistMan, "PairingSEQA", "sameevent-pairing"); // histograms for QA of the pairing + }; + if (fEnableBarrelMixingHistos) { + DefineHistograms(fHistMan, "PairingMEQA", "mixedevent-pairing"); // histograms for QA of the pairing + }; + dqhistograms::AddHistogramsFromJSON(fHistMan, fConfigAddJSONHistograms.value.c_str()); // ad-hoc histograms via JSON + VarManager::SetUseVars(fHistMan->GetUsedVars()); // provide the list of required variables so that VarManager knows what to fill fOutputList.setObject(fHistMan->GetMainHistogramList()); } } @@ -1741,6 +1770,7 @@ struct AnalysisSameEventPairing { constexpr bool eventHasQvectorCentr = ((TEventFillMap & VarManager::ObjTypes::CollisionQvect) > 0); constexpr bool trackHasCov = ((TTrackFillMap & VarManager::ObjTypes::TrackCov) > 0 || (TTrackFillMap & VarManager::ObjTypes::ReducedTrackBarrelCov) > 0); bool isSelectedBDT = false; + fNPairPerEvent = 0; for (auto& event : events) { if (!event.isEventSelected_bit(0)) { @@ -1788,6 +1818,7 @@ struct AnalysisSameEventPairing { twoTrackFilter |= (static_cast(1) << 31); } + fNPairPerEvent++; VarManager::FillPair(t1, t2); // compute quantities which depend on the associated collision, such as DCA if (fConfigOptions.propTrack) { @@ -2093,6 +2124,10 @@ struct AnalysisSameEventPairing { } } // end loop (cuts) } // end loop over pairs of track associations + VarManager::fgValues[VarManager::kNPairsPerEvent] = fNPairPerEvent; + if (fEnableBarrelHistos && fConfigQA) { + fHistMan->FillHistClass("PairingSEQA", VarManager::fgValues); + } } // end loop over events } @@ -2112,6 +2147,7 @@ struct AnalysisSameEventPairing { } auto t1 = a1.template reducedtrack_as(); auto t2 = a2.template reducedtrack_as(); + fNPairPerEvent++; VarManager::FillPairME(t1, t2); if constexpr ((TEventFillMap & VarManager::ObjTypes::ReducedEventQvector) > 0) { VarManager::FillPairVn(t1, t2); @@ -2277,102 +2313,12 @@ struct AnalysisSameEventPairing { auto assocs2 = assocs.sliceBy(preSlice, event2.globalIndex()); assocs2.bindExternalIndices(&events); + fNPairPerEvent = 0; runMixedPairing(assocs1, assocs2, tracks, tracks); - } // end event loop - } - - template - void runEmuMixedPairing(TAssoc1 const& assocs1, TAssoc2 const& assocs2, TTracks1 const& /*tracks1*/, TTracks2 const& /*tracks2*/) - { - const auto& histNames = fTrackMuonHistNames; - int sign1 = 0; - int sign2 = 0; - int nPairCuts = (fPairCuts.size() > 0) ? fPairCuts.size() : 1; - constexpr bool eventHasQvector = ((TEventFillMap & VarManager::ObjTypes::ReducedEventQvector) > 0); - constexpr bool eventHasQvectorCentr = ((TEventFillMap & VarManager::ObjTypes::CollisionQvect) > 0); - - for (auto& a1 : assocs1) { - if (!(a1.isBarrelSelected_raw() & fTrackFilterMask)) { - continue; + VarManager::fgValues[VarManager::kNPairsPerEvent] = fNPairPerEvent; + if (fEnableBarrelMixingHistos && fConfigQA) { + fHistMan->FillHistClass("PairingMEQA", VarManager::fgValues); } - for (auto& a2 : assocs2) { - if (!(a2.isMuonSelected_raw() & fMuonFilterMask)) { - continue; - } - - auto t1 = a1.template reducedtrack_as(); - auto t2 = a2.template reducedmuon_as(); - sign1 = t1.sign(); - sign2 = t2.sign(); - - VarManager::FillPairME(t1, t2); - if constexpr (eventHasQvector) { - VarManager::FillPairVn(t1, t2); - } - if constexpr (eventHasQvectorCentr) { - VarManager::FillPairVn(t1, t2); - } - - for (int iTrack = 0; iTrack < fNCutsBarrel; ++iTrack) { - if (!(a1.isBarrelSelected_raw() & (1u << iTrack))) { - continue; - } - for (int iMuon = 0; iMuon < fNCutsMuon; ++iMuon) { - if (!(a2.isMuonSelected_raw() & (1u << iMuon))) { - continue; - } - for (int iPairCut = 0; iPairCut < nPairCuts; ++iPairCut) { - if (!fPairCuts.empty()) { - AnalysisCompositeCut cut = fPairCuts.at(iPairCut); - if (!cut.IsSelected(VarManager::fgValues)) { - continue; - } - } - int index = iTrack * (fNCutsMuon * nPairCuts) + iMuon * nPairCuts + iPairCut; - auto itHist = histNames.find(index); - if (itHist == histNames.end() || itHist->second.size() < 6) { - continue; - } - if (sign1 * sign2 < 0) { - fHistMan->FillHistClass(itHist->second[3].Data(), VarManager::fgValues); - } else { - if (sign1 > 0) { - fHistMan->FillHistClass(itHist->second[4].Data(), VarManager::fgValues); - } else { - fHistMan->FillHistClass(itHist->second[5].Data(), VarManager::fgValues); - } - } - } // end pair cut loop - } // end muon cut loop - } // end barrel cut loop - } - } - } - - template - void runEmuSameSideMixing(TEvents& events, Preslice& preslice1, TTrackAssocs const& assocs1, TTracks const& tracks1, - Preslice& preslice2, TMuonAssocs const& assocs2, TMuons const& tracks2) - { - events.bindExternalIndices(&assocs1); - events.bindExternalIndices(&assocs2); - int mixingDepth = fConfigMixingDepth.value; - for (auto& [event1, event2] : selfCombinations(hashBin, mixingDepth, -1, events, events)) { - VarManager::ResetValues(0, VarManager::kNVars); - VarManager::FillEvent(event1, VarManager::fgValues); - - auto groupedAssocs1 = assocs1.sliceBy(preslice1, event1.globalIndex()); - groupedAssocs1.bindExternalIndices(&events); - if (groupedAssocs1.size() == 0) { - continue; - } - - auto groupedAssocs2 = assocs2.sliceBy(preslice2, event2.globalIndex()); - groupedAssocs2.bindExternalIndices(&events); - if (groupedAssocs2.size() == 0) { - continue; - } - - runEmuMixedPairing(groupedAssocs1, groupedAssocs2, tracks1, tracks2); } // end event loop } @@ -2493,6 +2439,101 @@ struct AnalysisSameEventPairing { } // end event loop } + template + void runEmuMixedPairing(TAssoc1 const& assocs1, TAssoc2 const& assocs2, TTracks1 const& /*tracks1*/, TTracks2 const& /*tracks2*/) + { + const auto& histNames = fTrackMuonHistNames; + int sign1 = 0; + int sign2 = 0; + int nPairCuts = (fPairCuts.size() > 0) ? fPairCuts.size() : 1; + constexpr bool eventHasQvector = ((TEventFillMap & VarManager::ObjTypes::ReducedEventQvector) > 0); + constexpr bool eventHasQvectorCentr = ((TEventFillMap & VarManager::ObjTypes::CollisionQvect) > 0); + + for (auto& a1 : assocs1) { + if (!(a1.isBarrelSelected_raw() & fTrackFilterMask)) { + continue; + } + for (auto& a2 : assocs2) { + if (!(a2.isMuonSelected_raw() & fMuonFilterMask)) { + continue; + } + + auto t1 = a1.template reducedtrack_as(); + auto t2 = a2.template reducedmuon_as(); + sign1 = t1.sign(); + sign2 = t2.sign(); + + VarManager::FillPairME(t1, t2); + if constexpr (eventHasQvector) { + VarManager::FillPairVn(t1, t2); + } + if constexpr (eventHasQvectorCentr) { + VarManager::FillPairVn(t1, t2); + } + + for (int iTrack = 0; iTrack < fNCutsBarrel; ++iTrack) { + if (!(a1.isBarrelSelected_raw() & (1u << iTrack))) { + continue; + } + for (int iMuon = 0; iMuon < fNCutsMuon; ++iMuon) { + if (!(a2.isMuonSelected_raw() & (1u << iMuon))) { + continue; + } + for (int iPairCut = 0; iPairCut < nPairCuts; ++iPairCut) { + if (!fPairCuts.empty()) { + AnalysisCompositeCut cut = fPairCuts.at(iPairCut); + if (!cut.IsSelected(VarManager::fgValues)) { + continue; + } + } + int index = iTrack * (fNCutsMuon * nPairCuts) + iMuon * nPairCuts + iPairCut; + auto itHist = histNames.find(index); + if (itHist == histNames.end() || itHist->second.size() < 6) { + continue; + } + if (sign1 * sign2 < 0) { + fHistMan->FillHistClass(itHist->second[3].Data(), VarManager::fgValues); + } else { + if (sign1 > 0) { + fHistMan->FillHistClass(itHist->second[4].Data(), VarManager::fgValues); + } else { + fHistMan->FillHistClass(itHist->second[5].Data(), VarManager::fgValues); + } + } + } // end pair cut loop + } // end muon cut loop + } // end barrel cut loop + } + } + } + + template + void runEmuSameSideMixing(TEvents& events, Preslice& preslice1, TTrackAssocs const& assocs1, TTracks const& tracks1, + Preslice& preslice2, TMuonAssocs const& assocs2, TMuons const& tracks2) + { + events.bindExternalIndices(&assocs1); + events.bindExternalIndices(&assocs2); + int mixingDepth = fConfigMixingDepth.value; + for (auto& [event1, event2] : selfCombinations(hashBin, mixingDepth, -1, events, events)) { + VarManager::ResetValues(0, VarManager::kNVars); + VarManager::FillEvent(event1, VarManager::fgValues); + + auto groupedAssocs1 = assocs1.sliceBy(preslice1, event1.globalIndex()); + groupedAssocs1.bindExternalIndices(&events); + if (groupedAssocs1.size() == 0) { + continue; + } + + auto groupedAssocs2 = assocs2.sliceBy(preslice2, event2.globalIndex()); + groupedAssocs2.bindExternalIndices(&events); + if (groupedAssocs2.size() == 0) { + continue; + } + + runEmuMixedPairing(groupedAssocs1, groupedAssocs2, tracks1, tracks2); + } // end event loop + } + void processAllSkimmed(MyEventsVtxCovSelected const& events, soa::Join const& barrelAssocs, MyBarrelTracksWithCovWithAmbiguities const& barrelTracks, soa::Join const& muonAssocs, MyMuonTracksWithCovWithAmbiguities const& muons) @@ -2617,8 +2658,12 @@ struct AnalysisSameEventPairing { runEmuSameSideMixing(events, trackEmuAssocsPerCollision, barrelAssocs, barrelTracks, muonAssocsPerCollision, muonAssocs, muons); } +<<<<<<< HEAD void processDummy(MyEvents&) >>>>>>> d32d3fbdc (Add electron-muon mixing) +======= + void processDummy(MyEventsBasic&) +>>>>>>> f947468c8 (Add e-mu mixing & solve conflicts) { // do nothing } @@ -3369,7 +3414,7 @@ struct AnalysisAsymmetricPairing { runThreeProng(events, trackAssocsPerCollision, barrelAssocs, barrelTracks, VarManager::kTripleCandidateToPKPi); } - void processDummy(MyEvents&) + void processDummy(MyEventsBasic&) { // do nothing } @@ -3416,6 +3461,10 @@ struct AnalysisDileptonTrack { Configurable> fConfigFitmassEC{"cfgTFitmassEC", std::vector{-0.541438, 2.8, 3.2}, "parameter from the fit fuction and fit range"}; Configurable> fConfigTransRange{"cfgTransRange", std::vector{0.333333, 0.666667}, "Transverse region for the energy correlstor analysis"}; + Configurable fConfigApplyEfficiency{"cfgApplyEfficiency", false, "If true, apply efficiency correction for the energy correlator study"}; + Configurable fConfigApplyEfficiencyME{"cfgApplyEfficiencyME", false, "If true, apply efficiency correction for the energy correlator study"}; + Configurable fConfigAccCCDBPath{"AccCCDBPath", "Users/y/yalin/pptest/test2", "Path of the efficiency corrections"}; + int fCurrentRun; // needed to detect if the run changed and trigger update of calibrations etc. int fNCuts; // number of dilepton leg cuts int fNLegCuts; @@ -3449,6 +3498,12 @@ struct AnalysisDileptonTrack { TF1* fMassBkg = nullptr; + TH2F* hAcceptance_rec; + TH2F* hAcceptance_gen; + TH1F* hEfficiency_dilepton; + TH1F* hEfficiency_hadron; + TH1F* hMasswindow; + void init(o2::framework::InitContext& context) { bool isBarrel = context.mOptions.get("processBarrelSkimmed"); @@ -3706,6 +3761,22 @@ struct AnalysisDileptonTrack { } } + void initAccFromCCDB(uint64_t timestamp) + { + TList* listAccs = fCCDB->getForTimeStamp(fConfigAccCCDBPath, timestamp); + if (!listAccs) { + LOG(fatal) << "Problem getting TList object with efficiencies!"; + } + hEfficiency_dilepton = static_cast(listAccs->FindObject("hEfficiency_dilepton")); + hEfficiency_hadron = static_cast(listAccs->FindObject("hEfficiency_hadron")); + hAcceptance_rec = static_cast(listAccs->FindObject("hAcceptance_rec")); + hAcceptance_gen = static_cast(listAccs->FindObject("hAcceptance_gen")); + hMasswindow = static_cast(listAccs->FindObject("hMasswindow")); + if (!hAcceptance_rec || !hAcceptance_gen || !hEfficiency_dilepton || !hEfficiency_hadron || !hMasswindow) { + LOG(fatal) << "Problem getting histograms from the TList object with efficiencies!"; + } + } + // Template function to run pair - hadron combinations template void runDileptonHadron(TEvent const& event, TTrackAssocs const& assocs, TTracks const& tracks, TDileptons const& dileptons) @@ -3725,7 +3796,7 @@ struct AnalysisDileptonTrack { } // dilepton rap cut float rap = dilepton.rap(); - if (fConfigUseRapcut && abs(rap) > fConfigDileptonRapCutAbs) + if (fConfigUseRapcut && std::abs(rap) > fConfigDileptonRapCutAbs) continue; VarManager::FillTrack(dilepton, fValuesDilepton); @@ -3782,8 +3853,21 @@ struct AnalysisDileptonTrack { VarManager::FillDileptonTrackVertexing(event, lepton1, lepton2, track, fValuesHadron); // for the energy correlator analysis + float Effweight_rec = 1.0f; + if (fConfigApplyEfficiency) { + float dilepton_eta = dilepton.eta(); + float dilepton_phi = dilepton.phi(); + float hadron_eta = track.eta(); + float hadron_phi = track.phi(); + float deltaphi = RecoDecay::constrainAngle(dilepton_phi - hadron_phi, -0.5 * o2::constants::math::PI); + Effweight_rec = hAcceptance_rec->Interpolate(dilepton_eta - hadron_eta, deltaphi); + float Effdilepton = hEfficiency_dilepton->Interpolate(dilepton.pt()); + float Masswindow = hMasswindow->Interpolate(dilepton.pt()); + float Effhadron = hEfficiency_hadron->Interpolate(track.pt()); + Effweight_rec = Effweight_rec * Effdilepton * Effhadron * Masswindow; + } std::vector fTransRange = fConfigTransRange; - VarManager::FillEnergyCorrelator(dilepton, track, fValuesHadron, fTransRange[0], fTransRange[1], fConfigApplyMassEC, fMassBkg->GetRandom()); + VarManager::FillEnergyCorrelatorTriple(lepton1, lepton2, track, fValuesHadron, fTransRange[0], fTransRange[1], fConfigApplyMassEC, fMassBkg->GetRandom(), 1. / Effweight_rec); // table to be written out for ML analysis BmesonsTable(event.runNumber(), event.globalIndex(), event.timestamp(), fValuesHadron[VarManager::kPairMass], dilepton.mass(), fValuesHadron[VarManager::kDeltaMass], fValuesHadron[VarManager::kPairPt], fValuesHadron[VarManager::kPairEta], fValuesHadron[VarManager::kPairPhi], fValuesHadron[VarManager::kPairRap], @@ -3888,6 +3972,9 @@ struct AnalysisDileptonTrack { } if (fCurrentRun != events.begin().runNumber()) { // start: runNumber initParamsFromCCDB(events.begin().timestamp()); + if (fConfigApplyEfficiency) { + initAccFromCCDB(events.begin().timestamp()); + } fCurrentRun = events.begin().runNumber(); } // end: runNumber for (auto& event : events) { @@ -3940,11 +4027,19 @@ struct AnalysisDileptonTrack { void processBarrelMixedEvent(soa::Filtered& events, soa::Filtered> const& assocs, - MyBarrelTracksWithCov const&, soa::Filtered const& dileptons) + MyBarrelTracksWithCov const& tracks, soa::Filtered const& dileptons) { if (events.size() == 0) { return; } + + if (fCurrentRun != events.begin().runNumber()) { // start: runNumber + if (fConfigApplyEfficiency) { + initAccFromCCDB(events.begin().timestamp()); + } + fCurrentRun = events.begin().runNumber(); + } // end: runNumber + events.bindExternalIndices(&dileptons); events.bindExternalIndices(&assocs); @@ -3976,18 +4071,41 @@ struct AnalysisDileptonTrack { // loop over dileptons for (auto dilepton : evDileptons) { - + // get full track info of tracks based on the index + auto lepton1 = tracks.rawIteratorAt(dilepton.index0Id()); + auto lepton2 = tracks.rawIteratorAt(dilepton.index1Id()); + // Check that the dilepton has zero charge + if (dilepton.sign() != 0) { + continue; + } // dilepton rap cut float rap = dilepton.rap(); - if (fConfigUseRapcut && abs(rap) > fConfigDileptonRapCutAbs) + if (fConfigUseRapcut && std::abs(rap) > fConfigDileptonRapCutAbs) continue; // compute dilepton - track quantities VarManager::FillDileptonHadron(dilepton, track, VarManager::fgValues); // for the energy correlator analysis + float Effweight_rec = 1.0f; + if (fConfigApplyEfficiency) { + float dilepton_eta = dilepton.eta(); + float dilepton_phi = dilepton.phi(); + float hadron_eta = track.eta(); + float hadron_phi = track.phi(); + float deltaphi = RecoDecay::constrainAngle(dilepton_phi - hadron_phi, -0.5 * o2::constants::math::PI); + Effweight_rec = hAcceptance_rec->Interpolate(dilepton_eta - hadron_eta, deltaphi); + float Effdilepton = hEfficiency_dilepton->Interpolate(dilepton.pt()); + float Masswindow = hMasswindow->Interpolate(dilepton.pt()); + float Effhadron = hEfficiency_hadron->Interpolate(track.pt()); + if (fConfigApplyEfficiencyME) { + Effweight_rec = Effdilepton * Effhadron * Masswindow; // for the moment, apply the efficiency correction also for the mixed event pairs, but this can be changed in case we want to apply it only for the same event pairs + } else { + Effweight_rec = Effweight_rec * Effdilepton * Effhadron * Masswindow; // apply acceptance and efficiency correction for the real pairs + } + } std::vector fTransRange = fConfigTransRange; - VarManager::FillEnergyCorrelator(dilepton, track, VarManager::fgValues, fTransRange[0], fTransRange[1], fConfigApplyMassEC, fMassBkg->GetRandom()); + VarManager::FillEnergyCorrelatorTriple(lepton1, lepton2, track, fValuesHadron, fTransRange[0], fTransRange[1], fConfigApplyMassEC, fMassBkg->GetRandom(), 1. / Effweight_rec); // loop over dilepton leg cuts and track cuts and fill histograms separately for each combination for (int icut = 0; icut < fNCuts; icut++) { @@ -4053,7 +4171,7 @@ struct AnalysisDileptonTrack { } // end event loop } - void processDummy(MyEvents&) + void processDummy(MyEventsBasic&) { // do nothing } @@ -4308,7 +4426,7 @@ struct AnalysisDileptonTrackTrack { } } - void processDummy(MyEvents&) + void processDummy(MyEventsBasic&) { // do nothing } @@ -4377,6 +4495,10 @@ void DefineHistograms(HistogramManager* histMan, TString histClasses, const char dqhistograms::DefineHistograms(histMan, objArray->At(iclass)->GetName(), "pair", histName); } + if (classStr.Contains("Pairing")) { + dqhistograms::DefineHistograms(histMan, objArray->At(iclass)->GetName(), "event", histName); + } + if (classStr.Contains("Triplets")) { dqhistograms::DefineHistograms(histMan, objArray->At(iclass)->GetName(), "pair", histName); } From d0553196de6e2c5d1ee8efbc525a43da72af4c8f Mon Sep 17 00:00:00 2001 From: Shunsuke-Kurita Date: Tue, 14 Apr 2026 16:20:13 +0900 Subject: [PATCH 3/4] Solve Formatting and MegaLinter errors --- PWGDQ/Tasks/tableReader_withAssoc.cxx | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/PWGDQ/Tasks/tableReader_withAssoc.cxx b/PWGDQ/Tasks/tableReader_withAssoc.cxx index 90c612bd44b..1c56dd850e3 100644 --- a/PWGDQ/Tasks/tableReader_withAssoc.cxx +++ b/PWGDQ/Tasks/tableReader_withAssoc.cxx @@ -1658,10 +1658,10 @@ struct AnalysisSameEventPairing { DefineHistograms(fHistMan, histNames.Data(), fConfigAddSEPHistogram.value.data()); // define all histograms if (fEnableBarrelHistos) { DefineHistograms(fHistMan, "PairingSEQA", "sameevent-pairing"); // histograms for QA of the pairing - }; + } if (fEnableBarrelMixingHistos) { DefineHistograms(fHistMan, "PairingMEQA", "mixedevent-pairing"); // histograms for QA of the pairing - }; + } dqhistograms::AddHistogramsFromJSON(fHistMan, fConfigAddJSONHistograms.value.c_str()); // ad-hoc histograms via JSON VarManager::SetUseVars(fHistMan->GetUsedVars()); // provide the list of required variables so that VarManager knows what to fill fOutputList.setObject(fHistMan->GetMainHistogramList()); @@ -2501,8 +2501,8 @@ struct AnalysisSameEventPairing { } } } // end pair cut loop - } // end muon cut loop - } // end barrel cut loop + } // end muon cut loop + } // end barrel cut loop } } } From c96fdbd5f65f6b4cc1a7e4a8225828414351fe4e Mon Sep 17 00:00:00 2001 From: Shunsuke-Kurita Date: Tue, 14 Apr 2026 19:09:52 +0900 Subject: [PATCH 4/4] Solve merged conflicts --- PWGDQ/Tasks/tableReader_withAssoc.cxx | 11 ----------- 1 file changed, 11 deletions(-) diff --git a/PWGDQ/Tasks/tableReader_withAssoc.cxx b/PWGDQ/Tasks/tableReader_withAssoc.cxx index 1c56dd850e3..8a55b51e97e 100644 --- a/PWGDQ/Tasks/tableReader_withAssoc.cxx +++ b/PWGDQ/Tasks/tableReader_withAssoc.cxx @@ -2642,15 +2642,12 @@ struct AnalysisSameEventPairing { runSameSideMixing(events, muonAssocs, muons, muonAssocsPerCollision); } -<<<<<<< HEAD void processMixingMuonSkimmedFlow(soa::Filtered& events, soa::Join const& muonAssocs, MyMuonTracksWithCovWithAmbiguities const& muons) { runSameSideMixing(events, muonAssocs, muons, muonAssocsPerCollision); } - void processDummy(MyEventsBasic&) -======= void processMixingElectronMuonSkimmed(soa::Filtered& events, soa::Join const& barrelAssocs, aod::ReducedTracks const& barrelTracks, soa::Join const& muonAssocs, MyMuonTracksWithCovWithAmbiguities const& muons) @@ -2658,12 +2655,7 @@ struct AnalysisSameEventPairing { runEmuSameSideMixing(events, trackEmuAssocsPerCollision, barrelAssocs, barrelTracks, muonAssocsPerCollision, muonAssocs, muons); } -<<<<<<< HEAD - void processDummy(MyEvents&) ->>>>>>> d32d3fbdc (Add electron-muon mixing) -======= void processDummy(MyEventsBasic&) ->>>>>>> f947468c8 (Add e-mu mixing & solve conflicts) { // do nothing } @@ -2684,11 +2676,8 @@ struct AnalysisSameEventPairing { PROCESS_SWITCH(AnalysisSameEventPairing, processMixingBarrelSkimmedFlow, "Run barrel type mixing pairing, with flow, with skimmed tracks", false); PROCESS_SWITCH(AnalysisSameEventPairing, processMixingBarrelWithQvectorCentrSkimmedNoCov, "Run barrel type mixing pairing, with skimmed tracks and with Qvector from central framework", false); PROCESS_SWITCH(AnalysisSameEventPairing, processMixingMuonSkimmed, "Run muon type mixing pairing, with skimmed muons", false); -<<<<<<< HEAD PROCESS_SWITCH(AnalysisSameEventPairing, processMixingMuonSkimmedFlow, "Run muon type mixing pairing, with skimmed muons and flow", false); -======= PROCESS_SWITCH(AnalysisSameEventPairing, processMixingElectronMuonSkimmed, "Run electron-muon mixing pairing, with skimmed tracks/muons", false); ->>>>>>> d32d3fbdc (Add electron-muon mixing) PROCESS_SWITCH(AnalysisSameEventPairing, processDummy, "Dummy function, enabled only if none of the others are enabled", true); };