1- #include < list>
21#include < vector>
32#include < limits>
43#include < string>
@@ -46,6 +45,7 @@ class DisplacedRegionSeedingVertexProducer : public edm::global::EDProducer<> {
4645 const double nearThreshold_;
4746 const double farThreshold_;
4847 const double discriminatorCut_;
48+ const unsigned int maxPseudoROIs_;
4949 const vector<string> input_names_;
5050 const vector<string> output_names_;
5151
@@ -63,6 +63,7 @@ DisplacedRegionSeedingVertexProducer::DisplacedRegionSeedingVertexProducer(const
6363 nearThreshold_(cfg.getParameter<double >(" nearThreshold" )),
6464 farThreshold_(cfg.getParameter<double >(" farThreshold" )),
6565 discriminatorCut_(cfg.getParameter<double >(" discriminatorCut" )),
66+ maxPseudoROIs_(cfg.getParameter<unsigned int >(" maxPseudoROIs" )),
6667 input_names_(cfg.getParameter<vector<string> >(" input_names" )),
6768 output_names_(cfg.getParameter<vector<string> >(" output_names" )),
6869 beamSpotToken_(consumes<reco::BeamSpot>(cfg.getParameter<edm::InputTag>(" beamSpot" ))),
@@ -93,16 +94,24 @@ void DisplacedRegionSeedingVertexProducer::produce(edm::StreamID streamID,
9394 const auto &trackClusters = event.get (trackClustersToken_);
9495
9596 // Initialize distances.
96- list<DisplacedVertexCluster> pseudoROIs;
97- list<Distance> distances;
97+ vector<DisplacedVertexCluster> pseudoROIs;
98+ pseudoROIs.reserve (std::min (maxPseudoROIs_, trackClusters.size ()));
99+ vector<Distance> distances;
98100 const double minTrackClusterRadius = minRadius_ - rParam_;
99101 for (unsigned i = 0 ; i < trackClusters.size (); i++) {
100102 const reco::VertexCompositeCandidate &trackCluster = trackClusters[i];
101103 const math::XYZVector x (trackCluster.vertex ());
102104 if (minRadius_ < 0.0 || minTrackClusterRadius < 0.0 || (x - bs).rho () > minTrackClusterRadius)
103105 pseudoROIs.emplace_back (&trackClusters.at (i), rParam_);
106+ if (pseudoROIs.size () == maxPseudoROIs_) {
107+ edm::LogWarning (" DisplacedRegionSeedingVertexProducer" )
108+ << " Truncated list of pseudoROIs at " << maxPseudoROIs_ << " out of " << trackClusters.size ()
109+ << " possible track clusters." ;
110+ break ;
111+ }
104112 }
105113 if (pseudoROIs.size () > 1 ) {
114+ distances.reserve (pseudoROIs.size () * std::max (1.0 , pseudoROIs.size () * 0.05 ));
106115 DisplacedVertexClusterItr secondToLast = pseudoROIs.end ();
107116 secondToLast--;
108117 for (DisplacedVertexClusterItr i = pseudoROIs.begin (); i != secondToLast; i++) {
@@ -119,11 +128,25 @@ void DisplacedRegionSeedingVertexProducer::produce(edm::StreamID streamID,
119128 }
120129 }
121130
131+ auto itBegin = distances.begin ();
132+ auto itLast = distances.end ();
133+
122134 // Do clustering.
123- while (!distances.empty ()) {
124- const auto comp = [](const Distance &a, const Distance &b) { return a.distance2 () <= b.distance2 (); };
125- distances.sort (comp);
126- DistanceItr dBest = distances.begin ();
135+ while (itBegin != itLast) {
136+ // find the lowest distance. Lots of repeatitive calculations done here
137+ // as from loop iteration to loop iteration only sqrt(distances.size()) distances
138+ // need to be recomputed (those involving best_i
139+ // but this is much better than sorting distances..
140+ DistanceItr dBest = itBegin;
141+ double distanceBest = dBest->distance2 ();
142+
143+ for (auto i = itBegin; i != itLast; i++) {
144+ if (distanceBest > i->distance2 ()) {
145+ dBest = i;
146+ distanceBest = i->distance2 ();
147+ }
148+ }
149+
127150 if (dBest->distance2 () > rParam_ * rParam_)
128151 break ;
129152
@@ -133,11 +156,12 @@ void DisplacedRegionSeedingVertexProducer::produce(edm::StreamID streamID,
133156 const auto distancePred = [](const Distance &a) {
134157 return (!a.entities ().first ->valid () || !a.entities ().second ->valid ());
135158 };
136- const auto pseudoROIPred = [](const DisplacedVertexCluster &a) { return !a.valid (); };
137- distances.remove_if (distancePred);
138- pseudoROIs.remove_if (pseudoROIPred);
159+ itLast = std::remove_if (itBegin, itLast, distancePred);
139160 }
140161
162+ const auto pseudoROIPred = [](const DisplacedVertexCluster &a) { return !a.valid (); };
163+ auto remove_invalid = std::remove_if (pseudoROIs.begin (), pseudoROIs.end (), pseudoROIPred);
164+
141165 // Remove invalid ROIs.
142166 const auto roiPred = [&](const DisplacedVertexCluster &roi) {
143167 if (!roi.valid ())
@@ -150,15 +174,16 @@ void DisplacedRegionSeedingVertexProducer::produce(edm::StreamID streamID,
150174 return true ;
151175 return false ;
152176 };
153- pseudoROIs.remove_if ( roiPred);
177+ auto remove_pred = std::remove_if ( pseudoROIs.begin (), remove_invalid, roiPred);
154178
155179 auto nearRegionsOfInterest = make_unique<vector<reco::Vertex> >();
156180 auto farRegionsOfInterest = make_unique<vector<reco::Vertex> >();
157181
158182 constexpr std::array<double , 6 > errorA{{1.0 , 0.0 , 1.0 , 0.0 , 0.0 , 1.0 }};
159183 static const reco::Vertex::Error errorRegion (errorA.begin (), errorA.end (), true , true );
160184
161- for (const auto &roi : pseudoROIs) {
185+ for (auto it = pseudoROIs.begin (); it != remove_pred; ++it) {
186+ auto const &roi = *it;
162187 const auto &x (roi.centerOfMass ());
163188 if ((x - bs).rho () < nearThreshold_)
164189 nearRegionsOfInterest->emplace_back (reco::Vertex::Point (roi.centerOfMass ()), errorRegion);
@@ -180,6 +205,7 @@ void DisplacedRegionSeedingVertexProducer::fillDescriptions(edm::ConfigurationDe
180205 desc.add <double >(" discriminatorCut" , -1.0 );
181206 desc.add <vector<string> >(" input_names" , {" phi_0" , " phi_1" });
182207 desc.add <vector<string> >(" output_names" , {" model_5/activation_10/Softmax" });
208+ desc.add <unsigned int >(" maxPseudoROIs" , 10000 );
183209 desc.addUntracked <unsigned >(" nThreads" , 1 );
184210 desc.add <edm::FileInPath>(
185211 " graph_path" ,
0 commit comments