Skip to content

Commit

Permalink
Followup to PR #107 (noise handling) (#114)
Browse files Browse the repository at this point in the history
* move TubeLayerModuleThetaCaloTool to RecFCCeeCalorimeter. Avoid hardcoding systemID in C++, pass it via configurables. Remove test based on legacy scripts with k4run

* add noise and noise-filtering to automatic test
  • Loading branch information
giovannimarchiori committed Sep 17, 2024
1 parent 26276b7 commit 5e8879e
Show file tree
Hide file tree
Showing 5 changed files with 76 additions and 658 deletions.
7 changes: 0 additions & 7 deletions RecFCCeeCalorimeter/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -50,10 +50,3 @@ add_test(NAME FCCeeLAr_benchmarkCorrection
WORKING_DIRECTORY ${PROJECT_SOURCE_DIR}
)
set_test_env(FCCeeLAr_benchmarkCorrection)

# test for ALLEGRO addNoise, to be moved ddsim+k4run
add_test(NAME FCCeeLAr_ecalNoise
COMMAND k4run RecFCCeeCalorimeter/tests/options/runEcalBarrel_addNoise_thetamodulemerged.py
WORKING_DIRECTORY ${PROJECT_SOURCE_DIR}
)
set_test_env(FCCeeLAr_ecalNoise)
Original file line number Diff line number Diff line change
@@ -1,20 +1,24 @@
#include "TubeLayerModuleThetaMergedCaloTool.h"
#include "TubeLayerModuleThetaCaloTool.h"

// segm
// dd4hep
#include "DD4hep/Detector.h"
#include "DD4hep/MultiSegmentation.h"

// k4geo
#include "detectorCommon/DetUtils_k4geo.h"

// k4FWCore
#include "k4Interface/IGeoSvc.h"

DECLARE_COMPONENT(TubeLayerModuleThetaMergedCaloTool)
DECLARE_COMPONENT(TubeLayerModuleThetaCaloTool)

