File indexing completed on 2021-10-01 22:41:03
0001 #include "FWCore/Framework/interface/stream/EDProducer.h"
0002 #include "FWCore/Framework/interface/ConsumesCollector.h"
0003 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0004 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0006 #include "FWCore/Utilities/interface/EDGetToken.h"
0007 #include "FWCore/Framework/interface/Event.h"
0008 #include "DataFormats/Common/interface/Handle.h"
0009 #include "FWCore/Utilities/interface/RunningAverage.h"
0010
0011 #include "RecoTracker/TkHitPairs/interface/IntermediateHitDoublets.h"
0012 #include "RecoTracker/TkHitPairs/interface/RegionsSeedingHitSets.h"
0013 #include "DataFormats/Common/interface/OwnVector.h"
0014 #include "RecoPixelVertexing/PixelTriplets/interface/LayerTriplets.h"
0015 #include "MultiHitGeneratorFromChi2.h"
0016
0017 class MultiHitFromChi2EDProducer : public edm::stream::EDProducer<> {
0018 public:
0019 MultiHitFromChi2EDProducer(const edm::ParameterSet& iConfig);
0020 ~MultiHitFromChi2EDProducer() override = default;
0021
0022 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0023
0024 void produce(edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0025
0026 private:
0027 edm::EDGetTokenT<IntermediateHitDoublets> doubletToken_;
0028
0029 edm::RunningAverage localRA_;
0030
0031 MultiHitGeneratorFromChi2 generator_;
0032 };
0033
0034 MultiHitFromChi2EDProducer::MultiHitFromChi2EDProducer(const edm::ParameterSet& iConfig)
0035 : doubletToken_(consumes<IntermediateHitDoublets>(iConfig.getParameter<edm::InputTag>("doublets"))),
0036 generator_(iConfig, consumesCollector()) {
0037 produces<RegionsSeedingHitSets>();
0038 produces<edm::OwnVector<BaseTrackerRecHit> >();
0039 }
0040
0041 void MultiHitFromChi2EDProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0042 edm::ParameterSetDescription desc;
0043
0044 desc.add<edm::InputTag>("doublets", edm::InputTag("hitPairEDProducer"));
0045
0046 MultiHitGeneratorFromChi2::fillDescriptions(desc);
0047
0048 auto label = MultiHitGeneratorFromChi2::fillDescriptionsLabel() + std::string("EDProducerDefault");
0049 descriptions.add(label, desc);
0050 }
0051
0052 void MultiHitFromChi2EDProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0053 edm::Handle<IntermediateHitDoublets> hdoublets;
0054 iEvent.getByToken(doubletToken_, hdoublets);
0055 const auto& regionDoublets = *hdoublets;
0056
0057 const SeedingLayerSetsHits& seedingLayerHits = regionDoublets.seedingLayerHits();
0058 if (seedingLayerHits.numberOfLayersInSet() < 3) {
0059 throw cms::Exception("LogicError")
0060 << "MultiHitFromChi2EDProducer expects SeedingLayerSetsHits::numberOfLayersInSet() to be >= 3, got "
0061 << seedingLayerHits.numberOfLayersInSet()
0062 << ". This is likely caused by a configuration error of this module, HitPairEDProducer, or "
0063 "SeedingLayersEDProducer.";
0064 }
0065
0066 auto seedingHitSets = std::make_unique<RegionsSeedingHitSets>();
0067 if (regionDoublets.empty()) {
0068 iEvent.put(std::move(seedingHitSets));
0069 return;
0070 }
0071 seedingHitSets->reserve(regionDoublets.regionSize(), localRA_.upper());
0072 generator_.initES(iSetup);
0073
0074
0075 std::vector<LayerTriplets::LayerSetAndLayers> trilayers = LayerTriplets::layers(seedingLayerHits);
0076
0077 OrderedMultiHits multihits;
0078 multihits.reserve(localRA_.upper());
0079 std::vector<std::unique_ptr<BaseTrackerRecHit> > refittedHitStorage;
0080 refittedHitStorage.reserve(localRA_.upper() * 2);
0081
0082 LogDebug("MultiHitFromChi2EDProducer") << "Creating multihits for " << regionDoublets.regionSize() << " regions, and "
0083 << trilayers.size() << " pair+3rd layers from "
0084 << regionDoublets.layerPairsSize() << " layer pairs";
0085
0086 LayerHitMapCache hitCache;
0087 for (const auto& regionLayerPairs : regionDoublets) {
0088 const TrackingRegion& region = regionLayerPairs.region();
0089
0090 auto seedingHitSetsFiller = RegionsSeedingHitSets::dummyFiller();
0091 seedingHitSetsFiller = seedingHitSets->beginRegion(®ion);
0092
0093 hitCache.clear();
0094 hitCache.extend(regionLayerPairs.layerHitMapCache());
0095
0096 LogTrace("MultiHitFromChi2EDProducer") << " starting region";
0097
0098 for (const auto& layerPair : regionLayerPairs) {
0099 LogTrace("MultiHitFromChi2EDProducer")
0100 << " starting layer pair " << layerPair.innerLayerIndex() << "," << layerPair.outerLayerIndex();
0101
0102 auto found = std::find_if(trilayers.begin(), trilayers.end(), [&](const LayerTriplets::LayerSetAndLayers& a) {
0103 return a.first[0].index() == layerPair.innerLayerIndex() && a.first[1].index() == layerPair.outerLayerIndex();
0104 });
0105 if (found == trilayers.end()) {
0106 auto exp = cms::Exception("LogicError") << "Did not find the layer pair from vector<pair+third layers>. This "
0107 "is a sign of some internal inconsistency\n";
0108 exp << "I was looking for layer pair " << layerPair.innerLayerIndex() << "," << layerPair.outerLayerIndex()
0109 << ". Triplets have the following pairs:\n";
0110 for (const auto& a : trilayers) {
0111 exp << " " << a.first[0].index() << "," << a.first[1].index() << ": 3rd layers";
0112 for (const auto& b : a.second) {
0113 exp << " " << b.index();
0114 }
0115 exp << "\n";
0116 }
0117 throw exp;
0118 }
0119 const auto& thirdLayers = found->second;
0120
0121 generator_.hitSets(region, multihits, layerPair.doublets(), thirdLayers, hitCache, refittedHitStorage);
0122
0123 #ifdef EDM_ML_DEBUG
0124 LogTrace("MultiHitFromChi2EDProducer")
0125 << " created " << multihits.size() << " multihits for layer pair " << layerPair.innerLayerIndex() << ","
0126 << layerPair.outerLayerIndex() << " and 3rd layers";
0127 for (const auto& l : thirdLayers) {
0128 LogTrace("MultiHitFromChi2EDProducer") << " " << l.index();
0129 }
0130 #endif
0131
0132 for (const SeedingHitSet& hitSet : multihits) {
0133 seedingHitSetsFiller.emplace_back(hitSet);
0134 }
0135 multihits.clear();
0136 }
0137 }
0138 localRA_.update(seedingHitSets->size());
0139
0140 auto storage = std::make_unique<edm::OwnVector<BaseTrackerRecHit> >();
0141 storage->reserve(refittedHitStorage.size());
0142 for (auto& ptr : refittedHitStorage)
0143 storage->push_back(ptr.release());
0144
0145 seedingHitSets->shrink_to_fit();
0146 storage->shrink_to_fit();
0147 iEvent.put(std::move(seedingHitSets));
0148 iEvent.put(std::move(storage));
0149 }
0150
0151 #include "FWCore/PluginManager/interface/ModuleDef.h"
0152 #include "FWCore/Framework/interface/MakerMacros.h"
0153 DEFINE_FWK_MODULE(MultiHitFromChi2EDProducer);