Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-04-04 01:27:03

0001 #include "FWCore/AbstractServices/interface/RandomNumberGenerator.h"
0002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0003 #include "FWCore/ServiceRegistry/interface/Service.h"
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 #include "DataFormats/Common/interface/Handle.h"
0007 #include "FWCore/Framework/interface/ConsumesCollector.h"
0008 #include "FWCore/Framework/interface/ProducesCollector.h"
0009 #include "SimGeneral/MixingModule/interface/PileUpEventPrincipal.h"
0010 
0011 //Data Formats
0012 #include "DataFormats/Common/interface/DetSetVector.h"
0013 #include "DataFormats/Common/interface/DetSet.h"
0014 #include "DataFormats/SiPixelDigi/interface/PixelDigi.h"
0015 #include "SimDataFormats/TrackerDigiSimLink/interface/PixelDigiSimLink.h"  // not really needed
0016 #include "SimDataFormats/PileupSummaryInfo/interface/PileupSummaryInfo.h"
0017 
0018 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0019 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0020 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0021 
0022 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0023 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0024 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0025 #include "Geometry/CommonTopologies/interface/PixelTopology.h"
0026 #include "Geometry/CommonDetUnit/interface/PixelGeomDetType.h"
0027 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0028 #include "CondFormats/SiPixelObjects/interface/PixelIndices.h"
0029 
0030 #include "CLHEP/Random/RandFlat.h"
0031 
0032 #include "SimGeneral/PreMixingModule/interface/PreMixingWorker.h"
0033 #include "SimGeneral/PreMixingModule/interface/PreMixingWorkerFactory.h"
0034 
0035 #include "SiPixelDigitizerAlgorithm.h"
0036 
0037 #include <map>
0038 #include <memory>
0039 
0040 class PreMixingSiPixelWorker : public PreMixingWorker {
0041 public:
0042   PreMixingSiPixelWorker(const edm::ParameterSet& ps, edm::ProducesCollector, edm::ConsumesCollector&& iC);
0043   ~PreMixingSiPixelWorker() override = default;
0044 
0045   void initializeEvent(edm::Event const& e, edm::EventSetup const& c) override;
0046   void addSignals(edm::Event const& e, edm::EventSetup const& es) override;
0047   void addPileups(PileUpEventPrincipal const&, edm::EventSetup const& es) override;
0048   void put(edm::Event& e, edm::EventSetup const& iSetup, std::vector<PileupSummaryInfo> const& ps, int bs) override;
0049 
0050 private:
0051   edm::InputTag pixeldigi_collectionSig_;   // secondary name given to collection of SiPixel digis
0052   edm::InputTag pixeldigi_collectionPile_;  // secondary name given to collection of SiPixel digis
0053   edm::InputTag pixeldigi_extraInfo_;       // secondary name given to collection of SiPixel digis
0054   edm::InputTag pixeldigi_extraInfoLite_;   // secondary name given to collection of SiPixel digis
0055   std::string PixelDigiCollectionDM_;       // secondary name to be given to new SiPixel digis
0056 
0057   edm::EDGetTokenT<edm::DetSetVector<PixelDigi>> PixelDigiToken_;   // Token to retrieve information
0058   edm::EDGetTokenT<edm::DetSetVector<PixelDigi>> PixelDigiPToken_;  // Token to retrieve information
0059   edm::EDGetTokenT<edm::DetSetVector<PixelSimHitExtraInfo>> PixelDigiPExtraToken_;
0060   edm::EDGetTokenT<edm::DetSetVector<PixelSimHitExtraInfoLite>> PixelDigiPExtraLiteToken_;  // Lite version
0061   const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> tTopoToken_;
0062   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> pDDToken_;
0063 
0064   SiPixelDigitizerAlgorithm digitizer_;
0065 
0066   //
0067   // Internal typedefs
0068 
0069   typedef int Amplitude;
0070   typedef std::map<int, Amplitude, std::less<int>> signal_map_type;  // from Digi.Skel.
0071   typedef std::map<uint32_t, signal_map_type> signalMaps;
0072 
0073   typedef std::multimap<int, PixelDigi>
0074       OneDetectorMap;  // maps by pixel ID for later combination - can have duplicate pixels
0075   typedef std::map<uint32_t, OneDetectorMap> SiGlobalIndex;  // map to all data for each detector ID
0076   typedef std::multimap<int, PixelSimHitExtraInfo> OneExtraInfoMap;
0077   typedef std::map<uint32_t, OneExtraInfoMap> SiPixelExtraInfo;
0078   typedef std::multimap<int, PixelSimHitExtraInfoLite> OneExtraInfoLiteMap;  // Lite version
0079   typedef std::map<uint32_t, OneExtraInfoLiteMap> SiPixelExtraInfoLite;      // Lite version
0080 
0081   SiGlobalIndex SiHitStorage_;
0082   SiPixelExtraInfo SiHitExtraStorage_;
0083   SiPixelExtraInfoLite SiHitExtraLiteStorage_;  // Lite version
0084 
0085   bool firstInitializeEvent_ = true;
0086   bool firstFinalizeEvent_ = true;
0087   bool applyLateReweighting_;
0088   bool usePixelExtraLiteFormat_;
0089 };
0090 
0091 // Constructor
0092 PreMixingSiPixelWorker::PreMixingSiPixelWorker(const edm::ParameterSet& ps,
0093                                                edm::ProducesCollector producesCollector,
0094                                                edm::ConsumesCollector&& iC)
0095     : tTopoToken_(iC.esConsumes()),
0096       pDDToken_(iC.esConsumes(edm::ESInputTag("", ps.getParameter<std::string>("PixGeometryType")))),
0097       digitizer_(ps, iC) {
0098   // declare the products to produce
0099 
0100   pixeldigi_collectionSig_ = ps.getParameter<edm::InputTag>("pixeldigiCollectionSig");
0101   pixeldigi_collectionPile_ = ps.getParameter<edm::InputTag>("pixeldigiCollectionPile");
0102   pixeldigi_extraInfo_ = ps.getParameter<edm::InputTag>("pixeldigiExtraCollectionPile");
0103   pixeldigi_extraInfoLite_ = ps.getParameter<edm::InputTag>("pixeldigiExtraLiteCollectionPile");  // Lite version
0104   PixelDigiCollectionDM_ = ps.getParameter<std::string>("PixelDigiCollectionDM");
0105   applyLateReweighting_ = ps.getParameter<bool>("applyLateReweighting");
0106   usePixelExtraLiteFormat_ = ps.getParameter<bool>("usePixelExtraLiteFormat");
0107   LogDebug("PreMixingSiPixelWorker") << "applyLateReweighting_ in PreMixingSiPixelWorker  " << applyLateReweighting_;
0108 
0109   LogDebug("PreMixingSiPixelWorker") << " applyLateReweighting_ in PreMixingSiPixelWorker " << applyLateReweighting_;
0110   LogDebug("PreMixingSiPixelWorker") << " usePixelExtraLiteFormat_ in PreMixingSiPixelWorker "
0111                                      << usePixelExtraLiteFormat_;
0112 
0113   PixelDigiToken_ = iC.consumes<edm::DetSetVector<PixelDigi>>(pixeldigi_collectionSig_);
0114   PixelDigiPToken_ = iC.consumes<edm::DetSetVector<PixelDigi>>(pixeldigi_collectionPile_);
0115   PixelDigiPExtraToken_ = iC.consumes<edm::DetSetVector<PixelSimHitExtraInfo>>(pixeldigi_extraInfo_);
0116   PixelDigiPExtraLiteToken_ =
0117       iC.consumes<edm::DetSetVector<PixelSimHitExtraInfoLite>>(pixeldigi_extraInfoLite_);  // Lite version
0118 
0119   producesCollector.produces<edm::DetSetVector<PixelDigi>>(PixelDigiCollectionDM_);
0120   producesCollector.produces<PixelFEDChannelCollection>(PixelDigiCollectionDM_);
0121 
0122   // clear local storage for this event
0123   SiHitStorage_.clear();
0124   SiHitExtraStorage_.clear();
0125   SiHitExtraLiteStorage_.clear();
0126 }
0127 
0128 // Need an event initialization
0129 
0130 void PreMixingSiPixelWorker::initializeEvent(edm::Event const& e, edm::EventSetup const& iSetup) {
0131   if (firstInitializeEvent_) {
0132     digitizer_.init(iSetup);
0133     firstInitializeEvent_ = false;
0134   }
0135   digitizer_.initializeEvent();
0136 }
0137 
0138 void PreMixingSiPixelWorker::addSignals(edm::Event const& e, edm::EventSetup const& es) {
0139   // fill in maps of hits
0140 
0141   LogDebug("PreMixingSiPixelWorker") << "===============> adding MC signals for " << e.id();
0142 
0143   edm::Handle<edm::DetSetVector<PixelDigi>> input;
0144 
0145   if (e.getByToken(PixelDigiToken_, input)) {
0146     //loop on all detsets (detectorIDs) inside the input collection
0147     edm::DetSetVector<PixelDigi>::const_iterator DSViter = input->begin();
0148     for (; DSViter != input->end(); DSViter++) {
0149 #ifdef DEBUG
0150       LogDebug("PreMixingSiPixelWorker") << "Processing DetID " << DSViter->id;
0151 #endif
0152 
0153       uint32_t detID = DSViter->id;
0154       edm::DetSet<PixelDigi>::const_iterator begin = (DSViter->data).begin();
0155       edm::DetSet<PixelDigi>::const_iterator end = (DSViter->data).end();
0156       edm::DetSet<PixelDigi>::const_iterator icopy;
0157 
0158       OneDetectorMap LocalMap;
0159 
0160       for (icopy = begin; icopy != end; icopy++) {
0161         LocalMap.insert(OneDetectorMap::value_type((icopy->channel()), *icopy));
0162       }
0163 
0164       SiHitStorage_.insert(SiGlobalIndex::value_type(detID, LocalMap));
0165     }
0166   }
0167 }  // end of addSiPixelSignals
0168 
0169 void PreMixingSiPixelWorker::addPileups(PileUpEventPrincipal const& pep, edm::EventSetup const& es) {
0170   LogDebug("PreMixingSiPixelWorker") << "\n===============> adding pileups from event  " << pep.principal().id()
0171                                      << " for bunchcrossing " << pep.bunchCrossing();
0172 
0173   // fill in maps of hits; same code as addSignals, except now applied to the pileup events
0174 
0175   edm::Handle<edm::DetSetVector<PixelDigi>> inputHandle;
0176   pep.getByLabel(pixeldigi_collectionPile_, inputHandle);
0177 
0178   // added for the Late CR
0179   edm::Handle<edm::DetSetVector<PixelSimHitExtraInfo>> pixelAddInfo;
0180   pep.getByLabel(pixeldigi_extraInfo_, pixelAddInfo);
0181   edm::Handle<edm::DetSetVector<PixelSimHitExtraInfoLite>> pixelAddInfoLite;  // Lite version
0182   pep.getByLabel(pixeldigi_extraInfoLite_, pixelAddInfoLite);                 // Lite version
0183   const TrackerTopology* tTopo = &es.getData(tTopoToken_);
0184   auto const& pDD = es.getData(pDDToken_);
0185   edm::Service<edm::RandomNumberGenerator> rng;
0186   CLHEP::HepRandomEngine* engine = &rng->getEngine(pep.principal().streamID());
0187 
0188   if (inputHandle.isValid()) {
0189     const auto& input = *inputHandle;
0190 
0191     bool loadExtraInformation = false;
0192 
0193     if (!usePixelExtraLiteFormat_) {
0194       if (pixelAddInfo.isValid() && applyLateReweighting_) {
0195         // access the extra information
0196         loadExtraInformation = true;
0197         // Iterate on detector units
0198         edm::DetSetVector<PixelSimHitExtraInfo>::const_iterator detIdIter;
0199         for (detIdIter = pixelAddInfo->begin(); detIdIter != pixelAddInfo->end(); detIdIter++) {
0200           uint32_t detid = detIdIter->id;  // = rawid
0201           OneExtraInfoMap LocalExtraMap;
0202           edm::DetSet<PixelSimHitExtraInfo>::const_iterator di;
0203           for (di = detIdIter->data.begin(); di != detIdIter->data.end(); di++) {
0204             LocalExtraMap.insert(OneExtraInfoMap::value_type((di->hitIndex()), *di));
0205           }
0206           SiHitExtraStorage_.insert(SiPixelExtraInfo::value_type(detid, LocalExtraMap));
0207         }  // end loop on detIdIter
0208       }  // end if applyLateReweighting_
0209       else if (!pixelAddInfo.isValid() && applyLateReweighting_) {
0210         edm::LogError("PreMixingSiPixelWorker") << " Problem in accessing the Extra Pixel SimHit Collection  !!!! ";
0211         edm::LogError("PreMixingSiPixelWorker") << " The Late Charge Reweighting can not be applied ";
0212         throw cms::Exception("PreMixingSiPixelWorker")
0213             << " Problem in accessing the Extra Pixel SimHit Collection for Late Charge Reweighting \n";
0214       }
0215     } else {  // Lite version
0216       if (pixelAddInfoLite.isValid() && applyLateReweighting_) {
0217         // access the extra information
0218         loadExtraInformation = true;
0219         // Iterate on detector units
0220         edm::DetSetVector<PixelSimHitExtraInfoLite>::const_iterator detIdIter;
0221         for (detIdIter = pixelAddInfoLite->begin(); detIdIter != pixelAddInfoLite->end(); detIdIter++) {
0222           uint32_t detid = detIdIter->id;  // = rawid
0223           OneExtraInfoLiteMap LocalExtraMap;
0224           edm::DetSet<PixelSimHitExtraInfoLite>::const_iterator di;
0225           for (di = detIdIter->data.begin(); di != detIdIter->data.end(); di++) {
0226             LocalExtraMap.insert(OneExtraInfoLiteMap::value_type((di->hitIndex()), *di));
0227           }
0228           SiHitExtraLiteStorage_.insert(SiPixelExtraInfoLite::value_type(detid, LocalExtraMap));
0229         }  // end loop on detIdIter
0230       }  // end if applyLateReweighting_
0231       else if (!pixelAddInfoLite.isValid() && applyLateReweighting_) {
0232         edm::LogError("PreMixingSiPixelWorker")
0233             << " Problem in accessing the Extra Pixel SimHit Lite Collection  !!!! ";
0234         edm::LogError("PreMixingSiPixelWorker") << " The Late Charge Reweighting can not be applied ";
0235         throw cms::Exception("PreMixingSiPixelWorker")
0236             << " Problem in accessing the Extra Pixel SimHit Lite Collection for Late Charge Reweighting \n";
0237       }
0238     }  // end if !usePixelExtraLiteFormat_
0239 
0240     //loop on all detsets (detectorIDs) inside the input collection
0241     edm::DetSetVector<PixelDigi>::const_iterator DSViter = input.begin();
0242     for (; DSViter != input.end(); DSViter++) {
0243 #ifdef DEBUG
0244       LogDebug("PreMixingSiPixelWorker") << "Pileups: Processing DetID " << DSViter->id;
0245 #endif
0246 
0247       uint32_t detID = DSViter->id;
0248       edm::DetSet<PixelDigi>::const_iterator begin = (DSViter->data).begin();
0249       edm::DetSet<PixelDigi>::const_iterator end = (DSViter->data).end();
0250       edm::DetSet<PixelDigi>::const_iterator icopy;
0251 
0252       // find correct local map (or new one) for this detector ID
0253 
0254       SiGlobalIndex::const_iterator itest;
0255 
0256       itest = SiHitStorage_.find(detID);
0257 
0258       std::vector<PixelDigi> TempDigis;
0259       for (icopy = begin; icopy != end; icopy++) {
0260         TempDigis.push_back(*icopy);
0261       }
0262       if (loadExtraInformation) {
0263         // apply the Late Charge Reweighthing on Pile-up digi
0264         if (!usePixelExtraLiteFormat_) {
0265           SiPixelExtraInfo::const_iterator jtest;
0266           jtest = SiHitExtraStorage_.find(detID);
0267           OneExtraInfoMap LocalSimHitExtraMap = jtest->second;
0268           std::vector<PixelSimHitExtraInfo> TempSimExtra;
0269           for (const auto& iLocal : LocalSimHitExtraMap) {
0270             TempSimExtra.push_back(iLocal.second);
0271           }
0272 
0273           for (const auto& iu : pDD.detUnits()) {
0274             if (iu->type().isTrackerPixel()) {
0275               uint32_t detIDinLoop = iu->geographicalId().rawId();
0276               if (detIDinLoop == detID) {
0277                 digitizer_.lateSignalReweight(
0278                     dynamic_cast<const PixelGeomDetUnit*>(iu), TempDigis, TempSimExtra, tTopo, engine);
0279                 break;
0280               }
0281             }
0282           }
0283         } else {  // Lite version
0284           SiPixelExtraInfoLite::const_iterator jtest;
0285           jtest = SiHitExtraLiteStorage_.find(detID);
0286           OneExtraInfoLiteMap LocalSimHitExtraMap = jtest->second;
0287           std::vector<PixelSimHitExtraInfoLite> TempSimExtra;
0288           for (const auto& iLocal : LocalSimHitExtraMap) {
0289             TempSimExtra.push_back(iLocal.second);
0290           }
0291 
0292           for (const auto& iu : pDD.detUnits()) {
0293             if (iu->type().isTrackerPixel()) {
0294               uint32_t detIDinLoop = iu->geographicalId().rawId();
0295               if (detIDinLoop == detID) {
0296                 digitizer_.lateSignalReweight(
0297                     dynamic_cast<const PixelGeomDetUnit*>(iu), TempDigis, TempSimExtra, tTopo, engine);
0298                 break;
0299               }
0300             }
0301           }
0302         }  // end if !usePixelExtraLiteFormat_
0303       }
0304 
0305       if (itest != SiHitStorage_.end()) {  // this detID already has hits, add to existing map
0306 
0307         OneDetectorMap LocalMap = itest->second;
0308 
0309         // fill in local map with extra channels
0310         for (unsigned int ij = 0; ij < TempDigis.size(); ij++) {
0311           LocalMap.insert(OneDetectorMap::value_type((TempDigis[ij].channel()), TempDigis[ij]));
0312         }
0313 
0314         SiHitStorage_[detID] = LocalMap;
0315 
0316       } else {  // fill local storage with this information, put in global collection
0317 
0318         OneDetectorMap LocalMap;
0319 
0320         for (unsigned int ij = 0; ij < TempDigis.size(); ij++) {
0321           LocalMap.insert(OneDetectorMap::value_type((TempDigis[ij].channel()), TempDigis[ij]));
0322         }
0323 
0324         SiHitStorage_.insert(SiGlobalIndex::value_type(detID, LocalMap));
0325       }
0326     }
0327   }
0328 }
0329 
0330 void PreMixingSiPixelWorker::put(edm::Event& e,
0331                                  edm::EventSetup const& iSetup,
0332                                  std::vector<PileupSummaryInfo> const& ps,
0333                                  int bs) {
0334   // collection of Digis to put in the event
0335 
0336   std::vector<edm::DetSet<PixelDigi>> vPixelDigi;
0337 
0338   // loop through our collection of detectors, merging hits and putting new ones in the output
0339   signalMaps signal;
0340 
0341   // big loop over Detector IDs:
0342 
0343   for (SiGlobalIndex::const_iterator IDet = SiHitStorage_.begin(); IDet != SiHitStorage_.end(); IDet++) {
0344     uint32_t detID = IDet->first;
0345 
0346     OneDetectorMap LocalMap = IDet->second;
0347 
0348     signal_map_type Signals;
0349     Signals.clear();
0350 
0351     //counter variables
0352     int formerPixel = -1;
0353     int currentPixel;
0354     int ADCSum = 0;
0355 
0356     OneDetectorMap::const_iterator iLocalchk;
0357 
0358     for (OneDetectorMap::const_iterator iLocal = LocalMap.begin(); iLocal != LocalMap.end(); ++iLocal) {
0359       currentPixel = iLocal->first;
0360 
0361       if (currentPixel == formerPixel) {  // we have to add these digis together
0362         ADCSum += (iLocal->second).adc();
0363       } else {
0364         if (formerPixel != -1) {  // ADC info stolen from SiStrips...
0365           if (ADCSum > 511)
0366             ADCSum = 255;
0367           else if (ADCSum > 253 && ADCSum < 512)
0368             ADCSum = 254;
0369 
0370           Signals.insert(std::make_pair(formerPixel, ADCSum));
0371         }
0372         // save pointers for next iteration
0373         formerPixel = currentPixel;
0374         ADCSum = (iLocal->second).adc();
0375       }
0376 
0377       iLocalchk = iLocal;
0378       if ((++iLocalchk) == LocalMap.end()) {  //make sure not to lose the last one
0379         if (ADCSum > 511)
0380           ADCSum = 255;
0381         else if (ADCSum > 253 && ADCSum < 512)
0382           ADCSum = 254;
0383         Signals.insert(std::make_pair(formerPixel, ADCSum));
0384       }
0385 
0386     }  // end of loop over one detector
0387 
0388     // stick this into the global vector of detector info
0389     signal.insert(std::make_pair(detID, Signals));
0390 
0391   }  // end of big loop over all detector IDs
0392 
0393   // put the collection of digis in the event
0394   edm::LogInfo("PreMixingSiPixelWorker") << "total # Merged Pixels: " << signal.size();
0395 
0396   std::vector<edm::DetSet<PixelDigi>> theDigiVector;
0397 
0398   // Load inefficiency constants (1st pass), set pileup information.
0399   if (firstFinalizeEvent_) {
0400     digitizer_.init_DynIneffDB(iSetup);
0401     firstFinalizeEvent_ = false;
0402   }
0403 
0404   digitizer_.calculateInstlumiFactor(ps, bs);
0405   digitizer_.setSimAccumulator(signal);
0406 
0407   edm::Service<edm::RandomNumberGenerator> rng;
0408   CLHEP::HepRandomEngine* engine = &rng->getEngine(e.streamID());
0409 
0410   auto const& pDD = iSetup.getData(pDDToken_);
0411   const TrackerTopology* tTopo = &iSetup.getData(tTopoToken_);
0412 
0413   if (digitizer_.killBadFEDChannels()) {
0414     std::unique_ptr<PixelFEDChannelCollection> PixelFEDChannelCollection_ = digitizer_.chooseScenario(ps, engine);
0415     if (PixelFEDChannelCollection_ == nullptr) {
0416       throw cms::Exception("NullPointerError") << "PixelFEDChannelCollection not set in chooseScenario function.\n";
0417     }
0418     e.put(std::move(PixelFEDChannelCollection_), PixelDigiCollectionDM_);
0419   }
0420 
0421   for (const auto& iu : pDD.detUnits()) {
0422     if (iu->type().isTrackerPixel()) {
0423       edm::DetSet<PixelDigi> collector(iu->geographicalId().rawId());
0424       edm::DetSet<PixelDigiSimLink> linkcollector(
0425           iu->geographicalId().rawId());                // ignored as DigiSimLinks are combined separately
0426       std::vector<PixelDigiAddTempInfo> tempcollector;  // to be ignored ?
0427 
0428       digitizer_.digitize(
0429           dynamic_cast<const PixelGeomDetUnit*>(iu), collector.data, linkcollector.data, tempcollector, tTopo, engine);
0430       if (!collector.data.empty()) {
0431         theDigiVector.push_back(std::move(collector));
0432       }
0433     }
0434   }
0435 
0436   e.put(std::make_unique<edm::DetSetVector<PixelDigi>>(theDigiVector), PixelDigiCollectionDM_);
0437 
0438   // clear local storage for this event
0439   SiHitStorage_.clear();
0440   SiHitExtraStorage_.clear();
0441   SiHitExtraLiteStorage_.clear();
0442 }
0443 
0444 DEFINE_PREMIXING_WORKER(PreMixingSiPixelWorker);