TubeLayerModuleThetaMergedCaloTool::TubeLayerModuleThetaMergedCaloTool(const std::string& type, const std::string& name,
TubeLayerModuleThetaCaloTool::TubeLayerModuleThetaCaloTool(const std::string& type, const std::string& name,
const IInterface* parent)
: AlgTool(type, name, parent), m_geoSvc("GeoSvc", name) {
declareInterface<ICalorimeterTool>(this);
}

StatusCode TubeLayerModuleThetaMergedCaloTool::initialize() {
StatusCode TubeLayerModuleThetaCaloTool::initialize() {
StatusCode sc = AlgTool::initialize();
if (sc.isFailure()) return sc;

Expand All @@ -24,7 +28,7 @@ StatusCode TubeLayerModuleThetaMergedCaloTool::initialize() {
return StatusCode::FAILURE;
}
if (m_readoutName != "") {
// Check if readouts exist
// Check if readout exists
info() << "Readout: " << m_readoutName << endmsg;
if (m_geoSvc->getDetector()->readouts().find(m_readoutName) == m_geoSvc->getDetector()->readouts().end()) {
error() << "Readout <<" << m_readoutName << ">> does not exist." << endmsg;
Expand All @@ -34,25 +38,17 @@ StatusCode TubeLayerModuleThetaMergedCaloTool::initialize() {
return sc;
}

StatusCode TubeLayerModuleThetaMergedCaloTool::finalize() { return AlgTool::finalize(); }
StatusCode TubeLayerModuleThetaCaloTool::finalize() { return AlgTool::finalize(); }

StatusCode TubeLayerModuleThetaMergedCaloTool::prepareEmptyCells(std::unordered_map<uint64_t, double>& aCells) {
StatusCode TubeLayerModuleThetaCaloTool::prepareEmptyCells(std::unordered_map<uint64_t, double>& aCells) {
// Get the total number of active volumes in the geometry
auto highestVol = gGeoManager->GetTopVolume();
unsigned int numLayers;
if (!m_activeVolumesNumber) {
numLayers = det::utils::countPlacedVolumes(highestVol, m_activeVolumeName);
} else {
// used when MergeLayers tool is used. To be removed once MergeLayer gets replaced by RedoSegmentation.
numLayers = m_activeVolumesNumber;
}
unsigned int numLayers = m_activeVolumesNumber;
info() << "Number of active layers " << numLayers << endmsg;

// get segmentation
dd4hep::DDSegmentation::Segmentation *aSegmentation = m_geoSvc->getDetector()->readout(m_readoutName).segmentation().segmentation();
std::string segmentationType = aSegmentation->type();
info() << "Segmentation type : " << segmentationType << endmsg;

dd4hep::DDSegmentation::FCCSWGridModuleThetaMerged_k4geo *moduleThetaSegmentation = nullptr;
if (segmentationType == "FCCSWGridModuleThetaMerged_k4geo") {
moduleThetaSegmentation = dynamic_cast<dd4hep::DDSegmentation::FCCSWGridModuleThetaMerged_k4geo *>(aSegmentation);
Expand All @@ -72,7 +68,9 @@ StatusCode TubeLayerModuleThetaMergedCaloTool::prepareEmptyCells(std::unordered_

for (unsigned int ilayer = 0; ilayer < numLayers; ilayer++) {
dd4hep::DDSegmentation::CellID volumeId = 0;
(*decoder)["system"].set(volumeId, 4);
for (unsigned int it = 0; it < m_fieldNames.size(); it++) {
decoder->set(volumeId, m_fieldNames[it], m_fieldValues[it]);
}
(*decoder)["layer"].set(volumeId, ilayer);
(*decoder)["theta"].set(volumeId, 0);
(*decoder)["module"].set(volumeId, 0);
Expand Down
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
#ifndef RECCALORIMETER_TUBELAYERMODULETHETAMERGEDCALOTOOL_H
#define RECCALORIMETER_TUBELAYERMODULETHETAMERGEDCALOTOOL_H
#ifndef RECFCCEECALORIMETER_TUBELAYERMODULETHETAMERGEDCALOTOOL_H
#define RECFCCEECALORIMETER_TUBELAYERMODULETHETAMERGEDCALOTOOL_H

// from Gaudi
#include "GaudiKernel/AlgTool.h"
Expand All @@ -9,8 +9,8 @@

class IGeoSvc;

/** @class TubeLayerModuleThetaMergedCaloTool Reconstruction/RecCalorimeter/src/components/TubeLayerModuleThetaMergedCaloTool.h
*TubeLayerModuleThetaMergedCaloTool.h
/** @class TubeLayerModuleThetaCaloTool k4RecCalorimeter/RecFCCeeCalorimeter/src/components/TubeLayerModuleThetaCaloTool.h
* TubeLayerModuleThetaCaloTool.h
*
* Tool for geometry-dependent settings of the digitisation.
* It assumes cylindrical geometry (layers) and phi-theta segmentation.
Expand All @@ -19,10 +19,10 @@ class IGeoSvc;
* @author Zhibo Wu
*/

class TubeLayerModuleThetaMergedCaloTool : public AlgTool, virtual public ICalorimeterTool {
class TubeLayerModuleThetaCaloTool : public AlgTool, virtual public ICalorimeterTool {
public:
TubeLayerModuleThetaMergedCaloTool(const std::string& type, const std::string& name, const IInterface* parent);
virtual ~TubeLayerModuleThetaMergedCaloTool() = default;
TubeLayerModuleThetaCaloTool(const std::string& type, const std::string& name, const IInterface* parent);
virtual ~TubeLayerModuleThetaCaloTool() = default;
virtual StatusCode initialize() final;
virtual StatusCode finalize() final;
/** Prepare a map of all existing cells in current geometry.
Expand Down Expand Up @@ -52,8 +52,8 @@ class TubeLayerModuleThetaMergedCaloTool : public AlgTool, virtual public ICalor
Gaudi::Property<std::vector<std::string>> m_fieldNames{this, "fieldNames"};
/// Values of the fields describing the segmented volume
Gaudi::Property<std::vector<int>> m_fieldValues{this, "fieldValues"};
/// Temporary: for use with MergeLayer tool
/// Number of layers
Gaudi::Property<unsigned int> m_activeVolumesNumber{this, "activeVolumesNumber", 0};
};

#endif /* RECCALORIMETER_TUBELAYERMODULETHETAMERGEDCALOTOOL_H */
#endif /* RECFCCEECALORIMETER_TUBELAYERMODULETHETAMERGEDCALOTOOL_H */
82 changes: 51 additions & 31 deletions RecFCCeeCalorimeter/tests/options/ALLEGRO_o1_v03_digi_reco.py
Original file line number Diff line number Diff line change
Expand Up @@ -53,7 +53,7 @@
# cluster energy corrections
# simple parametrisations of up/downstream losses for ECAL-only clusters
# not to be applied for ECAL+HCAL clustering
# superseded by MVA calibration, but turned on here for the purpose of testing that the code is not broken
# superseded by MVA calibration, but turned on here for the purpose of testing that the code is not broken - will end up in separate cluster collection
applyUpDownstreamCorrections = True

# BDT regression from total cluster energy and fraction of energy in each layer (after correction for sampling fraction)
Expand All @@ -68,6 +68,7 @@
# not run by default in production, but to be turned on here for the purpose of testing that the code is not broken
# currently off till we provide the onnx files
runPhotonIDTool = False
logEWeightInPhotonID = False

#
# ALGORITHMS AND SERVICES SETUP
Expand Down Expand Up @@ -147,6 +148,12 @@
readoutName=ecalBarrelReadoutName,
OutputLevel=INFO
)
# the noise tool needs the positioning tool, but if I reuse the previous one the code crashes..
cellPositionEcalBarrelToolForNoise = CellPositionsECalBarrelModuleThetaSegTool(
"CellPositionsECalBarrelForNoise",
readoutName=ecalBarrelReadoutName,
OutputLevel=INFO
)
if resegmentECalBarrel:
cellPositionEcalBarrelTool2 = CellPositionsECalBarrelModuleThetaSegTool(
"CellPositionsECalBarrel2",
Expand Down Expand Up @@ -273,50 +280,54 @@
ecalBarrelNoiseRMSHistName = "h_elecNoise_fcc_"
from Configurables import NoiseCaloCellsVsThetaFromFileTool
noiseBarrel = NoiseCaloCellsVsThetaFromFileTool("NoiseBarrel",
cellPositionsTool=cellPositionEcalBarrelTool,
cellPositionsTool=cellPositionEcalBarrelToolForNoise,
readoutName=ecalBarrelReadoutName,
noiseFileName=ecalBarrelNoisePath,
elecNoiseRMSHistoName=ecalBarrelNoiseRMSHistName,
setNoiseOffset=False,
activeFieldName="layer",
addPileup=False,
filterNoiseThreshold=0,
filterNoiseThreshold=1,
numRadialLayers=11,
scaleFactor=1 / 1000., # MeV to GeV
OutputLevel=INFO)

# needs to be migrated!
# from Configurables import TubeLayerPhiEtaCaloTool
# barrelGeometry = TubeLayerPhiEtaCaloTool("EcalBarrelGeo",
# readoutName=ecalBarrelReadoutNamePhiEta,
# activeVolumeName="LAr_sensitive",
# activeFieldName="layer",
# activeVolumesNumber=12,
# fieldNames=["system"],
# fieldValues=[4])
from Configurables import TubeLayerModuleThetaCaloTool
barrelGeometry = TubeLayerModuleThetaCaloTool("EcalBarrelGeo",
readoutName=ecalBarrelReadoutName,
activeVolumeName="LAr_sensitive",
activeFieldName="layer",
activeVolumesNumber=11,
fieldNames=["system"],
fieldValues=[4],
OutputLevel=INFO)

# cells with noise not filtered
# createEcalBarrelCellsNoise = CreateCaloCells("CreateECalBarrelCellsNoise",
# doCellCalibration=False,
# addCellNoise=True,
# filterCellNoise=False,
# OutputLevel=INFO,
# hits="ECalBarrelCellsStep2",
# noiseTool=noiseBarrel,
# geometryTool=barrelGeometry,
# cells=EcalBarrelCellsName)
createEcalBarrelCellsNoise = CreatePositionedCaloCells("CreatePositionedECalBarrelCellsWithNoise",
doCellCalibration=True,
positionsTool=cellPositionEcalBarrelTool,
calibTool=calibEcalBarrel,
addCellNoise=True,
filterCellNoise=False,
noiseTool=noiseBarrel,
geometryTool=barrelGeometry,
OutputLevel=INFO,
hits=ecalBarrelReadoutName, # uncalibrated & unpositioned cells without noise
cells=ecalBarrelPositionedCellsName + "WithNoise")

# cells with noise filtered
# createEcalBarrelCellsNoise = CreateCaloCells("CreateECalBarrelCellsNoise_filtered",
# doCellCalibration=False,
# addCellNoise=True,
# filterCellNoise=True,
# OutputLevel=INFO,
# hits="ECalBarrelCellsStep2",
# noiseTool=noiseBarrel,
# geometryTool=barrelGeometry,
# cells=EcalBarrelCellsName)

createEcalBarrelCellsNoiseFiltered = CreatePositionedCaloCells("CreateECalBarrelCellsWithNoiseFiltered",
doCellCalibration=True,
calibTool=calibEcalBarrel,
positionsTool=cellPositionEcalBarrelTool,
addCellNoise=True,
filterCellNoise=True,
noiseTool=noiseBarrel,
geometryTool=barrelGeometry,
OutputLevel=INFO,
hits=ecalBarrelReadoutName, # uncalibrated & unpositioned cells without noise
cells=ecalBarrelPositionedCellsName + "WithNoiseFiltered"
)

if runHCal:
# Apply calibration and positioning to cells in HCal barrel
Expand Down Expand Up @@ -526,6 +537,7 @@
layerFieldNames=["layer"],
thetaRecalcWeights=[ecalBarrelThetaWeights],
do_photon_shapeVar=runPhotonIDTool,
do_widthTheta_logE_weights=logEWeightInPhotonID,
OutputLevel=INFO
)

Expand Down Expand Up @@ -705,6 +717,7 @@
layerFieldNames=["layer"],
thetaRecalcWeights=[ecalBarrelThetaWeights],
do_photon_shapeVar=runPhotonIDTool,
do_widthTheta_logE_weights=logEWeightInPhotonID,
OutputLevel=INFO)

if applyMVAClusterEnergyCalibration:
Expand Down Expand Up @@ -882,6 +895,13 @@
]
createEcalBarrelCells.AuditExecute = True
createEcalEndcapCells.AuditExecute = True
if addNoise:
TopAlg += [
createEcalBarrelCellsNoise,
createEcalBarrelCellsNoiseFiltered
]
createEcalBarrelCellsNoise.AuditExecute = True
createEcalBarrelCellsNoiseFiltered.AuditExecute = True

if resegmentECalBarrel:
TopAlg += [
Expand Down
Loading

0 comments on commit 5e8879e

Please sign in to comment.