Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-06-11 03:00:22

0001 
0002 #include "SimCalorimetry/HGCalSimProducers/interface/HGCFEElectronics.h"
0003 #include "DataFormats/HGCDigi/interface/HGCDigiCollections.h"
0004 #include "FWCore/Utilities/interface/transform.h"
0005 
0006 #include "vdt/vdtMath.h"
0007 #include <cassert>
0008 
0009 using namespace hgc_digi;
0010 //#define EDM_ML_DEBUG
0011 
0012 //
0013 template <class DFr>
0014 HGCFEElectronics<DFr>::HGCFEElectronics(const edm::ParameterSet& ps)
0015     : fwVersion_{ps.getParameter<uint32_t>("fwVersion")},
0016       adcPulse_{},
0017       pulseAvgT_{},
0018       tdcForToAOnset_fC_{},
0019       adcSaturation_fC_{-1.0},
0020       adcLSB_fC_{},
0021       tdcLSB_fC_{},
0022       tdcSaturation_fC_{-1.0},
0023       adcThreshold_fC_{},
0024       tdcOnset_fC_{},
0025       toaLSB_ns_{},
0026       tdcResolutionInNs_{1e-9},  // set time resolution very small by default
0027       targetMIPvalue_ADC_{},
0028       jitterNoise_ns_{},
0029       jitterConstant_ns_{},
0030       eventTimeOffset_ns_{{0.02, 0.02, 0.02}},
0031       noise_fC_{},
0032       toaMode_(WEIGHTEDBYE) {
0033 #ifdef EDM_ML_DEBUG
0034   edm::LogVerbatim("HGCFE") << "[HGCFEElectronics] running with version " << fwVersion_;
0035 #endif
0036   if (ps.exists("adcPulse")) {
0037     auto temp = ps.getParameter<std::vector<double> >("adcPulse");
0038     for (unsigned i = 0; i < temp.size(); ++i) {
0039       adcPulse_[i] = (float)temp[i];
0040     }
0041     // normalize adc pulse
0042     for (unsigned i = 0; i < adcPulse_.size(); ++i) {
0043       adcPulse_[i] = adcPulse_[i] / adcPulse_[2];
0044     }
0045     temp = ps.getParameter<std::vector<double> >("pulseAvgT");
0046     for (unsigned i = 0; i < temp.size(); ++i) {
0047       pulseAvgT_[i] = (float)temp[i];
0048     }
0049   }
0050   if (ps.exists("adcNbits")) {
0051     uint32_t adcNbits = ps.getParameter<uint32_t>("adcNbits");
0052     adcSaturation_fC_ = ps.getParameter<double>("adcSaturation_fC");
0053     adcLSB_fC_ = adcSaturation_fC_ / pow(2., adcNbits);
0054 #ifdef EDM_ML_DEBUG
0055     edm::LogVerbatim("HGCFE") << "[HGCFEElectronics] " << adcNbits << " bit ADC defined"
0056                               << " with LSB=" << adcLSB_fC_ << " saturation to occur @ " << adcSaturation_fC_;
0057 #endif
0058   }
0059 
0060   if (ps.exists("tdcNbits")) {
0061     tdcNbits_ = ps.getParameter<uint32_t>("tdcNbits");
0062     setTDCfsc(ps.getParameter<double>("tdcSaturation_fC"));
0063 #ifdef EDM_ML_DEBUG
0064     edm::LogVerbatim("HGCFE") << "[HGCFEElectronics] " << tdcNbits_ << " bit TDC defined with LSB=" << tdcLSB_fC_
0065                               << " saturation to occur @ " << tdcSaturation_fC_
0066                               << " (NB lowered by 1 part in a million)";
0067 #endif
0068   }
0069   if (ps.exists("targetMIPvalue_ADC"))
0070     targetMIPvalue_ADC_ = ps.getParameter<uint32_t>("targetMIPvalue_ADC");
0071   if (ps.exists("adcThreshold_fC"))
0072     adcThreshold_fC_ = ps.getParameter<double>("adcThreshold_fC");
0073   if (ps.exists("tdcOnset_fC"))
0074     tdcOnset_fC_ = ps.getParameter<double>("tdcOnset_fC");
0075   if (ps.exists("tdcForToAOnset_fC")) {
0076     auto temp = ps.getParameter<std::vector<double> >("tdcForToAOnset_fC");
0077     if (temp.size() == tdcForToAOnset_fC_.size()) {
0078       std::copy_n(temp.begin(), temp.size(), tdcForToAOnset_fC_.begin());
0079     } else {
0080       throw cms::Exception("BadConfiguration") << " HGCFEElectronics wrong size for ToA thresholds ";
0081     }
0082   }
0083   if (ps.exists("toaLSB_ns"))
0084     toaLSB_ns_ = ps.getParameter<double>("toaLSB_ns");
0085   if (ps.exists("tdcChargeDrainParameterisation")) {
0086     for (auto val : ps.getParameter<std::vector<double> >("tdcChargeDrainParameterisation")) {
0087       tdcChargeDrainParameterisation_.push_back((float)val);
0088     }
0089   }
0090   if (ps.exists("tdcResolutionInPs"))
0091     tdcResolutionInNs_ = ps.getParameter<double>("tdcResolutionInPs") * 1e-3;  // convert to ns
0092   if (ps.exists("toaMode"))
0093     toaMode_ = ps.getParameter<uint32_t>("toaMode");
0094 
0095   if (ps.exists("jitterNoise_ns")) {
0096     auto temp = ps.getParameter<std::vector<double> >("jitterNoise_ns");
0097     if (temp.size() == jitterNoise_ns_.size()) {
0098       std::copy_n(temp.begin(), temp.size(), jitterNoise_ns_.begin());
0099     } else {
0100       throw cms::Exception("BadConfiguration") << " HGCFEElectronics wrong size for ToA jitterNoise ";
0101     }
0102   }
0103   if (ps.exists("jitterConstant_ns")) {
0104     auto temp = ps.getParameter<std::vector<double> >("jitterConstant_ns");
0105     if (temp.size() == jitterConstant_ns_.size()) {
0106       std::copy_n(temp.begin(), temp.size(), jitterConstant_ns_.begin());
0107     } else {
0108       throw cms::Exception("BadConfiguration") << " HGCFEElectronics wrong size for ToA jitterConstant ";
0109     }
0110   }
0111 }
0112 
0113 //
0114 template <class DFr>
0115 void HGCFEElectronics<DFr>::runTrivialShaper(
0116     DFr& dataFrame, HGCSimHitData& chargeColl, uint32_t thrADC, float lsbADC, uint32_t gainIdx, float maxADC) {
0117   bool debug(false);
0118 
0119 #ifdef EDM_ML_DEBUG
0120   for (int it = 0; it < (int)(chargeColl.size()); it++)
0121     debug |= (chargeColl[it] > adcThreshold_fC_);
0122 #endif
0123 
0124   if (debug)
0125     edm::LogVerbatim("HGCFE") << "[runTrivialShaper]";
0126 
0127   if (lsbADC < 0)
0128     lsbADC = adcLSB_fC_;
0129   if (maxADC < 0)
0130     // lower adcSaturation_fC_ by one part in a million
0131     // to ensure largest charge converted in bits is 0xfff==4095, not 0x1000
0132     // no effect on charges loewer than; no impact on cpu time, only done once
0133     maxADC = adcSaturation_fC_ * (1 - 1e-6);
0134   for (int it = 0; it < (int)(chargeColl.size()); it++) {
0135     //brute force saturation, maybe could to better with an exponential like saturation
0136     const uint32_t adc = std::floor(std::min(chargeColl[it], maxADC) / lsbADC);
0137     HGCSample newSample;
0138     newSample.set(adc > thrADC, false, gainIdx, 0, adc);
0139     dataFrame.setSample(it, newSample);
0140 
0141     if (debug)
0142       edm::LogVerbatim("HGCFE") << adc << " (" << chargeColl[it] << "/" << adcLSB_fC_ << ") ";
0143   }
0144 
0145   if (debug) {
0146     std::ostringstream msg;
0147     dataFrame.print(msg);
0148     edm::LogVerbatim("HGCFE") << msg.str();
0149   }
0150 }
0151 
0152 //
0153 template <class DFr>
0154 void HGCFEElectronics<DFr>::runSimpleShaper(DFr& dataFrame,
0155                                             HGCSimHitData& chargeColl,
0156                                             uint32_t thrADC,
0157                                             float lsbADC,
0158                                             uint32_t gainIdx,
0159                                             float maxADC,
0160                                             const hgc_digi::FEADCPulseShape& adcPulse) {
0161   //convolute with pulse shape to compute new ADCs
0162   newCharge_.fill(0.f);
0163   bool debug(false);
0164   for (int it = 0; it < (int)(chargeColl.size()); it++) {
0165     const float charge(chargeColl[it]);
0166     if (charge == 0.f)
0167       continue;
0168 
0169 #ifdef EDM_ML_DEBUG
0170     debug |= (charge > adcThreshold_fC_);
0171 #endif
0172 
0173     if (debug)
0174       edm::LogVerbatim("HGCFE") << "\t Redistributing SARS ADC" << charge << " @ " << it;
0175 
0176     for (int ipulse = -2; ipulse < (int)(adcPulse.size()) - 2; ipulse++) {
0177       if (it + ipulse < 0)
0178         continue;
0179       if (it + ipulse >= (int)(dataFrame.size()))
0180         continue;
0181       const float chargeLeak = charge * adcPulse[(ipulse + 2)];
0182       newCharge_[it + ipulse] += chargeLeak;
0183 
0184       if (debug)
0185         edm::LogVerbatim("HGCFE") << " | " << it + ipulse << " " << chargeLeak;
0186     }
0187 
0188     if (debug)
0189       edm::LogVerbatim("HGCFE") << std::endl;
0190   }
0191 
0192   for (int it = 0; it < (int)(newCharge_.size()); it++) {
0193     //brute force saturation, maybe could to better with an exponential like saturation
0194     const uint32_t adc = std::floor(std::min(newCharge_[it], maxADC) / lsbADC);
0195     HGCSample newSample;
0196     newSample.set(adc > thrADC, false, gainIdx, 0, adc);
0197     dataFrame.setSample(it, newSample);
0198 
0199     if (debug)
0200       edm::LogVerbatim("HGCFE") << adc << " (" << std::min(newCharge_[it], maxADC) << "/" << lsbADC << " ) ";
0201   }
0202 
0203   if (debug) {
0204     std::ostringstream msg;
0205     dataFrame.print(msg);
0206     edm::LogVerbatim("HGCFE") << msg.str();
0207   }
0208 }
0209 
0210 //
0211 template <class DFr>
0212 void HGCFEElectronics<DFr>::runShaperWithToT(DFr& dataFrame,
0213                                              HGCSimHitData& chargeColl,
0214                                              HGCSimHitData& toaColl,
0215                                              CLHEP::HepRandomEngine* engine,
0216                                              uint32_t thrADC,
0217                                              float lsbADC,
0218                                              uint32_t gainIdx,
0219                                              float maxADC,
0220                                              int thickness,
0221                                              float tdcOnsetAuto,
0222                                              float noiseWidth,
0223                                              const hgc_digi::FEADCPulseShape& adcPulse) {
0224   busyFlags_.fill(false);
0225   totFlags_.fill(false);
0226   toaFlags_.fill(false);
0227   newCharge_.fill(0.f);
0228   toaFromToT_.fill(0.f);
0229 
0230 #ifdef EDM_ML_DEBUG
0231   constexpr bool debug_state(true);
0232 #else
0233   constexpr bool debug_state(false);
0234 #endif
0235 
0236   bool debug = debug_state;
0237 
0238   float timeToA = 0.f;
0239 
0240   //configure the ADC <-> TDC transition depending on the value passed as argument
0241   double tdcOnset(maxADC);
0242   if (maxADC < 0) {
0243     maxADC = adcSaturation_fC_;
0244     tdcOnset = tdcOnset_fC_;
0245   }
0246 
0247   //configure the ADC LSB depending on the value passed as argument
0248   if (lsbADC < 0)
0249     lsbADC = adcLSB_fC_;
0250 
0251   //first look at time
0252   //for pileup look only at intime signals
0253   //ToA is in central BX if fired -- std::floor(BX/25.)+9;
0254   int fireBX = 9;
0255   //noise fluctuation on charge is added after ToA computation
0256   //to be done properly with realistic ToA shaper and jitter for the moment accounted in the smearing
0257   assert(static_cast<int>(tdcForToAOnset_fC_.size()) > thickness - 1);
0258   assert(static_cast<int>(noise_fC_.size()) > thickness - 1);
0259   assert(static_cast<int>(eventTimeOffset_ns_.size()) > thickness - 1);
0260   if (toaColl[fireBX] != 0.f && chargeColl[fireBX] > tdcForToAOnset_fC_[thickness - 1]) {
0261     timeToA = toaColl[fireBX];
0262     float sensor_noise = noiseWidth <= 0 ? noise_fC_[thickness - 1] : noiseWidth;
0263     float noise = jitterNoise_ns_[thickness - 1] * sensor_noise;
0264     float jitter = chargeColl[fireBX] == 0 ? 0 : (noise / chargeColl[fireBX]);
0265     if (jitter != 0)
0266       timeToA = CLHEP::RandGaussQ::shoot(engine, timeToA, jitter);
0267     else if (tdcResolutionInNs_ != 0)
0268       timeToA = CLHEP::RandGaussQ::shoot(engine, timeToA, tdcResolutionInNs_);
0269     timeToA += eventTimeOffset_ns_[thickness - 1];
0270     if (timeToA >= 0.f && timeToA <= 25.f)
0271       toaFlags_[fireBX] = true;
0272   }
0273 
0274   //now look at charge
0275   //first identify bunches which will trigger ToT
0276   //if(debug_state) edm::LogVerbatim("HGCFE") << "[runShaperWithToT]";
0277   assert(chargeColl.size() <= busyFlags_.size());
0278   assert(chargeColl.size() <= totFlags_.size());
0279   assert(chargeColl.size() <= newCharge_.size());
0280   for (int it = 0; it < (int)(chargeColl.size()); ++it) {
0281     debug = debug_state;
0282     //if already flagged as busy it can't be re-used to trigger the ToT
0283     if (busyFlags_[it])
0284       continue;
0285 
0286     if (tdcOnsetAuto < 0) {
0287       tdcOnsetAuto = tdcOnset_fC_;
0288     }
0289     //if below TDC onset will be handled by SARS ADC later
0290     float charge = chargeColl[it];
0291     if (charge < tdcOnset) {
0292       debug = false;
0293       continue;
0294     }
0295 
0296     //raise TDC mode for charge computation
0297     //ToA anyway fired independently will be sorted out with realistic ToA dedicated shaper
0298     float toa = timeToA;
0299     totFlags_[it] = true;
0300 
0301     if (debug)
0302       edm::LogVerbatim("HGCFE") << "\t q=" << charge << " fC with <toa>=" << toa << " ns, triggers ToT @ " << it;
0303 
0304     //compute total charge to be integrated and integration time
0305     //needs a loop as ToT will last as long as there is charge to dissipate
0306     int busyBxs(0);
0307     float totalCharge(charge), finalToA(toa), integTime(0);
0308     while (true) {
0309       //compute integration time in ns and # bunches
0310       //float newIntegTime(0);
0311       int poffset = 0;
0312       float charge_offset = 0.f;
0313       const float charge_kfC(totalCharge * 1e-3);
0314       if (charge_kfC < tdcChargeDrainParameterisation_[3]) {
0315         //newIntegTime=tdcChargeDrainParameterisation_[0]*pow(charge_kfC,2)+tdcChargeDrainParameterisation_[1]*charge_kfC+tdcChargeDrainParameterisation_[2];
0316       } else if (charge_kfC < tdcChargeDrainParameterisation_[7]) {
0317         poffset = 4;
0318         charge_offset = tdcChargeDrainParameterisation_[3];
0319         //newIntegTime=tdcChargeDrainParameterisation_[4]*pow(charge_kfC-tdcChargeDrainParameterisation_[3],2)+tdcChargeDrainParameterisation_[5]*(charge_kfC-tdcChargeDrainParameterisation_[3])+tdcChargeDrainParameterisation_[6];
0320       } else {
0321         poffset = 8;
0322         charge_offset = tdcChargeDrainParameterisation_[7];
0323         //newIntegTime=tdcChargeDrainParameterisation_[8]*pow(charge_kfC-tdcChargeDrainParameterisation_[7],2)+tdcChargeDrainParameterisation_[9]*(charge_kfC-tdcChargeDrainParameterisation_[7])+tdcChargeDrainParameterisation_[10];
0324       }
0325       const float charge_mod = charge_kfC - charge_offset;
0326       const float newIntegTime =
0327           ((tdcChargeDrainParameterisation_[poffset] * charge_mod + tdcChargeDrainParameterisation_[poffset + 1]) *
0328                charge_mod +
0329            tdcChargeDrainParameterisation_[poffset + 2]);
0330 
0331       const int newBusyBxs = std::floor(newIntegTime / 25.f) + 1;
0332 
0333       //if no update is needed regarding the number of bunches,
0334       //then the ToT integration time has converged
0335       integTime = newIntegTime;
0336       if (newBusyBxs == busyBxs)
0337         break;
0338 
0339       //update charge integrated during ToT
0340       if (debug) {
0341         if (busyBxs == 0)
0342           edm::LogVerbatim("HGCFE") << "\t Intial busy estimate=" << integTime << " ns = " << newBusyBxs << " bxs";
0343         else
0344           edm::LogVerbatim("HGCFE") << "\t ...integrated charge overflows initial busy estimate, interating again";
0345       }
0346 
0347       //update number of busy bunches
0348       busyBxs = newBusyBxs;
0349 
0350       //reset charge to be integrated
0351       totalCharge = charge;
0352       if (toaMode_ == WEIGHTEDBYE)
0353         finalToA = toa * charge;
0354 
0355       //add leakage from previous bunches in SARS ADC mode
0356       for (int jt = 0; jt < it; ++jt) {
0357         const unsigned int deltaT = (it - jt);
0358         if ((deltaT + 2) >= adcPulse.size() || chargeColl[jt] == 0.f || totFlags_[jt] || busyFlags_[jt])
0359           continue;
0360 
0361         const float leakCharge = chargeColl[jt] * adcPulse[deltaT + 2];
0362         totalCharge += leakCharge;
0363         if (toaMode_ == WEIGHTEDBYE)
0364           finalToA += leakCharge * pulseAvgT_[deltaT + 2];
0365 
0366         if (debug)
0367           edm::LogVerbatim("HGCFE") << "\t\t leaking " << chargeColl[jt] << " fC @ deltaT=-" << deltaT << " -> +"
0368                                     << leakCharge << " with avgT=" << pulseAvgT_[deltaT + 2];
0369       }
0370 
0371       //add contamination from posterior bunches
0372       for (int jt = it + 1; jt < it + busyBxs && jt < dataFrame.size(); ++jt) {
0373         //this charge will be integrated in TDC mode
0374         //disable for SARS ADC
0375         busyFlags_[jt] = true;
0376 
0377         const float extraCharge = chargeColl[jt];
0378         if (extraCharge == 0.f)
0379           continue;
0380         if (debug)
0381           edm::LogVerbatim("HGCFE") << "\t\t adding " << extraCharge << " fC @ deltaT=+" << (jt - it);
0382 
0383         totalCharge += extraCharge;
0384         if (toaMode_ == WEIGHTEDBYE)
0385           finalToA += extraCharge * toaColl[jt];
0386       }
0387 
0388       //finalize ToA contamination
0389       if (toaMode_ == WEIGHTEDBYE)
0390         finalToA /= totalCharge;
0391     }
0392 
0393     newCharge_[it] = (totalCharge - tdcOnset);
0394 
0395     if (debug)
0396       edm::LogVerbatim("HGCFE") << "\t Final busy estimate=" << integTime << " ns = " << busyBxs << " bxs" << std::endl
0397                                 << "\t Total integrated=" << totalCharge << " fC <toa>=" << toaFromToT_[it]
0398                                 << " (raw=" << finalToA << ") ns ";
0399 
0400     //last fC (tdcOnset) are dissipated trough pulse
0401     if (it + busyBxs < (int)(newCharge_.size())) {
0402       const float deltaT2nextBx((busyBxs * 25 - integTime));
0403       const float tdcOnsetLeakage(tdcOnset * vdt::fast_expf(-deltaT2nextBx / tdcChargeDrainParameterisation_[11]));
0404       if (debug)
0405         edm::LogVerbatim("HGCFE") << "\t Leaking remainder of TDC onset " << tdcOnset << " fC, to be dissipated in "
0406                                   << deltaT2nextBx << " DeltaT/tau=" << deltaT2nextBx << " / "
0407                                   << tdcChargeDrainParameterisation_[11] << " ns, adds " << tdcOnsetLeakage << " fC @ "
0408                                   << it + busyBxs << " bx (first free bx)";
0409       newCharge_[it + busyBxs] += tdcOnsetLeakage;
0410     }
0411   }
0412 
0413   //including the leakage from bunches in SARS ADC when not declared busy or in ToT
0414   auto runChargeSharing = [&]() {
0415     int ipulse = 0;
0416     for (int it = 0; it < (int)(chargeColl.size()); ++it) {
0417       //if busy, charge has been already integrated
0418       //if(debug) edm::LogVerbatim("HGCFE") << "\t SARS ADC pulse activated @ " << it << " : ";
0419       if (!totFlags_[it] && !busyFlags_[it]) {
0420         const int start = std::max(0, 2 - it);
0421         const int stop = std::min((int)adcPulse.size(), (int)newCharge_.size() - it + 2);
0422         for (ipulse = start; ipulse < stop; ++ipulse) {
0423           const int itoffset = it + ipulse - 2;
0424           //notice that if the channel is already busy,
0425           //it has already been affected by the leakage of the SARS ADC
0426           //if(totFlags_[itoffset] || busyFlags_[itoffset]) continue;
0427           if (!totFlags_[itoffset] && !busyFlags_[itoffset]) {
0428             newCharge_[itoffset] += chargeColl[it] * adcPulse[ipulse];
0429           }
0430           //if(debug) edm::LogVerbatim("HGCFE") << " | " << itoffset << " " << chargeColl[it]*adcPulse[ipulse] << "( " << chargeColl[it] << "->";
0431           //if(debug) edm::LogVerbatim("HGCFE") << newCharge_[itoffset] << ") ";
0432         }
0433       }
0434 
0435       if (debug)
0436         edm::LogVerbatim("HGCFE") << std::endl;
0437     }
0438   };
0439   runChargeSharing();
0440 
0441   //For the future need to understand how to deal with toa for out of time signals
0442   //and for that should keep track of the BX firing the ToA somewhere (also to restore the use of finalToA)
0443   /*
0444   float finalToA(0.);
0445   for(int it=0; it<(int)(newCharge_.size()); it++){
0446     if(toaFlags_[it]){
0447       finalToA = toaFromToT_[it];
0448       //to avoid +=25 for small negative time taken as 0
0449       while(finalToA < -1.e-5)  finalToA+=25.f;
0450       while(finalToA > 25.f) finalToA-=25.f;
0451       toaFromToT_[it] = finalToA;
0452     }
0453   }
0454   */
0455   //timeToA is already in 0-25ns range by construction
0456 
0457   //set new ADCs and ToA
0458   if (debug)
0459     edm::LogVerbatim("HGCFE") << "\t final result : ";
0460   assert(chargeColl.size() >= newCharge_.size());
0461   for (int it = 0; it < (int)(newCharge_.size()); it++) {
0462     if (debug)
0463       edm::LogVerbatim("HGCFE") << chargeColl[it] << " -> " << newCharge_[it] << " ";
0464 
0465     HGCSample newSample;
0466     if (totFlags_[it] || busyFlags_[it]) {
0467       if (totFlags_[it]) {
0468         //brute force saturation, maybe could to better with an exponential like saturation
0469         const float saturatedCharge(std::min(newCharge_[it], tdcSaturation_fC_));
0470         //working version for in-time PU and signal
0471         newSample.set(
0472             true, true, gainIdx, (uint16_t)(timeToA / toaLSB_ns_), (uint16_t)(std::floor(saturatedCharge / tdcLSB_fC_)));
0473         if (toaFlags_[it])
0474           newSample.setToAValid(true);
0475       } else {
0476         newSample.set(false, true, gainIdx, 0, 0);
0477       }
0478     } else {
0479       //brute force saturation, maybe could to better with an exponential like saturation
0480       const uint16_t adc = std::floor(std::min(newCharge_[it], maxADC) / lsbADC);
0481       //working version for in-time PU and signal
0482       newSample.set(adc > thrADC, false, gainIdx, (uint16_t)(timeToA / toaLSB_ns_), adc);
0483       if (toaFlags_[it])
0484         newSample.setToAValid(true);
0485     }
0486     dataFrame.setSample(it, newSample);
0487   }
0488 
0489   if (debug) {
0490     std::ostringstream msg;
0491     dataFrame.print(msg);
0492     edm::LogVerbatim("HGCFE") << msg.str();
0493   }
0494 }
0495 
0496 // cause the compiler to generate the appropriate code
0497 #include "DataFormats/HGCDigi/interface/HGCDigiCollections.h"
0498 template class HGCFEElectronics<HGCalDataFrame>;