CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
VertexProducer.cc
Go to the documentation of this file.
14 
15 using namespace l1tVertexFinder;
16 using namespace std;
17 
19  : l1TracksToken_(consumes<TTTrackRefCollectionType>(iConfig.getParameter<edm::InputTag>("l1TracksInputTag"))),
21  outputCollectionName_(iConfig.getParameter<std::string>("l1VertexCollectionName")),
22  settings_(AlgoSettings(iConfig)) {
23  // Get configuration parameters
24 
25  switch (settings_.vx_algo()) {
26  case Algorithm::fastHisto:
27  edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using the fastHisto binning algorithm";
28  break;
29  case Algorithm::fastHistoEmulation:
30  edm::LogInfo("VertexProducer")
31  << "VertexProducer::Finding vertices using the emulation version of the fastHisto binning algorithm";
32  break;
33  case Algorithm::fastHistoLooseAssociation:
34  edm::LogInfo("VertexProducer")
35  << "VertexProducer::Finding vertices using the fastHistoLooseAssociation binning algorithm";
36  break;
37  case Algorithm::GapClustering:
38  edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a gap clustering algorithm";
39  break;
40  case Algorithm::agglomerativeHierarchical:
41  edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a Simple Merge Clustering algorithm";
42  break;
43  case Algorithm::DBSCAN:
44  edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a DBSCAN algorithm";
45  break;
46  case Algorithm::PVR:
47  edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a PVR algorithm";
48  break;
49  case Algorithm::adaptiveVertexReconstruction:
50  edm::LogInfo("VertexProducer")
51  << "VertexProducer::Finding vertices using an AdaptiveVertexReconstruction algorithm";
52  break;
53  case Algorithm::HPV:
54  edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a Highest Pt Vertex algorithm";
55  break;
56  case Algorithm::Kmeans:
57  edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a kmeans algorithm";
58  break;
59  }
60 
61  //--- Define EDM output to be written to file (if required)
62  if (settings_.vx_algo() == Algorithm::fastHistoEmulation) {
63  produces<l1t::VertexWordCollection>(outputCollectionName_ + "Emulation");
64  } else {
65  produces<l1t::VertexCollection>(outputCollectionName_);
66  }
67 }
68 
71  iEvent.getByToken(l1TracksToken_, l1TracksHandle);
72 
73  std::vector<l1tVertexFinder::L1Track> l1Tracks;
74  l1Tracks.reserve(l1TracksHandle->size());
75  if (settings_.debug() > 1) {
76  edm::LogInfo("VertexProducer") << "produce::Processing " << l1TracksHandle->size() << " tracks";
77  }
78  for (const auto& track : *l1TracksHandle) {
79  auto l1track = L1Track(edm::refToPtr(track));
80  if (l1track.pt() >= settings_.vx_TrackMinPt()) {
81  l1Tracks.push_back(l1track);
82  } else {
83  if (settings_.debug() > 2) {
84  edm::LogInfo("VertexProducer") << "produce::Removing track with too low of a pt (" << l1track.pt() << ")\n"
85  << " word = " << l1track.getTTTrackPtr()->getTrackWord().to_string(2);
86  }
87  }
88  }
89 
90  if (settings_.debug() > 1) {
91  edm::LogInfo("VertexProducer") << "produce::Processing " << l1Tracks.size() << " tracks after minimum pt cut of"
92  << settings_.vx_TrackMinPt() << " GeV";
93  }
94 
96 
97  switch (settings_.vx_algo()) {
98  case Algorithm::fastHisto: {
99  const TrackerTopology& tTopo = iSetup.getData(tTopoToken);
100  vf.fastHisto(&tTopo);
101  break;
102  }
103  case Algorithm::fastHistoEmulation:
104  vf.fastHistoEmulation();
105  break;
106  case Algorithm::fastHistoLooseAssociation:
108  break;
109  case Algorithm::GapClustering:
110  vf.GapClustering();
111  break;
112  case Algorithm::agglomerativeHierarchical:
114  break;
115  case Algorithm::DBSCAN:
116  vf.DBSCAN();
117  break;
118  case Algorithm::PVR:
119  vf.PVR();
120  break;
121  case Algorithm::adaptiveVertexReconstruction:
123  break;
124  case Algorithm::HPV:
125  vf.HPV();
126  break;
127  case Algorithm::Kmeans:
128  vf.Kmeans();
129  break;
130  }
131 
132  vf.sortVerticesInPt();
133  vf.findPrimaryVertex();
134 
135  // //=== Store output EDM track and hardware stub collections.
136  if (settings_.vx_algo() == Algorithm::fastHistoEmulation) {
137  std::unique_ptr<l1t::VertexWordCollection> product_emulation =
138  std::make_unique<l1t::VertexWordCollection>(vf.verticesEmulation().begin(), vf.verticesEmulation().end());
139  iEvent.put(std::move(product_emulation), outputCollectionName_ + "Emulation");
140  } else {
141  std::unique_ptr<l1t::VertexCollection> product(new std::vector<l1t::Vertex>());
142  for (const auto& vtx : vf.vertices()) {
143  product->emplace_back(vtx.vertex());
144  }
145  iEvent.put(std::move(product), outputCollectionName_);
146  }
147 }
ESGetTokenH3DDVariant esConsumes(std::string const &Record, edm::ConsumesCollector &)
Definition: DeDxTools.cc:283
void fastHistoLooseAssociation()
High pT Vertex Algorithm.
Ptr< typename C::value_type > refToPtr(Ref< C, typename C::value_type, refhelper::FindUsingAdvance< C, typename C::value_type > > const &ref)
Definition: RefToPtr.h:18
T const & getData(const ESGetToken< T, R > &iToken) const noexcept(false)
Definition: EventSetup.h:119
void produce(edm::StreamID, edm::Event &, const edm::EventSetup &) const override
void GapClustering()
Gap Clustering Algorithm.
Definition: VertexFinder.cc:65
const edm::ESGetToken< TrackerTopology, TrackerTopologyRcd > tTopoToken
void agglomerativeHierarchicalClustering()
Simple Merge Algorithm.
const edm::EDGetTokenT< TTTrackRefCollectionType > l1TracksToken_
const std::string outputCollectionName_
unsigned int debug() const
Definition: AlgoSettings.h:77
void adaptiveVertexReconstruction()
Adaptive Vertex Reconstruction algorithm.
int iEvent
Definition: GenABIO.cc:224
void findPrimaryVertex()
Find the primary vertex.
void fastHisto(const TrackerTopology *tTopo)
Histogramming algorithm.
const std::vector< RecoVertex<> > & vertices() const
Returns the z positions of the reconstructed primary vertices.
Definition: VertexFinder.h:75
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
void sortVerticesInPt()
Sort vertices in pT.
VertexProducer(const edm::ParameterSet &)
Log< level::Info, false > LogInfo
void PVR()
Find maximum distance in two clusters of tracks.
size_type size() const
Size of the RefVector.
Definition: RefVector.h:102
void fastHistoEmulation()
Histogramming algorithm (emulation)
void HPV()
High pT Vertex Algorithm.
TTTrack< Ref_Phase2TrackerDigi_ > L1Track
Definition: L1Track.h:8
HLT enums.
void Kmeans()
Kmeans Algorithm.
const l1t::VertexWordCollection & verticesEmulation() const
Returns the emulation primary vertices.
Definition: VertexFinder.h:77
Algorithm vx_algo() const
Definition: AlgoSettings.h:34
void DBSCAN()
DBSCAN algorithm.
def move(src, dest)
Definition: eostools.py:511
l1tVertexFinder::AlgoSettings settings_