Skip to content
This repository was archived by the owner on Aug 31, 2026. It is now read-only.

Commit 35b66ff

Browse files
madbaronmadbaron
authored
Merge pull request #15 from MuonColliderSoft/missing-gaudi-algs
add BIB utils
2 parents e1d9b87 + b8dd1dc commit 35b66ff

40 files changed

Lines changed: 2677 additions & 1161 deletions

‎CMakeLists.txt‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -35,7 +35,7 @@ find_package(EDM4HEP REQUIRED)
3535
find_package(DD4hep REQUIRED)
3636
find_package(Gaudi REQUIRED)
3737
find_package(k4FWCore REQUIRED)
38-
find_package(ROOT REQUIRED COMPONENTS RIO Tree)
38+
find_package(ROOT REQUIRED COMPONENTS RIO Tree Hist)
3939
find_package(k4SimGeant4 REQUIRED)
4040
find_package(k4geo REQUIRED)
4141
find_package(CLHEP REQUIRED)

‎k4Reco/BIBUtils/README.md‎

Lines changed: 102 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,102 @@
1+
<!--
2+
Copyright (c) 2020-2024 Key4hep-Project.
3+
4+
This file is part of Key4hep.
5+
See https://key4hep.github.io/key4hep-doc/ for further info.
6+
7+
Licensed under the Apache License, Version 2.0 (the "License");
8+
you may not use this file except in compliance with the License.
9+
You may obtain a copy of the License at
10+
11+
http://www.apache.org/licenses/LICENSE-2.0
12+
13+
Unless required by applicable law or agreed to in writing, software
14+
distributed under the License is distributed on an "AS IS" BASIS,
15+
WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
16+
See the License for the specific language governing permissions and
17+
limitations under the License.
18+
-->
19+
# BIBUtils
20+
21+
Gaudi-native ports of the beam-induced-background (BIB) cleaning processors used
22+
in the Muon Collider reconstruction. They were originally Marlin processors:
23+
24+
| Gaudi algorithm | Original Marlin processor | Original package |
25+
|-------------------|---------------------------|------------------|
26+
| `FilterConeHits` | `FilterConeHits` | [MarlinTrkProcessors](https://github.com/MuonColliderSoft/MarlinTrkProcessors) |
27+
| `CaloConer` | `CaloConer` | [MyBIBUtils](https://github.com/madbaron/MyBIBUtils) |
28+
| `CaloHitSelector` | `CaloHitSelector` | [MyBIBUtils](https://github.com/madbaron/MyBIBUtils) |
29+
30+
All three are functional `k4FWCore::MultiTransformer`s. Selected hits are written
31+
to **subset** collections that reference the original hits, accompanied by a
32+
freshly built reco-to-sim link collection. The simulated hit of a given
33+
reconstructed hit is resolved through the input link collection rather than by
34+
relying on positional alignment between the hit and relation collections.
35+
36+
## FilterConeHits
37+
38+
Keeps the tracker hits that lie inside a cone around the trajectory of a
39+
generator-level MC particle. For each selected MC particle a helix is built from
40+
its production vertex, momentum and charge in the detector field (taken from the
41+
`GeoSvc`). A hit is kept when its angular distance to the helix is below
42+
`DeltaRCut` and/or its 3D distance to the helix is below `Dist3DCut`.
43+
44+
The helix math lives in the self-contained, header-only `TrackHelix` (see
45+
`include/TrackHelix.h`); it reproduces the point-to-helix distance and
46+
cylinder-crossing operations of MarlinUtil's `HelixClass` so the package carries
47+
no dependency on MarlinUtil or any other Marlin-era package.
48+
49+
Each instance handles a single tracker subdetector — configure one instance per
50+
collection (vertex/inner/outer × barrel/endcap), as in the original steering.
51+
52+
| Property | Default | Description |
53+
|---|---|---|
54+
| `MCParticleCollection` | `MCParticle` | input MC particles |
55+
| `TrackerHitInputCollections` | `VBTrackerHits` | input reco tracker hits |
56+
| `TrackerHitInputRelations` | `VBTrackerHitsRelations` | input reco→sim links |
57+
| `TrackerHitOutputCollections` | `VBTrackerHitsConed` | output reco hits (subset) |
58+
| `TrackerSimHitOutputCollections` | `VertexBarrelCollectionConed` | output sim hits (subset) |
59+
| `TrackerHitOutputRelations` | `VBTrackerHitsRelationsConed` | output reco→sim links |
60+
| `DeltaRCut` | `-1` | max angular distance to the helix [rad] (disabled if ≤ 0) |
61+
| `Dist3DCut` | `-1` | max 3D distance to the helix [mm] (disabled if ≤ 0) |
62+
| `ConeAroundStatus` | `[1]` | MC generator statuses to cone around |
63+
| `FillHistograms` | `false` | fill diagnostic histograms |
64+
| `TrackerOuterRadius` | `1500` | tracker barrel outer radius used to clip the helix [mm] |
65+
66+
## CaloConer
67+
68+
Keeps the calorimeter hits within a fixed angular cone (`ConeWidth`, in radians)
69+
around the direction of any generator-level (`generatorStatus == 1`) MC particle.
70+
71+
| Property | Default | Description |
72+
|---|---|---|
73+
| `MCParticleCollectionName` | `MCParticle` | input MC particles |
74+
| `CaloHitCollectionName` | `EcalBarrelCollectionRec` | input reco calo hits |
75+
| `CaloRelationCollectionName` | `EcalBarrelRelationsSimRec` | input reco→sim links |
76+
| `GoodHitCollection` | `EcalBarrelCollectionConed` | output reco hits (subset) |
77+
| `GoodRelationCollection` | `EcalBarrelRelationsSimConed` | output reco→sim links |
78+
| `ConeWidth` | `0.2` | half-opening angle of the cone [rad] |
79+
80+
## CaloHitSelector
81+
82+
Applies a per-cell energy threshold and a time-window selection to calorimeter
83+
hits. The threshold is read as a function of polar angle and layer from two ROOT
84+
histograms (`th_2dmode_sym`, `stddev_sym`) in `ThresholdsFilePath`:
85+
`threshold = mode + Nsigma * stddev`. A constant `FlatThreshold` (GeV) can be
86+
used instead. Surviving hits must fall inside `[TimeWindowMin, TimeWindowMax]`
87+
after a time-of-flight correction. The `GeoSvc` is used to decode the layer from
88+
the cellID.
89+
90+
| Property | Default | Description |
91+
|---|---|---|
92+
| `CaloHitCollectionName` | `EcalBarrelCollectionRec` | input reco calo hits |
93+
| `CaloRelationCollectionName` | `EcalBarrelRelationsSimRec` | input reco→sim links |
94+
| `GoodHitCollection` | `EcalBarrelCollectionSel` | output reco hits (subset) |
95+
| `GoodRelationCollection` | `EcalBarrelRelationsSimSel` | output reco→sim links |
96+
| `ThresholdsFilePath` | `""` | ROOT file with the threshold maps |
97+
| `Nsigma` | `3` | number of BIB-energy sigmas above the modal threshold |
98+
| `FlatThreshold` | `0` | constant threshold [GeV]; overrides the maps if > 0 |
99+
| `TimeWindowMin` / `TimeWindowMax` | `-0.5` / `10` | TOF-corrected time window [ns] |
100+
| `DoBIBsubtraction` | `false` | subtract the modal BIB energy from each cell |
101+
102+
See `options/` for runnable example steering files.

‎k4Reco/BIBUtils/include/BIBUtilsHelpers.h‎

Lines changed: 59 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,59 @@
1+
/*
2+
* Copyright (c) 2020-2024 Key4hep-Project.
3+
*
4+
* This file is part of Key4hep.
5+
* See https://key4hep.github.io/key4hep-doc/ for further info.
6+
*
7+
* Licensed under the Apache License, Version 2.0 (the "License");
8+
* you may not use this file except in compliance with the License.
9+
* You may obtain a copy of the License at
10+
*
11+
* http://www.apache.org/licenses/LICENSE-2.0
12+
*
13+
* Unless required by applicable law or agreed to in writing, software
14+
* distributed under the License is distributed on an "AS IS" BASIS,
15+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
16+
* See the License for the specific language governing permissions and
17+
* limitations under the License.
18+
*/
19+
#ifndef K4RECO_BIBUTILSHELPERS_H
20+
#define K4RECO_BIBUTILSHELPERS_H 1
21+
22+
#include <podio/ObjectID.h>
23+
24+
#include <cmath>
25+
#include <cstdint>
26+
27+
namespace k4reco::bibutils {
28+
29+
/// Pack a podio::ObjectID (collectionID, index) into a single 64-bit key so it
30+
/// can be used in unordered associative containers. This lets the BIB-cleaning
31+
/// algorithms map a reconstructed hit to its simulated hit via the input link
32+
/// collection, instead of relying on positional alignment between the hit and
33+
/// relation collections as the original Marlin processors did.
34+
inline std::uint64_t objectKey(const podio::ObjectID& id) {
35+
return (static_cast<std::uint64_t>(static_cast<std::uint32_t>(id.collectionID)) << 32) |
36+
static_cast<std::uint32_t>(id.index);
37+
}
38+
39+
/// Angle (in radians, in [0, pi]) between two 3-vectors. Reproduces the value
40+
/// returned by TVector3::Angle used in the original Marlin processors without
41+
/// pulling in the ROOT Physics library.
42+
inline double angleBetween(double ax, double ay, double az, double bx, double by, double bz) {
43+
const double na = std::sqrt(ax * ax + ay * ay + az * az);
44+
const double nb = std::sqrt(bx * bx + by * by + bz * bz);
45+
if (na == 0. || nb == 0.) {
46+
return 0.;
47+
}
48+
double cosAngle = (ax * bx + ay * by + az * bz) / (na * nb);
49+
if (cosAngle > 1.) {
50+
cosAngle = 1.;
51+
} else if (cosAngle < -1.) {
52+
cosAngle = -1.;
53+
}
54+
return std::acos(cosAngle);
55+
}
56+
57+
} // namespace k4reco::bibutils
58+
59+
#endif

‎k4Reco/BIBUtils/include/CaloConer.h‎

Lines changed: 60 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,60 @@
1+
/*
2+
* Copyright (c) 2020-2024 Key4hep-Project.
3+
*
4+
* This file is part of Key4hep.
5+
* See https://key4hep.github.io/key4hep-doc/ for further info.
6+
*
7+
* Licensed under the Apache License, Version 2.0 (the "License");
8+
* you may not use this file except in compliance with the License.
9+
* You may obtain a copy of the License at
10+
*
11+
* http://www.apache.org/licenses/LICENSE-2.0
12+
*
13+
* Unless required by applicable law or agreed to in writing, software
14+
* distributed under the License is distributed on an "AS IS" BASIS,
15+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
16+
* See the License for the specific language governing permissions and
17+
* limitations under the License.
18+
*/
19+
#ifndef K4RECO_CALOCONER_H
20+
#define K4RECO_CALOCONER_H 1
21+
22+
#include "Gaudi/Property.h"
23+
24+
#include <edm4hep/CaloHitSimCaloHitLinkCollection.h>
25+
#include <edm4hep/CalorimeterHitCollection.h>
26+
#include <edm4hep/MCParticleCollection.h>
27+
28+
#include <k4FWCore/Transformer.h>
29+
30+
#include <string>
31+
#include <tuple>
32+
33+
/** === CaloConer ===
34+
* Keeps only the calorimeter hits that fall within a fixed angular cone
35+
* (ConeWidth, in radians) around the direction of any generator-level
36+
* (generatorStatus == 1) MC particle. The selected hits are written to a
37+
* subset collection together with a freshly built reco-to-sim link collection.
38+
*
39+
* Gaudi-native port of the Marlin CaloConer processor from MyBIBUtils
40+
* (https://github.com/madbaron/MyBIBUtils), used to clean beam-induced
41+
* background (BIB) out of the calorimeters before particle flow.
42+
*
43+
* @author F. Meloni, DESY (original Marlin processor)
44+
*/
45+
struct CaloConer final : k4FWCore::MultiTransformer<
46+
std::tuple<edm4hep::CalorimeterHitCollection, edm4hep::CaloHitSimCaloHitLinkCollection>(
47+
const edm4hep::MCParticleCollection&, const edm4hep::CalorimeterHitCollection&,
48+
const edm4hep::CaloHitSimCaloHitLinkCollection&)> {
49+
CaloConer(const std::string& name, ISvcLocator* svcLoc);
50+
51+
std::tuple<edm4hep::CalorimeterHitCollection, edm4hep::CaloHitSimCaloHitLinkCollection>
52+
operator()(const edm4hep::MCParticleCollection& mcParticles, const edm4hep::CalorimeterHitCollection& caloHits,
53+
const edm4hep::CaloHitSimCaloHitLinkCollection& caloLinks) const override;
54+
55+
private:
56+
Gaudi::Property<double> m_coneSize{this, "ConeWidth", 0.2,
57+
"Half-opening angle of the cone around MC particles [rad]"};
58+
};
59+
60+
#endif

‎k4Reco/BIBUtils/include/CaloHitSelector.h‎

Lines changed: 88 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,88 @@
1+
/*
2+
* Copyright (c) 2020-2024 Key4hep-Project.
3+
*
4+
* This file is part of Key4hep.
5+
* See https://key4hep.github.io/key4hep-doc/ for further info.
6+
*
7+
* Licensed under the Apache License, Version 2.0 (the "License");
8+
* you may not use this file except in compliance with the License.
9+
* You may obtain a copy of the License at
10+
*
11+
* http://www.apache.org/licenses/LICENSE-2.0
12+
*
13+
* Unless required by applicable law or agreed to in writing, software
14+
* distributed under the License is distributed on an "AS IS" BASIS,
15+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
16+
* See the License for the specific language governing permissions and
17+
* limitations under the License.
18+
*/
19+
#ifndef K4RECO_CALOHITSELECTOR_H
20+
#define K4RECO_CALOHITSELECTOR_H 1
21+
22+
#include "Gaudi/Property.h"
23+
24+
#include <edm4hep/CaloHitSimCaloHitLinkCollection.h>
25+
#include <edm4hep/CalorimeterHitCollection.h>
26+
27+
#include <k4FWCore/Transformer.h>
28+
#include <k4Interface/IGeoSvc.h>
29+
30+
#include <GaudiKernel/StatusCode.h>
31+
32+
#include <memory>
33+
#include <string>
34+
#include <tuple>
35+
36+
class TH2D;
37+
38+
/** === CaloHitSelector ===
39+
* Applies per-cell energy and timing selections to calorimeter hits in order
40+
* to suppress beam-induced background (BIB). The energy threshold is read, as
41+
* a function of polar angle and calorimeter layer, from two ROOT histograms
42+
* (th_2dmode_sym, stddev_sym) stored in the file pointed to by
43+
* ThresholdsFilePath: threshold = mode + Nsigma * stddev. A constant
44+
* FlatThreshold (GeV) can be used instead. Hits passing the threshold are
45+
* further required to fall inside a time window, after correcting the hit time
46+
* for the time of flight from the origin. Selected hits are written to a
47+
* subset collection together with a freshly built reco-to-sim link collection.
48+
*
49+
* Gaudi-native port of the Marlin CaloHitSelector processor from MyBIBUtils
50+
* (https://github.com/madbaron/MyBIBUtils).
51+
*
52+
* @author F. Meloni, DESY (original Marlin processor)
53+
*/
54+
struct CaloHitSelector final
55+
: k4FWCore::MultiTransformer<
56+
std::tuple<edm4hep::CalorimeterHitCollection, edm4hep::CaloHitSimCaloHitLinkCollection>(
57+
const edm4hep::CalorimeterHitCollection&, const edm4hep::CaloHitSimCaloHitLinkCollection&)> {
58+
CaloHitSelector(const std::string& name, ISvcLocator* svcLoc);
59+
60+
StatusCode initialize() override;
61+
62+
std::tuple<edm4hep::CalorimeterHitCollection, edm4hep::CaloHitSimCaloHitLinkCollection>
63+
operator()(const edm4hep::CalorimeterHitCollection& caloHits,
64+
const edm4hep::CaloHitSimCaloHitLinkCollection& caloLinks) const override;
65+
66+
private:
67+
Gaudi::Property<std::string> m_thFile{this, "ThresholdsFilePath", "",
68+
"Path to the ROOT file holding the threshold maps (th_2dmode_sym, stddev_sym)"};
69+
Gaudi::Property<int> m_nSigma{this, "Nsigma", 3, "Number of BIB-energy sigmas added on top of the modal threshold"};
70+
Gaudi::Property<double> m_flatThreshold{this, "FlatThreshold", 0.,
71+
"Constant energy threshold [GeV]; if > 0 it overrides the map-based one"};
72+
Gaudi::Property<double> m_timeWindowMin{this, "TimeWindowMin", -0.5,
73+
"Lower edge of the time-of-flight corrected acceptance window [ns]"};
74+
Gaudi::Property<double> m_timeWindowMax{this, "TimeWindowMax", 10.,
75+
"Upper edge of the time-of-flight corrected acceptance window [ns]"};
76+
Gaudi::Property<bool> m_doBIBsubtraction{this, "DoBIBsubtraction", false,
77+
"Subtract the mean expected BIB energy (modal map) from each cell"};
78+
Gaudi::Property<std::string> m_encodingStringVariable{
79+
this, "EncodingStringParameterName", "GlobalCalorimeterReadoutID",
80+
"Name of the DD4hep constant holding the cellID encoding string for calorimeters"};
81+
82+
SmartIF<IGeoSvc> m_geoSvc;
83+
84+
std::unique_ptr<TH2D> m_thresholdMap;
85+
std::unique_ptr<TH2D> m_stddevMap;
86+
};
87+
88+
#endif

0 commit comments

Comments
 (0)