CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_6_1_2_SLHC4_patch1/src/SimCalorimetry/HcalSimProducers/src/HcalDigitizer.cc

Go to the documentation of this file.
00001 #include "SimCalorimetry/HcalSimProducers/interface/HcalDigitizer.h"
00002 #include "SimCalorimetry/HcalSimProducers/src/HcalTestHitGenerator.h"
00003 #include "SimDataFormats/CaloHit/interface/PCaloHitContainer.h"
00004 #include "SimCalorimetry/HcalSimAlgos/interface/HcalSimParameterMap.h"
00005 #include "SimCalorimetry/HcalSimAlgos/interface/HcalShapes.h"
00006 #include "SimCalorimetry/HcalSimAlgos/interface/HcalElectronicsSim.h"
00007 #include "SimCalorimetry/CaloSimAlgos/interface/CaloHitResponse.h"
00008 #include "SimCalorimetry/HcalSimAlgos/interface/HcalAmplifier.h"
00009 #include "SimCalorimetry/HcalSimAlgos/interface/HcalCoderFactory.h"
00010 #include "SimCalorimetry/HcalSimAlgos/interface/HcalHitCorrection.h"
00011 #include "SimCalorimetry/HcalSimAlgos/interface/HcalTimeSlewSim.h"
00012 #include "SimCalorimetry/HcalSimAlgos/interface/HcalSimParameterMap.h"
00013 #include "SimCalorimetry/HcalSimAlgos/interface/HcalSiPMHitResponse.h"
00014 #include "SimCalorimetry/HcalSimAlgos/interface/HPDIonFeedbackSim.h"
00015 #include "DataFormats/Common/interface/Handle.h"
00016 #include "DataFormats/Common/interface/Wrapper.h"
00017 #include "FWCore/Framework/interface/ESHandle.h"
00018 #include "FWCore/Framework/interface/Event.h"
00019 #include "SimGeneral/MixingModule/interface/PileUpEventPrincipal.h"
00020 #include "FWCore/Framework/interface/EventSetup.h"
00021 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00022 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00023 #include "SimCalorimetry/CaloSimAlgos/interface/CaloTDigitizer.h"
00024 #include "DataFormats/HcalDigi/interface/HcalDigiCollections.h"
00025 #include "Geometry/Records/interface/CaloGeometryRecord.h"
00026 #include "CalibFormats/HcalObjects/interface/HcalDbService.h"
00027 #include "CalibFormats/HcalObjects/interface/HcalDbRecord.h"
00028 #include "SimDataFormats/CrossingFrame/interface/CrossingFrame.h"
00029 #include "SimDataFormats/CrossingFrame/interface/MixCollection.h"
00030 #include "FWCore/Utilities/interface/RandomNumberGenerator.h"
00031 #include "FWCore/ServiceRegistry/interface/Service.h"
00032 #include "DataFormats/HcalDetId/interface/HcalZDCDetId.h"
00033 #include "SimCalorimetry/HcalSimAlgos/interface/HPDNoiseGenerator.h"
00034 #include "CondFormats/HcalObjects/interface/HcalCholeskyMatrix.h"
00035 #include "CondFormats/HcalObjects/interface/HcalCholeskyMatrices.h"
00036 #include <boost/foreach.hpp>
00037 #include "Geometry/CaloTopology/interface/HcalTopology.h"
00038 
00039 namespace HcalDigitizerImpl {
00040 
00041   template<typename SIPMDIGITIZER>
00042   void fillSiPMCells(std::vector<int> & siPMCells, SIPMDIGITIZER * siPMDigitizer)
00043   {
00044     std::vector<DetId> siPMDetIds;
00045     siPMDetIds.reserve(siPMCells.size());
00046     for(std::vector<int>::const_iterator idItr = siPMCells.begin();
00047         idItr != siPMCells.end(); ++idItr)
00048     {
00049       siPMDetIds.emplace_back(*idItr);
00050     }
00051     siPMDigitizer->setDetIds(siPMDetIds);
00052   }
00053 
00054   // if both exist, assume the SiPM one has cells filled, and
00055   // assign the rest to the HPD
00056   template<typename HPDDIGITIZER, typename SIPMDIGITIZER>
00057   void fillCells(std::vector<DetId>& allCells,
00058                  HPDDIGITIZER * hpdDigitizer,
00059                  SIPMDIGITIZER * siPMDigitizer)
00060   {
00061     // if both digitizers exist, split up the cells
00062     if(siPMDigitizer && hpdDigitizer)
00063     {
00064       std::vector<DetId> siPMDetIds = siPMDigitizer->detIds();
00065       std::sort(siPMDetIds.begin(), siPMDetIds.end());
00066       std::vector<DetId> sortedCells = allCells;
00067       std::sort(sortedCells.begin(), sortedCells.end());
00068       std::vector<DetId> hpdCells;
00069       std::set_difference(sortedCells.begin(), sortedCells.end(),
00070                           siPMDetIds.begin(), siPMDetIds.end(),
00071                           std::back_inserter(hpdCells) );
00072       hpdDigitizer->setDetIds(hpdCells);
00073     }
00074     else
00075     {
00076       if(siPMDigitizer) siPMDigitizer->setDetIds(allCells);
00077       if(hpdDigitizer) hpdDigitizer->setDetIds(allCells);
00078     }
00079   }
00080 } // namespace HcaiDigitizerImpl
00081 
00082 
00083 HcalDigitizer::HcalDigitizer(const edm::ParameterSet& ps) :
00084   theGeometry(0),
00085   theParameterMap(new HcalSimParameterMap(ps)),
00086   theShapes(new HcalShapes()),
00087   theHBHEResponse(0),
00088   theHBHESiPMResponse(0),
00089   theHOResponse(0),   
00090   theHOSiPMResponse(0),
00091   theHFResponse(new CaloHitResponse(theParameterMap, theShapes)),
00092   theZDCResponse(new CaloHitResponse(theParameterMap, theShapes)),
00093   theHBHEAmplifier(0),
00094   theHFAmplifier(0),
00095   theHOAmplifier(0),
00096   theZDCAmplifier(0),
00097   theIonFeedback(0),
00098   theCoderFactory(0),
00099   theUpgradeCoderFactory(0),
00100   theHBHEElectronicsSim(0),
00101   theHFElectronicsSim(0),
00102   theHOElectronicsSim(0),
00103   theZDCElectronicsSim(0),
00104   theUpgradeHBHEElectronicsSim(0),
00105   theUpgradeHFElectronicsSim(0),
00106   theHBHEHitFilter(),
00107   theHFHitFilter(ps.getParameter<bool>("doHFWindow")),
00108   theHOHitFilter(),
00109   theHOSiPMHitFilter(HcalOuter),
00110   theZDCHitFilter(),
00111   theHitCorrection(0),
00112   theNoiseGenerator(0),
00113   theNoiseHitGenerator(0),
00114   theHBHEDigitizer(0),
00115   theHBHESiPMDigitizer(0),
00116   theHODigitizer(0),
00117   theHOSiPMDigitizer(0),
00118   theHFDigitizer(0),
00119   theZDCDigitizer(0),
00120   theHBHEUpgradeDigitizer(0),
00121   theHFUpgradeDigitizer(0),
00122   theRelabeller(0),
00123   isZDC(true),
00124   isHCAL(true),
00125   zdcgeo(true),
00126   hbhegeo(true),
00127   hogeo(true),
00128   hfgeo(true),
00129   theHOSiPMCode(ps.getParameter<edm::ParameterSet>("ho").getParameter<int>("siPMCode"))
00130 {
00131   bool doNoise = ps.getParameter<bool>("doNoise");
00132   bool useOldNoiseHB = ps.getParameter<bool>("useOldHB");
00133   bool useOldNoiseHE = ps.getParameter<bool>("useOldHE");
00134   bool useOldNoiseHF = ps.getParameter<bool>("useOldHF");
00135   bool useOldNoiseHO = ps.getParameter<bool>("useOldHO");
00136   bool doEmpty = ps.getParameter<bool>("doEmpty");
00137   double HBtp = ps.getParameter<double>("HBTuningParameter");
00138   double HEtp = ps.getParameter<double>("HETuningParameter");
00139   double HFtp = ps.getParameter<double>("HFTuningParameter");
00140   double HOtp = ps.getParameter<double>("HOTuningParameter");
00141   bool doHBHEUpgrade = ps.getParameter<bool>("HBHEUpgradeQIE");
00142   bool doHFUpgrade   = ps.getParameter<bool>("HFUpgradeQIE");
00143 
00144   // need to make copies, because they might get different noise generators
00145   theHBHEAmplifier = new HcalAmplifier(theParameterMap, doNoise);
00146   theHFAmplifier = new HcalAmplifier(theParameterMap, doNoise);
00147   theHOAmplifier = new HcalAmplifier(theParameterMap, doNoise);
00148   theZDCAmplifier = new HcalAmplifier(theParameterMap, doNoise);
00149   theHBHEAmplifier->setHBtuningParameter(HBtp);
00150   theHBHEAmplifier->setHEtuningParameter(HEtp);
00151   theHFAmplifier->setHFtuningParameter(HFtp);
00152   theHOAmplifier->setHOtuningParameter(HOtp);
00153   theHBHEAmplifier->setUseOldHB(useOldNoiseHB);
00154   theHBHEAmplifier->setUseOldHE(useOldNoiseHE);
00155   theHFAmplifier->setUseOldHF(useOldNoiseHF);
00156   theHOAmplifier->setUseOldHO(useOldNoiseHO);
00157 
00158   theCoderFactory = new HcalCoderFactory(HcalCoderFactory::DB);
00159   theUpgradeCoderFactory = new HcalCoderFactory(HcalCoderFactory::UPGRADE);
00160 
00161   //  std::cout << "HcalDigitizer: theUpgradeCoderFactory created" << std::endl;
00162 
00163   theHBHEElectronicsSim = new HcalElectronicsSim(theHBHEAmplifier, theCoderFactory);
00164   theHFElectronicsSim = new HcalElectronicsSim(theHFAmplifier, theCoderFactory);
00165   theHOElectronicsSim = new HcalElectronicsSim(theHOAmplifier, theCoderFactory);
00166   theZDCElectronicsSim = new HcalElectronicsSim(theZDCAmplifier, theCoderFactory);
00167   theUpgradeHBHEElectronicsSim = new HcalElectronicsSim(theHBHEAmplifier, theUpgradeCoderFactory);
00168   theUpgradeHFElectronicsSim = new HcalElectronicsSim(theHFAmplifier, theUpgradeCoderFactory);
00169 
00170   //  std::cout << "HcalDigitizer: theUpgradeElectronicsSim created" <<  std::endl; 
00171 
00172   // a code of 1 means make all cells SiPM
00173   std::vector<int> hbSiPMCells(ps.getParameter<edm::ParameterSet>("hb").getParameter<std::vector<int> >("siPMCells"));
00174   //std::vector<int> hoSiPMCells(ps.getParameter<edm::ParameterSet>("ho").getParameter<std::vector<int> >("siPMCells"));
00175   // 0 means none, 1 means all, and 2 means use hardcoded
00176 
00177   //  std::cout << std::endl << " hbSiPMCells = " << hbSiPMCells[0] << std::endl;
00178 
00179   bool doHBHEHPD = hbSiPMCells.empty() || (hbSiPMCells[0] != 1);
00180   bool doHOHPD = (theHOSiPMCode != 1);
00181   bool doHBHESiPM = !hbSiPMCells.empty();
00182   bool doHOSiPM = (theHOSiPMCode != 0);
00183   if(doHBHEHPD) {
00184     theHBHEResponse = new CaloHitResponse(theParameterMap, theShapes);
00185     edm::LogInfo("HcalDigitizer") <<"Set scale for HB towers";
00186     theHBHEResponse->initHBHEScale();
00187 
00188     theHBHEResponse->setHitFilter(&theHBHEHitFilter);
00189     theHBHEDigitizer = new HBHEDigitizer(theHBHEResponse, theHBHEElectronicsSim, doEmpty);
00190     bool    changeResponse = ps.getParameter<bool>("ChangeResponse");
00191     edm::FileInPath fname  = ps.getParameter<edm::FileInPath>("CorrFactorFile");
00192     if (changeResponse) {
00193       std::string corrFileName = fname.fullPath();
00194       edm::LogInfo("HcalDigitizer") << "Set scale for HB towers from " << corrFileName;
00195       theHBHEResponse->setHBHEScale(corrFileName); //GMA
00196     }
00197   }
00198   if(doHOHPD) {
00199     theHOResponse = new CaloHitResponse(theParameterMap, theShapes);
00200     theHOResponse->setHitFilter(&theHOHitFilter);
00201     theHODigitizer = new HODigitizer(theHOResponse, theHOElectronicsSim, doEmpty);
00202   }
00203 
00204   if(doHBHESiPM) {
00205     theHBHESiPMResponse = new HcalSiPMHitResponse(theParameterMap, theShapes);
00206     theHBHESiPMResponse->setHitFilter(&theHBHEHitFilter);
00207     if (doHBHEUpgrade) {
00208       theHBHEUpgradeDigitizer = new UpgradeDigitizer(theHBHESiPMResponse, theUpgradeHBHEElectronicsSim, doEmpty);
00209 
00210       //      std::cout << "HcalDigitizer: theHBHEUpgradeDigitizer created" << std::endl;
00211 
00212     } else {
00213       theHBHESiPMDigitizer = new HBHEDigitizer(theHBHESiPMResponse, theHBHEElectronicsSim, doEmpty);
00214     }
00215   }
00216   if(doHOSiPM) {
00217     theHOSiPMResponse = new HcalSiPMHitResponse(theParameterMap, theShapes);
00218     theHOSiPMResponse->setHitFilter(&theHOSiPMHitFilter);
00219     theHOSiPMDigitizer = new HODigitizer(theHOSiPMResponse, theHOElectronicsSim, doEmpty);
00220   }
00221 
00222   // if both are present, fill the SiPM cells now
00223   if(doHBHEHPD && doHBHESiPM) {
00224 
00225     //    std::cout << "HcalDigitizer:  fill the SiPM cells now"  << std::endl;
00226 
00227     HcalDigitizerImpl::fillSiPMCells(hbSiPMCells, theHBHESiPMDigitizer);
00228   }
00229 
00230   theHFResponse->setHitFilter(&theHFHitFilter);
00231   theZDCResponse->setHitFilter(&theZDCHitFilter);
00232 
00233   bool doTimeSlew = ps.getParameter<bool>("doTimeSlew");
00234   if(doTimeSlew) {
00235     // no time slewing for HF
00236     theTimeSlewSim = new HcalTimeSlewSim(theParameterMap);
00237     theHBHEAmplifier->setTimeSlewSim(theTimeSlewSim);
00238     theHOAmplifier->setTimeSlewSim(theTimeSlewSim);
00239     theZDCAmplifier->setTimeSlewSim(theTimeSlewSim);
00240   }
00241 
00242   if (doHFUpgrade) {
00243     theHFUpgradeDigitizer = new UpgradeDigitizer(theHFResponse, theUpgradeHFElectronicsSim, doEmpty);
00244 
00245     //    std::cout << "HcalDigitizer: theHFUpgradeDigitizer created" << std::endl;
00246 
00247   } else {
00248     theHFDigitizer = new HFDigitizer(theHFResponse, theHFElectronicsSim, doEmpty);
00249   }
00250   theZDCDigitizer = new ZDCDigitizer(theZDCResponse, theZDCElectronicsSim, doEmpty);
00251 
00252   edm::ParameterSet ps0 = ps.getParameter<edm::ParameterSet>("HcalReLabel");
00253   relabel_ = ps0.getUntrackedParameter<bool>("RelabelHits");
00254   if (relabel_) {
00255     theRelabeller=new HcalHitRelabeller(ps0.getUntrackedParameter<edm::ParameterSet>("RelabelRules"));
00256   }     
00257 
00258   bool doHPDNoise = ps.getParameter<bool>("doHPDNoise");
00259   if(doHPDNoise) {
00260     //edm::ParameterSet hpdNoisePset = ps.getParameter<edm::ParameterSet>("HPDNoiseLibrary");
00261     theNoiseGenerator = new HPDNoiseGenerator(ps); 
00262     if(theHBHEDigitizer) theHBHEDigitizer->setNoiseSignalGenerator(theNoiseGenerator);
00263     if(theHBHESiPMDigitizer) theHBHESiPMDigitizer->setNoiseSignalGenerator(theNoiseGenerator);
00264   }
00265 
00266   if(ps.getParameter<bool>("doIonFeedback") && theHBHEResponse) {
00267     theIonFeedback = new HPDIonFeedbackSim(ps, theShapes);
00268     theHBHEResponse->setPECorrection(theIonFeedback);
00269     if(ps.getParameter<bool>("doThermalNoise")) {
00270       theHBHEAmplifier->setIonFeedbackSim(theIonFeedback);
00271     }
00272   }
00273 
00274   if(ps.getParameter<bool>("injectTestHits") ) {
00275     theNoiseHitGenerator = new HcalTestHitGenerator(ps);
00276     if(theHBHEDigitizer) theHBHEDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00277     if(theHBHESiPMDigitizer) theHBHESiPMDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00278     if(theHODigitizer) theHODigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00279     if(theHOSiPMDigitizer) theHOSiPMDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00280     if(theHBHEUpgradeDigitizer) {
00281       theHBHEUpgradeDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00282 
00283       //      std::cout << "HcalDigitizer: theHBHEUpgradeDigitizer setNoise" 
00284       //                << std::endl; 
00285     }
00286     if(theHFDigitizer) theHFDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00287     if(theHFUpgradeDigitizer) { 
00288       theHFUpgradeDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00289 
00290       //      std::cout << "HcalDigitizer: theHFUpgradeDigitizer setNoise" 
00291       //                << std::endl;
00292     }
00293     theZDCDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00294   }
00295 
00296   edm::Service<edm::RandomNumberGenerator> rng;
00297   if ( ! rng.isAvailable()) {
00298     throw cms::Exception("Configuration")
00299       << "HcalDigitizer requires the RandomNumberGeneratorService\n"
00300          "which is not present in the configuration file.  You must add the service\n"
00301          "in the configuration file or remove the modules that require it.";
00302   }
00303 
00304   CLHEP::HepRandomEngine& engine = rng->getEngine();
00305   if(theHBHEDigitizer) theHBHEDigitizer->setRandomEngine(engine);
00306   if(theHBHESiPMDigitizer) theHBHESiPMDigitizer->setRandomEngine(engine);
00307   if(theHODigitizer) theHODigitizer->setRandomEngine(engine);
00308   if(theHOSiPMDigitizer) theHOSiPMDigitizer->setRandomEngine(engine);
00309   if(theHBHEUpgradeDigitizer) theHBHEUpgradeDigitizer->setRandomEngine(engine);
00310   if(theIonFeedback) theIonFeedback->setRandomEngine(engine);
00311   if(theTimeSlewSim) theTimeSlewSim->setRandomEngine(engine);
00312   if(theHFUpgradeDigitizer) theHFUpgradeDigitizer->setRandomEngine(engine);
00313   if(theHFDigitizer) theHFDigitizer->setRandomEngine(engine);
00314   theZDCDigitizer->setRandomEngine(engine);
00315 
00316   if (theHitCorrection!=0) theHitCorrection->setRandomEngine(engine);
00317 
00318   hitsProducer_ = ps.getParameter<std::string>("hitsProducer");
00319 }
00320 
00321 
00322 HcalDigitizer::~HcalDigitizer() {
00323   if(theHBHEDigitizer)         delete theHBHEDigitizer;
00324   if(theHBHESiPMDigitizer)     delete theHBHESiPMDigitizer;
00325   if(theHODigitizer)           delete theHODigitizer;
00326   if(theHOSiPMDigitizer)       delete theHOSiPMDigitizer;
00327   if(theHFDigitizer)           delete theHFDigitizer;
00328   delete theZDCDigitizer;
00329   if(theHBHEUpgradeDigitizer)  delete theHBHEUpgradeDigitizer;
00330   if(theHFUpgradeDigitizer)    delete theHFUpgradeDigitizer;
00331   delete theParameterMap;
00332   delete theHBHEResponse;
00333   delete theHBHESiPMResponse;
00334   delete theHOResponse;
00335   delete theHOSiPMResponse;
00336   delete theHFResponse;
00337   delete theZDCResponse;
00338   delete theHBHEElectronicsSim;
00339   delete theHFElectronicsSim;
00340   delete theHOElectronicsSim;
00341   delete theZDCElectronicsSim;
00342   delete theUpgradeHBHEElectronicsSim;
00343   delete theUpgradeHFElectronicsSim;
00344   delete theHBHEAmplifier;
00345   delete theHFAmplifier;
00346   delete theHOAmplifier;
00347   delete theZDCAmplifier;
00348   delete theCoderFactory;
00349   delete theUpgradeCoderFactory;
00350   delete theHitCorrection;
00351   delete theNoiseGenerator;
00352   if (theRelabeller)           delete theRelabeller;
00353 }
00354 
00355 
00356 void HcalDigitizer::setHBHENoiseSignalGenerator(HcalBaseSignalGenerator * noiseGenerator) {
00357   noiseGenerator->setParameterMap(theParameterMap);
00358   noiseGenerator->setElectronicsSim(theHBHEElectronicsSim);
00359   if (theHBHEDigitizer) theHBHEDigitizer->setNoiseSignalGenerator(noiseGenerator);
00360   theHBHEAmplifier->setNoiseSignalGenerator(noiseGenerator);
00361 }
00362 
00363 void HcalDigitizer::setHFNoiseSignalGenerator(HcalBaseSignalGenerator * noiseGenerator) {
00364   noiseGenerator->setParameterMap(theParameterMap);
00365   noiseGenerator->setElectronicsSim(theHFElectronicsSim);
00366   if(theHFDigitizer) theHFDigitizer->setNoiseSignalGenerator(noiseGenerator);
00367   if(theHFUpgradeDigitizer) theHFUpgradeDigitizer->setNoiseSignalGenerator(noiseGenerator);
00368   theHFAmplifier->setNoiseSignalGenerator(noiseGenerator);
00369 }
00370 
00371 void HcalDigitizer::setHONoiseSignalGenerator(HcalBaseSignalGenerator * noiseGenerator) {
00372   noiseGenerator->setParameterMap(theParameterMap);
00373   noiseGenerator->setElectronicsSim(theHOElectronicsSim);
00374   theHODigitizer->setNoiseSignalGenerator(noiseGenerator);
00375   theHOAmplifier->setNoiseSignalGenerator(noiseGenerator);
00376 }
00377 
00378 void HcalDigitizer::setZDCNoiseSignalGenerator(HcalBaseSignalGenerator * noiseGenerator) {
00379   noiseGenerator->setParameterMap(theParameterMap);
00380   noiseGenerator->setElectronicsSim(theZDCElectronicsSim);
00381   theZDCDigitizer->setNoiseSignalGenerator(noiseGenerator);
00382   theZDCAmplifier->setNoiseSignalGenerator(noiseGenerator);
00383 }
00384 
00385 void HcalDigitizer::initializeEvent(edm::Event const& e, edm::EventSetup const& eventSetup) {
00386   // get the appropriate gains, noises, & widths for this event
00387   edm::ESHandle<HcalDbService> conditions;
00388   eventSetup.get<HcalDbRecord>().get(conditions);
00389   theHBHEAmplifier->setDbService(conditions.product());
00390   theHFAmplifier->setDbService(conditions.product());
00391   theHOAmplifier->setDbService(conditions.product());
00392   theZDCAmplifier->setDbService(conditions.product());
00393   theUpgradeHBHEElectronicsSim->setDbService(conditions.product());
00394   theUpgradeHFElectronicsSim->setDbService(conditions.product());
00395 
00396   theCoderFactory->setDbService(conditions.product());
00397   theUpgradeCoderFactory->setDbService(conditions.product());
00398   theParameterMap->setDbService(conditions.product());
00399 
00400   edm::ESHandle<HcalCholeskyMatrices> refCholesky;
00401   eventSetup.get<HcalCholeskyMatricesRcd>().get(refCholesky);
00402   const HcalCholeskyMatrices * myCholesky = refCholesky.product();
00403 
00404   edm::ESHandle<HcalPedestals> pedshandle;
00405   eventSetup.get<HcalPedestalsRcd>().get(pedshandle);
00406   const HcalPedestals *  myADCPedestals = pedshandle.product();
00407 
00408   theHBHEAmplifier->setCholesky(myCholesky);
00409   theHFAmplifier->setCholesky(myCholesky);
00410   theHOAmplifier->setCholesky(myCholesky);
00411   
00412   theHBHEAmplifier->setADCPeds(myADCPedestals);
00413   theHFAmplifier->setADCPeds(myADCPedestals);
00414   theHOAmplifier->setADCPeds(myADCPedestals);
00415 
00416   if(theHitCorrection != 0) {
00417     theHitCorrection->clear();
00418   }
00419 
00420   //initialize hits
00421   if(theHBHEDigitizer) theHBHEDigitizer->initializeHits();
00422   if(theHBHESiPMDigitizer) theHBHESiPMDigitizer->initializeHits();
00423   if(theHODigitizer) theHODigitizer->initializeHits();
00424   if(theHOSiPMDigitizer) theHOSiPMDigitizer->initializeHits();
00425   if(theHBHEUpgradeDigitizer) theHBHEUpgradeDigitizer->initializeHits();
00426   if(theHFUpgradeDigitizer) theHFUpgradeDigitizer->initializeHits();
00427   if(theHFDigitizer) theHFDigitizer->initializeHits();
00428   theZDCDigitizer->initializeHits();
00429 
00430 }
00431 
00432 void HcalDigitizer::accumulateCaloHits(edm::Handle<std::vector<PCaloHit> > const& hcalHandle, edm::Handle<std::vector<PCaloHit> > const& zdcHandle, int bunchCrossing) {
00433   // Step A: pass in inputs, and accumulate digirs
00434   if(isHCAL) {
00435     std::vector<PCaloHit> hcalHits = *hcalHandle.product();
00436     if (relabel_) {
00437       // Relabel PCaloHits
00438       edm::LogInfo("HcalDigitizer") << "Calling Relabller";
00439       theRelabeller->process(hcalHits);
00440     }
00441     if(theHitCorrection != 0) {
00442       theHitCorrection->fillChargeSums(hcalHits);
00443     }
00444     if(hbhegeo) {
00445       if(theHBHEDigitizer) theHBHEDigitizer->add(hcalHits, bunchCrossing);
00446       if(theHBHESiPMDigitizer) theHBHESiPMDigitizer->add(hcalHits, bunchCrossing);
00447       if(theHBHEUpgradeDigitizer) {
00448         //      std::cout << "HcalDigitizer::accumulateCaloHits  theHBHEUpgradeDigitizer->add" << std::endl;     
00449       theHBHEUpgradeDigitizer->add(hcalHits, bunchCrossing);
00450       }
00451     }
00452 
00453     if(hogeo) {
00454       if(theHODigitizer) theHODigitizer->add(hcalHits, bunchCrossing);
00455       if(theHOSiPMDigitizer) theHOSiPMDigitizer->add(hcalHits, bunchCrossing);
00456     }
00457 
00458     if(hfgeo) {
00459       if(theHFDigitizer) theHFDigitizer->add(hcalHits, bunchCrossing);
00460       if(theHFUpgradeDigitizer) theHFUpgradeDigitizer->add(hcalHits, bunchCrossing);
00461     } 
00462   } else {
00463     edm::LogInfo("HcalDigitizer") << "We don't have HCAL hit collection available ";
00464   }
00465 
00466   if(isZDC) {
00467     if(zdcgeo) {
00468       theZDCDigitizer->add(*zdcHandle.product(), bunchCrossing);
00469     } 
00470   } else {
00471     edm::LogInfo("HcalDigitizer") << "We don't have ZDC hit collection available ";
00472   }
00473 }
00474 
00475 void HcalDigitizer::accumulate(edm::Event const& e, edm::EventSetup const& eventSetup) {
00476   // Step A: Get Inputs
00477   edm::InputTag zdcTag(hitsProducer_, "ZDCHITS");
00478   edm::Handle<std::vector<PCaloHit> > zdcHandle;
00479   e.getByLabel(zdcTag, zdcHandle);
00480   isZDC = zdcHandle.isValid();
00481 
00482   edm::InputTag hcalTag(hitsProducer_, "HcalHits");
00483   edm::Handle<std::vector<PCaloHit> > hcalHandle;
00484   e.getByLabel(hcalTag, hcalHandle);
00485   isHCAL = hcalHandle.isValid();
00486 
00487   accumulateCaloHits(hcalHandle, zdcHandle, 0);
00488 }
00489 
00490 void HcalDigitizer::accumulate(PileUpEventPrincipal const& e, edm::EventSetup const& eventSetup) {
00491   // Step A: Get Inputs
00492   edm::InputTag zdcTag(hitsProducer_, "ZDCHITS");
00493   edm::Handle<std::vector<PCaloHit> > zdcHandle;
00494   e.getByLabel(zdcTag, zdcHandle);
00495   isZDC = zdcHandle.isValid();
00496 
00497   edm::InputTag hcalTag(hitsProducer_, "HcalHits");
00498   edm::Handle<std::vector<PCaloHit> > hcalHandle;
00499   e.getByLabel(hcalTag, hcalHandle);
00500   isHCAL = hcalHandle.isValid();
00501 
00502   accumulateCaloHits(hcalHandle, zdcHandle, e.bunchCrossing());
00503 }
00504 
00505 void HcalDigitizer::finalizeEvent(edm::Event& e, const edm::EventSetup& eventSetup) {
00506 
00507   // Step B: Create empty output
00508   std::auto_ptr<HBHEDigiCollection> hbheResult(new HBHEDigiCollection());
00509   std::auto_ptr<HODigiCollection> hoResult(new HODigiCollection());
00510   std::auto_ptr<HFDigiCollection> hfResult(new HFDigiCollection());
00511   std::auto_ptr<ZDCDigiCollection> zdcResult(new ZDCDigiCollection());
00512   std::auto_ptr<HBHEUpgradeDigiCollection> hbheupgradeResult(new HBHEUpgradeDigiCollection());
00513   std::auto_ptr<HFUpgradeDigiCollection> hfupgradeResult(new HFUpgradeDigiCollection());
00514 
00515   // Step C: Invoke the algorithm, getting back outputs.
00516   if(isHCAL&&hbhegeo){
00517     if(theHBHEDigitizer)        theHBHEDigitizer->run(*hbheResult);
00518     if(theHBHESiPMDigitizer)    theHBHESiPMDigitizer->run(*hbheResult);
00519     if(theHBHEUpgradeDigitizer) {
00520       theHBHEUpgradeDigitizer->run(*hbheupgradeResult);
00521 
00522       //      std::cout << "HcalDigitizer::finalizeEvent  theHBHEUpgradeDigitizer->run" << std::endl;     
00523     }
00524   }
00525   if(isHCAL&&hogeo) {
00526     if(theHODigitizer) theHODigitizer->run(*hoResult);
00527     if(theHOSiPMDigitizer) theHOSiPMDigitizer->run(*hoResult);
00528   }
00529   if(isHCAL&&hfgeo) {
00530     if(theHFDigitizer) theHFDigitizer->run(*hfResult);
00531     if(theHFUpgradeDigitizer) theHFUpgradeDigitizer->run(*hfupgradeResult);
00532   }
00533   if(isZDC&&zdcgeo) {
00534     theZDCDigitizer->run(*zdcResult);
00535   }
00536   
00537   edm::LogInfo("HcalDigitizer") << "HCAL HBHE digis : " << hbheResult->size();
00538   edm::LogInfo("HcalDigitizer") << "HCAL HO digis   : " << hoResult->size();
00539   edm::LogInfo("HcalDigitizer") << "HCAL HF digis   : " << hfResult->size();
00540   edm::LogInfo("HcalDigitizer") << "HCAL ZDC digis   : " << zdcResult->size();
00541   edm::LogInfo("HcalDigitizer") << "HCAL HBHE upgrade digis : " << hbheupgradeResult->size();
00542   edm::LogInfo("HcalDigitizer") << "HCAL HF upgrade digis : " << hfupgradeResult->size();
00543 
00544   /*
00545   std::cout << std::endl;
00546   std::cout << "HCAL HBHE digis : " << hbheResult->size() << std::endl;
00547   std::cout << "HCAL HO   digis : " << hoResult->size() << std::endl;
00548   std::cout << "HCAL HF   digis : " << hfResult->size() << std::endl;
00549  
00550   std::cout << "HCAL HBHE upgrade digis : " << hbheupgradeResult->size()
00551             << std::endl;
00552   std::cout << "HCAL HF   upgrade digis : " << hfupgradeResult->size()
00553             << std::endl;
00554   */
00555 
00556   // Step D: Put outputs into event
00557   e.put(hbheResult);
00558   e.put(hoResult);
00559   e.put(hfResult);
00560   e.put(zdcResult);
00561   e.put(hbheupgradeResult,"HBHEUpgradeDigiCollection");
00562   e.put(hfupgradeResult, "HFUpgradeDigiCollection");
00563 
00564   //  std::cout << std::endl << "========>  HcalDigitizer e.put " 
00565   //            << std::endl <<  std::endl;
00566 
00567   if(theHitCorrection) {
00568     theHitCorrection->clear();
00569   }
00570 }
00571 
00572 
00573 void HcalDigitizer::beginRun(const edm::EventSetup & es) {
00574   checkGeometry(es);
00575   theShapes->beginRun(es);
00576 }
00577 
00578 
00579 void HcalDigitizer::endRun() {
00580   theShapes->endRun();
00581 }
00582 
00583 
00584 void HcalDigitizer::checkGeometry(const edm::EventSetup & eventSetup) {
00585   // TODO find a way to avoid doing this every event
00586   edm::ESHandle<CaloGeometry> geometry;
00587   eventSetup.get<CaloGeometryRecord>().get(geometry);
00588   // See if it's been updated
00589   if(&*geometry != theGeometry)
00590   {
00591     theGeometry = &*geometry;
00592     updateGeometry(eventSetup);
00593   }
00594 }
00595 
00596 
00597 void  HcalDigitizer::updateGeometry(const edm::EventSetup & eventSetup) {
00598   if(theHBHEResponse) theHBHEResponse->setGeometry(theGeometry);
00599   if(theHBHESiPMResponse) theHBHESiPMResponse->setGeometry(theGeometry);
00600   if(theHOResponse) theHOResponse->setGeometry(theGeometry);
00601   if(theHOSiPMResponse) theHOSiPMResponse->setGeometry(theGeometry);
00602   theHFResponse->setGeometry(theGeometry);
00603   theZDCResponse->setGeometry(theGeometry);
00604   if(theRelabeller) theRelabeller->setGeometry(theGeometry);
00605 
00606   const std::vector<DetId>& hbCells = theGeometry->getValidDetIds(DetId::Hcal, HcalBarrel);
00607   const std::vector<DetId>& heCells = theGeometry->getValidDetIds(DetId::Hcal, HcalEndcap);
00608   const std::vector<DetId>& hoCells = theGeometry->getValidDetIds(DetId::Hcal, HcalOuter);
00609   const std::vector<DetId>& hfCells = theGeometry->getValidDetIds(DetId::Hcal, HcalForward);
00610   const std::vector<DetId>& zdcCells = theGeometry->getValidDetIds(DetId::Calo, HcalZDCDetId::SubdetectorId);
00611   //const std::vector<DetId>& hcalTrigCells = geometry->getValidDetIds(DetId::Hcal, HcalTriggerTower);
00612   //const std::vector<DetId>& hcalCalib = geometry->getValidDetIds(DetId::Calo, HcalCastorDetId::SubdetectorId);
00613   //std::cout<<"HcalDigitizer::CheckGeometry number of cells: "<<zdcCells.size()<<std::endl;
00614   if(zdcCells.empty()) zdcgeo = false;
00615   if(hbCells.empty() && heCells.empty()) hbhegeo = false;
00616   if(hoCells.empty()) hogeo = false;
00617   if(hfCells.empty()) hfgeo = false;
00618   // combine HB & HE
00619 
00620   theHBHEDetIds = hbCells;
00621   theHBHEDetIds.insert(theHBHEDetIds.end(), heCells.begin(), heCells.end());
00622 
00623   HcalDigitizerImpl::fillCells(theHBHEDetIds, theHBHEDigitizer, theHBHESiPMDigitizer);
00624   //HcalDigitizerImpl::fillCells(hoCells, theHODigitizer, theHOSiPMDigitizer);
00625   buildHOSiPMCells(hoCells, eventSetup);
00626   if(theHFDigitizer) theHFDigitizer->setDetIds(hfCells);
00627   if(theHFUpgradeDigitizer) theHFUpgradeDigitizer->setDetIds(hfCells);
00628   theZDCDigitizer->setDetIds(zdcCells); 
00629   if(theHBHEUpgradeDigitizer) {
00630     theHBHEUpgradeDigitizer->setDetIds(theHBHEDetIds);
00631 
00632     //    std::cout << " HcalDigitizer::updateGeometry  theHBHEUpgradeDigitizer->setDetIds(theHBHEDetIds)"<< std::endl;   
00633   }
00634 
00635 }
00636 
00637 
00638 void HcalDigitizer::buildHOSiPMCells(const std::vector<DetId>& allCells, const edm::EventSetup & eventSetup) {
00639   // all HPD
00640   if(theHOSiPMCode == 0) {
00641     theHODigitizer->setDetIds(allCells);
00642   } else if(theHOSiPMCode == 1) {
00643     theHOSiPMDigitizer->setDetIds(allCells);
00644     // FIXME pick Zecotek or hamamatsu?
00645   } else if(theHOSiPMCode == 2) {
00646     std::vector<HcalDetId> zecotekDetIds, hamamatsuDetIds;
00647     edm::ESHandle<HcalMCParams> p;
00648     eventSetup.get<HcalMCParamsRcd>().get(p);
00649     edm::ESHandle<HcalTopology> htopo;
00650     eventSetup.get<IdealGeometryRecord>().get(htopo);
00651    
00652     HcalMCParams mcParams(*p.product());
00653     if (mcParams.topo()==0) {
00654       mcParams.setTopo(htopo.product());
00655     }
00656 
00657     for(std::vector<DetId>::const_iterator detItr = allCells.begin();
00658         detItr != allCells.end(); ++detItr) {
00659       int shapeType = mcParams.getValues(*detItr)->signalShape();
00660       if(shapeType == HcalShapes::ZECOTEK) {
00661         zecotekDetIds.emplace_back(*detItr);
00662         theHOSiPMDetIds.push_back(*detItr);
00663       } else if(shapeType == HcalShapes::HAMAMATSU) {
00664         hamamatsuDetIds.emplace_back(*detItr);
00665         theHOSiPMDetIds.push_back(*detItr);
00666       } else {
00667         theHOHPDDetIds.push_back(*detItr);
00668       }
00669     }
00670 
00671     assert(theHODigitizer);
00672     assert(theHOSiPMDigitizer);
00673     theHODigitizer->setDetIds(theHOHPDDetIds);
00674     theHOSiPMDigitizer->setDetIds(theHOSiPMDetIds);
00675     theHOSiPMHitFilter.setDetIds(theHOSiPMDetIds);
00676     // FIXME not applying a HitFilter to the HPDs, for now
00677     theParameterMap->setHOZecotekDetIds(zecotekDetIds);
00678     theParameterMap->setHOHamamatsuDetIds(hamamatsuDetIds);
00679 
00680     // make sure we don't got through this exercise again
00681     theHOSiPMCode = -2;
00682   }
00683 }
00684 
00685       
00686     
00687