diff --git a/SimCalorimetry/HGCalSimProducers/interface/HGCDigitizerBase.h b/SimCalorimetry/HGCalSimProducers/interface/HGCDigitizerBase.h index 0a540bef0e402..cabd4d315a0c9 100644 --- a/SimCalorimetry/HGCalSimProducers/interface/HGCDigitizerBase.h +++ b/SimCalorimetry/HGCalSimProducers/interface/HGCDigitizerBase.h @@ -68,6 +68,10 @@ class HGCDigitizerBase { @short CTOR */ HGCDigitizerBase(const edm::ParameterSet& ps); + /** + @short Gaussian Noise Generation Member Function + */ + void GenerateGaussianNoise(CLHEP::HepRandomEngine* engine, const double NoiseMean, const double NoiseStd); /** @short steer digitization mode */ @@ -150,6 +154,16 @@ class HGCDigitizerBase { //if set to true, threshold will be computed based on the expected meap peak/2 bool thresholdFollowsMIP_; + + // New NoiseArray Parameters + + const double NoiseMean_, NoiseStd_; + static const size_t NoiseArrayLength_ = 200000; + static const size_t samplesize_ = 15; + std::array, NoiseArrayLength_> GaussianNoiseArray_; + bool RandNoiseGenerationFlag_; + // A parameter configurable from python configuration to decide which noise generation model to use + bool NoiseGeneration_Method_; }; #endif diff --git a/SimCalorimetry/HGCalSimProducers/interface/HGCFEElectronics.h b/SimCalorimetry/HGCalSimProducers/interface/HGCFEElectronics.h index 6c0d6511e30da..a5989c2bf7ce1 100644 --- a/SimCalorimetry/HGCalSimProducers/interface/HGCFEElectronics.h +++ b/SimCalorimetry/HGCalSimProducers/interface/HGCFEElectronics.h @@ -6,7 +6,7 @@ #include "CLHEP/Random/RandGauss.h" #include "CLHEP/Random/RandGaussQ.h" - +#include "CLHEP/Random/RandFlat.h" #include "SimCalorimetry/HGCalSimProducers/interface/HGCDigitizerTypes.h" /** diff --git a/SimCalorimetry/HGCalSimProducers/python/hgcalDigitizer_cfi.py b/SimCalorimetry/HGCalSimProducers/python/hgcalDigitizer_cfi.py index 593c7bdb36848..2a894ad326f86 100644 --- a/SimCalorimetry/HGCalSimProducers/python/hgcalDigitizer_cfi.py +++ b/SimCalorimetry/HGCalSimProducers/python/hgcalDigitizer_cfi.py @@ -66,6 +66,7 @@ accumulatorType = cms.string("HGCDigiProducer"), hitCollection = cms.string("HGCHitsEE"), digiCollection = cms.string("HGCDigisEE"), + NoiseGeneration_Method = cms.bool(True), maxSimHitsAccTime = cms.uint32(100), bxTime = cms.double(25), eVPerEleHolePair = cms.double(eV_per_eh_pair), @@ -135,6 +136,7 @@ accumulatorType = cms.string("HGCDigiProducer"), hitCollection = cms.string("HGCHitsHEfront"), digiCollection = cms.string("HGCDigisHEfront"), + NoiseGeneration_Method = cms.bool(True), maxSimHitsAccTime = cms.uint32(100), bxTime = cms.double(25), tofDelay = cms.double(5), @@ -203,6 +205,7 @@ accumulatorType = cms.string("HGCDigiProducer"), hitCollection = cms.string("HcalHits"), digiCollection = cms.string("HGCDigisHEback"), + NoiseGeneration_Method = cms.bool(True), maxSimHitsAccTime = cms.uint32(100), bxTime = cms.double(25), tofDelay = cms.double(1), @@ -246,6 +249,7 @@ accumulatorType = cms.string("HGCDigiProducer"), hitCollection = cms.string("HFNoseHits"), digiCollection = cms.string("HFNoseDigis"), + NoiseGeneration_Method = cms.bool(True), maxSimHitsAccTime = cms.uint32(100), bxTime = cms.double(25), eVPerEleHolePair = cms.double(eV_per_eh_pair), diff --git a/SimCalorimetry/HGCalSimProducers/src/HGCDigitizerBase.cc b/SimCalorimetry/HGCalSimProducers/src/HGCDigitizerBase.cc index 94623b346e439..4a40cdad1bf26 100644 --- a/SimCalorimetry/HGCalSimProducers/src/HGCDigitizerBase.cc +++ b/SimCalorimetry/HGCalSimProducers/src/HGCDigitizerBase.cc @@ -6,12 +6,13 @@ using namespace hgc_digi; using namespace hgc_digi_utils; template -HGCDigitizerBase::HGCDigitizerBase(const edm::ParameterSet& ps) : scaleByDose_(false) { +HGCDigitizerBase::HGCDigitizerBase(const edm::ParameterSet& ps) + : scaleByDose_(false), NoiseMean_(0.0), NoiseStd_(1.0) { bxTime_ = ps.getParameter("bxTime"); myCfg_ = ps.getParameter("digiCfg"); + NoiseGeneration_Method_ = ps.getParameter("NoiseGeneration_Method"); doTimeSamples_ = myCfg_.getParameter("doTimeSamples"); thresholdFollowsMIP_ = myCfg_.getParameter("thresholdFollowsMIP"); - if (myCfg_.exists("keV2fC")) keV2fC_ = myCfg_.getParameter("keV2fC"); else @@ -51,6 +52,18 @@ HGCDigitizerBase::HGCDigitizerBase(const edm::ParameterSet& ps) : scaleByDo edm::ParameterSet feCfg = myCfg_.getParameter("feCfg"); myFEelectronics_ = std::unique_ptr>(new HGCFEElectronics(feCfg)); myFEelectronics_->SetNoiseValues(noise_fC_); + RandNoiseGenerationFlag_ = 0; +} + +template +void HGCDigitizerBase::GenerateGaussianNoise(CLHEP::HepRandomEngine* engine, + const double NoiseMean, + const double NoiseStd) { + for (size_t i = 0; i < NoiseArrayLength_; i++) { + for (size_t j = 0; j < samplesize_; j++) { + GaussianNoiseArray_[i][j] = CLHEP::RandGaussQ::shoot(engine, NoiseMean, NoiseStd); + } + } } template @@ -62,6 +75,12 @@ void HGCDigitizerBase::run(std::unique_ptr& digiCo CLHEP::HepRandomEngine* engine) { if (scaleByDose_) scal_.setGeometry(theGeom); + if (NoiseGeneration_Method_ == true) { + if (RandNoiseGenerationFlag_ == false) { + GenerateGaussianNoise(engine, NoiseMean_, NoiseStd_); + RandNoiseGenerationFlag_ = true; + } + } if (digitizationType == 0) runSimple(digiColl, simData, theGeom, validIds, engine); else @@ -80,14 +99,20 @@ void HGCDigitizerBase::runSimple(std::unique_ptr& HGCCellInfo zeroData; zeroData.hit_info[0].fill(0.f); //accumulated energy zeroData.hit_info[1].fill(0.f); //time-of-flight - + std::array cellNoiseArray; + for (size_t i = 0; i < samplesize_; i++) + cellNoiseArray[i] = 0.0; for (const auto& id : validIds) { chargeColl.fill(0.f); toa.fill(0.f); HGCSimHitDataAccumulator::iterator it = simData.find(id); HGCCellInfo& cell = (simData.end() == it ? zeroData : it->second); addCellMetadata(cell, theGeom, id); + if (NoiseGeneration_Method_ == true) { + size_t hash_index = (CLHEP::RandFlat::shootInt(engine, (NoiseArrayLength_ - 1)) + id) % NoiseArrayLength_; + cellNoiseArray = GaussianNoiseArray_[hash_index]; + } //set the noise,cce, LSB and threshold to be used float cce(1.f), noiseWidth(0.f), lsbADC(-1.f), maxADC(-1.f); uint32_t thrADC(std::floor(myFEelectronics_->getTargetMipValue() / 2)); @@ -122,7 +147,11 @@ void HGCDigitizerBase::runSimple(std::unique_ptr& toa[i] = cell.hit_info[1][i] / rawCharge; //final charge estimation - float noise = CLHEP::RandGaussQ::shoot(engine, 0.0, noiseWidth); + float noise; + if (NoiseGeneration_Method_ == true) + noise = (float)cellNoiseArray[i] * noiseWidth; + else + noise = CLHEP::RandGaussQ::shoot(engine, cellNoiseArray[i], noiseWidth); float totalCharge(rawCharge * cce + noise); if (totalCharge < 0.f) totalCharge = 0.f; @@ -131,7 +160,8 @@ void HGCDigitizerBase::runSimple(std::unique_ptr& //run the shaper to create a new data frame DFr rawDataFrame(id); - myFEelectronics_->runShaper(rawDataFrame, chargeColl, toa, engine, thrADC, lsbADC, maxADC, cell.thickness); + int thickness = cell.thickness > 0 ? cell.thickness : 1; + myFEelectronics_->runShaper(rawDataFrame, chargeColl, toa, engine, thrADC, lsbADC, maxADC, thickness); //update the output according to the final shape updateOutput(coll, rawDataFrame);