From 0a558ae4f5f2c3f50201db01ce2de193bc23f551 Mon Sep 17 00:00:00 2001 From: akioogawa Date: Thu, 1 May 2025 15:00:31 -0400 Subject: [PATCH 1/6] Update FEMC for new geometry: Adding switch between Homogeneous and ScFi geometry implementations based on xml file loaded - For Homogenous, we keep enrgy smeasring as is. SF=1.0 - For ScFi, hits from fibers are summed to a tower and no enrgy smearing applied. SF=0.03. Added option to put SiPM saturation to CalorimeterHitDigi - Specify 2 new parameters : totalPixel and nPhotonPerGeV - default for totalPixel is 0, which case no attenuation is applied For FEMC for both geometry models, SiPM saturation is ON by default - Use "-PFEMC:SiPMSaturation=OFF" to turn it off --- .../calorimetry/CalorimeterHitDigi.cc | 12 +- .../calorimetry/CalorimeterHitDigiConfig.h | 3 + src/detectors/FEMC/FEMC.cc | 264 +++++++----------- .../calorimetry/CalorimeterHitDigi_factory.h | 2 + 4 files changed, 115 insertions(+), 166 deletions(-) diff --git a/src/algorithms/calorimetry/CalorimeterHitDigi.cc b/src/algorithms/calorimetry/CalorimeterHitDigi.cc index 95457082e9..f50a498b68 100644 --- a/src/algorithms/calorimetry/CalorimeterHitDigi.cc +++ b/src/algorithms/calorimetry/CalorimeterHitDigi.cc @@ -206,11 +206,19 @@ void CalorimeterHitDigi::process( ) : 0; double corrMeanScale_value = corrMeanScale(leading_hit); - + double ped = m_cfg.pedMeanADC + m_gaussian(m_generator) * m_cfg.pedSigmaADC; + //SiPM Saturation + double saturation = 1.0; + if(m_cfg.totalPixel>0){ + double nPhoton = edep * corrMeanScale_value * m_cfg.nPhotonPerGeV; + saturation = (1.0 - exp(-nPhoton/m_cfg.totalPixel))*m_cfg.totalPixel/nPhoton; + //printf("edep=%8.4f nPhoton=%8.2f totalPixel=%8.1f Saturation=%8.6f\n",edep,nPhoton,m_cfg.totalPixel,saturation); + } + // Note: both adc and tdc values must be positive numbers to avoid integer wraparound - unsigned long long adc = std::max(std::llround(ped + edep * corrMeanScale_value * (1.0 + eResRel) / m_cfg.dyRangeADC * m_cfg.capADC), 0LL); + unsigned long long adc = std::max(std::llround(ped + edep * corrMeanScale_value * saturation * (1.0 + eResRel) / m_cfg.dyRangeADC * m_cfg.capADC), 0LL); unsigned long long tdc = std::llround((time + m_gaussian(m_generator) * tRes) * stepTDC); if (edep> 1.e-3) trace("E sim {} \t adc: {} \t time: {}\t maxtime: {} \t tdc: {} \t corrMeanScale: {}", edep, adc, time, m_cfg.capTime, tdc, corrMeanScale_value); diff --git a/src/algorithms/calorimetry/CalorimeterHitDigiConfig.h b/src/algorithms/calorimetry/CalorimeterHitDigiConfig.h index 2b25d55c46..88f60d4898 100644 --- a/src/algorithms/calorimetry/CalorimeterHitDigiConfig.h +++ b/src/algorithms/calorimetry/CalorimeterHitDigiConfig.h @@ -29,6 +29,9 @@ namespace eicrecon { std::string readout{""}; std::vector fields{}; + //SiPM Saturation + double totalPixel{0}; + double nPhotonPerGeV{0}; }; } // eicrecon diff --git a/src/detectors/FEMC/FEMC.cc b/src/detectors/FEMC/FEMC.cc index 209a5fbd42..a25d9e54b6 100644 --- a/src/detectors/FEMC/FEMC.cc +++ b/src/detectors/FEMC/FEMC.cc @@ -1,3 +1,4 @@ + // SPDX-License-Identifier: LGPL-3.0-or-later // Copyright (C) 2021 - 2024, Chao Peng, Sylvester Joosten, Whitney Armstrong, David Lawrence, Friederike Bock, Wouter Deconinck, Kolja Kauder, Sebouh Paul @@ -8,6 +9,7 @@ #include #include "algorithms/calorimetry/CalorimeterHitDigiConfig.h" +#include "algorithms/calorimetry/CalorimeterHitRecoConfig.h" #include "extensions/jana/JOmniFactoryGeneratorT.h" #include "factories/calorimetry/CalorimeterClusterRecoCoG_factory.h" #include "factories/calorimetry/CalorimeterHitDigi_factory.h" @@ -23,55 +25,119 @@ extern "C" { using namespace eicrecon; InitJANAPlugin(app); + + auto log_service = app->GetService(); + auto mLog = log_service->logger("FEMC"); + + int SiPMSaturation=1; + std::string opt("ON"); + app->SetDefaultParameter("FEMC:SiPMSaturation",opt,"Turn ON(default) or OFF SiPM Saturation"); + if(opt.find("OFF") != std::string::npos || + opt.find("off") != std::string::npos || + opt.find("Off") != std::string::npos ){ + mLog->info("SiPM Saturation OFF"); + SiPMSaturation=0; + }else{ + mLog->info("SiPM Saturation ON"); + } + // Make sure digi and reco use the same value decltype(CalorimeterHitDigiConfig::capADC) EcalEndcapP_capADC = 16384; //16384, assuming 14 bits. For approximate HGCROC resolution use 65536 - decltype(CalorimeterHitDigiConfig::dyRangeADC) EcalEndcapP_dyRangeADC = 3 * dd4hep::GeV; + decltype(CalorimeterHitDigiConfig::dyRangeADC) EcalEndcapP_dyRangeADC = 100 * dd4hep::GeV; decltype(CalorimeterHitDigiConfig::pedMeanADC) EcalEndcapP_pedMeanADC = 200; decltype(CalorimeterHitDigiConfig::pedSigmaADC) EcalEndcapP_pedSigmaADC = 2.4576; decltype(CalorimeterHitDigiConfig::resolutionTDC) EcalEndcapP_resolutionTDC = 10 * dd4hep::picosecond; - app->Add(new JOmniFactoryGeneratorT( - "EcalEndcapPRawHits", - {"EcalEndcapPHits"}, + const double sampFrac=0.03; + decltype(CalorimeterHitDigiConfig::corrMeanScale) EcalEndcapP_corrMeanScale = Form("%f",1.0/sampFrac); + decltype(CalorimeterHitRecoConfig::sampFrac) EcalEndcapP_sampFrac = Form("%f",sampFrac); + decltype(CalorimeterHitDigiConfig::nPhotonPerGeV) EcalEndcapP_nPhotonPerGeV = 1500; + decltype(CalorimeterHitDigiConfig::totalPixel) EcalEndcapP_totalPixel = 4*159565; + if(SiPMSaturation==0) EcalEndcapP_totalPixel=0; + + int fEcalHomoScfi = 0; + try { + auto detector = app->GetService()->detector(); + fEcalHomoScfi = detector->constant("ForwardEcal_Homogeneous_Scfi"); + if(fEcalHomoScfi<=1) mLog->info("Homogeneous geometry loaded"); + else mLog->info("ScFi geometry loaded"); + } catch(...){}; + + if(fEcalHomoScfi<=1){ + app->Add(new JOmniFactoryGeneratorT( + "EcalEndcapPRawHits", + {"EcalEndcapPHits"}, #if EDM4EIC_VERSION_MAJOR >= 7 - {"EcalEndcapPRawHits", "EcalEndcapPRawHitAssociations"}, + {"EcalEndcapPRawHits", "EcalEndcapPRawHitAssociations"}, #else - {"EcalEndcapPRawHits"}, + {"EcalEndcapPRawHits"}, #endif - { - .eRes = {0.11333 * sqrt(dd4hep::GeV), 0.03, 0.0 * dd4hep::GeV}, // (11.333% / sqrt(E)) \oplus 3% - .tRes = 0.0, - .threshold = 0.0, - // .threshold = 15 * dd4hep::MeV for a single tower, applied on ADC level - .capADC = EcalEndcapP_capADC, - .capTime = 100, // given in ns, 4 samples in HGCROC - .dyRangeADC = EcalEndcapP_dyRangeADC, - .pedMeanADC = EcalEndcapP_pedMeanADC, - .pedSigmaADC = EcalEndcapP_pedSigmaADC, - .resolutionTDC = EcalEndcapP_resolutionTDC, - .corrMeanScale = "0.03", - .readout = "EcalEndcapPHits", - }, - app // TODO: Remove me once fixed - )); + { + .eRes = {0.11333 * sqrt(dd4hep::GeV), 0.03, 0.0 * dd4hep::GeV}, // (11.333% / sqrt(E)) \oplus 3% + .tRes = 0.0, + .threshold = 0.0, // 15MeV threshold for a single tower will be applied on ADC at Reco below + .capADC = EcalEndcapP_capADC, + .capTime = 100, // given in ns, 4 samples in HGCROC + .dyRangeADC = EcalEndcapP_dyRangeADC, + .pedMeanADC = EcalEndcapP_pedMeanADC, + .pedSigmaADC = EcalEndcapP_pedSigmaADC, + .resolutionTDC = EcalEndcapP_resolutionTDC, + .corrMeanScale = "1.0", + .readout = "EcalEndcapPHits", + .totalPixel = EcalEndcapP_totalPixel, + .nPhotonPerGeV = EcalEndcapP_nPhotonPerGeV, + }, + app // TODO: Remove me once fixed + )); + }else if(fEcalHomoScfi==2){ + app->Add(new JOmniFactoryGeneratorT( + "EcalEndcapPRawHits", + {"EcalEndcapPHits"}, +#if EDM4EIC_VERSION_MAJOR >= 7 + {"EcalEndcapPRawHits", "EcalEndcapPRawHitAssociations"}, +#else + {"EcalEndcapPRawHits"}, +#endif + { + .eRes = {0.0,0.0,0.0}, + .tRes = 0.0, + .threshold = 0.0, // 15MeV threshold for a single tower will be applied on ADC at Reco below + .capADC = EcalEndcapP_capADC, + .capTime = 100, // given in ns, 4 samples in HGCROC + .dyRangeADC = EcalEndcapP_dyRangeADC, + .pedMeanADC = EcalEndcapP_pedMeanADC, + .pedSigmaADC = EcalEndcapP_pedSigmaADC, + .resolutionTDC = EcalEndcapP_resolutionTDC, + .corrMeanScale = EcalEndcapP_corrMeanScale, + .readout = "EcalEndcapPHits", + .fields = {"fiberx","fibery","x","y"}, + .totalPixel = EcalEndcapP_totalPixel, + .nPhotonPerGeV = EcalEndcapP_nPhotonPerGeV, + }, + app // TODO: Remove me once fixed + )); + } + app->Add(new JOmniFactoryGeneratorT( "EcalEndcapPRecHits", {"EcalEndcapPRawHits"}, {"EcalEndcapPRecHits"}, { - .capADC = EcalEndcapP_capADC, - .dyRangeADC = EcalEndcapP_dyRangeADC, - .pedMeanADC = EcalEndcapP_pedMeanADC, - .pedSigmaADC = EcalEndcapP_pedSigmaADC, - .resolutionTDC = EcalEndcapP_resolutionTDC, - .thresholdFactor = 0.0, - .thresholdValue = 2, // The ADC of a 15 MeV particle is adc = 200 + 15 * 0.03 * ( 1.0 + 0) / 3000 * 16384 = 200 + 2.4576 - .sampFrac = "0.03", - .readout = "EcalEndcapPHits", - }, - app // TODO: Remove me once fixed - )); + .capADC = EcalEndcapP_capADC, + .dyRangeADC = EcalEndcapP_dyRangeADC, + .pedMeanADC = EcalEndcapP_pedMeanADC, + .pedSigmaADC = EcalEndcapP_pedSigmaADC, + .resolutionTDC = EcalEndcapP_resolutionTDC, + .thresholdFactor = 0.0, + .thresholdValue = 2, // The ADC of a 15 MeV particle is adc = 200 + 15 * 0.03 * ( 1.0 + 0) / 3000 * 16384 = 200 + 2.4576 + .sampFrac = "1.0", // alerady taken care in DIGI code above + .readout = "EcalEndcapPHits", + }, + app // TODO: Remove me once fixed )); + )); + app->Add(new JOmniFactoryGeneratorT( "EcalEndcapPTruthProtoClusters", {"EcalEndcapPRecHits", "EcalEndcapPHits"}, {"EcalEndcapPTruthProtoClusters"}, app // TODO: Remove me once fixed )); + app->Add(new JOmniFactoryGeneratorT( "EcalEndcapPIslandProtoClusters", {"EcalEndcapPRecHits"}, {"EcalEndcapPIslandProtoClusters"}, { @@ -208,135 +274,5 @@ extern "C" { app ) ); - - // Insert is identical to regular Ecal - app->Add(new JOmniFactoryGeneratorT( - "EcalEndcapPInsertRawHits", - {"EcalEndcapPInsertHits"}, -#if EDM4EIC_VERSION_MAJOR >= 7 - {"EcalEndcapPInsertRawHits", "EcalEndcapPInsertRawHitAssociations"}, -#else - {"EcalEndcapPInsertRawHits"}, -#endif - { - .eRes = {0.11333 * sqrt(dd4hep::GeV), 0.03, 0.0 * dd4hep::GeV}, // (11.333% / sqrt(E)) \oplus 3% - .tRes = 0.0, - .threshold = 0.0, - // .threshold = 15 * dd4hep::MeV for a single tower, applied on ADC level - .capADC = EcalEndcapP_capADC, - .capTime = 100, // given in ns, 4 samples in HGCROC - .dyRangeADC = EcalEndcapP_dyRangeADC, - .pedMeanADC = EcalEndcapP_pedMeanADC, - .pedSigmaADC = EcalEndcapP_pedSigmaADC, - .resolutionTDC = EcalEndcapP_resolutionTDC, - .corrMeanScale = "0.03", - .readout = "EcalEndcapPInsertHits", - }, - app // TODO: Remove me once fixed - )); - app->Add(new JOmniFactoryGeneratorT( - "EcalEndcapPInsertRecHits", {"EcalEndcapPInsertRawHits"}, {"EcalEndcapPInsertRecHits"}, - { - .capADC = EcalEndcapP_capADC, - .dyRangeADC = EcalEndcapP_dyRangeADC, - .pedMeanADC = EcalEndcapP_pedMeanADC, - .pedSigmaADC = EcalEndcapP_pedSigmaADC, - .resolutionTDC = EcalEndcapP_resolutionTDC, - .thresholdFactor = 0.0, - .thresholdValue = 2, // The ADC of a 15 MeV particle is adc = 200 + 15 * 0.03 * ( 1.0 + 0) / 3000 * 16384 = 200 + 2.4576 - .sampFrac = "0.03", - .readout = "EcalEndcapPInsertHits", - }, - app // TODO: Remove me once fixed - )); - app->Add(new JOmniFactoryGeneratorT( - "EcalEndcapPInsertTruthProtoClusters", {"EcalEndcapPInsertRecHits", "EcalEndcapPInsertHits"}, {"EcalEndcapPInsertTruthProtoClusters"}, - app // TODO: Remove me once fixed - )); - app->Add(new JOmniFactoryGeneratorT( - "EcalEndcapPInsertIslandProtoClusters", {"EcalEndcapPInsertRecHits"}, {"EcalEndcapPInsertIslandProtoClusters"}, - { - .sectorDist = 5.0 * dd4hep::cm, - .dimScaledLocalDistXY = {1.5,1.5}, - .splitCluster = false, - .minClusterHitEdep = 0.0 * dd4hep::MeV, - .minClusterCenterEdep = 60.0 * dd4hep::MeV, - .transverseEnergyProfileMetric = "dimScaledLocalDistXY", - .transverseEnergyProfileScale = 1., - }, - app // TODO: Remove me once fixed - )); - - app->Add( - new JOmniFactoryGeneratorT( - "EcalEndcapPInsertTruthClustersWithoutShapes", - {"EcalEndcapPInsertTruthProtoClusters", // edm4eic::ProtoClusterCollection -#if EDM4EIC_VERSION_MAJOR >= 7 - "EcalEndcapPInsertRawHitAssociations"}, // edm4eic::MCRecoCalorimeterHitCollection -#else - "EcalEndcapPInsertHits"}, // edm4hep::SimCalorimeterHitCollection -#endif - {"EcalEndcapPInsertTruthClustersWithoutShapes", // edm4eic::Cluster - "EcalEndcapPInsertTruthClusterAssociationsWithoutShapes"}, // edm4eic::MCRecoClusterParticleAssociation - { - .energyWeight = "log", - .sampFrac = 1.0, - .logWeightBase = 6.2, - .enableEtaBounds = true - }, - app // TODO: Remove me once fixed - ) - ); - - app->Add( - new JOmniFactoryGeneratorT( - "EcalEndcapPInsertTruthClusters", - {"EcalEndcapPInsertTruthClustersWithoutShapes", - "EcalEndcapPInsertTruthClusterAssociationsWithoutShapes"}, - {"EcalEndcapPInsertTruthClusters", - "EcalEndcapPInsertTruthClusterAssociations"}, - { - .energyWeight = "log", - .logWeightBase = 6.2 - }, - app - ) - ); - - app->Add( - new JOmniFactoryGeneratorT( - "EcalEndcapPInsertClustersWithoutShapes", - {"EcalEndcapPInsertIslandProtoClusters", // edm4eic::ProtoClusterCollection -#if EDM4EIC_VERSION_MAJOR >= 7 - "EcalEndcapPInsertRawHitAssociations"}, // edm4eic::MCRecoCalorimeterHitCollection -#else - "EcalEndcapPInsertHits"}, // edm4hep::SimCalorimeterHitCollection -#endif - {"EcalEndcapPInsertClustersWithoutShapes", // edm4eic::Cluster - "EcalEndcapPInsertClusterAssociationsWithoutShapes"}, // edm4eic::MCRecoClusterParticleAssociation - { - .energyWeight = "log", - .sampFrac = 1.0, - .logWeightBase = 3.6, - .enableEtaBounds = false, - }, - app // TODO: Remove me once fixed - ) - ); - - app->Add( - new JOmniFactoryGeneratorT( - "EcalEndcapPInsertClusters", - {"EcalEndcapPInsertClustersWithoutShapes", - "EcalEndcapPInsertClusterAssociationsWithoutShapes"}, - {"EcalEndcapPInsertClusters", - "EcalEndcapPInsertClusterAssociations"}, - { - .energyWeight = "log", - .logWeightBase = 3.6 - }, - app - ) - ); } } diff --git a/src/factories/calorimetry/CalorimeterHitDigi_factory.h b/src/factories/calorimetry/CalorimeterHitDigi_factory.h index 1fc28bc2da..b55dbc830a 100644 --- a/src/factories/calorimetry/CalorimeterHitDigi_factory.h +++ b/src/factories/calorimetry/CalorimeterHitDigi_factory.h @@ -34,6 +34,8 @@ class CalorimeterHitDigi_factory : public JOmniFactory m_corrMeanScale {this, "scaleResponse", config().corrMeanScale}; ParameterRef> m_fields {this, "signalSumFields", config().fields}; ParameterRef m_readout {this, "readoutClass", config().readout}; + ParameterRef m_nPhotonPerGeV {this, "numberOfPhotonPerGeV", config().nPhotonPerGeV}; + ParameterRef m_totalPixel {this, "totalNnumberOfPixel", config().totalPixel}; Service m_algorithmsInit {this}; From 93febe8ae7d5c6be20d23fbe8d8401558f1d0055 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 2 May 2025 17:42:38 +0000 Subject: [PATCH 2/6] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- .../calorimetry/CalorimeterHitDigi.cc | 18 +- .../calorimetry/CalorimeterHitDigiConfig.h | 4 +- src/detectors/FEMC/FEMC.cc | 173 +++++++++--------- .../calorimetry/CalorimeterHitDigi_factory.h | 4 +- 4 files changed, 103 insertions(+), 96 deletions(-) diff --git a/src/algorithms/calorimetry/CalorimeterHitDigi.cc b/src/algorithms/calorimetry/CalorimeterHitDigi.cc index a46b13dc0e..9813a77aaf 100644 --- a/src/algorithms/calorimetry/CalorimeterHitDigi.cc +++ b/src/algorithms/calorimetry/CalorimeterHitDigi.cc @@ -194,7 +194,7 @@ void CalorimeterHitDigi::process(const CalorimeterHitDigi::Input& input, #endif } if (time > m_cfg.capTime) - continue; + continue; // safety check const double eResRel = @@ -207,14 +207,14 @@ void CalorimeterHitDigi::process(const CalorimeterHitDigi::Input& input, double ped = m_cfg.pedMeanADC + m_gaussian(m_generator) * m_cfg.pedSigmaADC; - //SiPM Saturation - double saturation = 1.0; - if(m_cfg.totalPixel>0){ - double nPhoton = edep * corrMeanScale_value * (1.0 + eResRel) * m_cfg.nPhotonPerGeV; - saturation = (1.0 - exp(-nPhoton/m_cfg.totalPixel))*m_cfg.totalPixel/nPhoton; - //printf("edep=%8.4f nPhoton=%8.2f totalPixel=%8.1f Saturation=%8.6f\n",edep,nPhoton,m_cfg.totalPixel,saturation); - } - + //SiPM Saturation + double saturation = 1.0; + if (m_cfg.totalPixel > 0) { + double nPhoton = edep * corrMeanScale_value * (1.0 + eResRel) * m_cfg.nPhotonPerGeV; + saturation = (1.0 - exp(-nPhoton / m_cfg.totalPixel)) * m_cfg.totalPixel / nPhoton; + //printf("edep=%8.4f nPhoton=%8.2f totalPixel=%8.1f Saturation=%8.6f\n",edep,nPhoton,m_cfg.totalPixel,saturation); + } + // Note: both adc and tdc values must be positive numbers to avoid integer wraparound unsigned long long adc = std::max(std::llround(ped + edep * corrMeanScale_value * (1.0 + eResRel) * saturation / diff --git a/src/algorithms/calorimetry/CalorimeterHitDigiConfig.h b/src/algorithms/calorimetry/CalorimeterHitDigiConfig.h index c3d850e31b..b033892f7b 100644 --- a/src/algorithms/calorimetry/CalorimeterHitDigiConfig.h +++ b/src/algorithms/calorimetry/CalorimeterHitDigiConfig.h @@ -30,7 +30,7 @@ struct CalorimeterHitDigiConfig { std::vector fields{}; //SiPM Saturation - double totalPixel{0}; - double nPhotonPerGeV{0}; + double totalPixel{0}; + double nPhotonPerGeV{0}; }; } // namespace eicrecon diff --git a/src/detectors/FEMC/FEMC.cc b/src/detectors/FEMC/FEMC.cc index 56c2f5e081..f501dc2686 100644 --- a/src/detectors/FEMC/FEMC.cc +++ b/src/detectors/FEMC/FEMC.cc @@ -23,100 +23,107 @@ extern "C" { void InitPlugin(JApplication* app) { using namespace eicrecon; - + InitJANAPlugin(app); auto log_service = app->GetService(); - auto mLog = log_service->logger("FEMC"); - - int SiPMSaturation=1; - std::string opt("ON"); - app->SetDefaultParameter("FEMC:SiPMSaturation",opt,"Turn ON(default) or OFF SiPM Saturation"); - if(opt.find("OFF") != std::string::npos || - opt.find("off") != std::string::npos || - opt.find("Off") != std::string::npos ){ - mLog->info("SiPM Saturation OFF"); - SiPMSaturation=0; - }else{ - mLog->info("SiPM Saturation ON"); - } + auto mLog = log_service->logger("FEMC"); + + int SiPMSaturation = 1; + std::string opt("ON"); + app->SetDefaultParameter("FEMC:SiPMSaturation", opt, "Turn ON(default) or OFF SiPM Saturation"); + if (opt.find("OFF") != std::string::npos || opt.find("off") != std::string::npos || + opt.find("Off") != std::string::npos) { + mLog->info("SiPM Saturation OFF"); + SiPMSaturation = 0; + } else { + mLog->info("SiPM Saturation ON"); + } // Make sure digi and reco use the same value - decltype(CalorimeterHitDigiConfig::capADC) EcalEndcapP_capADC = 16384; //16384, assuming 14 bits. For approximate HGCROC resolution use 65536 - decltype(CalorimeterHitDigiConfig::dyRangeADC) EcalEndcapP_dyRangeADC = 100 * dd4hep::GeV; - decltype(CalorimeterHitDigiConfig::pedMeanADC) EcalEndcapP_pedMeanADC = 200; - decltype(CalorimeterHitDigiConfig::pedSigmaADC) EcalEndcapP_pedSigmaADC = 2.4576; - decltype(CalorimeterHitDigiConfig::resolutionTDC) EcalEndcapP_resolutionTDC = 10 * dd4hep::picosecond; - const double sampFrac=0.03; - decltype(CalorimeterHitDigiConfig::corrMeanScale) EcalEndcapP_corrMeanScale = Form("%f",1.0/sampFrac); - decltype(CalorimeterHitRecoConfig::sampFrac) EcalEndcapP_sampFrac = Form("%f",sampFrac); + decltype(CalorimeterHitDigiConfig::capADC) EcalEndcapP_capADC = + 16384; //16384, assuming 14 bits. For approximate HGCROC resolution use 65536 + decltype(CalorimeterHitDigiConfig::dyRangeADC) EcalEndcapP_dyRangeADC = 100 * dd4hep::GeV; + decltype(CalorimeterHitDigiConfig::pedMeanADC) EcalEndcapP_pedMeanADC = 200; + decltype(CalorimeterHitDigiConfig::pedSigmaADC) EcalEndcapP_pedSigmaADC = 2.4576; + decltype(CalorimeterHitDigiConfig::resolutionTDC) EcalEndcapP_resolutionTDC = + 10 * dd4hep::picosecond; + const double sampFrac = 0.03; + decltype(CalorimeterHitDigiConfig::corrMeanScale) EcalEndcapP_corrMeanScale = + Form("%f", 1.0 / sampFrac); + decltype(CalorimeterHitRecoConfig::sampFrac) EcalEndcapP_sampFrac = Form("%f", sampFrac); decltype(CalorimeterHitDigiConfig::nPhotonPerGeV) EcalEndcapP_nPhotonPerGeV = 1500; - decltype(CalorimeterHitDigiConfig::totalPixel) EcalEndcapP_totalPixel = 4*159565; - if(SiPMSaturation==0) EcalEndcapP_totalPixel=0; - - int fEcalHomoScfi = 0; + decltype(CalorimeterHitDigiConfig::totalPixel) EcalEndcapP_totalPixel = 4 * 159565; + if (SiPMSaturation == 0) + EcalEndcapP_totalPixel = 0; + + int fEcalHomoScfi = 0; try { auto detector = app->GetService()->detector(); - fEcalHomoScfi = detector->constant("ForwardEcal_Homogeneous_Scfi"); - if(fEcalHomoScfi<=1) mLog->info("Homogeneous geometry loaded"); - else mLog->info("ScFi geometry loaded"); - } catch(...){}; - - if(fEcalHomoScfi<=1){ - app->Add(new JOmniFactoryGeneratorT( - "EcalEndcapPRawHits", - {"EcalEndcapPHits"}, + fEcalHomoScfi = detector->constant("ForwardEcal_Homogeneous_Scfi"); + if (fEcalHomoScfi <= 1) + mLog->info("Homogeneous geometry loaded"); + else + mLog->info("ScFi geometry loaded"); + } catch (...) { + }; + + if (fEcalHomoScfi <= 1) { + app->Add(new JOmniFactoryGeneratorT( + "EcalEndcapPRawHits", {"EcalEndcapPHits"}, #if EDM4EIC_VERSION_MAJOR >= 7 - {"EcalEndcapPRawHits", "EcalEndcapPRawHitAssociations"}, + {"EcalEndcapPRawHits", "EcalEndcapPRawHitAssociations"}, #else - {"EcalEndcapPRawHits"}, + {"EcalEndcapPRawHits"}, #endif - { - .eRes = {0.11333 * sqrt(dd4hep::GeV), 0.03, 0.0 * dd4hep::GeV}, // (11.333% / sqrt(E)) \oplus 3% - .tRes = 0.0, - .threshold = 0.0, // 15MeV threshold for a single tower will be applied on ADC at Reco below - .capADC = EcalEndcapP_capADC, - .capTime = 100, // given in ns, 4 samples in HGCROC - .dyRangeADC = EcalEndcapP_dyRangeADC, - .pedMeanADC = EcalEndcapP_pedMeanADC, - .pedSigmaADC = EcalEndcapP_pedSigmaADC, - .resolutionTDC = EcalEndcapP_resolutionTDC, - .corrMeanScale = "1.0", - .readout = "EcalEndcapPHits", - .totalPixel = EcalEndcapP_totalPixel, - .nPhotonPerGeV = EcalEndcapP_nPhotonPerGeV, - }, - app // TODO: Remove me once fixed - )); - }else if(fEcalHomoScfi==2){ - app->Add(new JOmniFactoryGeneratorT( - "EcalEndcapPRawHits", - {"EcalEndcapPHits"}, + { + .eRes = {0.11333 * sqrt(dd4hep::GeV), 0.03, + 0.0 * dd4hep::GeV}, // (11.333% / sqrt(E)) \oplus 3% + .tRes = 0.0, + .threshold = + 0.0, // 15MeV threshold for a single tower will be applied on ADC at Reco below + .capADC = EcalEndcapP_capADC, + .capTime = 100, // given in ns, 4 samples in HGCROC + .dyRangeADC = EcalEndcapP_dyRangeADC, + .pedMeanADC = EcalEndcapP_pedMeanADC, + .pedSigmaADC = EcalEndcapP_pedSigmaADC, + .resolutionTDC = EcalEndcapP_resolutionTDC, + .corrMeanScale = "1.0", + .readout = "EcalEndcapPHits", + .totalPixel = EcalEndcapP_totalPixel, + .nPhotonPerGeV = EcalEndcapP_nPhotonPerGeV, + }, + app // TODO: Remove me once fixed + )); + } else if (fEcalHomoScfi == 2) { + app->Add(new JOmniFactoryGeneratorT( + "EcalEndcapPRawHits", {"EcalEndcapPHits"}, #if EDM4EIC_VERSION_MAJOR >= 7 - {"EcalEndcapPRawHits", "EcalEndcapPRawHitAssociations"}, + {"EcalEndcapPRawHits", "EcalEndcapPRawHitAssociations"}, #else - {"EcalEndcapPRawHits"}, + {"EcalEndcapPRawHits"}, #endif - { - .eRes = {0.0,0.0,0.0}, - .tRes = 0.0, - .threshold = 0.0, // 15MeV threshold for a single tower will be applied on ADC at Reco below - .capADC = EcalEndcapP_capADC, - .capTime = 100, // given in ns, 4 samples in HGCROC - .dyRangeADC = EcalEndcapP_dyRangeADC, - .pedMeanADC = EcalEndcapP_pedMeanADC, - .pedSigmaADC = EcalEndcapP_pedSigmaADC, - .resolutionTDC = EcalEndcapP_resolutionTDC, - .corrMeanScale = EcalEndcapP_corrMeanScale, - .readout = "EcalEndcapPHits", - .fields = {"fiberx","fibery","x","y"}, - .totalPixel = EcalEndcapP_totalPixel, - .nPhotonPerGeV = EcalEndcapP_nPhotonPerGeV, - }, - app // TODO: Remove me once fixed - )); - } - + { + .eRes = {0.0, 0.0, 0.0}, + .tRes = 0.0, + .threshold = + 0.0, // 15MeV threshold for a single tower will be applied on ADC at Reco below + .capADC = EcalEndcapP_capADC, + .capTime = 100, // given in ns, 4 samples in HGCROC + .dyRangeADC = EcalEndcapP_dyRangeADC, + .pedMeanADC = EcalEndcapP_pedMeanADC, + .pedSigmaADC = EcalEndcapP_pedSigmaADC, + .resolutionTDC = EcalEndcapP_resolutionTDC, + .corrMeanScale = EcalEndcapP_corrMeanScale, + .readout = "EcalEndcapPHits", + .fields = {"fiberx", "fibery", "x", "y"}, + .totalPixel = EcalEndcapP_totalPixel, + .nPhotonPerGeV = EcalEndcapP_nPhotonPerGeV, + }, + app // TODO: Remove me once fixed + )); + } + app->Add(new JOmniFactoryGeneratorT( "EcalEndcapPRecHits", {"EcalEndcapPRawHits"}, {"EcalEndcapPRecHits"}, { @@ -126,8 +133,9 @@ void InitPlugin(JApplication* app) { .pedSigmaADC = EcalEndcapP_pedSigmaADC, .resolutionTDC = EcalEndcapP_resolutionTDC, .thresholdFactor = 0.0, - .thresholdValue = 2, // The ADC of a 15 MeV particle is adc = 200 + 15 * 0.03 * ( 1.0 + 0) / 3000 * 16384 = 200 + 2.4576 - .sampFrac = "1.00", // alerady taken care in DIGI code above + .thresholdValue = + 2, // The ADC of a 15 MeV particle is adc = 200 + 15 * 0.03 * ( 1.0 + 0) / 3000 * 16384 = 200 + 2.4576 + .sampFrac = "1.00", // alerady taken care in DIGI code above .readout = "EcalEndcapPHits", }, app // TODO: Remove me once fixed @@ -239,5 +247,4 @@ void InitPlugin(JApplication* app) { "EcalEndcapPSplitMergeClusterAssociationsWithoutShapes"}, {"EcalEndcapPSplitMergeClusters", "EcalEndcapPSplitMergeClusterAssociations"}, {.energyWeight = "log", .logWeightBase = 3.6}, app)); - } diff --git a/src/factories/calorimetry/CalorimeterHitDigi_factory.h b/src/factories/calorimetry/CalorimeterHitDigi_factory.h index 839d487c2f..533069547d 100644 --- a/src/factories/calorimetry/CalorimeterHitDigi_factory.h +++ b/src/factories/calorimetry/CalorimeterHitDigi_factory.h @@ -36,8 +36,8 @@ class CalorimeterHitDigi_factory ParameterRef m_corrMeanScale{this, "scaleResponse", config().corrMeanScale}; ParameterRef> m_fields{this, "signalSumFields", config().fields}; ParameterRef m_readout{this, "readoutClass", config().readout}; - ParameterRef m_nPhotonPerGeV {this, "numberOfPhotonPerGeV", config().nPhotonPerGeV}; - ParameterRef m_totalPixel {this, "totalNnumberOfPixel", config().totalPixel}; + ParameterRef m_nPhotonPerGeV{this, "numberOfPhotonPerGeV", config().nPhotonPerGeV}; + ParameterRef m_totalPixel{this, "totalNnumberOfPixel", config().totalPixel}; Service m_algorithmsInit{this}; From 581c4b447a390273ab66a58415c1774594df9490 Mon Sep 17 00:00:00 2001 From: akioogawa Date: Fri, 2 May 2025 14:13:07 -0400 Subject: [PATCH 3/6] typo --- src/detectors/FEMC/FEMC.cc | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/detectors/FEMC/FEMC.cc b/src/detectors/FEMC/FEMC.cc index f501dc2686..7fe801513d 100644 --- a/src/detectors/FEMC/FEMC.cc +++ b/src/detectors/FEMC/FEMC.cc @@ -135,7 +135,7 @@ void InitPlugin(JApplication* app) { .thresholdFactor = 0.0, .thresholdValue = 2, // The ADC of a 15 MeV particle is adc = 200 + 15 * 0.03 * ( 1.0 + 0) / 3000 * 16384 = 200 + 2.4576 - .sampFrac = "1.00", // alerady taken care in DIGI code above + .sampFrac = "1.00", // already taken care in DIGI code above .readout = "EcalEndcapPHits", }, app // TODO: Remove me once fixed @@ -248,3 +248,4 @@ void InitPlugin(JApplication* app) { {"EcalEndcapPSplitMergeClusters", "EcalEndcapPSplitMergeClusterAssociations"}, {.energyWeight = "log", .logWeightBase = 3.6}, app)); } +} From 3096e5e2988d9a54d1cb20446302fcf2e65a0d29 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 8 May 2025 13:06:10 +0000 Subject: [PATCH 4/6] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- src/detectors/FEMC/FEMC.cc | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/detectors/FEMC/FEMC.cc b/src/detectors/FEMC/FEMC.cc index 6c0461190c..18cab9a680 100644 --- a/src/detectors/FEMC/FEMC.cc +++ b/src/detectors/FEMC/FEMC.cc @@ -254,7 +254,7 @@ void InitPlugin(JApplication* app) { "EcalEndcapPSplitMergeClusterAssociationsWithoutShapes"}, {"EcalEndcapPSplitMergeClusters", "EcalEndcapPSplitMergeClusterAssociations"}, {.energyWeight = "log", .logWeightBase = 3.6}, app)); - + // Insert is identical to regular Ecal app->Add(new JOmniFactoryGeneratorT( "EcalEndcapPInsertRawHits", {"EcalEndcapPInsertHits"}, From 10bd33890f4f04a0ab7e3f578235a4f3535c50cc Mon Sep 17 00:00:00 2001 From: akioogawa Date: Tue, 13 May 2025 11:52:24 -0400 Subject: [PATCH 5/6] - Update ScFi model's Sampling Faction from 3% to 2.9% based on ratio to Homogeneous model (see https://www.star.bnl.gov/~akio/epic/reco/index.html) - Adding ScFi model's resolution smearing of 2.2% constant term based on high energy end (see https://www.star.bnl.gov/~akio/epic/reco/index.html) - Changing ADC threshold from 2 (~9MeV which was too low) to 3 (~15MeV as intended) - Rename an unfortunate variable name --- src/detectors/FEMC/FEMC.cc | 19 +++++++++++-------- 1 file changed, 11 insertions(+), 8 deletions(-) diff --git a/src/detectors/FEMC/FEMC.cc b/src/detectors/FEMC/FEMC.cc index 18cab9a680..3ee3749efd 100644 --- a/src/detectors/FEMC/FEMC.cc +++ b/src/detectors/FEMC/FEMC.cc @@ -48,7 +48,8 @@ void InitPlugin(JApplication* app) { decltype(CalorimeterHitDigiConfig::pedSigmaADC) EcalEndcapP_pedSigmaADC = 2.4576; decltype(CalorimeterHitDigiConfig::resolutionTDC) EcalEndcapP_resolutionTDC = 10 * dd4hep::picosecond; - const double sampFrac = 0.03; + //const double sampFrac = 0.03; + const double sampFrac = 0.029043; //updated with ratio to ScFi model decltype(CalorimeterHitDigiConfig::corrMeanScale) EcalEndcapP_corrMeanScale = Form("%f", 1.0 / sampFrac); decltype(CalorimeterHitRecoConfig::sampFrac) EcalEndcapP_sampFrac = Form("%f", sampFrac); @@ -57,18 +58,18 @@ void InitPlugin(JApplication* app) { if (SiPMSaturation == 0) EcalEndcapP_totalPixel = 0; - int fEcalHomoScfi = 0; + int FEMCHomoScfi = 0; try { auto detector = app->GetService()->detector(); - fEcalHomoScfi = detector->constant("ForwardEcal_Homogeneous_Scfi"); - if (fEcalHomoScfi <= 1) + FEMCHomoScfi = detector->constant("ForwardEcal_Homogeneous_Scfi"); + if (FEMCHomoScfi <= 1) mLog->info("Homogeneous geometry loaded"); else mLog->info("ScFi geometry loaded"); } catch (...) { }; - if (fEcalHomoScfi <= 1) { + if (FEMCHomoScfi <= 1) { app->Add(new JOmniFactoryGeneratorT( "EcalEndcapPRawHits", {"EcalEndcapPHits"}, #if EDM4EIC_VERSION_MAJOR >= 7 @@ -95,7 +96,7 @@ void InitPlugin(JApplication* app) { }, app // TODO: Remove me once fixed )); - } else if (fEcalHomoScfi == 2) { + } else if (FEMCHomoScfi == 2) { app->Add(new JOmniFactoryGeneratorT( "EcalEndcapPRawHits", {"EcalEndcapPHits"}, #if EDM4EIC_VERSION_MAJOR >= 7 @@ -104,7 +105,8 @@ void InitPlugin(JApplication* app) { {"EcalEndcapPRawHits"}, #endif { - .eRes = {0.0, 0.0, 0.0}, + //.eRes = {0.0, 0.0, 0.0}, // No smearing for ScFi + .eRes = {0.0, 0.022, 0.0}, // just constant term 2.2% based on MC data comparison .tRes = 0.0, .threshold = 0.0, // 15MeV threshold for a single tower will be applied on ADC at Reco below @@ -134,7 +136,8 @@ void InitPlugin(JApplication* app) { .resolutionTDC = EcalEndcapP_resolutionTDC, .thresholdFactor = 0.0, .thresholdValue = - 2, // The ADC of a 15 MeV particle is adc = 200 + 15 * 0.03 * ( 1.0 + 0) / 3000 * 16384 = 200 + 2.4576 + // 2, // The ADC of a 15 MeV particle is adc = 200 + 15 * 0.03 * ( 1.0 + 0) / 3000 * 16384 = 200 + 2.4576 + 3, // 15 MeV = 2.4576, but adc=llround(dE) and cut off is "<". So 3 here = 15.25MeV .sampFrac = "1.00", // already taken care in DIGI code above .readout = "EcalEndcapPHits", }, From 1b7f3b70af67fde285ea1ba22b7b0a3bbcc82e1e Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Tue, 13 May 2025 15:59:45 +0000 Subject: [PATCH 6/6] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- src/detectors/FEMC/FEMC.cc | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/src/detectors/FEMC/FEMC.cc b/src/detectors/FEMC/FEMC.cc index 3ee3749efd..2c6b5d3cd0 100644 --- a/src/detectors/FEMC/FEMC.cc +++ b/src/detectors/FEMC/FEMC.cc @@ -48,8 +48,8 @@ void InitPlugin(JApplication* app) { decltype(CalorimeterHitDigiConfig::pedSigmaADC) EcalEndcapP_pedSigmaADC = 2.4576; decltype(CalorimeterHitDigiConfig::resolutionTDC) EcalEndcapP_resolutionTDC = 10 * dd4hep::picosecond; - //const double sampFrac = 0.03; - const double sampFrac = 0.029043; //updated with ratio to ScFi model + //const double sampFrac = 0.03; + const double sampFrac = 0.029043; //updated with ratio to ScFi model decltype(CalorimeterHitDigiConfig::corrMeanScale) EcalEndcapP_corrMeanScale = Form("%f", 1.0 / sampFrac); decltype(CalorimeterHitRecoConfig::sampFrac) EcalEndcapP_sampFrac = Form("%f", sampFrac); @@ -61,7 +61,7 @@ void InitPlugin(JApplication* app) { int FEMCHomoScfi = 0; try { auto detector = app->GetService()->detector(); - FEMCHomoScfi = detector->constant("ForwardEcal_Homogeneous_Scfi"); + FEMCHomoScfi = detector->constant("ForwardEcal_Homogeneous_Scfi"); if (FEMCHomoScfi <= 1) mLog->info("Homogeneous geometry loaded"); else @@ -105,8 +105,8 @@ void InitPlugin(JApplication* app) { {"EcalEndcapPRawHits"}, #endif { - //.eRes = {0.0, 0.0, 0.0}, // No smearing for ScFi - .eRes = {0.0, 0.022, 0.0}, // just constant term 2.2% based on MC data comparison + //.eRes = {0.0, 0.0, 0.0}, // No smearing for ScFi + .eRes = {0.0, 0.022, 0.0}, // just constant term 2.2% based on MC data comparison .tRes = 0.0, .threshold = 0.0, // 15MeV threshold for a single tower will be applied on ADC at Reco below @@ -136,8 +136,8 @@ void InitPlugin(JApplication* app) { .resolutionTDC = EcalEndcapP_resolutionTDC, .thresholdFactor = 0.0, .thresholdValue = - // 2, // The ADC of a 15 MeV particle is adc = 200 + 15 * 0.03 * ( 1.0 + 0) / 3000 * 16384 = 200 + 2.4576 - 3, // 15 MeV = 2.4576, but adc=llround(dE) and cut off is "<". So 3 here = 15.25MeV + // 2, // The ADC of a 15 MeV particle is adc = 200 + 15 * 0.03 * ( 1.0 + 0) / 3000 * 16384 = 200 + 2.4576 + 3, // 15 MeV = 2.4576, but adc=llround(dE) and cut off is "<". So 3 here = 15.25MeV .sampFrac = "1.00", // already taken care in DIGI code above .readout = "EcalEndcapPHits", },