Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 10:03:01

0001 #include "FWCore/Framework/interface/stream/EDProducer.h"
0002 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0003 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0004 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0005 #include "FWCore/Framework/interface/ConsumesCollector.h"
0006 #include "FWCore/Framework/interface/ProducesCollector.h"
0007 #include "FWCore/Utilities/interface/EDGetToken.h"
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "DataFormats/Common/interface/Handle.h"
0010 #include "FWCore/Utilities/interface/RunningAverage.h"
0011 
0012 #include "RecoTracker/TkTrackingRegions/interface/TrackingRegion.h"
0013 #include "DataFormats/Common/interface/OwnVector.h"
0014 #include "TrackingTools/TransientTrackingRecHit/interface/SeedingLayerSetsHits.h"
0015 #include "RecoTracker/TkTrackingRegions/interface/TrackingRegionsSeedingLayerSets.h"
0016 #include "RecoTracker/TkHitPairs/interface/LayerHitMapCache.h"
0017 #include "RecoTracker/TkHitPairs/interface/HitPairGeneratorFromLayerPair.h"
0018 #include "RecoTracker/TkHitPairs/interface/IntermediateHitDoublets.h"
0019 #include "RecoTracker/TkHitPairs/interface/RegionsSeedingHitSets.h"
0020 
0021 namespace {
0022   class ImplBase;
0023 }
0024 
0025 class HitPairEDProducer : public edm::stream::EDProducer<> {
0026 public:
0027   HitPairEDProducer(const edm::ParameterSet& iConfig);
0028   ~HitPairEDProducer() override = default;
0029 
0030   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0031 
0032   void produce(edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0033 
0034 private:
0035   edm::EDGetTokenT<bool> clusterCheckToken_;
0036 
0037   std::unique_ptr<::ImplBase> impl_;
0038 };
0039 
0040 namespace {
0041   class ImplBase {
0042   public:
0043     ImplBase(const edm::ParameterSet& iConfig, edm::ConsumesCollector iC);
0044     virtual ~ImplBase() = default;
0045 
0046     virtual void produces(edm::ProducesCollector) const = 0;
0047 
0048     virtual void produce(const bool clusterCheckOk, edm::Event& iEvent, const edm::EventSetup& iSetup) = 0;
0049 
0050   protected:
0051     edm::RunningAverage localRA_;
0052     const unsigned int maxElement_;
0053     const unsigned int maxElementTotal_;
0054     const bool putEmptyIfMaxElementReached_;
0055 
0056     HitPairGeneratorFromLayerPair generator_;
0057     std::vector<unsigned> layerPairBegins_;
0058   };
0059   ImplBase::ImplBase(const edm::ParameterSet& iConfig, edm::ConsumesCollector iC)
0060       : maxElement_(iConfig.getParameter<unsigned int>("maxElement")),
0061         maxElementTotal_(iConfig.getParameter<unsigned int>("maxElementTotal")),
0062         putEmptyIfMaxElementReached_(iConfig.getParameter<bool>("putEmptyIfMaxElementReached")),
0063         generator_(
0064             iC, 0, 1, nullptr, maxElement_),  // these indices are dummy, TODO: cleanup HitPairGeneratorFromLayerPair
0065         layerPairBegins_(iConfig.getParameter<std::vector<unsigned>>("layerPairs")) {
0066     if (layerPairBegins_.empty())
0067       throw cms::Exception("Configuration")
0068           << "HitPairEDProducer requires at least index for layer pairs (layerPairs parameter), none was given";
0069   }
0070 
0071   /////
0072   template <typename T_SeedingHitSets, typename T_IntermediateHitDoublets, typename T_RegionLayers>
0073   class Impl : public ImplBase {
0074   public:
0075     template <typename... Args>
0076     Impl(const edm::ParameterSet& iConfig, edm::ConsumesCollector iC, Args&&... args)
0077         : ImplBase(iConfig, iC), regionsLayers_(&layerPairBegins_, std::forward<Args>(args)..., iC) {}
0078     ~Impl() override = default;
0079 
0080     void produces(edm::ProducesCollector producesCollector) const override {
0081       T_SeedingHitSets::produces(producesCollector);
0082       T_IntermediateHitDoublets::produces(producesCollector);
0083     }
0084 
0085     void produce(const bool clusterCheckOk, edm::Event& iEvent, const edm::EventSetup& iSetup) override {
0086       auto regionsLayers = regionsLayers_.beginEvent(iEvent);
0087 
0088       auto seedingHitSetsProducer = T_SeedingHitSets(&localRA_);
0089       auto intermediateHitDoubletsProducer = T_IntermediateHitDoublets(regionsLayers.seedingLayerSetsHitsPtr());
0090 
0091       if (!clusterCheckOk) {
0092         seedingHitSetsProducer.putEmpty(iEvent);
0093         intermediateHitDoubletsProducer.putEmpty(iEvent);
0094         return;
0095       }
0096 
0097       seedingHitSetsProducer.reserve(regionsLayers.regionsSize());
0098       intermediateHitDoubletsProducer.reserve(regionsLayers.regionsSize());
0099 
0100       unsigned int nDoublets = 0;
0101 
0102       for (const auto& regionLayers : regionsLayers) {
0103         const TrackingRegion& region = regionLayers.region();
0104         auto hitCachePtr_filler_shs = seedingHitSetsProducer.beginRegion(&region, nullptr);
0105         auto hitCachePtr_filler_ihd =
0106             intermediateHitDoubletsProducer.beginRegion(&region, std::get<0>(hitCachePtr_filler_shs));
0107         auto hitCachePtr = std::get<0>(hitCachePtr_filler_ihd);
0108 
0109         for (SeedingLayerSetsHits::SeedingLayerSet layerSet : regionLayers.layerPairs()) {
0110           auto doubletsOpt = generator_.doublets(region, iEvent, iSetup, layerSet, *hitCachePtr);
0111           if (not doubletsOpt) {
0112             if (putEmptyIfMaxElementReached_) {
0113               putEmpty(iEvent, regionsLayers);
0114               return;
0115             } else {
0116               continue;
0117             }
0118           }
0119           auto& doublets = *doubletsOpt;
0120           LogTrace("HitPairEDProducer") << " created " << doublets.size() << " doublets for layers "
0121                                         << layerSet[0].index() << "," << layerSet[1].index();
0122           if (doublets.empty())
0123             continue;  // don't bother if no pairs from these layers
0124           nDoublets += doublets.size();
0125           if (nDoublets >= maxElementTotal_) {
0126             edm::LogError("TooManyPairs") << "number of total pairs exceed maximum, no pairs produced";
0127             putEmpty(iEvent, regionsLayers);
0128             return;
0129           }
0130           seedingHitSetsProducer.fill(std::get<1>(hitCachePtr_filler_shs), doublets);
0131           intermediateHitDoubletsProducer.fill(std::get<1>(hitCachePtr_filler_ihd), layerSet, std::move(doublets));
0132         }
0133       }
0134 
0135       seedingHitSetsProducer.put(iEvent);
0136       intermediateHitDoubletsProducer.put(iEvent);
0137     }
0138 
0139   private:
0140     template <typename T>
0141     void putEmpty(edm::Event& iEvent, T& regionsLayers) {
0142       auto seedingHitSetsProducerDummy = T_SeedingHitSets(&localRA_);
0143       auto intermediateHitDoubletsProducerDummy = T_IntermediateHitDoublets(regionsLayers.seedingLayerSetsHitsPtr());
0144       seedingHitSetsProducerDummy.putEmpty(iEvent);
0145       intermediateHitDoubletsProducerDummy.putEmpty(iEvent);
0146     }
0147 
0148     T_RegionLayers regionsLayers_;
0149   };
0150 
0151   /////
0152   class DoNothing {
0153   public:
0154     DoNothing(const SeedingLayerSetsHits*) {}
0155     DoNothing(edm::RunningAverage*) {}
0156 
0157     static void produces(edm::ProducesCollector){};
0158 
0159     void reserve(size_t) {}
0160 
0161     auto beginRegion(const TrackingRegion*, LayerHitMapCache* ptr) { return std::make_tuple(ptr, 0); }
0162 
0163     void fill(int, const HitDoublets&) {}
0164     void fill(int, const SeedingLayerSetsHits::SeedingLayerSet&, HitDoublets&&) {}
0165 
0166     void put(edm::Event&) {}
0167     void putEmpty(edm::Event&) {}
0168   };
0169 
0170   /////
0171   class ImplSeedingHitSets {
0172   public:
0173     ImplSeedingHitSets(edm::RunningAverage* localRA)
0174         : seedingHitSets_(std::make_unique<RegionsSeedingHitSets>()), localRA_(localRA) {}
0175 
0176     static void produces(edm::ProducesCollector producesCollector) {
0177       producesCollector.produces<RegionsSeedingHitSets>();
0178     }
0179 
0180     void reserve(size_t regionsSize) { seedingHitSets_->reserve(regionsSize, localRA_->upper()); }
0181 
0182     auto beginRegion(const TrackingRegion* region, LayerHitMapCache*) {
0183       hitCacheTmp_.clear();
0184       return std::make_tuple(&hitCacheTmp_, seedingHitSets_->beginRegion(region));
0185     }
0186 
0187     void fill(RegionsSeedingHitSets::RegionFiller& filler, const HitDoublets& doublets) {
0188       for (size_t i = 0, size = doublets.size(); i < size; ++i) {
0189         filler.emplace_back(doublets.hit(i, HitDoublets::inner), doublets.hit(i, HitDoublets::outer));
0190       }
0191     }
0192 
0193     void put(edm::Event& iEvent) {
0194       seedingHitSets_->shrink_to_fit();
0195       localRA_->update(seedingHitSets_->size());
0196       putEmpty(iEvent);
0197     }
0198 
0199     void putEmpty(edm::Event& iEvent) { iEvent.put(std::move(seedingHitSets_)); }
0200 
0201   private:
0202     std::unique_ptr<RegionsSeedingHitSets> seedingHitSets_;
0203     edm::RunningAverage* localRA_;
0204     LayerHitMapCache hitCacheTmp_;  // used if !produceIntermediateHitDoublets
0205   };
0206 
0207   /////
0208   class ImplIntermediateHitDoublets {
0209   public:
0210     ImplIntermediateHitDoublets(const SeedingLayerSetsHits* layers)
0211         : intermediateHitDoublets_(std::make_unique<IntermediateHitDoublets>(layers)), layers_(layers) {}
0212 
0213     static void produces(edm::ProducesCollector producesCollector) {
0214       producesCollector.produces<IntermediateHitDoublets>();
0215     }
0216 
0217     void reserve(size_t regionsSize) { intermediateHitDoublets_->reserve(regionsSize, layers_->size()); }
0218 
0219     auto beginRegion(const TrackingRegion* region, LayerHitMapCache*) {
0220       auto filler = intermediateHitDoublets_->beginRegion(region);
0221       return std::make_tuple(&(filler.layerHitMapCache()), filler);
0222     }
0223 
0224     void fill(IntermediateHitDoublets::RegionFiller& filler,
0225               const SeedingLayerSetsHits::SeedingLayerSet& layerSet,
0226               HitDoublets&& doublets) {
0227       filler.addDoublets(layerSet, std::move(doublets));
0228     }
0229 
0230     void put(edm::Event& iEvent) {
0231       intermediateHitDoublets_->shrink_to_fit();
0232       putEmpty(iEvent);
0233     }
0234 
0235     void putEmpty(edm::Event& iEvent) { iEvent.put(std::move(intermediateHitDoublets_)); }
0236 
0237   private:
0238     std::unique_ptr<IntermediateHitDoublets> intermediateHitDoublets_;
0239     const SeedingLayerSetsHits* layers_;
0240   };
0241 
0242   /////
0243   // For the usual case that TrackingRegions and seeding layers are read separately
0244   class RegionsLayersSeparate {
0245   public:
0246     class RegionLayers {
0247     public:
0248       RegionLayers(const TrackingRegion* region, const std::vector<SeedingLayerSetsHits::SeedingLayerSet>* layerPairs)
0249           : region_(region), layerPairs_(layerPairs) {}
0250 
0251       const TrackingRegion& region() const { return *region_; }
0252       const std::vector<SeedingLayerSetsHits::SeedingLayerSet>& layerPairs() const { return *layerPairs_; }
0253 
0254     private:
0255       const TrackingRegion* region_;
0256       const std::vector<SeedingLayerSetsHits::SeedingLayerSet>* layerPairs_;
0257     };
0258 
0259     class EventTmp {
0260     public:
0261       class const_iterator {
0262       public:
0263         using internal_iterator_type = edm::OwnVector<TrackingRegion>::const_iterator;
0264         using value_type = RegionLayers;
0265         using difference_type = internal_iterator_type::difference_type;
0266 
0267         const_iterator(internal_iterator_type iter,
0268                        const std::vector<SeedingLayerSetsHits::SeedingLayerSet>* layerPairs)
0269             : iter_(iter), layerPairs_(layerPairs) {}
0270 
0271         value_type operator*() const { return value_type(&(*iter_), layerPairs_); }
0272 
0273         const_iterator& operator++() {
0274           ++iter_;
0275           return *this;
0276         }
0277         const_iterator operator++(int) {
0278           const_iterator clone(*this);
0279           ++(*this);
0280           return clone;
0281         }
0282 
0283         bool operator==(const const_iterator& other) const { return iter_ == other.iter_; }
0284         bool operator!=(const const_iterator& other) const { return !operator==(other); }
0285 
0286       private:
0287         internal_iterator_type iter_;
0288         const std::vector<SeedingLayerSetsHits::SeedingLayerSet>* layerPairs_;
0289       };
0290 
0291       EventTmp(const SeedingLayerSetsHits* seedingLayerSetsHits,
0292                const edm::OwnVector<TrackingRegion>* regions,
0293                const std::vector<unsigned>& layerPairBegins)
0294           : seedingLayerSetsHits_(seedingLayerSetsHits), regions_(regions) {
0295         // construct the pairs from the sets
0296         if (seedingLayerSetsHits_->numberOfLayersInSet() > 2) {
0297           for (const auto& layerSet : *seedingLayerSetsHits_) {
0298             for (const auto pairBeginIndex : layerPairBegins) {
0299               if (pairBeginIndex + 1 >= seedingLayerSetsHits->numberOfLayersInSet()) {
0300                 throw cms::Exception("LogicError")
0301                     << "Layer pair index " << pairBeginIndex
0302                     << " is out of bounds, input SeedingLayerSetsHits has only "
0303                     << seedingLayerSetsHits->numberOfLayersInSet()
0304                     << " layers per set, and the index+1 must be < than the number of layers in set";
0305               }
0306 
0307               // Take only the requested pair of the set
0308               SeedingLayerSetsHits::SeedingLayerSet pairCandidate = layerSet.slice(pairBeginIndex, pairBeginIndex + 1);
0309 
0310               // it would be trivial to use 128-bit bitfield for O(1) check
0311               // if a layer pair has been inserted, but let's test first how
0312               // a "straightforward" solution works
0313               auto found = std::find_if(
0314                   layerPairs.begin(), layerPairs.end(), [&](const SeedingLayerSetsHits::SeedingLayerSet& pair) {
0315                     return pair[0].index() == pairCandidate[0].index() && pair[1].index() == pairCandidate[1].index();
0316                   });
0317               if (found != layerPairs.end())
0318                 continue;
0319 
0320               layerPairs.push_back(pairCandidate);
0321             }
0322           }
0323         } else {
0324           if (layerPairBegins.size() != 1) {
0325             throw cms::Exception("LogicError")
0326                 << "With pairs of input layers, it doesn't make sense to specify more than one input layer pair, got "
0327                 << layerPairBegins.size();
0328           }
0329           if (layerPairBegins[0] != 0) {
0330             throw cms::Exception("LogicError")
0331                 << "With pairs of input layers, it doesn't make sense to specify other input layer pair than 0; got "
0332                 << layerPairBegins[0];
0333           }
0334 
0335           layerPairs.reserve(seedingLayerSetsHits->size());
0336           for (const auto& set : *seedingLayerSetsHits)
0337             layerPairs.push_back(set);
0338         }
0339       }
0340 
0341       const SeedingLayerSetsHits* seedingLayerSetsHitsPtr() const { return seedingLayerSetsHits_; }
0342 
0343       size_t regionsSize() const { return regions_->size(); }
0344 
0345       const_iterator begin() const { return const_iterator(regions_->begin(), &layerPairs); }
0346       const_iterator cbegin() const { return begin(); }
0347       const_iterator end() const { return const_iterator(regions_->end(), &layerPairs); }
0348       const_iterator cend() const { return end(); }
0349 
0350     private:
0351       const SeedingLayerSetsHits* seedingLayerSetsHits_;
0352       const edm::OwnVector<TrackingRegion>* regions_;
0353       std::vector<SeedingLayerSetsHits::SeedingLayerSet> layerPairs;
0354     };
0355 
0356     RegionsLayersSeparate(const std::vector<unsigned>* layerPairBegins,
0357                           const edm::InputTag& seedingLayerTag,
0358                           const edm::InputTag& regionTag,
0359                           edm::ConsumesCollector iC)
0360         : layerPairBegins_(layerPairBegins),
0361           seedingLayerToken_(iC.consumes<SeedingLayerSetsHits>(seedingLayerTag)),
0362           regionToken_(iC.consumes<edm::OwnVector<TrackingRegion>>(regionTag)) {}
0363 
0364     EventTmp beginEvent(const edm::Event& iEvent) const {
0365       edm::Handle<SeedingLayerSetsHits> hlayers;
0366       iEvent.getByToken(seedingLayerToken_, hlayers);
0367       const auto* layers = hlayers.product();
0368       if (layers->numberOfLayersInSet() < 2)
0369         throw cms::Exception("LogicError")
0370             << "HitPairEDProducer expects SeedingLayerSetsHits::numberOfLayersInSet() to be >= 2, got "
0371             << layers->numberOfLayersInSet()
0372             << ". This is likely caused by a configuration error of this module, or SeedingLayersEDProducer.";
0373       edm::Handle<edm::OwnVector<TrackingRegion>> hregions;
0374       iEvent.getByToken(regionToken_, hregions);
0375 
0376       return EventTmp(layers, hregions.product(), *layerPairBegins_);
0377     }
0378 
0379   private:
0380     const std::vector<unsigned>* layerPairBegins_;
0381     edm::EDGetTokenT<SeedingLayerSetsHits> seedingLayerToken_;
0382     edm::EDGetTokenT<edm::OwnVector<TrackingRegion>> regionToken_;
0383   };
0384 
0385   /////
0386   // For the case of automated pixel inactive region recovery where
0387   // TrackingRegions and seeding layers become together
0388   class RegionsLayersTogether {
0389   public:
0390     class EventTmp {
0391     public:
0392       using const_iterator = TrackingRegionsSeedingLayerSets::const_iterator;
0393 
0394       explicit EventTmp(const TrackingRegionsSeedingLayerSets* regionsLayers) : regionsLayers_(regionsLayers) {
0395         if (regionsLayers->seedingLayerSetsHits().numberOfLayersInSet() != 2) {
0396           throw cms::Exception("LogicError")
0397               << "With trackingRegionsSeedingLayers input, the seeding layer sets may only contain layer pairs, now "
0398                  "got sets with "
0399               << regionsLayers->seedingLayerSetsHits().numberOfLayersInSet() << " layers";
0400         }
0401       }
0402 
0403       const SeedingLayerSetsHits* seedingLayerSetsHitsPtr() const { return &(regionsLayers_->seedingLayerSetsHits()); }
0404 
0405       size_t regionsSize() const { return regionsLayers_->regionsSize(); }
0406 
0407       const_iterator begin() const { return regionsLayers_->begin(); }
0408       const_iterator cbegin() const { return begin(); }
0409       const_iterator end() const { return regionsLayers_->end(); }
0410       const_iterator cend() const { return end(); }
0411 
0412     private:
0413       const TrackingRegionsSeedingLayerSets* regionsLayers_;
0414     };
0415 
0416     RegionsLayersTogether(const std::vector<unsigned>* layerPairBegins,
0417                           const edm::InputTag& regionLayerTag,
0418                           edm::ConsumesCollector iC)
0419         : regionLayerToken_(iC.consumes<TrackingRegionsSeedingLayerSets>(regionLayerTag)) {
0420       if (layerPairBegins->size() != 1) {
0421         throw cms::Exception("LogicError") << "With trackingRegionsSeedingLayers mode, it doesn't make sense to "
0422                                               "specify more than one input layer pair, got "
0423                                            << layerPairBegins->size();
0424       }
0425       if ((*layerPairBegins)[0] != 0) {
0426         throw cms::Exception("LogicError") << "With trackingRegionsSeedingLayers mode, it doesn't make sense to "
0427                                               "specify other input layer pair than 0; got "
0428                                            << (*layerPairBegins)[0];
0429       }
0430     }
0431 
0432     EventTmp beginEvent(const edm::Event& iEvent) const {
0433       edm::Handle<TrackingRegionsSeedingLayerSets> hregions;
0434       iEvent.getByToken(regionLayerToken_, hregions);
0435       return EventTmp(hregions.product());
0436     }
0437 
0438   private:
0439     edm::EDGetTokenT<TrackingRegionsSeedingLayerSets> regionLayerToken_;
0440   };
0441 }  // namespace
0442 
0443 HitPairEDProducer::HitPairEDProducer(const edm::ParameterSet& iConfig) {
0444   auto layersTag = iConfig.getParameter<edm::InputTag>("seedingLayers");
0445   auto regionTag = iConfig.getParameter<edm::InputTag>("trackingRegions");
0446   auto regionLayerTag = iConfig.getParameter<edm::InputTag>("trackingRegionsSeedingLayers");
0447   const bool useRegionLayers = !regionLayerTag.label().empty();
0448   if (useRegionLayers) {
0449     if (!regionTag.label().empty()) {
0450       throw cms::Exception("Configuration")
0451           << "HitPairEDProducer requires either trackingRegions or trackingRegionsSeedingLayers to be set, now both "
0452              "are set to non-empty value. Set the unneeded parameter to empty value.";
0453     }
0454     if (!layersTag.label().empty()) {
0455       throw cms::Exception("Configuration")
0456           << "With non-empty trackingRegionsSeedingLayers, please set also seedingLayers to empty value to reduce "
0457              "confusion, because the parameter is not used";
0458     }
0459   }
0460   if (regionTag.label().empty() && regionLayerTag.label().empty()) {
0461     throw cms::Exception("Configuration")
0462         << "HitPairEDProducer requires either trackingRegions or trackingRegionsSeedingLayers to be set, now both are "
0463            "set to empty value. Set the needed parameter to a non-empty value.";
0464   }
0465 
0466   const bool produceSeedingHitSets = iConfig.getParameter<bool>("produceSeedingHitSets");
0467   const bool produceIntermediateHitDoublets = iConfig.getParameter<bool>("produceIntermediateHitDoublets");
0468 
0469   if (produceSeedingHitSets && produceIntermediateHitDoublets) {
0470     if (useRegionLayers)
0471       throw cms::Exception("Configuration")
0472           << "Mode 'trackingRegionsSeedingLayers' makes sense only with 'produceSeedingHitsSets', now also "
0473              "'produceIntermediateHitDoublets is active";
0474     impl_ = std::make_unique<::Impl<::ImplSeedingHitSets, ::ImplIntermediateHitDoublets, ::RegionsLayersSeparate>>(
0475         iConfig, consumesCollector(), layersTag, regionTag);
0476   } else if (produceSeedingHitSets) {
0477     if (useRegionLayers) {
0478       impl_ = std::make_unique<::Impl<::ImplSeedingHitSets, ::DoNothing, ::RegionsLayersTogether>>(
0479           iConfig, consumesCollector(), regionLayerTag);
0480     } else {
0481       impl_ = std::make_unique<::Impl<::ImplSeedingHitSets, ::DoNothing, ::RegionsLayersSeparate>>(
0482           iConfig, consumesCollector(), layersTag, regionTag);
0483     }
0484   } else if (produceIntermediateHitDoublets) {
0485     if (useRegionLayers)
0486       throw cms::Exception("Configuration")
0487           << "Mode 'trackingRegionsSeedingLayers' makes sense only with 'produceSeedingHitsSets', now "
0488              "'produceIntermediateHitDoublets is active instead";
0489     impl_ = std::make_unique<::Impl<::DoNothing, ::ImplIntermediateHitDoublets, ::RegionsLayersSeparate>>(
0490         iConfig, consumesCollector(), layersTag, regionTag);
0491   } else
0492     throw cms::Exception("Configuration")
0493         << "HitPairEDProducer requires either produceIntermediateHitDoublets or produceSeedingHitSets to be True. If "
0494            "neither are needed, just remove this module from your sequence/path as it doesn't do anything useful";
0495 
0496   auto clusterCheckTag = iConfig.getParameter<edm::InputTag>("clusterCheck");
0497   if (!clusterCheckTag.label().empty())
0498     clusterCheckToken_ = consumes<bool>(clusterCheckTag);
0499 
0500   impl_->produces(producesCollector());
0501 }
0502 
0503 void HitPairEDProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0504   edm::ParameterSetDescription desc;
0505 
0506   desc.add<edm::InputTag>("seedingLayers", edm::InputTag("seedingLayersEDProducer"))
0507       ->setComment("Set this empty if 'trackingRegionsSeedingLayers' is non-empty");
0508   desc.add<edm::InputTag>("trackingRegions", edm::InputTag("globalTrackingRegionFromBeamSpot"))
0509       ->setComment(
0510           "Input tracking regions when using all layer sets in 'seedingLayers' (conflicts with "
0511           "'trackingRegionsSeedingLayers', set this empty to use the other)");
0512   desc.add<edm::InputTag>("trackingRegionsSeedingLayers", edm::InputTag(""))
0513       ->setComment(
0514           "Input tracking regions and corresponding layer sets in case of dynamically limiting the seeding layers "
0515           "(conflicts with 'trackingRegions', set this empty to use the other; if using this set also 'seedingLayers' "
0516           "to empty)");
0517   desc.add<edm::InputTag>("clusterCheck", edm::InputTag("trackerClusterCheck"));
0518   desc.add<bool>("produceSeedingHitSets", false);
0519   desc.add<bool>("produceIntermediateHitDoublets", false);
0520   desc.add<unsigned int>("maxElement", 1000000);
0521   desc.add<unsigned int>("maxElementTotal", 50000000);
0522   desc.add<bool>("putEmptyIfMaxElementReached", false)
0523       ->setComment(
0524           "If set to true (default is 'false'), abort processing and put empty data products also if any layer pair "
0525           "yields at least maxElement doublets, in addition to aborting processing if the sum of doublets from all "
0526           "layer pairs reaches maxElementTotal.");
0527   desc.add<std::vector<unsigned>>("layerPairs", std::vector<unsigned>{0})
0528       ->setComment("Indices to the pairs of consecutive layers, i.e. 0 means (0,1), 1 (1,2) etc.");
0529 
0530   descriptions.add("hitPairEDProducerDefault", desc);
0531 }
0532 
0533 void HitPairEDProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0534   bool clusterCheckOk = true;
0535   if (!clusterCheckToken_.isUninitialized()) {
0536     edm::Handle<bool> hclusterCheck;
0537     iEvent.getByToken(clusterCheckToken_, hclusterCheck);
0538     clusterCheckOk = *hclusterCheck;
0539   }
0540 
0541   impl_->produce(clusterCheckOk, iEvent, iSetup);
0542 }
0543 
0544 #include "FWCore/PluginManager/interface/ModuleDef.h"
0545 #include "FWCore/Framework/interface/MakerMacros.h"
0546 DEFINE_FWK_MODULE(HitPairEDProducer);