CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_6_1_2_SLHC2_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   theUpgradeElectronicsSim(0),
00105   theHBHEHitFilter(),
00106   theHFHitFilter(ps.getParameter<bool>("doHFWindow")),
00107   theHOHitFilter(),
00108   theHOSiPMHitFilter(HcalOuter),
00109   theZDCHitFilter(),
00110   theHitCorrection(0),
00111   theNoiseGenerator(0),
00112   theNoiseHitGenerator(0),
00113   theHBHEDigitizer(0),
00114   theHBHESiPMDigitizer(0),
00115   theHODigitizer(0),
00116   theHOSiPMDigitizer(0),
00117   theHFDigitizer(0),
00118   theZDCDigitizer(0),
00119   theHBHEUpgradeDigitizer(0),
00120   theHFUpgradeDigitizer(0),
00121   theRelabeller(0),
00122   isZDC(true),
00123   isHCAL(true),
00124   zdcgeo(true),
00125   hbhegeo(true),
00126   hogeo(true),
00127   hfgeo(true),
00128   theHOSiPMCode(ps.getParameter<edm::ParameterSet>("ho").getParameter<int>("siPMCode"))
00129 {
00130   bool doNoise = ps.getParameter<bool>("doNoise");
00131   bool useOldNoiseHB = ps.getParameter<bool>("useOldHB");
00132   bool useOldNoiseHE = ps.getParameter<bool>("useOldHE");
00133   bool useOldNoiseHF = ps.getParameter<bool>("useOldHF");
00134   bool useOldNoiseHO = ps.getParameter<bool>("useOldHO");
00135   bool doEmpty = ps.getParameter<bool>("doEmpty");
00136   double HBtp = ps.getParameter<double>("HBTuningParameter");
00137   double HEtp = ps.getParameter<double>("HETuningParameter");
00138   double HFtp = ps.getParameter<double>("HFTuningParameter");
00139   double HOtp = ps.getParameter<double>("HOTuningParameter");
00140   bool doHBHEUpgrade = ps.getParameter<bool>("HBHEUpgradeQIE");
00141   bool doHFUpgrade   = ps.getParameter<bool>("HFUpgradeQIE");
00142 
00143   // need to make copies, because they might get different noise generators
00144   theHBHEAmplifier = new HcalAmplifier(theParameterMap, doNoise);
00145   theHFAmplifier = new HcalAmplifier(theParameterMap, doNoise);
00146   theHOAmplifier = new HcalAmplifier(theParameterMap, doNoise);
00147   theZDCAmplifier = new HcalAmplifier(theParameterMap, doNoise);
00148   theHBHEAmplifier->setHBtuningParameter(HBtp);
00149   theHBHEAmplifier->setHEtuningParameter(HEtp);
00150   theHFAmplifier->setHFtuningParameter(HFtp);
00151   theHOAmplifier->setHOtuningParameter(HOtp);
00152   theHBHEAmplifier->setUseOldHB(useOldNoiseHB);
00153   theHBHEAmplifier->setUseOldHE(useOldNoiseHE);
00154   theHFAmplifier->setUseOldHF(useOldNoiseHF);
00155   theHOAmplifier->setUseOldHO(useOldNoiseHO);
00156 
00157   theCoderFactory = new HcalCoderFactory(HcalCoderFactory::DB);
00158   theUpgradeCoderFactory = new HcalCoderFactory(HcalCoderFactory::UPGRADE);
00159   theHBHEElectronicsSim = new HcalElectronicsSim(theHBHEAmplifier, theCoderFactory);
00160   theHFElectronicsSim = new HcalElectronicsSim(theHFAmplifier, theCoderFactory);
00161   theHOElectronicsSim = new HcalElectronicsSim(theHOAmplifier, theCoderFactory);
00162   theZDCElectronicsSim = new HcalElectronicsSim(theZDCAmplifier, theCoderFactory);
00163   theUpgradeElectronicsSim = new HcalElectronicsSim(theHBHEAmplifier, theUpgradeCoderFactory);
00164 
00165   // a code of 1 means make all cells SiPM
00166   std::vector<int> hbSiPMCells(ps.getParameter<edm::ParameterSet>("hb").getParameter<std::vector<int> >("siPMCells"));
00167   //std::vector<int> hoSiPMCells(ps.getParameter<edm::ParameterSet>("ho").getParameter<std::vector<int> >("siPMCells"));
00168   // 0 means none, 1 means all, and 2 means use hardcoded
00169 
00170   bool doHBHEHPD = hbSiPMCells.empty() || (hbSiPMCells[0] != 1);
00171   bool doHOHPD = (theHOSiPMCode != 1);
00172   bool doHBHESiPM = !hbSiPMCells.empty();
00173   bool doHOSiPM = (theHOSiPMCode != 0);
00174   if(doHBHEHPD) {
00175     theHBHEResponse = new CaloHitResponse(theParameterMap, theShapes);
00176     edm::LogInfo("HcalDigitizer") <<"Set scale for HB towers";
00177     theHBHEResponse->initHBHEScale();
00178 
00179     theHBHEResponse->setHitFilter(&theHBHEHitFilter);
00180     theHBHEDigitizer = new HBHEDigitizer(theHBHEResponse, theHBHEElectronicsSim, doEmpty);
00181     bool    changeResponse = ps.getParameter<bool>("ChangeResponse");
00182     edm::FileInPath fname  = ps.getParameter<edm::FileInPath>("CorrFactorFile");
00183     if (changeResponse) {
00184       std::string corrFileName = fname.fullPath();
00185       edm::LogInfo("HcalDigitizer") << "Set scale for HB towers from " << corrFileName;
00186       theHBHEResponse->setHBHEScale(corrFileName); //GMA
00187     }
00188   }
00189   if(doHOHPD) {
00190     theHOResponse = new CaloHitResponse(theParameterMap, theShapes);
00191     theHOResponse->setHitFilter(&theHOHitFilter);
00192     theHODigitizer = new HODigitizer(theHOResponse, theHOElectronicsSim, doEmpty);
00193   }
00194 
00195   if(doHBHESiPM) {
00196     theHBHESiPMResponse = new HcalSiPMHitResponse(theParameterMap, theShapes);
00197     theHBHESiPMResponse->setHitFilter(&theHBHEHitFilter);
00198     if (doHBHEUpgrade) {
00199       theHBHEUpgradeDigitizer = new UpgradeDigitizer(theHBHESiPMResponse, theUpgradeElectronicsSim, doEmpty);
00200     } else {
00201       theHBHESiPMDigitizer = new HBHEDigitizer(theHBHESiPMResponse, theHBHEElectronicsSim, doEmpty);
00202     }
00203   }
00204   if(doHOSiPM) {
00205     theHOSiPMResponse = new HcalSiPMHitResponse(theParameterMap, theShapes);
00206     theHOSiPMResponse->setHitFilter(&theHOSiPMHitFilter);
00207     theHOSiPMDigitizer = new HODigitizer(theHOSiPMResponse, theHOElectronicsSim, doEmpty);
00208   }
00209 
00210   // if both are present, fill the SiPM cells now
00211   if(doHBHEHPD && doHBHESiPM) {
00212     HcalDigitizerImpl::fillSiPMCells(hbSiPMCells, theHBHESiPMDigitizer);
00213   }
00214 
00215   theHFResponse->setHitFilter(&theHFHitFilter);
00216   theZDCResponse->setHitFilter(&theZDCHitFilter);
00217 
00218   bool doTimeSlew = ps.getParameter<bool>("doTimeSlew");
00219   if(doTimeSlew) {
00220     // no time slewing for HF
00221     theTimeSlewSim = new HcalTimeSlewSim(theParameterMap);
00222     theHBHEAmplifier->setTimeSlewSim(theTimeSlewSim);
00223     theHOAmplifier->setTimeSlewSim(theTimeSlewSim);
00224     theZDCAmplifier->setTimeSlewSim(theTimeSlewSim);
00225   }
00226 
00227   if (doHFUpgrade) {
00228     theHFUpgradeDigitizer = new UpgradeDigitizer(theHFResponse, theUpgradeElectronicsSim, doEmpty);
00229   } else {
00230     theHFDigitizer = new HFDigitizer(theHFResponse, theHFElectronicsSim, doEmpty);
00231   }
00232   theZDCDigitizer = new ZDCDigitizer(theZDCResponse, theZDCElectronicsSim, doEmpty);
00233 
00234   edm::ParameterSet ps0 = ps.getParameter<edm::ParameterSet>("HcalReLabel");
00235   relabel_ = ps0.getUntrackedParameter<bool>("RelabelHits");
00236   if (relabel_) {
00237     theRelabeller=new HcalHitRelabeller(ps0.getUntrackedParameter<edm::ParameterSet>("RelabelRules"));
00238   }     
00239 
00240   bool doHPDNoise = ps.getParameter<bool>("doHPDNoise");
00241   if(doHPDNoise) {
00242     //edm::ParameterSet hpdNoisePset = ps.getParameter<edm::ParameterSet>("HPDNoiseLibrary");
00243     theNoiseGenerator = new HPDNoiseGenerator(ps); 
00244     if(theHBHEDigitizer) theHBHEDigitizer->setNoiseSignalGenerator(theNoiseGenerator);
00245     if(theHBHESiPMDigitizer) theHBHESiPMDigitizer->setNoiseSignalGenerator(theNoiseGenerator);
00246   }
00247 
00248   if(ps.getParameter<bool>("doIonFeedback") && theHBHEResponse) {
00249     theIonFeedback = new HPDIonFeedbackSim(ps, theShapes);
00250     theHBHEResponse->setPECorrection(theIonFeedback);
00251     if(ps.getParameter<bool>("doThermalNoise")) {
00252       theHBHEAmplifier->setIonFeedbackSim(theIonFeedback);
00253     }
00254   }
00255 
00256   if(ps.getParameter<bool>("injectTestHits") ) {
00257     theNoiseHitGenerator = new HcalTestHitGenerator(ps);
00258     if(theHBHEDigitizer) theHBHEDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00259     if(theHBHESiPMDigitizer) theHBHESiPMDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00260     if(theHODigitizer) theHODigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00261     if(theHOSiPMDigitizer) theHOSiPMDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00262     if(theHBHEUpgradeDigitizer) theHBHEUpgradeDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00263     if(theHFDigitizer) theHFDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00264     if(theHFUpgradeDigitizer) theHFUpgradeDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00265     theZDCDigitizer->setNoiseHitGenerator(theNoiseHitGenerator);
00266   }
00267 
00268   edm::Service<edm::RandomNumberGenerator> rng;
00269   if ( ! rng.isAvailable()) {
00270     throw cms::Exception("Configuration")
00271       << "HcalDigitizer requires the RandomNumberGeneratorService\n"
00272          "which is not present in the configuration file.  You must add the service\n"
00273          "in the configuration file or remove the modules that require it.";
00274   }
00275 
00276   CLHEP::HepRandomEngine& engine = rng->getEngine();
00277   if(theHBHEDigitizer) theHBHEDigitizer->setRandomEngine(engine);
00278   if(theHBHESiPMDigitizer) theHBHESiPMDigitizer->setRandomEngine(engine);
00279   if(theHODigitizer) theHODigitizer->setRandomEngine(engine);
00280   if(theHOSiPMDigitizer) theHOSiPMDigitizer->setRandomEngine(engine);
00281   if(theHBHEUpgradeDigitizer) theHBHEUpgradeDigitizer->setRandomEngine(engine);
00282   if(theIonFeedback) theIonFeedback->setRandomEngine(engine);
00283   if(theTimeSlewSim) theTimeSlewSim->setRandomEngine(engine);
00284   if(theHFUpgradeDigitizer) theHFUpgradeDigitizer->setRandomEngine(engine);
00285   if(theHFDigitizer) theHFDigitizer->setRandomEngine(engine);
00286   theZDCDigitizer->setRandomEngine(engine);
00287 
00288   if (theHitCorrection!=0) theHitCorrection->setRandomEngine(engine);
00289 
00290   hitsProducer_ = ps.getParameter<std::string>("hitsProducer");
00291 }
00292 
00293 
00294 HcalDigitizer::~HcalDigitizer() {
00295   if(theHBHEDigitizer)         delete theHBHEDigitizer;
00296   if(theHBHESiPMDigitizer)     delete theHBHESiPMDigitizer;
00297   if(theHODigitizer)           delete theHODigitizer;
00298   if(theHOSiPMDigitizer)       delete theHOSiPMDigitizer;
00299   if(theHFDigitizer)           delete theHFDigitizer;
00300   delete theZDCDigitizer;
00301   if(theHBHEUpgradeDigitizer)  delete theHBHEUpgradeDigitizer;
00302   if(theHFUpgradeDigitizer)    delete theHFUpgradeDigitizer;
00303   delete theParameterMap;
00304   delete theHBHEResponse;
00305   delete theHBHESiPMResponse;
00306   delete theHOResponse;
00307   delete theHOSiPMResponse;
00308   delete theHFResponse;
00309   delete theZDCResponse;
00310   delete theHBHEElectronicsSim;
00311   delete theHFElectronicsSim;
00312   delete theHOElectronicsSim;
00313   delete theZDCElectronicsSim;
00314   delete theUpgradeElectronicsSim;
00315   delete theHBHEAmplifier;
00316   delete theHFAmplifier;
00317   delete theHOAmplifier;
00318   delete theZDCAmplifier;
00319   delete theCoderFactory;
00320   delete theUpgradeCoderFactory;
00321   delete theHitCorrection;
00322   delete theNoiseGenerator;
00323   if (theRelabeller)           delete theRelabeller;
00324 }
00325 
00326 
00327 void HcalDigitizer::setHBHENoiseSignalGenerator(HcalBaseSignalGenerator * noiseGenerator) {
00328   noiseGenerator->setParameterMap(theParameterMap);
00329   noiseGenerator->setElectronicsSim(theHBHEElectronicsSim);
00330   if (theHBHEDigitizer) theHBHEDigitizer->setNoiseSignalGenerator(noiseGenerator);
00331   theHBHEAmplifier->setNoiseSignalGenerator(noiseGenerator);
00332 }
00333 
00334 void HcalDigitizer::setHFNoiseSignalGenerator(HcalBaseSignalGenerator * noiseGenerator) {
00335   noiseGenerator->setParameterMap(theParameterMap);
00336   noiseGenerator->setElectronicsSim(theHFElectronicsSim);
00337   if(theHFDigitizer) theHFDigitizer->setNoiseSignalGenerator(noiseGenerator);
00338   if(theHFUpgradeDigitizer) theHFUpgradeDigitizer->setNoiseSignalGenerator(noiseGenerator);
00339   theHFAmplifier->setNoiseSignalGenerator(noiseGenerator);
00340 }
00341 
00342 void HcalDigitizer::setHONoiseSignalGenerator(HcalBaseSignalGenerator * noiseGenerator) {
00343   noiseGenerator->setParameterMap(theParameterMap);
00344   noiseGenerator->setElectronicsSim(theHOElectronicsSim);
00345   theHODigitizer->setNoiseSignalGenerator(noiseGenerator);
00346   theHOAmplifier->setNoiseSignalGenerator(noiseGenerator);
00347 }
00348 
00349 void HcalDigitizer::setZDCNoiseSignalGenerator(HcalBaseSignalGenerator * noiseGenerator) {
00350   noiseGenerator->setParameterMap(theParameterMap);
00351   noiseGenerator->setElectronicsSim(theZDCElectronicsSim);
00352   theZDCDigitizer->setNoiseSignalGenerator(noiseGenerator);
00353   theZDCAmplifier->setNoiseSignalGenerator(noiseGenerator);
00354 }
00355 
00356 void HcalDigitizer::initializeEvent(edm::Event const& e, edm::EventSetup const& eventSetup) {
00357   // get the appropriate gains, noises, & widths for this event
00358   edm::ESHandle<HcalDbService> conditions;
00359   eventSetup.get<HcalDbRecord>().get(conditions);
00360   theHBHEAmplifier->setDbService(conditions.product());
00361   theHFAmplifier->setDbService(conditions.product());
00362   theHOAmplifier->setDbService(conditions.product());
00363   theZDCAmplifier->setDbService(conditions.product());
00364   theUpgradeElectronicsSim->setDbService(conditions.product());
00365 
00366   theCoderFactory->setDbService(conditions.product());
00367   theUpgradeCoderFactory->setDbService(conditions.product());
00368   theParameterMap->setDbService(conditions.product());
00369 
00370   edm::ESHandle<HcalCholeskyMatrices> refCholesky;
00371   eventSetup.get<HcalCholeskyMatricesRcd>().get(refCholesky);
00372   const HcalCholeskyMatrices * myCholesky = refCholesky.product();
00373 
00374   edm::ESHandle<HcalPedestals> pedshandle;
00375   eventSetup.get<HcalPedestalsRcd>().get(pedshandle);
00376   const HcalPedestals *  myADCPedestals = pedshandle.product();
00377 
00378   theHBHEAmplifier->setCholesky(myCholesky);
00379   theHFAmplifier->setCholesky(myCholesky);
00380   theHOAmplifier->setCholesky(myCholesky);
00381   
00382   theHBHEAmplifier->setADCPeds(myADCPedestals);
00383   theHFAmplifier->setADCPeds(myADCPedestals);
00384   theHOAmplifier->setADCPeds(myADCPedestals);
00385 
00386   if(theHitCorrection != 0) {
00387     theHitCorrection->clear();
00388   }
00389 }
00390 
00391 void HcalDigitizer::accumulateCaloHits(edm::Handle<std::vector<PCaloHit> > const& hcalHandle, edm::Handle<std::vector<PCaloHit> > const& zdcHandle, int bunchCrossing) {
00392   // Step A: pass in inputs, and accumulate digirs
00393   if(isHCAL) {
00394     std::vector<PCaloHit> hcalHits = *hcalHandle.product();
00395     if (relabel_) {
00396       // Relabel PCaloHits
00397       edm::LogInfo("HcalDigitizer") << "Calling Relabller";
00398       theRelabeller->process(hcalHits);
00399     }
00400     if(theHitCorrection != 0) {
00401       theHitCorrection->fillChargeSums(hcalHits);
00402     }
00403     if(hbhegeo) {
00404       if(theHBHEDigitizer) theHBHEDigitizer->add(hcalHits, bunchCrossing);
00405       if(theHBHESiPMDigitizer) theHBHESiPMDigitizer->add(hcalHits, bunchCrossing);
00406       if(theHBHEUpgradeDigitizer) theHBHEUpgradeDigitizer->add(hcalHits, bunchCrossing);
00407     }
00408 
00409     if(hogeo) {
00410       if(theHODigitizer) theHODigitizer->add(hcalHits, bunchCrossing);
00411       if(theHOSiPMDigitizer) theHOSiPMDigitizer->add(hcalHits, bunchCrossing);
00412     }
00413 
00414     if(hfgeo) {
00415       if(theHFDigitizer) theHFDigitizer->add(hcalHits, bunchCrossing);
00416       if(theHFUpgradeDigitizer) theHFUpgradeDigitizer->add(hcalHits, bunchCrossing);
00417     } 
00418   } else {
00419     edm::LogInfo("HcalDigitizer") << "We don't have HCAL hit collection available ";
00420   }
00421 
00422   if(isZDC) {
00423     if(zdcgeo) {
00424       theZDCDigitizer->add(*zdcHandle.product(), bunchCrossing);
00425     } 
00426   } else {
00427     edm::LogInfo("HcalDigitizer") << "We don't have ZDC hit collection available ";
00428   }
00429 }
00430 
00431 void HcalDigitizer::accumulate(edm::Event const& e, edm::EventSetup const& eventSetup) {
00432   // Step A: Get Inputs
00433   edm::InputTag zdcTag(hitsProducer_, "ZDCHITS");
00434   edm::Handle<std::vector<PCaloHit> > zdcHandle;
00435   e.getByLabel(zdcTag, zdcHandle);
00436   isZDC = zdcHandle.isValid();
00437 
00438   edm::InputTag hcalTag(hitsProducer_, "HcalHits");
00439   edm::Handle<std::vector<PCaloHit> > hcalHandle;
00440   e.getByLabel(hcalTag, hcalHandle);
00441   isHCAL = hcalHandle.isValid();
00442 
00443   accumulateCaloHits(hcalHandle, zdcHandle, 0);
00444 }
00445 
00446 void HcalDigitizer::accumulate(PileUpEventPrincipal const& e, edm::EventSetup const& eventSetup) {
00447   // Step A: Get Inputs
00448   edm::InputTag zdcTag(hitsProducer_, "ZDCHITS");
00449   edm::Handle<std::vector<PCaloHit> > zdcHandle;
00450   e.getByLabel(zdcTag, zdcHandle);
00451   isZDC = zdcHandle.isValid();
00452 
00453   edm::InputTag hcalTag(hitsProducer_, "HcalHits");
00454   edm::Handle<std::vector<PCaloHit> > hcalHandle;
00455   e.getByLabel(hcalTag, hcalHandle);
00456   isHCAL = hcalHandle.isValid();
00457 
00458   accumulateCaloHits(hcalHandle, zdcHandle, e.bunchCrossing());
00459 }
00460 
00461 void HcalDigitizer::finalizeEvent(edm::Event& e, const edm::EventSetup& eventSetup) {
00462 
00463   // Step B: Create empty output
00464   std::auto_ptr<HBHEDigiCollection> hbheResult(new HBHEDigiCollection());
00465   std::auto_ptr<HODigiCollection> hoResult(new HODigiCollection());
00466   std::auto_ptr<HFDigiCollection> hfResult(new HFDigiCollection());
00467   std::auto_ptr<ZDCDigiCollection> zdcResult(new ZDCDigiCollection());
00468   std::auto_ptr<HBHEUpgradeDigiCollection> hbheupgradeResult(new HBHEUpgradeDigiCollection());
00469   std::auto_ptr<HFUpgradeDigiCollection> hfupgradeResult(new HFUpgradeDigiCollection());
00470 
00471   // Step C: Invoke the algorithm, getting back outputs.
00472   if(isHCAL&&hbhegeo){
00473     if(theHBHEDigitizer)        theHBHEDigitizer->run(*hbheResult);
00474     if(theHBHESiPMDigitizer)    theHBHESiPMDigitizer->run(*hbheResult);
00475     if(theHBHEUpgradeDigitizer) theHBHEUpgradeDigitizer->run(*hbheupgradeResult);
00476   }
00477   if(isHCAL&&hogeo) {
00478     if(theHODigitizer) theHODigitizer->run(*hoResult);
00479     if(theHOSiPMDigitizer) theHOSiPMDigitizer->run(*hoResult);
00480   }
00481   if(isHCAL&&hfgeo) {
00482     if(theHFDigitizer) theHFDigitizer->run(*hfResult);
00483     if(theHFUpgradeDigitizer) theHFUpgradeDigitizer->run(*hfupgradeResult);
00484   }
00485   if(isZDC&&zdcgeo) {
00486     theZDCDigitizer->run(*zdcResult);
00487   }
00488   
00489   edm::LogInfo("HcalDigitizer") << "HCAL HBHE digis : " << hbheResult->size();
00490   edm::LogInfo("HcalDigitizer") << "HCAL HO digis   : " << hoResult->size();
00491   edm::LogInfo("HcalDigitizer") << "HCAL HF digis   : " << hfResult->size();
00492   edm::LogInfo("HcalDigitizer") << "HCAL ZDC digis   : " << zdcResult->size();
00493   edm::LogInfo("HcalDigitizer") << "HCAL HBHE upgrade digis : " << hbheupgradeResult->size();
00494   edm::LogInfo("HcalDigitizer") << "HCAL HF upgrade digis : " << hfupgradeResult->size();
00495 
00496   // Step D: Put outputs into event
00497   e.put(hbheResult);
00498   e.put(hoResult);
00499   e.put(hfResult);
00500   e.put(zdcResult);
00501   e.put(hbheupgradeResult,"HBHEUpgradeDigiCollection");
00502   e.put(hfupgradeResult, "HFUpgradeDigiCollection");
00503 
00504   if(theHitCorrection) {
00505     theHitCorrection->clear();
00506   }
00507 }
00508 
00509 
00510 void HcalDigitizer::beginRun(const edm::EventSetup & es) {
00511   checkGeometry(es);
00512   theShapes->beginRun(es);
00513 }
00514 
00515 
00516 void HcalDigitizer::endRun() {
00517   theShapes->endRun();
00518 }
00519 
00520 
00521 void HcalDigitizer::checkGeometry(const edm::EventSetup & eventSetup) {
00522   // TODO find a way to avoid doing this every event
00523   edm::ESHandle<CaloGeometry> geometry;
00524   eventSetup.get<CaloGeometryRecord>().get(geometry);
00525   // See if it's been updated
00526   if(&*geometry != theGeometry)
00527   {
00528     theGeometry = &*geometry;
00529     updateGeometry(eventSetup);
00530   }
00531 }
00532 
00533 
00534 void  HcalDigitizer::updateGeometry(const edm::EventSetup & eventSetup) {
00535   if(theHBHEResponse) theHBHEResponse->setGeometry(theGeometry);
00536   if(theHBHESiPMResponse) theHBHESiPMResponse->setGeometry(theGeometry);
00537   if(theHOResponse) theHOResponse->setGeometry(theGeometry);
00538   if(theHOSiPMResponse) theHOSiPMResponse->setGeometry(theGeometry);
00539   theHFResponse->setGeometry(theGeometry);
00540   theZDCResponse->setGeometry(theGeometry);
00541   if(theRelabeller) theRelabeller->setGeometry(theGeometry);
00542 
00543   const std::vector<DetId>& hbCells = theGeometry->getValidDetIds(DetId::Hcal, HcalBarrel);
00544   const std::vector<DetId>& heCells = theGeometry->getValidDetIds(DetId::Hcal, HcalEndcap);
00545   const std::vector<DetId>& hoCells = theGeometry->getValidDetIds(DetId::Hcal, HcalOuter);
00546   const std::vector<DetId>& hfCells = theGeometry->getValidDetIds(DetId::Hcal, HcalForward);
00547   const std::vector<DetId>& zdcCells = theGeometry->getValidDetIds(DetId::Calo, HcalZDCDetId::SubdetectorId);
00548   //const std::vector<DetId>& hcalTrigCells = geometry->getValidDetIds(DetId::Hcal, HcalTriggerTower);
00549   //const std::vector<DetId>& hcalCalib = geometry->getValidDetIds(DetId::Calo, HcalCastorDetId::SubdetectorId);
00550   //std::cout<<"HcalDigitizer::CheckGeometry number of cells: "<<zdcCells.size()<<std::endl;
00551   if(zdcCells.empty()) zdcgeo = false;
00552   if(hbCells.empty() && heCells.empty()) hbhegeo = false;
00553   if(hoCells.empty()) hogeo = false;
00554   if(hfCells.empty()) hfgeo = false;
00555   // combine HB & HE
00556 
00557   theHBHEDetIds = hbCells;
00558   theHBHEDetIds.insert(theHBHEDetIds.end(), heCells.begin(), heCells.end());
00559 
00560   HcalDigitizerImpl::fillCells(theHBHEDetIds, theHBHEDigitizer, theHBHESiPMDigitizer);
00561   //HcalDigitizerImpl::fillCells(hoCells, theHODigitizer, theHOSiPMDigitizer);
00562   buildHOSiPMCells(hoCells, eventSetup);
00563   if(theHFDigitizer) theHFDigitizer->setDetIds(hfCells);
00564   if(theHFUpgradeDigitizer) theHFUpgradeDigitizer->setDetIds(hfCells);
00565   theZDCDigitizer->setDetIds(zdcCells); 
00566   if(theHBHEUpgradeDigitizer) theHBHEUpgradeDigitizer->setDetIds(theHBHEDetIds);
00567 }
00568 
00569 
00570 void HcalDigitizer::buildHOSiPMCells(const std::vector<DetId>& allCells, const edm::EventSetup & eventSetup) {
00571   // all HPD
00572   if(theHOSiPMCode == 0) {
00573     theHODigitizer->setDetIds(allCells);
00574   } else if(theHOSiPMCode == 1) {
00575     theHOSiPMDigitizer->setDetIds(allCells);
00576     // FIXME pick Zecotek or hamamatsu?
00577   } else if(theHOSiPMCode == 2) {
00578     std::vector<HcalDetId> zecotekDetIds, hamamatsuDetIds;
00579     edm::ESHandle<HcalMCParams> p;
00580     eventSetup.get<HcalMCParamsRcd>().get(p);
00581     edm::ESHandle<HcalTopology> htopo;
00582     eventSetup.get<IdealGeometryRecord>().get(htopo);
00583    
00584     HcalMCParams mcParams(*p.product());
00585     if (mcParams.topo()==0) {
00586       mcParams.setTopo(htopo.product());
00587     }
00588 
00589     for(std::vector<DetId>::const_iterator detItr = allCells.begin();
00590         detItr != allCells.end(); ++detItr) {
00591       int shapeType = mcParams.getValues(*detItr)->signalShape();
00592       if(shapeType == HcalShapes::ZECOTEK) {
00593         zecotekDetIds.emplace_back(*detItr);
00594         theHOSiPMDetIds.push_back(*detItr);
00595       } else if(shapeType == HcalShapes::HAMAMATSU) {
00596         hamamatsuDetIds.emplace_back(*detItr);
00597         theHOSiPMDetIds.push_back(*detItr);
00598       } else {
00599         theHOHPDDetIds.push_back(*detItr);
00600       }
00601     }
00602 
00603     assert(theHODigitizer);
00604     assert(theHOSiPMDigitizer);
00605     theHODigitizer->setDetIds(theHOHPDDetIds);
00606     theHOSiPMDigitizer->setDetIds(theHOSiPMDetIds);
00607     theHOSiPMHitFilter.setDetIds(theHOSiPMDetIds);
00608     // FIXME not applying a HitFilter to the HPDs, for now
00609     theParameterMap->setHOZecotekDetIds(zecotekDetIds);
00610     theParameterMap->setHOHamamatsuDetIds(hamamatsuDetIds);
00611 
00612     // make sure we don't got through this exercise again
00613     theHOSiPMCode = -2;
00614   }
00615 }
00616 
00617       
00618     
00619