CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_6_1_2_SLHC2_patch1/src/SimTracker/SiStripDigitizer/plugins/SiStripDigitizer.cc

Go to the documentation of this file.
00001 // File: SiStripDigitizerAlgorithm.cc
00002 // Description:  Class for digitization.
00003 
00004 // system include files
00005 #include <memory>
00006 
00007 #include "SimTracker/Common/interface/SimHitSelectorFromDB.h"
00008 
00009 #include "SiStripDigitizer.h"
00010 #include "SiStripDigitizerAlgorithm.h"
00011 #include "SimGeneral/MixingModule/interface/PileUpEventPrincipal.h"
00012 
00013 #include "DataFormats/Common/interface/DetSetVector.h"
00014 #include "DataFormats/SiStripDigi/interface/SiStripDigi.h"
00015 #include "DataFormats/SiStripDigi/interface/SiStripRawDigi.h"
00016 
00017 // user include files
00018 #include "FWCore/Framework/interface/Frameworkfwd.h"
00019 #include "FWCore/Framework/interface/EDProducer.h"
00020 
00021 #include "FWCore/Framework/interface/Event.h"
00022 #include "FWCore/Framework/interface/MakerMacros.h"
00023 
00024 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00025 #include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
00026 #include "SimDataFormats/TrackingHit/interface/PSimHit.h"
00027 
00028 //needed for the geometry:
00029 #include "FWCore/Framework/interface/EventSetup.h"
00030 #include "FWCore/Framework/interface/ESHandle.h"
00031 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
00032 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
00033 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
00034 #include "Geometry/TrackerNumberingBuilder/interface/GeometricDet.h"
00035 #include "Geometry/CommonTopologies/interface/StripTopology.h"
00036 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetType.h"
00037 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
00038 //needed for the magnetic field:
00039 #include "MagneticField/Engine/interface/MagneticField.h"
00040 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
00041 
00042 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
00043 
00044 //Data Base infromations
00045 #include "CalibFormats/SiStripObjects/interface/SiStripDetCabling.h"
00046 #include "CalibTracker/Records/interface/SiStripDependentRecords.h"
00047 #include "CondFormats/DataRecord/interface/SiStripCondDataRecords.h"
00048 
00049 //Random Number
00050 #include "FWCore/ServiceRegistry/interface/Service.h"
00051 #include "FWCore/Utilities/interface/RandomNumberGenerator.h"
00052 #include "FWCore/Utilities/interface/Exception.h"
00053 #include "CLHEP/Random/RandomEngine.h"
00054 
00055 SiStripDigitizer::SiStripDigitizer(const edm::ParameterSet& conf, edm::EDProducer& mixMod) : 
00056   gainLabel(conf.getParameter<std::string>("Gain")),
00057   hitsProducer(conf.getParameter<std::string>("hitsProducer")),
00058   trackerContainers(conf.getParameter<std::vector<std::string> >("ROUList")),
00059   ZSDigi(conf.getParameter<edm::ParameterSet>("DigiModeList").getParameter<std::string>("ZSDigi")),
00060   SCDigi(conf.getParameter<edm::ParameterSet>("DigiModeList").getParameter<std::string>("SCDigi")),
00061   VRDigi(conf.getParameter<edm::ParameterSet>("DigiModeList").getParameter<std::string>("VRDigi")),
00062   PRDigi(conf.getParameter<edm::ParameterSet>("DigiModeList").getParameter<std::string>("PRDigi")),
00063   geometryType(conf.getParameter<std::string>("GeometryType")),
00064   useConfFromDB(conf.getParameter<bool>("TrackerConfigurationFromDB")),
00065   zeroSuppression(conf.getParameter<bool>("ZeroSuppression"))
00066 { 
00067   const std::string alias("simSiStripDigis");
00068   
00069   mixMod.produces<edm::DetSetVector<SiStripDigi> >(ZSDigi).setBranchAlias(ZSDigi);
00070   mixMod.produces<edm::DetSetVector<SiStripRawDigi> >(SCDigi).setBranchAlias(alias + SCDigi);
00071   mixMod.produces<edm::DetSetVector<SiStripRawDigi> >(VRDigi).setBranchAlias(alias + VRDigi);
00072   mixMod.produces<edm::DetSetVector<SiStripRawDigi> >(PRDigi).setBranchAlias(alias + PRDigi);
00073   edm::Service<edm::RandomNumberGenerator> rng;
00074   if ( ! rng.isAvailable()) {
00075     throw cms::Exception("Configuration")
00076       << "SiStripDigitizer requires the RandomNumberGeneratorService\n"
00077       "which is not present in the configuration file.  You must add the service\n"
00078       "in the configuration file or remove the modules that require it.";
00079   }
00080   
00081   rndEngine = &(rng->getEngine());
00082   theDigiAlgo.reset(new SiStripDigitizerAlgorithm(conf,(*rndEngine)));
00083 
00084 }
00085 
00086 // Virtual destructor needed.
00087 SiStripDigitizer::~SiStripDigitizer() { 
00088 }  
00089 
00090 void SiStripDigitizer::accumulateStripHits(edm::Handle<std::vector<PSimHit> > hSimHits) {
00091   if(hSimHits.isValid()) {
00092     std::set<unsigned int> detIds;
00093     std::vector<PSimHit> const& simHits = *hSimHits.product();
00094     for(std::vector<PSimHit>::const_iterator it = simHits.begin(), itEnd = simHits.end(); it != itEnd; ++it) {
00095       unsigned int detId = (*it).detUnitId();
00096       if(detIds.insert(detId).second) {
00097         // The insert succeeded, so this detector element has not yet been processed.
00098         StripGeomDetUnit* stripdet = detectorUnits[detId];
00099         //access to magnetic field in global coordinates
00100         GlobalVector bfield = pSetup->inTesla(stripdet->surface().position());
00101         LogDebug ("Digitizer ") << "B-field(T) at " << stripdet->surface().position() << "(cm): "
00102                                 << pSetup->inTesla(stripdet->surface().position());
00103         theDigiAlgo->accumulateSimHits(it, itEnd, stripdet, bfield);
00104       }
00105     }
00106   }
00107 }
00108 
00109 // Functions that gets called by framework every event
00110   void
00111   SiStripDigitizer::accumulate(edm::Event const& iEvent, edm::EventSetup const& iSetup) {
00112     // Step A: Get Inputs
00113     for(vstring::const_iterator i = trackerContainers.begin(), iEnd = trackerContainers.end(); i != iEnd; ++i) {
00114       edm::Handle<std::vector<PSimHit> > simHits;
00115       edm::InputTag tag(hitsProducer, *i);
00116 
00117       iEvent.getByLabel(tag, simHits);
00118       accumulateStripHits(simHits);
00119     }
00120   }
00121 
00122   void
00123   SiStripDigitizer::accumulate(PileUpEventPrincipal const& iEvent, edm::EventSetup const& iSetup) {
00124     // Step A: Get Inputs
00125     for(vstring::const_iterator i = trackerContainers.begin(), iEnd = trackerContainers.end(); i != iEnd; ++i) {
00126       edm::Handle<std::vector<PSimHit> > simHits;
00127       edm::InputTag tag(hitsProducer, *i);
00128 
00129       iEvent.getByLabel(tag, simHits);
00130       accumulateStripHits(simHits);
00131     }
00132   }
00133 
00134 
00135 void SiStripDigitizer::initializeEvent(edm::Event const& iEvent, edm::EventSetup const& iSetup) {
00136   // Step A: Get Inputs
00137 
00138   if(useConfFromDB){
00139     edm::ESHandle<SiStripDetCabling> detCabling;
00140     iSetup.get<SiStripDetCablingRcd>().get(detCabling);
00141     detCabling->addConnected(theDetIdList);
00142   }
00143 
00144   theDigiAlgo->initializeEvent(iSetup);
00145 
00146   iSetup.get<TrackerDigiGeometryRecord>().get(geometryType,pDD);
00147   iSetup.get<IdealMagneticFieldRecord>().get(pSetup);
00148 
00149   // FIX THIS! We only need to clear and (re)fill detectorUnits when the geometry type IOV changes.  Use ESWatcher to determine this.
00150   bool changes = true;
00151   if(changes) { // Replace with ESWatcher
00152     detectorUnits.clear();
00153   }
00154   for(TrackingGeometry::DetUnitContainer::const_iterator iu = pDD->detUnits().begin(); iu != pDD->detUnits().end(); ++iu) {
00155     unsigned int detId = (*iu)->geographicalId().rawId();
00156     DetId idet=DetId(detId);
00157     unsigned int isub=idet.subdetId();
00158     if((isub == StripSubdetector::TIB) ||
00159        (isub == StripSubdetector::TID) ||
00160        (isub == StripSubdetector::TOB) ||
00161        (isub == StripSubdetector::TEC)) {
00162       StripGeomDetUnit* stripdet = dynamic_cast<StripGeomDetUnit*>((*iu));
00163       assert(stripdet != 0);
00164       if(changes) { // Replace with ESWatcher
00165         detectorUnits.insert(std::make_pair(detId, stripdet));
00166       }
00167       theDigiAlgo->initializeDetUnit(stripdet, iSetup);
00168     }
00169   }
00170 }
00171 
00172 void SiStripDigitizer::finalizeEvent(edm::Event& iEvent, edm::EventSetup const& iSetup) {
00173   edm::ESHandle<SiStripGain> gainHandle;
00174   edm::ESHandle<SiStripNoises> noiseHandle;
00175   edm::ESHandle<SiStripThreshold> thresholdHandle;
00176   edm::ESHandle<SiStripPedestals> pedestalHandle;
00177   iSetup.get<SiStripGainSimRcd>().get(gainLabel,gainHandle);
00178   iSetup.get<SiStripNoisesRcd>().get(noiseHandle);
00179   iSetup.get<SiStripThresholdRcd>().get(thresholdHandle);
00180   iSetup.get<SiStripPedestalsRcd>().get(pedestalHandle);
00181 
00182   std::vector<edm::DetSet<SiStripDigi> > theDigiVector;
00183   std::vector<edm::DetSet<SiStripRawDigi> > theRawDigiVector;
00184 
00185   
00186   // Step B: LOOP on StripGeomDetUnit
00187   theDigiVector.reserve(10000);
00188   theDigiVector.clear();
00189 
00190   for(TrackingGeometry::DetUnitContainer::const_iterator iu = pDD->detUnits().begin(); iu != pDD->detUnits().end(); iu ++){
00191     if(useConfFromDB){
00192       //apply the cable map _before_ digitization: consider only the detis that are connected 
00193       if(theDetIdList.find((*iu)->geographicalId().rawId())==theDetIdList.end())
00194         continue;
00195     }
00196     StripGeomDetUnit* sgd = dynamic_cast<StripGeomDetUnit*>((*iu));
00197     if (sgd != 0){
00198       edm::DetSet<SiStripDigi> collectorZS((*iu)->geographicalId().rawId());
00199       edm::DetSet<SiStripRawDigi> collectorRaw((*iu)->geographicalId().rawId());
00200       theDigiAlgo->digitize(collectorZS,collectorRaw,sgd,
00201                        gainHandle,thresholdHandle,noiseHandle,pedestalHandle);
00202       if(zeroSuppression){
00203         if(collectorZS.data.size()>0){
00204           theDigiVector.push_back(collectorZS);
00205         }
00206       }else{
00207         if(collectorRaw.data.size()>0){
00208           theRawDigiVector.push_back(collectorRaw);
00209         }
00210       }
00211     }
00212   }
00213 
00214   if(zeroSuppression){
00215     // Step C: create output collection
00216     std::auto_ptr<edm::DetSetVector<SiStripRawDigi> > output_virginraw(new edm::DetSetVector<SiStripRawDigi>());
00217     std::auto_ptr<edm::DetSetVector<SiStripRawDigi> > output_scopemode(new edm::DetSetVector<SiStripRawDigi>());
00218     std::auto_ptr<edm::DetSetVector<SiStripRawDigi> > output_processedraw(new edm::DetSetVector<SiStripRawDigi>());
00219     std::auto_ptr<edm::DetSetVector<SiStripDigi> > output(new edm::DetSetVector<SiStripDigi>(theDigiVector) );
00220     // Step D: write output to file
00221     iEvent.put(output, ZSDigi);
00222     iEvent.put(output_scopemode, SCDigi);
00223     iEvent.put(output_virginraw, VRDigi);
00224     iEvent.put(output_processedraw, PRDigi);
00225   }else{
00226     // Step C: create output collection
00227     std::auto_ptr<edm::DetSetVector<SiStripRawDigi> > output_virginraw(new edm::DetSetVector<SiStripRawDigi>(theRawDigiVector));
00228     std::auto_ptr<edm::DetSetVector<SiStripRawDigi> > output_scopemode(new edm::DetSetVector<SiStripRawDigi>());
00229     std::auto_ptr<edm::DetSetVector<SiStripRawDigi> > output_processedraw(new edm::DetSetVector<SiStripRawDigi>());
00230     std::auto_ptr<edm::DetSetVector<SiStripDigi> > output(new edm::DetSetVector<SiStripDigi>() );
00231     // Step D: write output to file
00232     iEvent.put(output, ZSDigi);
00233     iEvent.put(output_scopemode, SCDigi);
00234     iEvent.put(output_virginraw, VRDigi);
00235     iEvent.put(output_processedraw, PRDigi);
00236   }
00237 }