Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-11-17 23:30:01

0001 #include "SimCalorimetry/HcalSimProducers/interface/HcalDigitizer.h"
0002 #include "CalibFormats/HcalObjects/interface/HcalDbService.h"
0003 #include "DataFormats/Common/interface/Handle.h"
0004 #include "DataFormats/Common/interface/Wrapper.h"
0005 #include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
0006 #include "DataFormats/HcalDetId/interface/HcalTestNumbering.h"
0007 #include "DataFormats/HcalDetId/interface/HcalZDCDetId.h"
0008 #include "DataFormats/HcalDigi/interface/HcalDigiCollections.h"
0009 #include "DataFormats/HcalDigi/interface/HcalQIENum.h"
0010 #include "FWCore/Framework/interface/ConsumesCollector.h"
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/EventSetup.h"
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "Geometry/CaloTopology/interface/HcalTopology.h"
0016 #include "Geometry/ForwardGeometry/interface/ZdcTopology.h"
0017 #include "SimCalorimetry/CaloSimAlgos/interface/CaloHitResponse.h"
0018 #include "SimCalorimetry/CaloSimAlgos/interface/CaloTDigitizer.h"
0019 #include "SimCalorimetry/HcalSimAlgos/interface/HPDIonFeedbackSim.h"
0020 #include "SimCalorimetry/HcalSimAlgos/interface/HcalAmplifier.h"
0021 #include "SimCalorimetry/HcalSimAlgos/interface/HcalBaseSignalGenerator.h"
0022 #include "SimCalorimetry/HcalSimAlgos/interface/HcalCoderFactory.h"
0023 #include "SimCalorimetry/HcalSimAlgos/interface/HcalElectronicsSim.h"
0024 #include "SimCalorimetry/HcalSimAlgos/interface/HcalSiPMHitResponse.h"
0025 #include "SimCalorimetry/HcalSimAlgos/interface/HcalSimParameterMap.h"
0026 #include "SimCalorimetry/HcalSimAlgos/interface/HcalTimeSlewSim.h"
0027 #include "SimDataFormats/CrossingFrame/interface/CrossingFrame.h"
0028 #include "SimDataFormats/CrossingFrame/interface/MixCollection.h"
0029 #include "SimGeneral/MixingModule/interface/PileUpEventPrincipal.h"
0030 #include <boost/foreach.hpp>
0031 #include <memory>
0032 
0033 // #define EDM_ML_DEBUG
0034 
0035 HcalDigitizer::HcalDigitizer(const edm::ParameterSet &ps, edm::ConsumesCollector &iC)
0036     : conditionsToken_(iC.esConsumes()),
0037       topoToken_(iC.esConsumes()),
0038       topoZToken_(iC.esConsumes()),
0039       hcalTimeSlew_delay_token_(iC.esConsumes(edm::ESInputTag("", "HBHE"))),
0040       theGeometryToken(iC.esConsumes()),
0041       theRecNumberToken(iC.esConsumes()),
0042       qieTypesToken_(iC.esConsumes()),
0043       theGeometry(nullptr),
0044       theRecNumber(nullptr),
0045       theParameterMap(ps),
0046       theShapes(),
0047       theHBHEResponse(std::make_unique<CaloHitResponse>(&theParameterMap, &theShapes)),
0048       theHBHESiPMResponse(std::make_unique<HcalSiPMHitResponse>(
0049           &theParameterMap, &theShapes, ps.getParameter<bool>("HcalPreMixStage1"), true)),
0050       theHOResponse(std::make_unique<CaloHitResponse>(&theParameterMap, &theShapes)),
0051       theHOSiPMResponse(std::make_unique<HcalSiPMHitResponse>(
0052           &theParameterMap, &theShapes, ps.getParameter<bool>("HcalPreMixStage1"), false)),
0053       theHFResponse(std::make_unique<CaloHitResponse>(&theParameterMap, &theShapes)),
0054       theHFQIE10Response(std::make_unique<CaloHitResponse>(
0055           &theParameterMap, &theShapes, ps.getParameter<bool>("HcalPreMixStage1"), true)),
0056       theZDCResponse(std::make_unique<CaloHitResponse>(
0057           &theParameterMap, &theShapes, ps.getParameter<bool>("HcalPreMixStage1"), false)),
0058       theHBHEAmplifier(nullptr),
0059       theHFAmplifier(nullptr),
0060       theHOAmplifier(nullptr),
0061       theZDCAmplifier(nullptr),
0062       theHFQIE10Amplifier(nullptr),
0063       theHBHEQIE11Amplifier(nullptr),
0064       theIonFeedback(nullptr),
0065       theCoderFactory(nullptr),
0066       theHBHEElectronicsSim(nullptr),
0067       theHFElectronicsSim(nullptr),
0068       theHOElectronicsSim(nullptr),
0069       theZDCElectronicsSim(nullptr),
0070       theHFQIE10ElectronicsSim(nullptr),
0071       theHBHEQIE11ElectronicsSim(nullptr),
0072       theHBHEHitFilter(),
0073       theHBHEQIE11HitFilter(),
0074       theHFHitFilter(),
0075       theHFQIE10HitFilter(),
0076       theHOHitFilter(),
0077       theHOSiPMHitFilter(),
0078       theZDCHitFilter(),
0079       theHBHEDigitizer(nullptr),
0080       theHODigitizer(nullptr),
0081       theHOSiPMDigitizer(nullptr),
0082       theHFDigitizer(nullptr),
0083       theZDCDigitizer(nullptr),
0084       theHFQIE10Digitizer(nullptr),
0085       theHBHEQIE11Digitizer(nullptr),
0086       theRelabeller(nullptr),
0087       isZDC(ps.getParameter<bool>("doZDCDigi")),
0088       isHCAL(true),
0089       zdcgeo(true),
0090       hbhegeo(true),
0091       hogeo(true),
0092       hfgeo(true),
0093       doHFWindow_(ps.getParameter<bool>("doHFWindow")),
0094       killHE_(ps.getParameter<bool>("killHE")),
0095       debugCS_(ps.getParameter<bool>("debugCaloSamples")),
0096       ignoreTime_(ps.getParameter<bool>("ignoreGeantTime")),
0097       injectTestHits_(ps.getParameter<bool>("injectTestHits")),
0098       hitsProducer_(ps.getParameter<std::string>("hitsProducer")),
0099       theHOSiPMCode(ps.getParameter<edm::ParameterSet>("ho").getParameter<int>("siPMCode")),
0100       deliveredLumi(0.),
0101       agingFlagHB(ps.getParameter<bool>("HBDarkening")),
0102       agingFlagHE(ps.getParameter<bool>("HEDarkening")),
0103       zdcToken_(iC.consumes(edm::InputTag(hitsProducer_, "ZDCHITS"))),
0104       hcalToken_(iC.consumes(edm::InputTag(hitsProducer_, "HcalHits"))),
0105       m_HBDarkening(nullptr),
0106       m_HEDarkening(nullptr),
0107       m_HFRecalibration(nullptr),
0108       injectedHitsEnergy_(ps.getParameter<std::vector<double>>("injectTestHitsEnergy")),
0109       injectedHitsTime_(ps.getParameter<std::vector<double>>("injectTestHitsTime")),
0110       injectedHitsCells_(ps.getParameter<std::vector<int>>("injectTestHitsCells")) {
0111   if (agingFlagHB) {
0112     m_HBDarkeningToken = iC.esConsumes(edm::ESInputTag("", "HB"));
0113   }
0114   if (agingFlagHE) {
0115     m_HEDarkeningToken = iC.esConsumes(edm::ESInputTag("", "HE"));
0116   }
0117   if (theHOSiPMCode == 2) {
0118     mcParamsToken_ = iC.esConsumes();
0119   }
0120 
0121   bool doNoise = ps.getParameter<bool>("doNoise");
0122 
0123   bool PreMix1 = ps.getParameter<bool>("HcalPreMixStage1");  // special threshold/pedestal treatment
0124   bool PreMix2 = ps.getParameter<bool>("HcalPreMixStage2");  // special threshold/pedestal treatment
0125   bool doEmpty = ps.getParameter<bool>("doEmpty");
0126   deliveredLumi = ps.getParameter<double>("DelivLuminosity");
0127   bool agingFlagHF = ps.getParameter<bool>("HFDarkening");
0128   double minFCToDelay = ps.getParameter<double>("minFCToDelay");
0129 
0130   if (PreMix1 && PreMix2) {
0131     throw cms::Exception("Configuration") << "HcalDigitizer cannot operate in PreMixing digitization and "
0132                                              "PreMixing\n"
0133                                              "digi combination modes at the same time.  Please set one mode to "
0134                                              "False\n"
0135                                              "in the configuration file.";
0136   }
0137 
0138   // need to make copies, because they might get different noise generators
0139   theHBHEAmplifier = std::make_unique<HcalAmplifier>(&theParameterMap, doNoise, PreMix1, PreMix2);
0140   theHFAmplifier = std::make_unique<HcalAmplifier>(&theParameterMap, doNoise, PreMix1, PreMix2);
0141   theHOAmplifier = std::make_unique<HcalAmplifier>(&theParameterMap, doNoise, PreMix1, PreMix2);
0142   theZDCAmplifier = std::make_unique<HcalAmplifier>(&theParameterMap, doNoise, PreMix1, PreMix2);
0143   theHFQIE10Amplifier = std::make_unique<HcalAmplifier>(&theParameterMap, doNoise, PreMix1, PreMix2);
0144   theHBHEQIE11Amplifier = std::make_unique<HcalAmplifier>(&theParameterMap, doNoise, PreMix1, PreMix2);
0145 
0146   theCoderFactory = std::make_unique<HcalCoderFactory>(HcalCoderFactory::DB);
0147 
0148   theHBHEElectronicsSim =
0149       std::make_unique<HcalElectronicsSim>(&theParameterMap, theHBHEAmplifier.get(), theCoderFactory.get(), PreMix1);
0150   theHFElectronicsSim =
0151       std::make_unique<HcalElectronicsSim>(&theParameterMap, theHFAmplifier.get(), theCoderFactory.get(), PreMix1);
0152   theHOElectronicsSim =
0153       std::make_unique<HcalElectronicsSim>(&theParameterMap, theHOAmplifier.get(), theCoderFactory.get(), PreMix1);
0154   theZDCElectronicsSim =
0155       std::make_unique<HcalElectronicsSim>(&theParameterMap, theZDCAmplifier.get(), theCoderFactory.get(), PreMix1);
0156   theHFQIE10ElectronicsSim =
0157       std::make_unique<HcalElectronicsSim>(&theParameterMap,
0158                                            theHFQIE10Amplifier.get(),
0159                                            theCoderFactory.get(),
0160                                            PreMix1);  // should this use a different coder factory?
0161   theHBHEQIE11ElectronicsSim =
0162       std::make_unique<HcalElectronicsSim>(&theParameterMap,
0163                                            theHBHEQIE11Amplifier.get(),
0164                                            theCoderFactory.get(),
0165                                            PreMix1);  // should this use a different coder factory?
0166 
0167   bool doHOHPD = (theHOSiPMCode != 1);
0168   bool doHOSiPM = (theHOSiPMCode != 0);
0169   if (doHOHPD) {
0170     theHOResponse = std::make_unique<CaloHitResponse>(&theParameterMap, &theShapes);
0171     theHOResponse->setHitFilter(&theHOHitFilter);
0172     theHODigitizer = std::make_unique<HODigitizer>(theHOResponse.get(), theHOElectronicsSim.get(), doEmpty);
0173   }
0174   if (doHOSiPM) {
0175     theHOSiPMResponse->setHitFilter(&theHOSiPMHitFilter);
0176     theHOSiPMDigitizer = std::make_unique<HODigitizer>(theHOSiPMResponse.get(), theHOElectronicsSim.get(), doEmpty);
0177   }
0178 
0179   theHBHEResponse->setHitFilter(&theHBHEHitFilter);
0180   theHBHESiPMResponse->setHitFilter(&theHBHEQIE11HitFilter);
0181 
0182   // QIE8 and QIE11 can coexist in HBHE
0183   theHBHEQIE11Digitizer =
0184       std::make_unique<QIE11Digitizer>(theHBHESiPMResponse.get(), theHBHEQIE11ElectronicsSim.get(), doEmpty);
0185   theHBHEDigitizer = std::make_unique<HBHEDigitizer>(theHBHEResponse.get(), theHBHEElectronicsSim.get(), doEmpty);
0186 
0187   bool doTimeSlew = ps.getParameter<bool>("doTimeSlew");
0188   // initialize: they won't be called later if flag is set
0189   hcalTimeSlew_delay_ = nullptr;
0190   theTimeSlewSim.reset(nullptr);
0191   if (doTimeSlew) {
0192     // no time slewing for HF
0193     theTimeSlewSim = std::make_unique<HcalTimeSlewSim>(&theParameterMap, minFCToDelay);
0194     theHBHEAmplifier->setTimeSlewSim(theTimeSlewSim.get());
0195     theHBHEQIE11Amplifier->setTimeSlewSim(theTimeSlewSim.get());
0196     theHOAmplifier->setTimeSlewSim(theTimeSlewSim.get());
0197     theZDCAmplifier->setTimeSlewSim(theTimeSlewSim.get());
0198   }
0199 
0200   theHFResponse->setHitFilter(&theHFHitFilter);
0201   theHFQIE10Response->setHitFilter(&theHFQIE10HitFilter);
0202   theZDCResponse->setHitFilter(&theZDCHitFilter);
0203 
0204   // QIE8 and QIE10 can coexist in HF
0205   theHFQIE10Digitizer =
0206       std::make_unique<QIE10Digitizer>(theHFQIE10Response.get(), theHFQIE10ElectronicsSim.get(), doEmpty);
0207   theHFDigitizer = std::make_unique<HFDigitizer>(theHFResponse.get(), theHFElectronicsSim.get(), doEmpty);
0208 
0209   // temporary move until Run3 ZDC will be added
0210   if (isZDC) {
0211     theZDCDigitizer = std::make_unique<ZDCDigitizer>(theZDCResponse.get(), theZDCElectronicsSim.get(), doEmpty);
0212   } else {
0213     theZDCDigitizer = nullptr;
0214     edm::LogVerbatim("HcalDigitizer") << "Inhibit ZDC Digitization";
0215   }
0216   edm::LogVerbatim("HcalDigitizer") << "iszDC: " << isZDC << " theZDCDigitizer: " << theZDCDigitizer;
0217 
0218   testNumbering_ = ps.getParameter<bool>("TestNumbering");
0219 #ifdef EDM_ML_DEBUG
0220   edm::LogVerbatim("HcalSim") << "Flag to see if Hit Relabeller to be initiated " << testNumbering_;
0221 #endif
0222   if (testNumbering_)
0223     theRelabeller = std::make_unique<HcalHitRelabeller>(ps.getParameter<bool>("doNeutralDensityFilter"));
0224 
0225   if (ps.getParameter<bool>("doIonFeedback") && theHBHEResponse) {
0226     theIonFeedback = std::make_unique<HPDIonFeedbackSim>(ps, &theShapes);
0227     theHBHEResponse->setPECorrection(theIonFeedback.get());
0228     if (ps.getParameter<bool>("doThermalNoise")) {
0229       theHBHEAmplifier->setIonFeedbackSim(theIonFeedback.get());
0230     }
0231   }
0232 
0233   // option to save CaloSamples as event product for debugging
0234   if (debugCS_) {
0235     if (theHBHEDigitizer)
0236       theHBHEDigitizer->setDebugCaloSamples(true);
0237     if (theHBHEQIE11Digitizer)
0238       theHBHEQIE11Digitizer->setDebugCaloSamples(true);
0239     if (theHODigitizer)
0240       theHODigitizer->setDebugCaloSamples(true);
0241     if (theHOSiPMDigitizer)
0242       theHOSiPMDigitizer->setDebugCaloSamples(true);
0243     if (theHFDigitizer)
0244       theHFDigitizer->setDebugCaloSamples(true);
0245     if (theHFQIE10Digitizer)
0246       theHFQIE10Digitizer->setDebugCaloSamples(true);
0247     if (theZDCDigitizer)
0248       theZDCDigitizer->setDebugCaloSamples(true);
0249   }
0250 
0251   // option to ignore Geant time distribution in SimHits, for debugging
0252   if (ignoreTime_) {
0253     theHBHEResponse->setIgnoreGeantTime(ignoreTime_);
0254     theHBHESiPMResponse->setIgnoreGeantTime(ignoreTime_);
0255     theHOResponse->setIgnoreGeantTime(ignoreTime_);
0256     theHOSiPMResponse->setIgnoreGeantTime(ignoreTime_);
0257     theHFResponse->setIgnoreGeantTime(ignoreTime_);
0258     theHFQIE10Response->setIgnoreGeantTime(ignoreTime_);
0259     theZDCResponse->setIgnoreGeantTime(ignoreTime_);
0260   }
0261 
0262   if (agingFlagHF)
0263     m_HFRecalibration = std::make_unique<HFRecalibration>(ps.getParameter<edm::ParameterSet>("HFRecalParameterBlock"));
0264 }
0265 
0266 HcalDigitizer::~HcalDigitizer() {}
0267 
0268 void HcalDigitizer::setHBHENoiseSignalGenerator(HcalBaseSignalGenerator *noiseGenerator) {
0269   noiseGenerator->setParameterMap(&theParameterMap);
0270   noiseGenerator->setElectronicsSim(theHBHEElectronicsSim.get());
0271   if (theHBHEDigitizer)
0272     theHBHEDigitizer->setNoiseSignalGenerator(noiseGenerator);
0273   theHBHEAmplifier->setNoiseSignalGenerator(noiseGenerator);
0274 }
0275 
0276 void HcalDigitizer::setQIE11NoiseSignalGenerator(HcalBaseSignalGenerator *noiseGenerator) {
0277   noiseGenerator->setParameterMap(&theParameterMap);
0278   noiseGenerator->setElectronicsSim(theHBHEQIE11ElectronicsSim.get());
0279   if (theHBHEQIE11Digitizer)
0280     theHBHEQIE11Digitizer->setNoiseSignalGenerator(noiseGenerator);
0281   theHBHEQIE11Amplifier->setNoiseSignalGenerator(noiseGenerator);
0282 }
0283 
0284 void HcalDigitizer::setHFNoiseSignalGenerator(HcalBaseSignalGenerator *noiseGenerator) {
0285   noiseGenerator->setParameterMap(&theParameterMap);
0286   noiseGenerator->setElectronicsSim(theHFElectronicsSim.get());
0287   if (theHFDigitizer)
0288     theHFDigitizer->setNoiseSignalGenerator(noiseGenerator);
0289   theHFAmplifier->setNoiseSignalGenerator(noiseGenerator);
0290 }
0291 
0292 void HcalDigitizer::setQIE10NoiseSignalGenerator(HcalBaseSignalGenerator *noiseGenerator) {
0293   noiseGenerator->setParameterMap(&theParameterMap);
0294   noiseGenerator->setElectronicsSim(theHFQIE10ElectronicsSim.get());
0295   if (theHFQIE10Digitizer)
0296     theHFQIE10Digitizer->setNoiseSignalGenerator(noiseGenerator);
0297   theHFQIE10Amplifier->setNoiseSignalGenerator(noiseGenerator);
0298 }
0299 
0300 void HcalDigitizer::setHONoiseSignalGenerator(HcalBaseSignalGenerator *noiseGenerator) {
0301   noiseGenerator->setParameterMap(&theParameterMap);
0302   noiseGenerator->setElectronicsSim(theHOElectronicsSim.get());
0303   if (theHODigitizer)
0304     theHODigitizer->setNoiseSignalGenerator(noiseGenerator);
0305   if (theHOSiPMDigitizer)
0306     theHOSiPMDigitizer->setNoiseSignalGenerator(noiseGenerator);
0307   theHOAmplifier->setNoiseSignalGenerator(noiseGenerator);
0308 }
0309 
0310 void HcalDigitizer::setZDCNoiseSignalGenerator(HcalBaseSignalGenerator *noiseGenerator) {
0311   noiseGenerator->setParameterMap(&theParameterMap);
0312   noiseGenerator->setElectronicsSim(theZDCElectronicsSim.get());
0313   if (theZDCDigitizer) {
0314     theZDCDigitizer->setNoiseSignalGenerator(noiseGenerator);
0315     theZDCAmplifier->setNoiseSignalGenerator(noiseGenerator);
0316   }
0317 }
0318 
0319 void HcalDigitizer::initializeEvent(edm::Event const &e, edm::EventSetup const &eventSetup) {
0320   setup(eventSetup);
0321 
0322   // get the appropriate gains, noises, & widths for this event
0323   const HcalDbService *conditions = &eventSetup.getData(conditionsToken_);
0324 
0325   theShapes.setDbService(conditions);
0326 
0327   theHBHEAmplifier->setDbService(conditions);
0328   theHFAmplifier->setDbService(conditions);
0329   theHOAmplifier->setDbService(conditions);
0330   theZDCAmplifier->setDbService(conditions);
0331   theHFQIE10Amplifier->setDbService(conditions);
0332   theHBHEQIE11Amplifier->setDbService(conditions);
0333 
0334   theHFQIE10ElectronicsSim->setDbService(conditions);
0335   theHBHEQIE11ElectronicsSim->setDbService(conditions);
0336 
0337   theCoderFactory->setDbService(conditions);
0338   theParameterMap.setDbService(conditions);
0339 
0340   // initialize hits
0341   if (theHBHEDigitizer)
0342     theHBHEDigitizer->initializeHits();
0343   if (theHBHEQIE11Digitizer)
0344     theHBHEQIE11Digitizer->initializeHits();
0345   if (theHODigitizer)
0346     theHODigitizer->initializeHits();
0347   if (theHOSiPMDigitizer)
0348     theHOSiPMDigitizer->initializeHits();
0349   if (theHFQIE10Digitizer)
0350     theHFQIE10Digitizer->initializeHits();
0351   if (theHFDigitizer)
0352     theHFDigitizer->initializeHits();
0353   if (theZDCDigitizer)
0354     theZDCDigitizer->initializeHits();
0355 }
0356 
0357 void HcalDigitizer::accumulateCaloHits(edm::Handle<std::vector<PCaloHit>> const &hcalHandle,
0358                                        edm::Handle<std::vector<PCaloHit>> const &zdcHandle,
0359                                        int bunchCrossing,
0360                                        CLHEP::HepRandomEngine *engine,
0361                                        const HcalTopology *htopoP,
0362                                        const ZdcTopology *ztopoP) {
0363   // Step A: pass in inputs, and accumulate digis
0364   if (isHCAL) {
0365     std::vector<PCaloHit> hcalHitsOrig = *hcalHandle.product();
0366     if (injectTestHits_)
0367       hcalHitsOrig = injectedHits_;
0368     std::vector<PCaloHit> hcalHits;
0369     hcalHits.reserve(hcalHitsOrig.size());
0370 
0371     // evaluate darkening before relabeling
0372     if (testNumbering_) {
0373       if (m_HBDarkening || m_HEDarkening || m_HFRecalibration) {
0374         darkening(hcalHitsOrig);
0375       }
0376       // Relabel PCaloHits if necessary
0377       edm::LogVerbatim("HcalDigitizer") << "Calling Relabeller";
0378       theRelabeller->process(hcalHitsOrig);
0379     }
0380 
0381     // eliminate bad hits
0382     for (unsigned int i = 0; i < hcalHitsOrig.size(); i++) {
0383       DetId id(hcalHitsOrig[i].id());
0384       HcalDetId hid(id);
0385       if (!htopoP->validHcal(hid)) {
0386         edm::LogError("HcalDigitizer") << "bad hcal id found in digitizer. Skipping " << id.rawId() << " " << hid;
0387         continue;
0388       } else if (hid.subdet() == HcalForward && !doHFWindow_ && hcalHitsOrig[i].depth() != 0) {
0389         // skip HF window hits unless desired
0390         continue;
0391       } else if (killHE_ && hid.subdet() == HcalEndcap) {
0392         // remove HE hits if asked for (phase 2)
0393         continue;
0394       } else {
0395 #ifdef EDM_ML_DEBUG
0396         edm::LogVerbatim("HcalSim") << "HcalDigitizer format " << hid.oldFormat() << " for " << hid;
0397 #endif
0398         DetId newid = DetId(hid.newForm());
0399 #ifdef EDM_ML_DEBUG
0400         edm::LogVerbatim("HcalSim") << "Hit " << i << " out of " << hcalHits.size() << " " << std::hex << id.rawId()
0401                                     << " --> " << newid.rawId() << std::dec << " " << HcalDetId(newid.rawId()) << '\n';
0402 #endif
0403         hcalHitsOrig[i].setID(newid.rawId());
0404         hcalHits.push_back(hcalHitsOrig[i]);
0405       }
0406     }
0407 
0408     if (hbhegeo) {
0409       if (theHBHEDigitizer)
0410         theHBHEDigitizer->add(hcalHits, bunchCrossing, engine);
0411       if (theHBHEQIE11Digitizer)
0412         theHBHEQIE11Digitizer->add(hcalHits, bunchCrossing, engine);
0413     }
0414 
0415     if (hogeo) {
0416       if (theHODigitizer)
0417         theHODigitizer->add(hcalHits, bunchCrossing, engine);
0418       if (theHOSiPMDigitizer)
0419         theHOSiPMDigitizer->add(hcalHits, bunchCrossing, engine);
0420     }
0421 
0422     if (hfgeo) {
0423       if (theHFDigitizer)
0424         theHFDigitizer->add(hcalHits, bunchCrossing, engine);
0425       if (theHFQIE10Digitizer)
0426         theHFQIE10Digitizer->add(hcalHits, bunchCrossing, engine);
0427     }
0428   } else {
0429     edm::LogVerbatim("HcalDigitizer") << "We don't have HCAL hit collection available ";
0430   }
0431 
0432   if (isZDC) {
0433     if (zdcgeo && theZDCDigitizer) {
0434       std::vector<PCaloHit> zdcHitsOrig = *zdcHandle.product();
0435       std::vector<PCaloHit> zdcHits;
0436       zdcHits.reserve(zdcHitsOrig.size());
0437       // eliminate bad hits
0438       for (unsigned int i = 0; i < zdcHitsOrig.size(); i++) {
0439         DetId id(zdcHitsOrig[i].id());
0440         HcalZDCDetId hid(id);
0441         if (!ztopoP->valid(hid)) {
0442           edm::LogError("HcalDigitizer") << "bad zdc id found in digitizer. Skipping " << std::hex << id.rawId()
0443                                          << std::dec << " " << hid;
0444           continue;
0445         }
0446         zdcHits.push_back(zdcHitsOrig[i]);
0447 #ifdef EDM_ML_DEBUG
0448         edm::LogVerbatim("HcalSim") << "Hit " << i << " out of " << zdcHitsOrig.size() << " " << std::hex << id.rawId()
0449                                     << " " << hid;
0450 #endif
0451       }
0452       theZDCDigitizer->add(zdcHits, bunchCrossing, engine);
0453     }
0454   } else {
0455     edm::LogVerbatim("HcalDigitizer") << "We don't have ZDC hit collection available ";
0456   }
0457 }
0458 
0459 void HcalDigitizer::accumulate(edm::Event const &e, edm::EventSetup const &eventSetup, CLHEP::HepRandomEngine *engine) {
0460   // Step A: Get Inputs
0461   const edm::Handle<std::vector<PCaloHit>> &zdcHandle = e.getHandle(zdcToken_);
0462   isZDC = zdcHandle.isValid();
0463 
0464   const edm::Handle<std::vector<PCaloHit>> &hcalHandle = e.getHandle(hcalToken_);
0465   isHCAL = hcalHandle.isValid() or injectTestHits_;
0466 
0467   const HcalTopology *htopoP = &eventSetup.getData(topoToken_);
0468   const ZdcTopology *ztopoP = &eventSetup.getData(topoZToken_);
0469 
0470   accumulateCaloHits(hcalHandle, zdcHandle, 0, engine, htopoP, ztopoP);
0471 }
0472 
0473 void HcalDigitizer::accumulate(PileUpEventPrincipal const &e,
0474                                edm::EventSetup const &eventSetup,
0475                                CLHEP::HepRandomEngine *engine) {
0476   // Step A: Get Inputs
0477   edm::InputTag zdcTag(hitsProducer_, "ZDCHITS");
0478   edm::Handle<std::vector<PCaloHit>> zdcHandle;
0479   e.getByLabel(zdcTag, zdcHandle);
0480   isZDC = zdcHandle.isValid();
0481 
0482   edm::InputTag hcalTag(hitsProducer_, "HcalHits");
0483   edm::Handle<std::vector<PCaloHit>> hcalHandle;
0484   e.getByLabel(hcalTag, hcalHandle);
0485   isHCAL = hcalHandle.isValid();
0486 
0487   const HcalTopology *htopoP = &eventSetup.getData(topoToken_);
0488   const ZdcTopology *ztopoP = &eventSetup.getData(topoZToken_);
0489 
0490   accumulateCaloHits(hcalHandle, zdcHandle, e.bunchCrossing(), engine, htopoP, ztopoP);
0491 }
0492 
0493 void HcalDigitizer::finalizeEvent(edm::Event &e, const edm::EventSetup &eventSetup, CLHEP::HepRandomEngine *engine) {
0494   // Step B: Create empty output
0495   std::unique_ptr<HBHEDigiCollection> hbheResult(new HBHEDigiCollection());
0496   std::unique_ptr<HODigiCollection> hoResult(new HODigiCollection());
0497   std::unique_ptr<HFDigiCollection> hfResult(new HFDigiCollection());
0498   std::unique_ptr<ZDCDigiCollection> zdcResult(new ZDCDigiCollection());
0499   std::unique_ptr<QIE10DigiCollection> hfQIE10Result(new QIE10DigiCollection(
0500       !theHFQIE10DetIds.empty() ? theHFQIE10Response.get()->getReadoutFrameSize(theHFQIE10DetIds[0])
0501                                 : QIE10DigiCollection::MAXSAMPLES));
0502   std::unique_ptr<QIE11DigiCollection> hbheQIE11Result(new QIE11DigiCollection(
0503       !theHBHEQIE11DetIds.empty() ? theHBHESiPMResponse.get()->getReadoutFrameSize(theHBHEQIE11DetIds[0]) :
0504                                   //      theParameterMap->simParameters(theHBHEQIE11DetIds[0]).readoutFrameSize()
0505           //      :
0506           QIE11DigiCollection::MAXSAMPLES));
0507 
0508   // Step C: Invoke the algorithm, getting back outputs.
0509   if (isHCAL && hbhegeo) {
0510     if (theHBHEDigitizer)
0511       theHBHEDigitizer->run(*hbheResult, engine);
0512     if (theHBHEQIE11Digitizer)
0513       theHBHEQIE11Digitizer->run(*hbheQIE11Result, engine);
0514   }
0515   if (isHCAL && hogeo) {
0516     if (theHODigitizer)
0517       theHODigitizer->run(*hoResult, engine);
0518     if (theHOSiPMDigitizer)
0519       theHOSiPMDigitizer->run(*hoResult, engine);
0520   }
0521   if (isHCAL && hfgeo) {
0522     if (theHFDigitizer)
0523       theHFDigitizer->run(*hfResult, engine);
0524     if (theHFQIE10Digitizer)
0525       theHFQIE10Digitizer->run(*hfQIE10Result, engine);
0526   }
0527   if (isZDC && zdcgeo && theZDCDigitizer) {
0528     theZDCDigitizer->run(*zdcResult, engine);
0529   }
0530 
0531   edm::LogVerbatim("HcalDigitizer") << "HCAL HBHE digis : " << hbheResult->size();
0532   edm::LogVerbatim("HcalDigitizer") << "HCAL HO digis   : " << hoResult->size();
0533   edm::LogVerbatim("HcalDigitizer") << "HCAL HF digis   : " << hfResult->size();
0534   edm::LogVerbatim("HcalDigitizer") << "HCAL ZDC digis  : " << zdcResult->size();
0535   edm::LogVerbatim("HcalDigitizer") << "HCAL HF QIE10 digis : " << hfQIE10Result->size();
0536   edm::LogVerbatim("HcalDigitizer") << "HCAL HBHE QIE11 digis : " << hbheQIE11Result->size();
0537 
0538 #ifdef EDM_ML_DEBUG
0539   edm::LogVerbatim("HcalSim") << "\nHCAL HBHE digis : " << hbheResult->size();
0540   edm::LogVerbatim("HcalSim") << "HCAL HO   digis : " << hoResult->size();
0541   edm::LogVerbatim("HcalSim") << "HCAL HF   digis : " << hfResult->size();
0542   edm::LogVerbatim("HcalSim") << "HCAL ZDC  digis : " << zdcResult->size();
0543   edm::LogVerbatim("HcalSim") << "HCAL HF QIE10 digis : " << hfQIE10Result->size();
0544   edm::LogVerbatim("HcalSim") << "HCAL HBHE QIE11 digis : " << hbheQIE11Result->size();
0545 #endif
0546 
0547   // Step D: Put outputs into event
0548   e.put(std::move(hbheResult));
0549   e.put(std::move(hoResult));
0550   e.put(std::move(hfResult));
0551   e.put(std::move(zdcResult));
0552   e.put(std::move(hfQIE10Result), "HFQIE10DigiCollection");
0553   e.put(std::move(hbheQIE11Result), "HBHEQIE11DigiCollection");
0554 
0555   if (debugCS_) {
0556     std::unique_ptr<CaloSamplesCollection> csResult(new CaloSamplesCollection());
0557     // smush together all the results
0558     if (theHBHEDigitizer)
0559       csResult->insert(
0560           csResult->end(), theHBHEDigitizer->getCaloSamples().begin(), theHBHEDigitizer->getCaloSamples().end());
0561     if (theHBHEQIE11Digitizer)
0562       csResult->insert(csResult->end(),
0563                        theHBHEQIE11Digitizer->getCaloSamples().begin(),
0564                        theHBHEQIE11Digitizer->getCaloSamples().end());
0565     if (theHODigitizer)
0566       csResult->insert(
0567           csResult->end(), theHODigitizer->getCaloSamples().begin(), theHODigitizer->getCaloSamples().end());
0568     if (theHOSiPMDigitizer)
0569       csResult->insert(
0570           csResult->end(), theHOSiPMDigitizer->getCaloSamples().begin(), theHOSiPMDigitizer->getCaloSamples().end());
0571     if (theHFDigitizer)
0572       csResult->insert(
0573           csResult->end(), theHFDigitizer->getCaloSamples().begin(), theHFDigitizer->getCaloSamples().end());
0574     if (theHFQIE10Digitizer)
0575       csResult->insert(
0576           csResult->end(), theHFQIE10Digitizer->getCaloSamples().begin(), theHFQIE10Digitizer->getCaloSamples().end());
0577     if (theZDCDigitizer)
0578       csResult->insert(
0579           csResult->end(), theZDCDigitizer->getCaloSamples().begin(), theZDCDigitizer->getCaloSamples().end());
0580     e.put(std::move(csResult), "HcalSamples");
0581   }
0582 
0583   if (injectTestHits_) {
0584     std::unique_ptr<edm::PCaloHitContainer> pcResult(new edm::PCaloHitContainer());
0585     pcResult->insert(pcResult->end(), injectedHits_.begin(), injectedHits_.end());
0586     e.put(std::move(pcResult), "HcalHits");
0587   }
0588 
0589 #ifdef EDM_ML_DEBUG
0590   edm::LogVerbatim("HcalSim") << "\n========>  HcalDigitizer e.put\n";
0591 #endif
0592 }
0593 
0594 void HcalDigitizer::setup(const edm::EventSetup &es) {
0595   checkGeometry(es);
0596 
0597   if (agingFlagHB) {
0598     m_HBDarkening = &es.getData(m_HBDarkeningToken);
0599   }
0600   if (agingFlagHE) {
0601     m_HEDarkening = &es.getData(m_HEDarkeningToken);
0602   }
0603 
0604   hcalTimeSlew_delay_ = &es.getData(hcalTimeSlew_delay_token_);
0605 
0606   theHBHEAmplifier->setTimeSlew(hcalTimeSlew_delay_);
0607   theHBHEQIE11Amplifier->setTimeSlew(hcalTimeSlew_delay_);
0608   theHOAmplifier->setTimeSlew(hcalTimeSlew_delay_);
0609   theZDCAmplifier->setTimeSlew(hcalTimeSlew_delay_);
0610 }
0611 
0612 void HcalDigitizer::checkGeometry(const edm::EventSetup &eventSetup) {
0613   theGeometry = &eventSetup.getData(theGeometryToken);
0614   theRecNumber = &eventSetup.getData(theRecNumberToken);
0615 
0616   if (theHBHEResponse)
0617     theHBHEResponse->setGeometry(theGeometry);
0618   if (theHBHESiPMResponse)
0619     theHBHESiPMResponse->setGeometry(theGeometry);
0620   if (theHOResponse)
0621     theHOResponse->setGeometry(theGeometry);
0622   if (theHOSiPMResponse)
0623     theHOSiPMResponse->setGeometry(theGeometry);
0624   theHFResponse->setGeometry(theGeometry);
0625   theHFQIE10Response->setGeometry(theGeometry);
0626   theZDCResponse->setGeometry(theGeometry);
0627   if (theRelabeller)
0628     theRelabeller->setGeometry(theRecNumber);
0629 
0630   // See if it's been updated
0631   bool check1 = theGeometryWatcher_.check(eventSetup);
0632   bool check2 = theRecNumberWatcher_.check(eventSetup);
0633   if (check1 or check2) {
0634     updateGeometry(eventSetup);
0635   }
0636 }
0637 
0638 void HcalDigitizer::updateGeometry(const edm::EventSetup &eventSetup) {
0639   const std::vector<DetId> &hbCells = theGeometry->getValidDetIds(DetId::Hcal, HcalBarrel);
0640   const std::vector<DetId> &heCells = theGeometry->getValidDetIds(DetId::Hcal, HcalEndcap);
0641   const std::vector<DetId> &hoCells = theGeometry->getValidDetIds(DetId::Hcal, HcalOuter);
0642   const std::vector<DetId> &hfCells = theGeometry->getValidDetIds(DetId::Hcal, HcalForward);
0643   const std::vector<DetId> &zdcCells = theGeometry->getValidDetIds(DetId::Calo, HcalZDCDetId::SubdetectorId);
0644   // const std::vector<DetId>& hcalTrigCells =
0645   // geometry->getValidDetIds(DetId::Hcal, HcalTriggerTower); const
0646   // std::vector<DetId>& hcalCalib = geometry->getValidDetIds(DetId::Calo,
0647   // HcalCastorDetId::SubdetectorId);
0648 #ifdef EDM_ML_DEBUG
0649   edm::LogVerbatim("HcalSim") << "HcalDigitizer::CheckGeometry number of cells: " << zdcCells.size();
0650 #endif
0651   if (zdcCells.empty())
0652     zdcgeo = false;
0653   if (hbCells.empty() && heCells.empty())
0654     hbhegeo = false;
0655   if (hoCells.empty())
0656     hogeo = false;
0657   if (hfCells.empty())
0658     hfgeo = false;
0659   // combine HB & HE
0660 
0661   hbheCells = hbCells;
0662   if (!killHE_) {
0663     hbheCells.insert(hbheCells.end(), heCells.begin(), heCells.end());
0664   }
0665   // handle mixed QIE8/11 scenario in HBHE
0666   buildHBHEQIECells(hbheCells, eventSetup);
0667   if (theHBHESiPMResponse)
0668     theHBHESiPMResponse->setDetIds(theHBHEQIE11DetIds);
0669 
0670   if (theHOSiPMDigitizer) {
0671     buildHOSiPMCells(hoCells, eventSetup);
0672     if (theHOSiPMResponse)
0673       theHOSiPMResponse->setDetIds(hoCells);
0674   }
0675 
0676   // handle mixed QIE8/10 scenario in HF
0677   buildHFQIECells(hfCells, eventSetup);
0678 
0679   if (theZDCDigitizer)
0680     theZDCDigitizer->setDetIds(zdcCells);
0681 
0682   // fill test hits collection if desired and empty
0683   if (injectTestHits_ && injectedHits_.empty() && !injectedHitsCells_.empty() && !injectedHitsEnergy_.empty()) {
0684     // make list of specified cells if desired
0685     std::vector<DetId> testCells;
0686     if (injectedHitsCells_.size() >= 4) {
0687       testCells.reserve(injectedHitsCells_.size() / 4);
0688       for (unsigned ic = 0; ic < injectedHitsCells_.size(); ic += 4) {
0689         if (ic + 4 > injectedHitsCells_.size())
0690           break;
0691         testCells.push_back(HcalDetId((HcalSubdetector)injectedHitsCells_[ic],
0692                                       injectedHitsCells_[ic + 1],
0693                                       injectedHitsCells_[ic + 2],
0694                                       injectedHitsCells_[ic + 3]));
0695       }
0696     } else {
0697       int testSubdet = injectedHitsCells_[0];
0698       if (testSubdet == HcalBarrel)
0699         testCells = hbCells;
0700       else if (testSubdet == HcalEndcap)
0701         testCells = heCells;
0702       else if (testSubdet == HcalForward)
0703         testCells = hfCells;
0704       else if (testSubdet == HcalOuter)
0705         testCells = hoCells;
0706       else
0707         throw cms::Exception("Configuration") << "Unknown subdet " << testSubdet << " for HCAL test hit injection";
0708     }
0709     bool useHitTimes = (injectedHitsTime_.size() == injectedHitsEnergy_.size());
0710     injectedHits_.reserve(testCells.size() * injectedHitsEnergy_.size());
0711     for (unsigned ih = 0; ih < injectedHitsEnergy_.size(); ++ih) {
0712       double tmp = useHitTimes ? injectedHitsTime_[ih] : 0.;
0713       for (auto &aCell : testCells) {
0714         injectedHits_.emplace_back(aCell, injectedHitsEnergy_[ih], tmp);
0715       }
0716     }
0717   }
0718 }
0719 
0720 void HcalDigitizer::buildHFQIECells(const std::vector<DetId> &allCells, const edm::EventSetup &eventSetup) {
0721   // if results are already cached, no need to look again
0722   if (!theHFQIE8DetIds.empty() || !theHFQIE10DetIds.empty())
0723     return;
0724 
0725   // get the QIETypes
0726   // intentional copy
0727   HcalQIETypes qieTypes = eventSetup.getData(qieTypesToken_);
0728   if (qieTypes.topo() == nullptr) {
0729     qieTypes.setTopo(&eventSetup.getData(topoToken_));
0730   }
0731 
0732   for (std::vector<DetId>::const_iterator detItr = allCells.begin(); detItr != allCells.end(); ++detItr) {
0733     HcalQIENum qieType = HcalQIENum(qieTypes.getValues(*detItr)->getValue());
0734     if (qieType == QIE8) {
0735       theHFQIE8DetIds.push_back(*detItr);
0736     } else if (qieType == QIE10) {
0737       theHFQIE10DetIds.push_back(*detItr);
0738     } else {  // default is QIE8
0739       theHFQIE8DetIds.push_back(*detItr);
0740     }
0741   }
0742 
0743   if (!theHFQIE8DetIds.empty())
0744     theHFDigitizer->setDetIds(theHFQIE8DetIds);
0745   else {
0746     theHFDigitizer.reset();
0747   }
0748 
0749   if (!theHFQIE10DetIds.empty())
0750     theHFQIE10Digitizer->setDetIds(theHFQIE10DetIds);
0751   else {
0752     theHFQIE10Digitizer.reset();
0753   }
0754 }
0755 
0756 void HcalDigitizer::buildHBHEQIECells(const std::vector<DetId> &allCells, const edm::EventSetup &eventSetup) {
0757   // if results are already cached, no need to look again
0758   if (!theHBHEQIE8DetIds.empty() || !theHBHEQIE11DetIds.empty())
0759     return;
0760 
0761   // get the QIETypes
0762   // intentional copy
0763   HcalQIETypes qieTypes = eventSetup.getData(qieTypesToken_);
0764   if (qieTypes.topo() == nullptr) {
0765     qieTypes.setTopo(&eventSetup.getData(topoToken_));
0766   }
0767 
0768   for (std::vector<DetId>::const_iterator detItr = allCells.begin(); detItr != allCells.end(); ++detItr) {
0769     HcalQIENum qieType = HcalQIENum(qieTypes.getValues(*detItr)->getValue());
0770     if (qieType == QIE8) {
0771       theHBHEQIE8DetIds.push_back(*detItr);
0772     } else if (qieType == QIE11) {
0773       theHBHEQIE11DetIds.push_back(*detItr);
0774     } else {  // default is QIE8
0775       theHBHEQIE8DetIds.push_back(*detItr);
0776     }
0777   }
0778 
0779   if (!theHBHEQIE8DetIds.empty())
0780     theHBHEDigitizer->setDetIds(theHBHEQIE8DetIds);
0781   else {
0782     theHBHEDigitizer.reset();
0783   }
0784 
0785   if (!theHBHEQIE11DetIds.empty())
0786     theHBHEQIE11Digitizer->setDetIds(theHBHEQIE11DetIds);
0787   else {
0788     theHBHEQIE11Digitizer.reset();
0789   }
0790 
0791   if (!theHBHEQIE8DetIds.empty() && !theHBHEQIE11DetIds.empty()) {
0792     theHBHEHitFilter.setDetIds(theHBHEQIE8DetIds);
0793     theHBHEQIE11HitFilter.setDetIds(theHBHEQIE11DetIds);
0794   }
0795 }
0796 
0797 void HcalDigitizer::buildHOSiPMCells(const std::vector<DetId> &allCells, const edm::EventSetup &eventSetup) {
0798   // all HPD
0799 
0800   if (theHOSiPMCode == 0) {
0801     theHODigitizer->setDetIds(allCells);
0802   } else if (theHOSiPMCode == 1) {
0803     theHOSiPMDigitizer->setDetIds(allCells);
0804     // FIXME pick Zecotek or hamamatsu?
0805   } else if (theHOSiPMCode == 2) {
0806     std::vector<HcalDetId> zecotekDetIds, hamamatsuDetIds;
0807 
0808     // intentional copy
0809     HcalMCParams mcParams = eventSetup.getData(mcParamsToken_);
0810     if (mcParams.topo() == nullptr) {
0811       mcParams.setTopo(&eventSetup.getData(topoToken_));
0812     }
0813 
0814     for (std::vector<DetId>::const_iterator detItr = allCells.begin(); detItr != allCells.end(); ++detItr) {
0815       int shapeType = mcParams.getValues(*detItr)->signalShape();
0816       if (shapeType == HcalShapes::ZECOTEK) {
0817         zecotekDetIds.emplace_back(*detItr);
0818         theHOSiPMDetIds.push_back(*detItr);
0819       } else if (shapeType == HcalShapes::HAMAMATSU) {
0820         hamamatsuDetIds.emplace_back(*detItr);
0821         theHOSiPMDetIds.push_back(*detItr);
0822       } else {
0823         theHOHPDDetIds.push_back(*detItr);
0824       }
0825     }
0826 
0827     if (!theHOHPDDetIds.empty())
0828       theHODigitizer->setDetIds(theHOHPDDetIds);
0829     else {
0830       theHODigitizer.reset();
0831     }
0832 
0833     if (!theHOSiPMDetIds.empty())
0834       theHOSiPMDigitizer->setDetIds(theHOSiPMDetIds);
0835     else {
0836       theHOSiPMDigitizer.reset();
0837     }
0838 
0839     if (!theHOHPDDetIds.empty() && !theHOSiPMDetIds.empty()) {
0840       theHOSiPMHitFilter.setDetIds(theHOSiPMDetIds);
0841       theHOHitFilter.setDetIds(theHOHPDDetIds);
0842     }
0843 
0844     theParameterMap.setHOZecotekDetIds(zecotekDetIds);
0845     theParameterMap.setHOHamamatsuDetIds(hamamatsuDetIds);
0846 
0847     // make sure we don't got through this exercise again
0848     theHOSiPMCode = -2;
0849   }
0850 }
0851 
0852 void HcalDigitizer::darkening(std::vector<PCaloHit> &hcalHits) {
0853   for (unsigned int ii = 0; ii < hcalHits.size(); ++ii) {
0854     uint32_t tmpId = hcalHits[ii].id();
0855     int det, z, depth, ieta, phi, lay;
0856     HcalTestNumbering::unpackHcalIndex(tmpId, det, z, depth, ieta, phi, lay);
0857 
0858     bool darkened = false;
0859     float dweight = 1.;
0860 
0861     if (det == int(HcalBarrel) && m_HBDarkening) {
0862       // HB darkening
0863       dweight = m_HBDarkening->degradation(deliveredLumi, ieta, lay);
0864       darkened = true;
0865     } else if (det == int(HcalEndcap) && m_HEDarkening) {
0866       // HE darkening
0867       dweight = m_HEDarkening->degradation(deliveredLumi, ieta, lay);
0868       darkened = true;
0869     } else if (det == int(HcalForward) && m_HFRecalibration) {
0870       // HF darkening - approximate: invert recalibration factor
0871       dweight = 1.0 / m_HFRecalibration->getCorr(ieta, depth, deliveredLumi);
0872       darkened = true;
0873     }
0874 
0875     // reset hit energy
0876     if (darkened)
0877       hcalHits[ii].setEnergy(hcalHits[ii].energy() * dweight);
0878   }
0879 }