CMS 3D CMS Logo

/afs/cern.ch/work/c/cnuttens/private/SiStripDev/DocumentationProduction/CMSSW_6_1_0/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                                            const TrackerTopology *tTopo) {
00092   if(hSimHits.isValid()) {
00093     std::set<unsigned int> detIds;
00094     std::vector<PSimHit> const& simHits = *hSimHits.product();
00095     for(std::vector<PSimHit>::const_iterator it = simHits.begin(), itEnd = simHits.end(); it != itEnd; ++it) {
00096       unsigned int detId = (*it).detUnitId();
00097       if(detIds.insert(detId).second) {
00098         // The insert succeeded, so this detector element has not yet been processed.
00099         StripGeomDetUnit* stripdet = detectorUnits[detId];
00100         //access to magnetic field in global coordinates
00101         GlobalVector bfield = pSetup->inTesla(stripdet->surface().position());
00102         LogDebug ("Digitizer ") << "B-field(T) at " << stripdet->surface().position() << "(cm): "
00103                                 << pSetup->inTesla(stripdet->surface().position());
00104         theDigiAlgo->accumulateSimHits(it, itEnd, stripdet, bfield, tTopo);
00105       }
00106     }
00107   }
00108 }
00109 
00110 // Functions that gets called by framework every event
00111   void
00112   SiStripDigitizer::accumulate(edm::Event const& iEvent, edm::EventSetup const& iSetup) {
00113     //Retrieve tracker topology from geometry
00114     edm::ESHandle<TrackerTopology> tTopoHand;
00115     iSetup.get<IdealGeometryRecord>().get(tTopoHand);
00116     const TrackerTopology *tTopo=tTopoHand.product();
00117 
00118     // Step A: Get Inputs
00119     for(vstring::const_iterator i = trackerContainers.begin(), iEnd = trackerContainers.end(); i != iEnd; ++i) {
00120       edm::Handle<std::vector<PSimHit> > simHits;
00121       edm::InputTag tag(hitsProducer, *i);
00122 
00123       iEvent.getByLabel(tag, simHits);
00124       accumulateStripHits(simHits,tTopo);
00125     }
00126   }
00127 
00128   void
00129   SiStripDigitizer::accumulate(PileUpEventPrincipal const& iEvent, edm::EventSetup const& iSetup) {
00130 
00131     edm::ESHandle<TrackerTopology> tTopoHand;
00132     iSetup.get<IdealGeometryRecord>().get(tTopoHand);
00133     const TrackerTopology *tTopo=tTopoHand.product();
00134 
00135     // Step A: Get Inputs
00136     for(vstring::const_iterator i = trackerContainers.begin(), iEnd = trackerContainers.end(); i != iEnd; ++i) {
00137       edm::Handle<std::vector<PSimHit> > simHits;
00138       edm::InputTag tag(hitsProducer, *i);
00139 
00140       iEvent.getByLabel(tag, simHits);
00141       accumulateStripHits(simHits,tTopo);
00142     }
00143   }
00144 
00145 
00146 void SiStripDigitizer::initializeEvent(edm::Event const& iEvent, edm::EventSetup const& iSetup) {
00147   // Step A: Get Inputs
00148 
00149   if(useConfFromDB){
00150     edm::ESHandle<SiStripDetCabling> detCabling;
00151     iSetup.get<SiStripDetCablingRcd>().get(detCabling);
00152     detCabling->addConnected(theDetIdList);
00153   }
00154 
00155   theDigiAlgo->initializeEvent(iSetup);
00156 
00157   iSetup.get<TrackerDigiGeometryRecord>().get(geometryType,pDD);
00158   iSetup.get<IdealMagneticFieldRecord>().get(pSetup);
00159 
00160   // FIX THIS! We only need to clear and (re)fill detectorUnits when the geometry type IOV changes.  Use ESWatcher to determine this.
00161   bool changes = true;
00162   if(changes) { // Replace with ESWatcher
00163     detectorUnits.clear();
00164   }
00165   for(TrackingGeometry::DetUnitContainer::const_iterator iu = pDD->detUnits().begin(); iu != pDD->detUnits().end(); ++iu) {
00166     unsigned int detId = (*iu)->geographicalId().rawId();
00167     DetId idet=DetId(detId);
00168     unsigned int isub=idet.subdetId();
00169     if((isub == StripSubdetector::TIB) ||
00170        (isub == StripSubdetector::TID) ||
00171        (isub == StripSubdetector::TOB) ||
00172        (isub == StripSubdetector::TEC)) {
00173       StripGeomDetUnit* stripdet = dynamic_cast<StripGeomDetUnit*>((*iu));
00174       assert(stripdet != 0);
00175       if(changes) { // Replace with ESWatcher
00176         detectorUnits.insert(std::make_pair(detId, stripdet));
00177       }
00178       theDigiAlgo->initializeDetUnit(stripdet, iSetup);
00179     }
00180   }
00181 }
00182 
00183 void SiStripDigitizer::finalizeEvent(edm::Event& iEvent, edm::EventSetup const& iSetup) {
00184   edm::ESHandle<SiStripGain> gainHandle;
00185   edm::ESHandle<SiStripNoises> noiseHandle;
00186   edm::ESHandle<SiStripThreshold> thresholdHandle;
00187   edm::ESHandle<SiStripPedestals> pedestalHandle;
00188   iSetup.get<SiStripGainSimRcd>().get(gainLabel,gainHandle);
00189   iSetup.get<SiStripNoisesRcd>().get(noiseHandle);
00190   iSetup.get<SiStripThresholdRcd>().get(thresholdHandle);
00191   iSetup.get<SiStripPedestalsRcd>().get(pedestalHandle);
00192 
00193   std::vector<edm::DetSet<SiStripDigi> > theDigiVector;
00194   std::vector<edm::DetSet<SiStripRawDigi> > theRawDigiVector;
00195 
00196   
00197   // Step B: LOOP on StripGeomDetUnit
00198   theDigiVector.reserve(10000);
00199   theDigiVector.clear();
00200 
00201   for(TrackingGeometry::DetUnitContainer::const_iterator iu = pDD->detUnits().begin(); iu != pDD->detUnits().end(); iu ++){
00202     if(useConfFromDB){
00203       //apply the cable map _before_ digitization: consider only the detis that are connected 
00204       if(theDetIdList.find((*iu)->geographicalId().rawId())==theDetIdList.end())
00205         continue;
00206     }
00207     StripGeomDetUnit* sgd = dynamic_cast<StripGeomDetUnit*>((*iu));
00208     if (sgd != 0){
00209       edm::DetSet<SiStripDigi> collectorZS((*iu)->geographicalId().rawId());
00210       edm::DetSet<SiStripRawDigi> collectorRaw((*iu)->geographicalId().rawId());
00211       theDigiAlgo->digitize(collectorZS,collectorRaw,sgd,
00212                        gainHandle,thresholdHandle,noiseHandle,pedestalHandle);
00213       if(zeroSuppression){
00214         if(collectorZS.data.size()>0){
00215           theDigiVector.push_back(collectorZS);
00216         }
00217       }else{
00218         if(collectorRaw.data.size()>0){
00219           theRawDigiVector.push_back(collectorRaw);
00220         }
00221       }
00222     }
00223   }
00224 
00225   if(zeroSuppression){
00226     // Step C: create output collection
00227     std::auto_ptr<edm::DetSetVector<SiStripRawDigi> > output_virginraw(new edm::DetSetVector<SiStripRawDigi>());
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>(theDigiVector) );
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   }else{
00237     // Step C: create output collection
00238     std::auto_ptr<edm::DetSetVector<SiStripRawDigi> > output_virginraw(new edm::DetSetVector<SiStripRawDigi>(theRawDigiVector));
00239     std::auto_ptr<edm::DetSetVector<SiStripRawDigi> > output_scopemode(new edm::DetSetVector<SiStripRawDigi>());
00240     std::auto_ptr<edm::DetSetVector<SiStripRawDigi> > output_processedraw(new edm::DetSetVector<SiStripRawDigi>());
00241     std::auto_ptr<edm::DetSetVector<SiStripDigi> > output(new edm::DetSetVector<SiStripDigi>() );
00242     // Step D: write output to file
00243     iEvent.put(output, ZSDigi);
00244     iEvent.put(output_scopemode, SCDigi);
00245     iEvent.put(output_virginraw, VRDigi);
00246     iEvent.put(output_processedraw, PRDigi);
00247   }
00248 }