File indexing completed on 2024-04-06 12:22:12
0001 #include "L1Trigger/VertexFinder/interface/VertexProducer.h"
0002 #include "DataFormats/L1TrackTrigger/interface/TTTypes.h"
0003 #include "DataFormats/L1Trigger/interface/Vertex.h"
0004 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "FWCore/Framework/interface/EventSetup.h"
0007 #include "FWCore/Framework/interface/MakerMacros.h"
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0010 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0011 #include "L1Trigger/VertexFinder/interface/AlgoSettings.h"
0012 #include "L1Trigger/VertexFinder/interface/RecoVertex.h"
0013 #include "L1Trigger/VertexFinder/interface/VertexFinder.h"
0014
0015 using namespace l1tVertexFinder;
0016 using namespace std;
0017
0018 VertexProducer::VertexProducer(const edm::ParameterSet& iConfig)
0019 : l1TracksToken_(consumes<TTTrackRefCollectionType>(iConfig.getParameter<edm::InputTag>("l1TracksInputTag"))),
0020 tTopoToken(esConsumes<TrackerTopology, TrackerTopologyRcd>()),
0021 outputCollectionName_(iConfig.getParameter<std::string>("l1VertexCollectionName")),
0022 settings_(AlgoSettings(iConfig)) {
0023
0024
0025 switch (settings_.vx_algo()) {
0026 case Algorithm::fastHisto:
0027 edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using the fastHisto binning algorithm";
0028 break;
0029 case Algorithm::fastHistoEmulation:
0030 edm::LogInfo("VertexProducer")
0031 << "VertexProducer::Finding vertices using the emulation version of the fastHisto binning algorithm";
0032 break;
0033 case Algorithm::fastHistoLooseAssociation:
0034 edm::LogInfo("VertexProducer")
0035 << "VertexProducer::Finding vertices using the fastHistoLooseAssociation binning algorithm";
0036 break;
0037 case Algorithm::GapClustering:
0038 edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a gap clustering algorithm";
0039 break;
0040 case Algorithm::agglomerativeHierarchical:
0041 edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a Simple Merge Clustering algorithm";
0042 break;
0043 case Algorithm::DBSCAN:
0044 edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a DBSCAN algorithm";
0045 break;
0046 case Algorithm::PVR:
0047 edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a PVR algorithm";
0048 break;
0049 case Algorithm::adaptiveVertexReconstruction:
0050 edm::LogInfo("VertexProducer")
0051 << "VertexProducer::Finding vertices using an AdaptiveVertexReconstruction algorithm";
0052 break;
0053 case Algorithm::HPV:
0054 edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a Highest Pt Vertex algorithm";
0055 break;
0056 case Algorithm::Kmeans:
0057 edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a kmeans algorithm";
0058 break;
0059 case Algorithm::NNEmulation:
0060 edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using the Neural Network Emulation";
0061 break;
0062 }
0063
0064
0065 if (settings_.vx_algo() == Algorithm::fastHistoEmulation || settings_.vx_algo() == Algorithm::NNEmulation) {
0066 produces<l1t::VertexWordCollection>(outputCollectionName_ + "Emulation");
0067 } else {
0068 produces<l1t::VertexCollection>(outputCollectionName_);
0069 }
0070
0071 if (settings_.vx_algo() == Algorithm::NNEmulation) {
0072
0073 if (settings_.debug() > 1) {
0074 edm::LogInfo("VertexProducer") << "loading TrkWeight graph from " << settings_.vx_trkw_graph() << std::endl;
0075 edm::LogInfo("VertexProducer") << "loading PatternRec graph from " << settings_.vx_pattrec_graph() << std::endl;
0076 }
0077 TrkWGraph_ = tensorflow::loadGraphDef(settings_.vx_trkw_graph());
0078 TrkWSesh_ = tensorflow::createSession(TrkWGraph_);
0079 PattRecGraph_ = tensorflow::loadGraphDef(settings_.vx_pattrec_graph());
0080 PattRecSesh_ = tensorflow::createSession(PattRecGraph_);
0081 }
0082 }
0083
0084 void VertexProducer::produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0085 edm::Handle<TTTrackRefCollectionType> l1TracksHandle;
0086 iEvent.getByToken(l1TracksToken_, l1TracksHandle);
0087
0088 std::vector<l1tVertexFinder::L1Track> l1Tracks;
0089 l1Tracks.reserve(l1TracksHandle->size());
0090 if (settings_.debug() > 1) {
0091 edm::LogInfo("VertexProducer") << "produce::Processing " << l1TracksHandle->size() << " tracks";
0092 }
0093 for (const auto& track : *l1TracksHandle) {
0094 auto l1track = L1Track(edm::refToPtr(track));
0095 if (l1track.pt() >= settings_.vx_TrackMinPt()) {
0096 l1Tracks.push_back(l1track);
0097 } else {
0098 if (settings_.debug() > 2) {
0099 edm::LogInfo("VertexProducer") << "produce::Removing track with too low of a pt (" << l1track.pt() << ")\n"
0100 << " word = " << l1track.getTTTrackPtr()->getTrackWord().to_string(2);
0101 }
0102 }
0103 }
0104
0105 if (settings_.debug() > 1) {
0106 edm::LogInfo("VertexProducer") << "produce::Processing " << l1Tracks.size() << " tracks after minimum pt cut of"
0107 << settings_.vx_TrackMinPt() << " GeV";
0108 }
0109
0110 VertexFinder vf(l1Tracks, settings_);
0111
0112 switch (settings_.vx_algo()) {
0113 case Algorithm::fastHisto: {
0114 const TrackerTopology& tTopo = iSetup.getData(tTopoToken);
0115 vf.fastHisto(&tTopo);
0116 break;
0117 }
0118 case Algorithm::fastHistoEmulation:
0119 vf.fastHistoEmulation();
0120 break;
0121 case Algorithm::fastHistoLooseAssociation:
0122 vf.fastHistoLooseAssociation();
0123 break;
0124 case Algorithm::GapClustering:
0125 vf.GapClustering();
0126 break;
0127 case Algorithm::agglomerativeHierarchical:
0128 vf.agglomerativeHierarchicalClustering();
0129 break;
0130 case Algorithm::DBSCAN:
0131 vf.DBSCAN();
0132 break;
0133 case Algorithm::PVR:
0134 vf.PVR();
0135 break;
0136 case Algorithm::adaptiveVertexReconstruction:
0137 vf.adaptiveVertexReconstruction();
0138 break;
0139 case Algorithm::HPV:
0140 vf.HPV();
0141 break;
0142 case Algorithm::Kmeans:
0143 vf.Kmeans();
0144 break;
0145 case Algorithm::NNEmulation:
0146 vf.NNVtxEmulation(TrkWSesh_, PattRecSesh_);
0147 break;
0148 }
0149
0150 vf.sortVerticesInPt();
0151 vf.findPrimaryVertex();
0152
0153
0154 if (settings_.vx_algo() == Algorithm::fastHistoEmulation || settings_.vx_algo() == Algorithm::NNEmulation) {
0155 std::unique_ptr<l1t::VertexWordCollection> product_emulation =
0156 std::make_unique<l1t::VertexWordCollection>(vf.verticesEmulation().begin(), vf.verticesEmulation().end());
0157 iEvent.put(std::move(product_emulation), outputCollectionName_ + "Emulation");
0158 } else {
0159 std::unique_ptr<l1t::VertexCollection> product(new std::vector<l1t::Vertex>());
0160 for (const auto& vtx : vf.vertices()) {
0161 product->emplace_back(vtx.vertex());
0162 }
0163 iEvent.put(std::move(product), outputCollectionName_);
0164 }
0165 }
0166 DEFINE_FWK_MODULE(VertexProducer);