From af4c4fe9ac12e33bc9fe85e0a87e95b0130abee8 Mon Sep 17 00:00:00 2001
From: Volker Friese <v.friese@gsi.de>
Date: Fri, 4 Dec 2020 20:02:26 +0100
Subject: [PATCH] Enable global deactivation of random channels in the STS
simulation.
---
core/detectors/sts/CbmStsParAsic.cxx | 58 ++++++++++++++++++-----
core/detectors/sts/CbmStsParAsic.h | 24 ++++++++--
core/detectors/sts/CbmStsParModule.cxx | 15 +++++-
core/detectors/sts/CbmStsParModule.h | 7 +++
core/detectors/sts/CbmStsParSetModule.cxx | 26 +++++++++-
core/detectors/sts/CbmStsParSetModule.h | 15 +++++-
sim/detectors/sts/CbmStsDigitize.cxx | 37 +++++++++------
sim/detectors/sts/CbmStsDigitize.h | 16 ++++++-
8 files changed, 165 insertions(+), 33 deletions(-)
diff --git a/core/detectors/sts/CbmStsParAsic.cxx b/core/detectors/sts/CbmStsParAsic.cxx
index b45d072377..75e07c4945 100644
--- a/core/detectors/sts/CbmStsParAsic.cxx
+++ b/core/detectors/sts/CbmStsParAsic.cxx
@@ -6,28 +6,31 @@
#include <TF1.h> // for TF1
#include <TMath.h> // for Exp
+#include <TRandom.h>
#include <assert.h> // for assert
#include <sstream> // for operator<<, basic_ostream, stringstream
ClassImp(CbmStsParAsic)
- // ----- Constructor ---------------------------------------------------
- CbmStsParAsic::CbmStsParAsic(UShort_t nAdc,
- Double_t dynRange,
- Double_t threshold,
- Double_t timeResol,
- Double_t deadTime,
- Double_t noise,
- Double_t znr) {
- Set(nAdc, dynRange, threshold, timeResol, deadTime, noise, znr);
+// ----- Constructor ---------------------------------------------------
+CbmStsParAsic::CbmStsParAsic(UShort_t nChannels,
+ UShort_t nAdc,
+ Double_t dynRange,
+ Double_t threshold,
+ Double_t timeResol,
+ Double_t deadTime,
+ Double_t noise,
+ Double_t znr) {
+ Set(nChannels, nAdc, dynRange, threshold, timeResol, deadTime, noise, znr);
}
// -------------------------------------------------------------------------
// ----- Copy constructor ----------------------------------------------
CbmStsParAsic::CbmStsParAsic(const CbmStsParAsic& other) {
- Set(other.fNofAdc,
+ Set(other.fNofChannels,
+ other.fNofAdc,
other.fDynRange,
other.fThreshold,
other.fTimeResolution,
@@ -40,7 +43,8 @@ CbmStsParAsic::CbmStsParAsic(const CbmStsParAsic& other) {
// ----- Copy assignment operator --------------------------------------
CbmStsParAsic& CbmStsParAsic::operator=(const CbmStsParAsic& other) {
- Set(other.fNofAdc,
+ Set(other.fNofChannels,
+ other.fNofAdc,
other.fDynRange,
other.fThreshold,
other.fTimeResolution,
@@ -68,6 +72,34 @@ Short_t CbmStsParAsic::ChargeToAdc(Double_t charge) const {
// -------------------------------------------------------------------------
+// ----- Deactivate channels -------------------------------------------
+UInt_t CbmStsParAsic::DeactivateRandomChannels(Double_t fraction) {
+
+ if ( fraction < 0. ) return 0;
+
+ // --- Average number of dead channels
+ Double_t meanDead = fraction * Double_t(fNofChannels);
+ if ( meanDead > fNofChannels) meanDead = fNofChannels;
+
+ // --- Sample actual number of dead channels from Poissonian
+ Int_t nDead = gRandom->Poisson(meanDead);
+
+ // --- Deactivate the given number of channels
+ Int_t nDeactivated = 0;
+ while(nDeactivated < nDead) {
+ Int_t channel = Int_t(gRandom->Uniform(0, fNofChannels));
+ if ( IsChannelActive(channel) ) {
+ fDeadChannels.insert(channel);
+ nDeactivated++;
+ } //? Channel was active
+ } //# Deactivated channels
+
+ assert(nDeactivated == nDead);
+ return nDead;
+}
+// -------------------------------------------------------------------------
+
+
// ----- Single-channel noise rate -------------------------------------
Double_t CbmStsParAsic::GetNoiseRate() const {
if (fNoise == 0.) return 0.;
@@ -97,7 +129,8 @@ void CbmStsParAsic::Init() {
// ----- Set the parameters ---------------------------------------------
-void CbmStsParAsic::Set(UShort_t nAdc,
+void CbmStsParAsic::Set(UShort_t nChannels,
+ UShort_t nAdc,
Double_t dynRange,
Double_t threshold,
Double_t timeResol,
@@ -114,6 +147,7 @@ void CbmStsParAsic::Set(UShort_t nAdc,
assert(noise >= 0.);
assert(zeroNoiseRate >= 0.);
+ fNofChannels = nChannels;
fNofAdc = nAdc;
fDynRange = dynRange;
fThreshold = threshold;
diff --git a/core/detectors/sts/CbmStsParAsic.h b/core/detectors/sts/CbmStsParAsic.h
index 5678bb3df2..79977e6eb5 100644
--- a/core/detectors/sts/CbmStsParAsic.h
+++ b/core/detectors/sts/CbmStsParAsic.h
@@ -30,6 +30,7 @@ public:
/** @brief Constructor with parameters
+ ** @param nChannels Number of readout channels
** @param nAdc Number of ADC channels
** @param dynRange Dynamic range of ADC [e]
** @param threshold ADC threshold [e]
@@ -38,7 +39,8 @@ public:
** @param noise Noise RMS [e]
** @param znr Zero-crossing noise rate [1/ns]
**/
- CbmStsParAsic(UShort_t nAdc,
+ CbmStsParAsic(UShort_t nChannels,
+ UShort_t nAdc,
Double_t dynRange,
Double_t threshold,
Double_t timeResol,
@@ -76,6 +78,13 @@ public:
}
+ /** @brief Randomly deactivate a fraction of the channels
+ ** @param fraction Fraction of channels to deactivate
+ ** @return Number of deactivated channels
+ **/
+ UInt_t DeactivateRandomChannels(Double_t fraction);
+
+
/** @brief ADC channel for a given charge
** @param charge Charge [e]
** @return ADC channel number
@@ -103,6 +112,12 @@ public:
UShort_t GetNofAdc() const { return fNofAdc; }
+ /** @brief Number of readout channels
+ ** @return Number of readout channels
+ **/
+ UShort_t GetNofChannels() const { return fNofChannels; }
+
+
/** @brief Electronic noise RMS
** @return Noise RMS [e]
**/
@@ -160,6 +175,7 @@ public:
/** @brief Set parameters
+ ** @param nChannels Number of readout channels
** @param nAdc Number of ADC channels
** @param dynRange Dynamic range [e]
** @param threshold Threshold [e]
@@ -169,7 +185,8 @@ public:
** @param zeroNoiseRate Zero-crossing noise rate
** @param deadChannels Set of dead channels
**/
- void Set(UShort_t nAdc,
+ void Set(UShort_t nChannels,
+ UShort_t nAdc,
Double_t dynRange,
Double_t threshold,
Double_t timeResol,
@@ -184,6 +201,7 @@ public:
private:
+ UShort_t fNofChannels = 0; ///< Number of readout channels
UShort_t fNofAdc = 0; ///< Number of ADC channels
Double_t fDynRange = 0.; ///< Dynamic range [e]
Double_t fThreshold = 0.; ///< Threshold [e]
@@ -199,7 +217,7 @@ private:
** method in order to avoid frequent re-calculation. **/
TF1* fNoiseCharge = nullptr; //!
- ClassDefNV(CbmStsParAsic, 2);
+ ClassDefNV(CbmStsParAsic, 3);
};
#endif /* CBMSTSPARASIC_H */
diff --git a/core/detectors/sts/CbmStsParModule.cxx b/core/detectors/sts/CbmStsParModule.cxx
index f15b92cc99..03340beebf 100644
--- a/core/detectors/sts/CbmStsParModule.cxx
+++ b/core/detectors/sts/CbmStsParModule.cxx
@@ -25,6 +25,18 @@ ClassImp(CbmStsParModule)
// -------------------------------------------------------------------------
+// ----- Randomly deactivate channels ----------------------------------
+UInt_t CbmStsParModule::DeactivateRandomChannels(Double_t fraction) {
+ if ( fraction <= 0. ) return 0;
+ UInt_t nDeactivated = 0;
+ for ( auto& asic : fAsicPars ) {
+ nDeactivated += asic.DeactivateRandomChannels(fraction);
+ }
+ return nDeactivated;
+}
+// -------------------------------------------------------------------------
+
+
// ----- Get ASIC parameters -------------------------------------------
const CbmStsParAsic& CbmStsParModule::GetParAsic(UInt_t channel) const {
assert(!fAsicPars.empty());
@@ -36,7 +48,7 @@ const CbmStsParAsic& CbmStsParModule::GetParAsic(UInt_t channel) const {
// -------------------------------------------------------------------------
-// ----- Check for a chennel being active ------------------------------
+// ----- Check for a channel being active ------------------------------
Bool_t CbmStsParModule::IsChannelActive(UInt_t channel) const {
const CbmStsParAsic& parAsic = GetParAsic(channel);
UShort_t asicChannel = channel % fNofAsicChannels;
@@ -47,6 +59,7 @@ Bool_t CbmStsParModule::IsChannelActive(UInt_t channel) const {
// ----- Set parameters for all ASICs ----------------------------------
void CbmStsParModule::SetAllAsics(const CbmStsParAsic& asicPar) {
+ assert(asicPar.GetNofChannels() == fNofAsicChannels);
for (auto& par : fAsicPars) {
par = asicPar;
par.Init();
diff --git a/core/detectors/sts/CbmStsParModule.h b/core/detectors/sts/CbmStsParModule.h
index 850a159acb..29a666e1f1 100644
--- a/core/detectors/sts/CbmStsParModule.h
+++ b/core/detectors/sts/CbmStsParModule.h
@@ -43,6 +43,13 @@ public:
~CbmStsParModule() {};
+ /** @brief Randomly deactivate a fraction of the channels
+ ** @param fraction Fraction of channels to deactivate
+ ** @return Number of deactivated channels
+ **/
+ UInt_t DeactivateRandomChannels(Double_t fraction);
+
+
/** @brief ASIC parameters for a given channel
** @param channel Channel number
** @return ASIC parameters
diff --git a/core/detectors/sts/CbmStsParSetModule.cxx b/core/detectors/sts/CbmStsParSetModule.cxx
index f4d1289e29..d370dfd31e 100644
--- a/core/detectors/sts/CbmStsParSetModule.cxx
+++ b/core/detectors/sts/CbmStsParSetModule.cxx
@@ -12,6 +12,7 @@
#include <sstream> // for operator<<, basic_ostream, stringstream
#include <string> // for char_traits
+
ClassImp(CbmStsParSetModule)
// ----- Constructor ----------------------------------------------------
@@ -37,6 +38,18 @@ void CbmStsParSetModule::clear() {
// --------------------------------------------------------------------------
+// ----- Randomly deactivate channels -----------------------------------
+UInt_t CbmStsParSetModule::DeactivateRandomChannels(Double_t fraction) {
+ if ( fraction <= 0. ) return 0;
+ UInt_t nDeactivated = 0;
+ for ( auto& entry : fParams ) {
+ nDeactivated += entry.second.DeactivateRandomChannels(fraction);
+ }
+ return nDeactivated;
+}
+// --------------------------------------------------------------------------
+
+
// ----- Read parameters from ASCII file --------------------------------
Bool_t CbmStsParSetModule::getParams(FairParamList*) {
LOG(fatal) << GetName() << ": ASCII input is not defined!";
@@ -45,7 +58,7 @@ Bool_t CbmStsParSetModule::getParams(FairParamList*) {
// --------------------------------------------------------------------------
-// ----- Get condition parameters of a sensor ---------------------------
+// ----- Get condition parameters of a module ---------------------------
const CbmStsParModule& CbmStsParSetModule::GetParModule(UInt_t address) {
if (fUseGlobal) return fGlobalParams;
assert(fParams.count(address));
@@ -61,6 +74,17 @@ void CbmStsParSetModule::putParams(FairParamList*) {
// --------------------------------------------------------------------------
+// ----- Set module parameters ------------------------------------------
+void CbmStsParSetModule::SetParModule(UInt_t address,
+ const CbmStsParModule& par) {
+ if (fParams.count(address))
+ LOG(fatal) << GetName() << ": Replacing parameters for sensor address "
+ << address;
+ fParams[address] = par;
+}
+// --------------------------------------------------------------------------
+
+
// ----- Info to string ------------------------------------------------
std::string CbmStsParSetModule::ToString() const {
std::stringstream ss;
diff --git a/core/detectors/sts/CbmStsParSetModule.h b/core/detectors/sts/CbmStsParSetModule.h
index 987f228377..47054ec4f0 100644
--- a/core/detectors/sts/CbmStsParSetModule.h
+++ b/core/detectors/sts/CbmStsParSetModule.h
@@ -56,6 +56,13 @@ public:
virtual void clear();
+ /** @brief Randomly deactivate a fraction of the channels
+ ** @param fraction Fraction of channels to deactivate
+ ** @return Number of deactivated channels
+ **/
+ UInt_t DeactivateRandomChannels(Double_t fraction);
+
+
/** @brief Reading parameters from ASCII. Abstract in base class.
**
** An ASCII I/O is not implemented. The method throws an error.
@@ -102,6 +109,12 @@ public:
fUseGlobal = kTRUE;
}
+ /** @brief Set the parameters for a module
+ ** @param address Module address
+ ** @param par Module parameter object
+ **/
+ void SetParModule(UInt_t address, const CbmStsParModule& par);
+
/** @brief Info to string **/
std::string ToString() const;
@@ -114,7 +127,7 @@ private:
/** @brief Global parameters, used for all modules **/
CbmStsParModule fGlobalParams {};
- /** @brief Map of parameters. Key is sensor address. **/
+ /** @brief Map of parameters. Key is module address. **/
std::map<UInt_t, CbmStsParModule> fParams {};
diff --git a/sim/detectors/sts/CbmStsDigitize.cxx b/sim/detectors/sts/CbmStsDigitize.cxx
index d3956d2c10..7554c2bf33 100644
--- a/sim/detectors/sts/CbmStsDigitize.cxx
+++ b/sim/detectors/sts/CbmStsDigitize.cxx
@@ -438,12 +438,20 @@ void CbmStsDigitize::InitParams() {
assert(fUserParModule);
assert(fUserParAsic);
fUserParModule->SetAllAsics(*fUserParAsic);
- fParSetModule->SetGlobalPar(*fUserParModule);
+ for (Int_t iModule = 0; iModule < fSetup->GetNofModules(); iModule++) {
+ UInt_t address = fSetup->GetModule(iModule)->GetAddress();
+ fParSetModule->SetParModule(address, *fUserParModule);
+ }
+ UInt_t deactivated = 0;
+ if ( fUserFracDeadChan > 0. ) {
+ deactivated = fParSetModule->DeactivateRandomChannels(fUserFracDeadChan);
+ }
fParSetModule->setChanged();
fParSetModule->setInputVersion(-2, 1);
LOG(info) << "--- Using global ASIC parameters: \n "
<< fUserParAsic->ToString();
LOG(info) << "--- Module parameters: " << fParSetModule->ToString();
+ LOG(info) << "--- Deactive channels: " << deactivated << " " << fUserFracDeadChan;
// --- Sensor parameters
// TODO: The code above is highly unsatisfactory. A better implementation
@@ -714,24 +722,23 @@ void CbmStsDigitize::SetGlobalDefaults() {
// through the base class CbmDigitizeBase, such that the settings here
// will be overwritten later (in Init).
+ // --- Module parameters
+ if (fUserParModule) delete fUserParModule;
+ UInt_t nChannels = 2048; // Number of module readout channels
+ UInt_t nAsicChannels = 128; // Number of readout channels per ASIC
+ fUserParModule = new CbmStsParModule(nChannels, nAsicChannels);
+
// --- ASIC parameters
if (fUserParAsic) delete fUserParAsic;
- UInt_t nAdc = 32; // Number of ADC channels (5 bit)
+ UShort_t nAdc = 32; // Number of ADC channels (5 bit)
Double_t dynRange = 75000.; // Dynamic range [e]
Double_t threshold = 3000.; // Threshold [e]
Double_t timeResol = 5.; // Time resolution [ns]
Double_t deadTime = 800.; // Channel dead time [ns]
Double_t noiseRms = 1000.; // RMS of noise [e]
Double_t znr = 3.9789e-3; // Zero-crossing noise rate [1/ns]
- fUserParAsic = new CbmStsParAsic(
- nAdc, dynRange, threshold, timeResol, deadTime, noiseRms, znr);
-
- // --- Module parameters
- if (fUserParModule) delete fUserParModule;
- UInt_t nChannels = 2048; // Number of module readout channels
- UInt_t nAsicChannels = 128; // Number of readout channels per ASIC
- fUserParModule = new CbmStsParModule(nChannels, nAsicChannels);
-
+ fUserParAsic = new CbmStsParAsic(nAsicChannels, nAdc, dynRange, threshold,
+ timeResol, deadTime, noiseRms, znr);
// --- Sensor parameters
// --- Here, only the default pitch and stereo angles are defined. The
// --- other parameters are extracted from the geometry.
@@ -762,7 +769,8 @@ void CbmStsDigitize::SetGlobalDefaults() {
// ----- Set the global module parameters ------------------------------
-void CbmStsDigitize::SetGlobalAsicParams(UInt_t nAdc,
+void CbmStsDigitize::SetGlobalAsicParams(UShort_t nChannels,
+ UShort_t nAdc,
Double_t dynRange,
Double_t threshold,
Double_t timeResolution,
@@ -772,8 +780,9 @@ void CbmStsDigitize::SetGlobalAsicParams(UInt_t nAdc,
assert(!fIsInitialised);
assert(nAdc > 0);
if (fUserParAsic) delete fUserParAsic;
- fUserParAsic = new CbmStsParAsic(
- nAdc, dynRange, threshold, timeResolution, deadTime, noise, zeroNoiseRate);
+ fUserParAsic = new CbmStsParAsic(nChannels, nAdc, dynRange, threshold,
+ timeResolution, deadTime, noise,
+ zeroNoiseRate);
}
// -------------------------------------------------------------------------
diff --git a/sim/detectors/sts/CbmStsDigitize.h b/sim/detectors/sts/CbmStsDigitize.h
index f4b96e90a3..f8ce957de1 100644
--- a/sim/detectors/sts/CbmStsDigitize.h
+++ b/sim/detectors/sts/CbmStsDigitize.h
@@ -116,6 +116,7 @@ public:
/** @brief Set the global ASIC parameters
+ ** @param nChannels Number of readout channels
** @param nAdc Number of ADC channels
** @param dynRange Dynamic range [e]
** @param threshold Threshold [e]
@@ -126,7 +127,8 @@ public:
**
** These parameters will be applied to all ASICS in all modules.
**/
- void SetGlobalAsicParams(UInt_t nAdc,
+ void SetGlobalAsicParams(UShort_t nChannels,
+ UShort_t nAdc,
Double_t dynRange,
Double_t threshold,
Double_t timeResolution,
@@ -135,6 +137,17 @@ public:
Double_t zeroNoiseRate);
+ /** @brief Set global fraction of dead channels
+ ** @param fraction Fraction of dead channels
+ **
+ ** If this number is different from zero, in each ASIC a number of
+ ** channels corresponding to this fraction are deactivated.
+ **/
+ void SetGlobalFracDeadChannels(Double_t fraction) {
+ fUserFracDeadChan = fraction;
+ }
+
+
/** @brief Set the global module parameters
** @param nChannels Number of readout channels
** @param nAsicChannels Number of readout channels per ASIC
@@ -253,6 +266,7 @@ private:
CbmStsParSensor* fUserParSensor = nullptr; ///< User defined, global
CbmStsParSensorCond* fUserParCond = nullptr; ///< User defined, global
Double_t fUserDinactive = 0.; ///< Size of inactive sensor border [cm]
+ Double_t fUserFracDeadChan = 0.; ///< Fraction of inactive ASIC channels
// --- Module and sensor parameters for runtime DB output
CbmStsParSim* fParSim = nullptr; ///< Simulation settings
--
GitLab