RPDigiProducer/plugins/RPDigiProducer.cc

0001 // user include files
0002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0003 #include "FWCore/Framework/interface/Frameworkfwd.h"
0004 #include "FWCore/Framework/interface/EDProducer.h"
0005 #include "FWCore/Framework/interface/EventSetup.h"
0006 #include "FWCore/Framework/interface/Event.h"
0007 #include "FWCore/Framework/interface/MakerMacros.h"
0008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0009 #include "FWCore/Utilities/interface/Exception.h"
0010 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0012
0013 #include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
0014 #include "SimDataFormats/TrackingHit/interface/PSimHit.h"
0015 #include "SimDataFormats/CrossingFrame/interface/CrossingFrame.h"
0016 #include "SimDataFormats/CrossingFrame/interface/MixCollection.h"
0017
0018 #include "DataFormats/CTPPSDigi/interface/TotemRPDigi.h"
0019 #include "DataFormats/CTPPSDigi/interface/TotemRPDigi.h"
0020 #include "DataFormats/Common/interface/DetSetVector.h"
0021 #include "DataFormats/CTPPSDetId/interface/TotemRPDetId.h"
0022 #include "DataFormats/Common/interface/DetSet.h"
0023
0024 //Random Number
0025 #include "FWCore/ServiceRegistry/interface/Service.h"
0026 #include "FWCore/Utilities/interface/RandomNumberGenerator.h"
0027
0028 #include "CondFormats/PPSObjects/interface/TotemAnalysisMask.h"
0029 #include "CondFormats/DataRecord/interface/TotemReadoutRcd.h"
0030 #include "SimPPS/RPDigiProducer/interface/RPSimTypes.h"
0031 #include "SimPPS/RPDigiProducer/plugins/RPDetDigitizer.h"
0032 #include "SimPPS/RPDigiProducer/plugins/DeadChannelsManager.h"
0033
0034 // system include files
0035 #include <memory>
0036 #include <vector>
0037 #include <map>
0038 #include <string>
0039 #include <iostream>
0040 #include <cstdlib>  // I need it for random numbers
0041
0042 // user include files
0043
0044 //
0045 // class decleration
0046 //
0047
0048 namespace CLHEP {
0049   class HepRandomEngine;
0050 }
0051
0052 class RPDigiProducer : public edm::EDProducer {
0053 public:
0054   explicit RPDigiProducer(const edm::ParameterSet&);
0055   ~RPDigiProducer() override = default;
0056
0057   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0058
0059 private:
0060   void beginRun(const edm::Run&, const edm::EventSetup&) override;
0061   void produce(edm::Event&, const edm::EventSetup&) override;
0062
0063   edm::DetSet<TotemRPDigi> convertRPStripDetSet(const edm::DetSet<TotemRPDigi>&);
0064
0065   // ----------member data ---------------------------
0066   std::vector<std::string> RP_hit_containers_;
0067   typedef std::map<unsigned int, std::vector<PSimHit>> simhit_map;
0068   typedef simhit_map::iterator simhit_map_iterator;
0069
0070   edm::ParameterSet conf_;
0071   std::map<RPDetId, std::unique_ptr<RPDetDigitizer>> theAlgoMap;
0072
0073   CLHEP::HepRandomEngine* rndEngine_ = nullptr;
0074   int verbosity_;
0075
0076   /**
0077        * this variable answers the question whether given channel is dead or not
0078        */
0079   DeadChannelsManager deadChannelsManager;
0080   /**
0081        * this variable indicates whether we take into account dead channels or simulate as if all
0082        * channels work ok (by default we do not simulate dead channels)
0083        */
0084   bool simulateDeadChannels;
0085
0086   edm::EDGetTokenT<CrossingFrame<PSimHit>> tokenCrossingFrameTotemRP;
0087   edm::ESGetToken<TotemAnalysisMask, TotemReadoutRcd> tokenAnalysisMask;
0088 };
0089
0090 RPDigiProducer::RPDigiProducer(const edm::ParameterSet& conf) : conf_(conf) {
0091   //now do what ever other initialization is needed
0092   produces<edm::DetSetVector<TotemRPDigi>>();
0093
0094   // register data to consume
0095   tokenCrossingFrameTotemRP = consumes<CrossingFrame<PSimHit>>(edm::InputTag("mix", "g4SimHitsTotemHitsRP", ""));
0096
0097   RP_hit_containers_ = conf.getParameter<std::vector<std::string>>("ROUList");
0098   verbosity_ = conf.getParameter<int>("RPVerbosity");
0099
0100   simulateDeadChannels = false;
0101   if (conf.exists(
0102           "simulateDeadChannels")) {  //check if "simulateDeadChannels" variable is defined in configuration file
0103     simulateDeadChannels = conf.getParameter<bool>("simulateDeadChannels");
0104   }
0105   if (simulateDeadChannels) {
0106     tokenAnalysisMask = esConsumes();
0107   }
0108 }
0109
0110 //
0111 // member functions
0112 //
0113
0114 // ------------ method called to produce the data  ------------
0115 void RPDigiProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0116   using namespace edm;
0117
0118   // initialize random engine
0119   if (!rndEngine_) {
0120     Service<RandomNumberGenerator> rng;
0121     if (!rng.isAvailable()) {
0122       throw cms::Exception("Configuration")
0123           << "This class requires the RandomNumberGeneratorService\n"
0124              "which is not present in the configuration file.  You must add the service\n"
0125              "in the configuration file or remove the modules that require it.";
0126     }
0127     rndEngine_ = &(rng->getEngine(iEvent.streamID()));
0128   }
0129
0130   // Step A: Get Inputs
0131   edm::Handle<CrossingFrame<PSimHit>> cf;
0132   iEvent.getByLabel("mix", "g4SimHitsTotemHitsRP", cf);
0133
0134   if (verbosity_) {
0135     edm::LogInfo("RPDigiProducer") << "\n\n=================== Starting SimHit access"
0136                                    << "  ==================="
0137                                    << "\n";
0138
0139     MixCollection<PSimHit> col{cf.product(), std::pair(-0, 0)};
0140     MixCollection<PSimHit>::iterator cfi;
0141     int count = 0;
0142     for (cfi = col.begin(); cfi != col.end(); cfi++) {
0143       edm::LogInfo("RPDigiProducer") << " Hit " << count << " has tof " << cfi->timeOfFlight() << " trackid "
0144                                      << cfi->trackId() << " bunchcr " << cfi.bunch() << " trigger " << cfi.getTrigger()
0145                                      << ", from EncodedEventId: " << cfi->eventId().bunchCrossing() << " "
0146                                      << cfi->eventId().event() << " bcr from MixCol " << cfi.bunch() << "\n";
0147       edm::LogInfo("RPDigiProducer") << " Hit: " << (*cfi) << "\n";
0148       count++;
0149     }
0150   }
0151
0152   MixCollection<PSimHit> allRPHits{cf.product(), std::pair(0, 0)};
0153
0154   if (verbosity_)
0155     edm::LogInfo("RPDigiProducer") << "Input MixCollection size = " << allRPHits.size() << "\n";
0156
0157   //Loop on PSimHit
0158   simhit_map simHitMap_;
0159   simHitMap_.clear();
0160
0161   MixCollection<PSimHit>::iterator isim;
0162   for (isim = allRPHits.begin(); isim != allRPHits.end(); ++isim) {
0163     simHitMap_[(*isim).detUnitId()].push_back((*isim));
0164   }
0165
0166   // Step B: LOOP on hits in event
0167   std::vector<edm::DetSet<TotemRPDigi>> DigiVector;
0168   DigiVector.reserve(400);
0169   DigiVector.clear();
0170
0171   for (simhit_map_iterator it = simHitMap_.begin(); it != simHitMap_.end(); ++it) {
0172     edm::DetSet<TotemRPDigi> digi_collector(it->first);
0173
0174     if (theAlgoMap.find(it->first) == theAlgoMap.end()) {
0175       theAlgoMap[it->first] = std::make_unique<RPDetDigitizer>(conf_, *rndEngine_, it->first, iSetup);
0176     }
0177
0178     std::vector<int> input_links;
0179     simromanpot::DigiPrimaryMapType output_digi_links;
0180
0181     (theAlgoMap.find(it->first)->second)
0182         ->run(simHitMap_[it->first], input_links, digi_collector.data, output_digi_links);
0183
0184     if (!digi_collector.data.empty()) {
0185       DigiVector.push_back(convertRPStripDetSet(digi_collector));
0186     }
0187   }
0188
0189   // Step C: create empty output collection
0190   std::unique_ptr<edm::DetSetVector<TotemRPDigi>> digi_output(new edm::DetSetVector<TotemRPDigi>(DigiVector));
0191
0192   if (verbosity_) {
0193     edm::LogInfo("RPDigiProducer") << "digi_output->size()=" << digi_output->size() << "\n";
0194   }
0195   // Step D: write output to file
0196   iEvent.put(std::move(digi_output));
0197 }
0198
0199 // ------------ method called once each job just before starting event loop  ------------
0200 void RPDigiProducer::beginRun(const edm::Run& beginrun, const edm::EventSetup& es) {
0201   // get analysis mask to mask channels
0202   if (simulateDeadChannels) {
0203     //set analysisMask in deadChannelsManager
0204     deadChannelsManager = DeadChannelsManager(&es.getData(tokenAnalysisMask));
0205   }
0206 }
0207
0208 edm::DetSet<TotemRPDigi> RPDigiProducer::convertRPStripDetSet(const edm::DetSet<TotemRPDigi>& rpstrip_detset) {
0209   edm::DetSet<TotemRPDigi> rpdigi_detset(rpstrip_detset.detId());
0210   rpdigi_detset.reserve(rpstrip_detset.size());
0211
0212   for (std::vector<TotemRPDigi>::const_iterator stripIterator = rpstrip_detset.data.begin();
0213        stripIterator < rpstrip_detset.data.end();
0214        ++stripIterator) {
0215     rpdigi_detset.push_back(TotemRPDigi(stripIterator->stripNumber()));
0216   }
0217
0218   return rpdigi_detset;
0219 }
0220
0221 void RPDigiProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0222   //RPSiDetDigitizer
0223   //all distances in [mm]
0224   edm::ParameterSetDescription desc;
0225   desc.add<bool>("RPLandauFluctuations", true);
0226   desc.add<bool>("RPDisplacementOn", false);
0227   desc.add<int>("RPVerbosity", 0);
0228   desc.add<double>("RPVFATThreshold", 9000.0);
0229   desc.add<double>("RPTopEdgePosition", 1.5);
0230   desc.add<double>("RPActiveEdgeSmearing", 0.013);
0231   desc.add<double>("RPEquivalentNoiseCharge300um", 1000.0);
0232   desc.add<int>("RPVFATTriggerMode", 2);
0233   desc.add<std::vector<double>>("RPInterStripSmearing",
0234                                 {
0235                                     0.011,
0236                                 });
0237   desc.add<double>("RPSharingSigmas", 5.0);  //how many sigmas taken into account for the edges and inter strips
0238   desc.add<double>("RPGeVPerElectron", 3.61e-09);
0239   desc.add<double>("RPActiveEdgePosition", 0.034);  //from the physical edge
0240   desc.add<bool>("RPDeadStripSimulationOn", false);
0241   desc.add<std::vector<std::string>>("ROUList",
0242                                      {
0243                                          "TotemHitsRP",
0244                                      });
0245   desc.add<bool>("RPNoNoise", false);
0246   desc.add<bool>("RPDigiSimHitRelationsPresistence", false);  //save links betweend digi, clusters and OSCAR/Geant4 hits
0247   desc.add<std::string>("mixLabel", "mix");
0248   desc.add<int>("RPChargeDivisionsPerThickness", 5);
0249   desc.add<double>("RPDeltaProductionCut", 0.120425);  //[MeV]
0250   desc.add<double>("RPBottomEdgePosition", 1.5);
0251   desc.add<double>("RPBottomEdgeSmearing", 0.011);
0252   desc.add<double>("RPTopEdgeSmearing", 0.011);
0253   desc.add<std::string>("InputCollection", "g4SimHitsTotemHitsRP");
0254   desc.add<double>("RPInterStripCoupling",
0255                    1.0);  //fraction of charge going to the strip, the missing part is taken by its neighbours
0256   desc.add<double>("RPDeadStripProbability", 0.001);
0257   desc.add<int>("RPChargeDivisionsPerStrip", 15);
0258   descriptions.add("RPSiDetDigitizer", desc);
0259 }
0260
0261 DEFINE_FWK_MODULE(RPDigiProducer);