Skip to content

Refactor far-forward neutrals reconstruction #2736

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
baptistefraisse wants to merge 3 commits into
base: main
Choose a base branch
from
Open

Refactor far-forward neutrals reconstruction #2736

Show file tree
Hide file tree
Changes from 2 commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
1127c0c
Refactor far-forward neutral reconstruction for multiple calorimeters
baptistefraisse Jul 9, 2026
33ae708
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Jul 9, 2026
fb18593
Merge branch 'main' into reco/ff-lambda-multicalo
veprbl Jul 23, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
97 changes: 32 additions & 65 deletions src/algorithms/reco/FarForwardNeutralsReconstruction.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,7 @@
#include <stdexcept>
#include <tuple>
#include <vector>
#include <string>

#include "FarForwardNeutralsReconstruction.h"

Expand Down Expand Up @@ -75,10 +76,8 @@ int FarForwardNeutralsReconstruction::processNeutralCalo(
edm4eic::ReconstructedParticleCollection* out_neutrals,
const std::vector<double>& gammaScaleCoeff, const std::vector<double>& neutronScaleCoeff,
bool canDetectGammas, bool canDetectNeutrons, const CorrFunc& gammaCorr,
const CorrFunc& neutronCorr, GammaMode gammaMode, double gammaLeaderFracMin, double clusterEmin,
NeutronMode neutronMode, bool associateAllClustersToNeutron) const {

(void)gammaLeaderFracMin;
const CorrFunc& neutronCorr, GammaMode gammaMode, double clusterEmin, NeutronMode neutronMode,
bool associateAllClustersToNeutron) const {

const double m_neutron = m_particleSvc.particle(2112).mass;

Expand Down Expand Up @@ -248,67 +247,35 @@ void FarForwardNeutralsReconstruction::process(
const FarForwardNeutralsReconstruction::Output& output) const {

// Unpacking
const auto [clustersHcal, clustersB0, clustersEcalEndcapP, clustersLFHCAL] = input;
auto [out_neutralsHcal, out_neutralsB0, out_neutralsEcalEndcapP, out_neutralsLFHCAL] = output;

// Global
int n_neutrons = 0;

// ZDC-Hcal
n_neutrons += processNeutralCalo(clustersHcal, out_neutralsHcal,
/*gammaScaleCoeff=*/m_cfg.gammaScaleCorrCoeffHcalZDC,
/*neutronScaleCoeff=*/m_cfg.neutronScaleCorrCoeffHcalZDC,
/*canDetectGammas=*/true,
/*canDetectNeutrons=*/true,
/*gammaCorr=*/corrPower,
/*neutronCorr=*/corrPower,
/*gammaMode=*/GammaMode::AllPassing,
/*gammaLeaderFracMin=*/0.0,
/*clusterEmin=*/m_cfg.clusterEminHcalZDC,
/*neutronMode=*/NeutronMode::SumAll,
/*associateAllClustersToNeutron=*/true);

// B0-Ecal
n_neutrons += processNeutralCalo(clustersB0, out_neutralsB0,
/*gammaScaleCoeff=*/m_cfg.gammaScaleCorrCoeffB0Ecal,
/*neutronScaleCoeff=*/m_cfg.neutronScaleCorrCoeffB0Ecal,
/*canDetectGammas=*/true,
/*canDetectNeutrons=*/false,
/*gammaCorr=*/corrPower,
/*neutronCorr=*/corrPower,
/*gammaMode=*/GammaMode::LeaderOnly,
/*gammaLeaderFracMin=*/0.0,
/*clusterEmin=*/m_cfg.clusterEminB0Ecal,
/*neutronMode=*/NeutronMode::None,
/*associateAllClustersToNeutron=*/false);

// EndcapP-Ecal
n_neutrons += processNeutralCalo(clustersEcalEndcapP, out_neutralsEcalEndcapP,
/*gammaScaleCoeff=*/m_cfg.gammaScaleCorrCoeffEcalEndcapP,
/*neutronScaleCoeff=*/m_cfg.neutronScaleCorrCoeffEcalEndcapP,
/*canDetectGammas=*/true,
/*canDetectNeutrons=*/false,
/*gammaCorr=*/corrPower,
/*neutronCorr=*/corrPower,
/*gammaMode=*/GammaMode::LeaderOnly,
/*gammaLeaderFracMin=*/0.0,
/*clusterEmin=*/m_cfg.clusterEminEcalEndcapP,
/*neutronMode=*/NeutronMode::None,
/*associateAllClustersToNeutron=*/false);

// LFHCAL
n_neutrons += processNeutralCalo(clustersLFHCAL, out_neutralsLFHCAL,
/*gammaScaleCoeff=*/m_cfg.gammaScaleCorrCoeffLFHCAL,
/*neutronScaleCoeff=*/m_cfg.neutronScaleCorrCoeffLFHCAL,
/*canDetectGammas=*/false,
/*canDetectNeutrons=*/true,
/*gammaCorr=*/corrPower,
/*neutronCorr=*/corrPower,
/*gammaMode=*/GammaMode::None,
/*gammaLeaderFracMin=*/0.0,
/*clusterEmin=*/m_cfg.clusterEminLFHCAL,
/*neutronMode=*/NeutronMode::LeaderOnly,
/*associateAllClustersToNeutron=*/false);
const auto [clusters] = input;
auto [out_neutrals] = output;

GammaMode gammaMode = GammaMode::None;
if (m_cfg.gammaMode == "LeaderOnly") {
gammaMode = GammaMode::LeaderOnly;
} else if (m_cfg.gammaMode == "AllPassing") {
gammaMode = GammaMode::AllPassing;
}

NeutronMode neutronMode = NeutronMode::None;
if (m_cfg.neutronMode == "SumAll") {
neutronMode = NeutronMode::SumAll;
} else if (m_cfg.neutronMode == "LeaderOnly") {
neutronMode = NeutronMode::LeaderOnly;
}

const int n_neutrons =
processNeutralCalo(clusters, out_neutrals,
/*gammaScaleCoeff=*/m_cfg.gammaScaleCorrCoeff,
/*neutronScaleCoeff=*/m_cfg.neutronScaleCorrCoeff,
/*canDetectGammas=*/m_cfg.canDetectGammas,
/*canDetectNeutrons=*/m_cfg.canDetectNeutrons,
/*gammaCorr=*/corrPower,
/*neutronCorr=*/corrPower,
/*gammaMode=*/gammaMode,
/*clusterEmin=*/m_cfg.clusterEmin,
/*neutronMode=*/neutronMode,
/*associateAllClustersToNeutron=*/m_cfg.associateAllClustersToNeutron);

debug("Found {} neutron candidates", n_neutrons);
}
Expand Down
27 changes: 6 additions & 21 deletions src/algorithms/reco/FarForwardNeutralsReconstruction.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -21,18 +21,8 @@
namespace eicrecon {

using FarForwardNeutralsReconstructionAlgorithm =

algorithms::Algorithm<
algorithms::Input<const edm4eic::ClusterCollection, // clusters ZDC-Hcal
const edm4eic::ClusterCollection, // clusters B0-Ecal
const edm4eic::ClusterCollection, // clusters EndcapP-Ecal
const edm4eic::ClusterCollection>, // clusters LFHCAL

algorithms::Output<
edm4eic::ReconstructedParticleCollection, // neutrons/gamma in ZDC-Hcal
edm4eic::ReconstructedParticleCollection, // neutrons/gamma in B0-Ecal
edm4eic::ReconstructedParticleCollection, // neutrons/gamma in EndcapP-Ecal
edm4eic::ReconstructedParticleCollection>>; // neutrons/gamma in LFHCAL
algorithms::Algorithm<algorithms::Input<const edm4eic::ClusterCollection>,
algorithms::Output<edm4eic::ReconstructedParticleCollection>>;
/**
* Reconstructs far-forward neutral candidates from multiple calorimeter cluster collections.
*
Expand All @@ -59,12 +49,8 @@ class FarForwardNeutralsReconstruction
FarForwardNeutralsReconstruction(std::string_view name)
: FarForwardNeutralsReconstructionAlgorithm{
name,

{"clustersHcal", "clustersB0", "clustersEcalEndCapP", "clustersLFHCAL"},

{"outputNeutralsHcal", "outputNeutralsB0", "outputNeutralsEcalEndCapP",
"outputNeutralsLFHCAL"},

{"clusters"},
{"outputNeutrals"},
"Convert EMCal and HCal clusters into neutron or photon candidates"} {}

void init() final;
Expand All @@ -89,9 +75,8 @@ class FarForwardNeutralsReconstruction
const std::vector<double>& gammaScaleCoeff,
const std::vector<double>& neutronScaleCoeff, bool canDetectGammas,
bool canDetectNeutrons, const CorrFunc& gammaCorr,
const CorrFunc& neutronCorr, GammaMode gammaMode,
double gammaLeaderFracMin, double clusterEmin, NeutronMode neutronMode,
bool associateAllClustersToNeutron) const;
const CorrFunc& neutronCorr, GammaMode gammaMode, double clusterEmin,
NeutronMode neutronMode, bool associateAllClustersToNeutron) const;
};

} // namespace eicrecon
47 changes: 24 additions & 23 deletions src/algorithms/reco/FarForwardNeutralsReconstructionConfig.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3,36 +3,37 @@
#pragma once
#include <float.h>
#include <DD4hep/Detector.h>
#include <string>
#include <vector>

namespace eicrecon {

struct FarForwardNeutralsReconstructionConfig {
/** detector constant describing distance to the ZDC */
/** detector constant describing distance reference position */
std::string offsetPositionName = "HcalFarForwardZDC_SiPMonTile_r_pos";
/** Correction factors for neutrons in the Hcal (ZDC) */
std::vector<double> neutronScaleCorrCoeffHcalZDC = {2.4, 0.89};
/** Correction factors for gammas in the Hcal (ZDC) */
std::vector<double> gammaScaleCorrCoeffHcalZDC = {1.1, 0.98};
/** Correction factors for neutrons in the LFHCAL */
std::vector<double> neutronScaleCorrCoeffLFHCAL = {2.55, 0.95};
/** Correction factors for gammas in the LFHCAL */
std::vector<double> gammaScaleCorrCoeffLFHCAL = {0., 0.};
/** Correction factors for neutrons in the B0-Ecal */
std::vector<double> neutronScaleCorrCoeffB0Ecal = {0., 0.};
/** Correction factors for gammas in the B0-Ecal */
std::vector<double> gammaScaleCorrCoeffB0Ecal = {0.99, 1.14};
/** Correction factors for neutrons in the Endcap-Ecal */
std::vector<double> neutronScaleCorrCoeffEcalEndcapP = {0., 0.};
/** Correction factors for gammas in the Endcap-Ecal */
std::vector<double> gammaScaleCorrCoeffEcalEndcapP = {1.05, 1.01};
/** Cluster thresholds */
double clusterEminHcalZDC = 0.0; // GeV
double clusterEminB0Ecal = 1.0; // GeV
double clusterEminEcalEndcapP = 1.0; // GeV
double clusterEminLFHCAL = 7.0; // GeV

/** Correction factors */
std::vector<double> neutronScaleCorrCoeff = {0.0, 0.0};
std::vector<double> gammaScaleCorrCoeff = {0.0, 0.0};

/** Detector capabilities */
bool canDetectGammas = true;
bool canDetectNeutrons = true;

/** Reconstruction modes */
std::string gammaMode = "None";
std::string neutronMode = "None";

/** Cluster threshold */
double clusterEmin = 0.0; // GeV

/** Whether all non-gamma clusters are summed into one neutron candidate */
bool associateAllClustersToNeutron = false;

/** rotation from global to local coordinates */
double globalToProtonRotation = -0.025;
/** Neutron-photon separation in HcalFarForwardZDC */

/** Neutron-photon separation used for ZDC */
double gammaZMaxOffset = 400;
double gammaMaxLength = 100;
double gammaMaxWidth = 12;
Expand Down
41 changes: 21 additions & 20 deletions src/algorithms/reco/LambdaReconstruction.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -373,10 +373,6 @@ void LambdaReconstruction::process(const LambdaReconstruction::Input& input,
const double mL = std::sqrt(m2L);
const double dL = mL - m_lambda;

if (std::abs(dL) > m_cfg.lambdaMassWindow * m_lambda) {
continue;
}

const double pi0_term = pp.dmpi0_cand / (m_cfg.pi0Window * m_pi0);
const double lambda_term = dL / (m_cfg.lambdaMassWindow * m_lambda);

Expand Down Expand Up @@ -417,36 +413,41 @@ void LambdaReconstruction::process(const LambdaReconstruction::Input& input,
// --------------------------------------------------------------------------

auto better = [&](const LambdaCandidate& a, const LambdaCandidate& b) -> bool {
// Primary criterion: physical compatibility.
if (a.chi2 != b.chi2) {
return a.chi2 < b.chi2;
}
// Tie-breakers only: prefer ZDC information when candidates are similarly good.
if (a.n_cat != b.n_cat) {
return static_cast<int>(a.n_cat) < static_cast<int>(b.n_cat);
return a.n_cat == NeutronCategory::ZDC;
}
if (a.g_cat != b.g_cat) {
return static_cast<int>(a.g_cat) < static_cast<int>(b.g_cat);
}
if (a.chi2 != b.chi2) {
return a.chi2 < b.chi2;
}
if (a.pz != b.pz) {
return a.pz > b.pz;
}
return a.E > b.E;
};

int best_k = 0;
for (int k = 1; k < static_cast<int>(cands.size()); ++k) {
if (better(cands[k], cands[best_k])) {
best_k = k;
}
}
std::sort(cands.begin(), cands.end(), better);

@github-actions github-actions Bot Jul 9, 2026

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

⚠️ modernize-use-ranges ⚠️
use a ranges version of this algorithm

Suggested change
std::sort(cands.begin(), cands.end(), better);
std::ranges::sort(cands,, better);


const int max_trials = cands.size();

const auto& best = cands[best_k];
for (int k = 0; k < max_trials; ++k) {
const auto& best = cands[k];

const auto& g1 = gamma_pool[best.g_i];
const auto& g2 = gamma_pool[best.g_j];
const auto& n =
(best.n_cat == NeutronCategory::ZDC) ? neutrons_zdc[best.n_idx] : neutrons_other[best.n_idx];
const auto& g1 = gamma_pool[best.g_i];
const auto& g2 = gamma_pool[best.g_j];
const auto& n = (best.n_cat == NeutronCategory::ZDC) ? neutrons_zdc[best.n_idx]
: neutrons_other[best.n_idx];

reconstruct_from_triplet(n, g1, g2, out_lambdas, out_decay_products);
if (!reconstruct_from_triplet(n, g1, g2, out_lambdas, out_decay_products)) {
continue;
}

break;
}
}

} // namespace eicrecon
58 changes: 15 additions & 43 deletions src/factories/reco/FarForwardNeutralsReconstruction_factory.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,52 +20,30 @@ class FarForwardNeutralsReconstruction_factory
private:
std::unique_ptr<AlgoT> m_algo;

PodioInput<edm4eic::Cluster> m_clusters_hcal_input{this};
PodioInput<edm4eic::Cluster> m_clusters_b0_input{this};
PodioInput<edm4eic::Cluster> m_clusters_ecalendcapp_input{this};
PodioInput<edm4eic::Cluster> m_clusters_lfhcal_input{this};

PodioOutput<edm4eic::ReconstructedParticle> m_hcal_neutrals_output{this};
PodioOutput<edm4eic::ReconstructedParticle> m_b0_neutrals_output{this};
PodioOutput<edm4eic::ReconstructedParticle> m_ecalendcapp_neutrals_output{this};
PodioOutput<edm4eic::ReconstructedParticle> m_lfhcal_neutrals_output{this};
PodioInput<edm4eic::Cluster> m_clusters_input{this};
PodioOutput<edm4eic::ReconstructedParticle> m_neutrals_output{this};

ParameterRef<std::string> m_offset_position_name{this, "offsetPositionName",
config().offsetPositionName};

ParameterRef<std::vector<double>> m_n_scale_corr_coeff_hcal_zdc{
this, "neutronScaleCorrCoeffHcalZDC", config().neutronScaleCorrCoeffHcalZDC};

ParameterRef<std::vector<double>> m_gamma_scale_corr_coeff_hcal_zdc{
this, "gammaScaleCorrCoeffHcalZDC", config().gammaScaleCorrCoeffHcalZDC};

ParameterRef<std::vector<double>> m_n_scale_corr_coeff_b0ecal{
this, "neutronScaleCorrCoeffB0Ecal", config().neutronScaleCorrCoeffB0Ecal};

ParameterRef<std::vector<double>> m_gamma_scale_corr_coeff_b0ecal{
this, "gammaScaleCorrCoeffB0Ecal", config().gammaScaleCorrCoeffB0Ecal};

ParameterRef<std::vector<double>> m_n_scale_corr_coeff_ecalendcapp{
this, "neutronScaleCorrCoeffEcalEndcapP", config().neutronScaleCorrCoeffEcalEndcapP};
ParameterRef<std::vector<double>> m_neutron_scale_corr_coeff{this, "neutronScaleCorrCoeff",
config().neutronScaleCorrCoeff};

ParameterRef<std::vector<double>> m_gamma_scale_corr_coeff_ecalendcapp{
this, "gammaScaleCorrCoeffEcalEndcapP", config().gammaScaleCorrCoeffEcalEndcapP};
ParameterRef<std::vector<double>> m_gamma_scale_corr_coeff{this, "gammaScaleCorrCoeff",
config().gammaScaleCorrCoeff};

ParameterRef<std::vector<double>> m_n_scale_corr_coeff_lfhcal{
this, "neutronScaleCorrCoeffLFHCAL", config().neutronScaleCorrCoeffLFHCAL};
ParameterRef<bool> m_can_detect_gammas{this, "canDetectGammas", config().canDetectGammas};

ParameterRef<std::vector<double>> m_gamma_scale_corr_coeff_lfhcal{
this, "gammaScaleCorrCoeffLFHCAL", config().gammaScaleCorrCoeffLFHCAL};
ParameterRef<bool> m_can_detect_neutrons{this, "canDetectNeutrons", config().canDetectNeutrons};

ParameterRef<double> m_cluster_emin_hcal_zdc{this, "clusterEminHcalZDC",
config().clusterEminHcalZDC};
ParameterRef<std::string> m_gamma_mode{this, "gammaMode", config().gammaMode};

ParameterRef<double> m_cluster_emin_b0ecal{this, "clusterEminB0Ecal", config().clusterEminB0Ecal};
ParameterRef<std::string> m_neutron_mode{this, "neutronMode", config().neutronMode};

ParameterRef<double> m_cluster_emin_ecalendcapp{this, "clusterEminEcalEndcapP",
config().clusterEminEcalEndcapP};
ParameterRef<double> m_cluster_emin{this, "clusterEmin", config().clusterEmin};

ParameterRef<double> m_cluster_emin_lfhcal{this, "clusterEminLFHCAL", config().clusterEminLFHCAL};
ParameterRef<bool> m_associate_all_clusters_to_neutron{this, "associateAllClustersToNeutron",
config().associateAllClustersToNeutron};

ParameterRef<double> m_global_to_proton_rotation{this, "globalToProtonRotation",
config().globalToProtonRotation};
Expand All @@ -90,17 +68,11 @@ class FarForwardNeutralsReconstruction_factory
void Process(int32_t /* run_number */, uint64_t /* event_number */) {
m_algo->process(
{
m_clusters_hcal_input(),
m_clusters_b0_input(),
m_clusters_ecalendcapp_input(),
m_clusters_lfhcal_input(),
m_clusters_input(),

},
{
m_hcal_neutrals_output().get(),
m_b0_neutrals_output().get(),
m_ecalendcapp_neutrals_output().get(),
m_lfhcal_neutrals_output().get(),
m_neutrals_output().get(),

});
}
Expand Down
Loading
Loading