CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
VertexAssociatorByTracks.cc
Go to the documentation of this file.
6 
13 
17 
20 
22 
23 
24 /* Constructor */
26 {
27  R2SMatchedSimRatio_ = config.getParameter<double>("R2SMatchedSimRatio");
28  R2SMatchedRecoRatio_ = config.getParameter<double>("R2SMatchedRecoRatio");
29  S2RMatchedSimRatio_ = config.getParameter<double>("S2RMatchedSimRatio");
30  S2RMatchedRecoRatio_ = config.getParameter<double>("S2RMatchedRecoRatio");
31 
32  std::string trackQualityType = config.getParameter<std::string>("trackQuality");
34 
35  edm::ParameterSet param = config.getParameter<edm::ParameterSet>("trackingParticleSelector");
36 
38  param.getParameter<double>("ptMinTP"),
39  param.getParameter<double>("minRapidityTP"),
40  param.getParameter<double>("maxRapidityTP"),
41  param.getParameter<double>("tipTP"),
42  param.getParameter<double>("lipTP"),
43  param.getParameter<int>("minHitTP"),
44  param.getParameter<bool>("signalOnlyTP"),
45  param.getParameter<bool>("chargedOnlyTP"),
46  param.getParameter<bool>("stableOnlyTP"),
47  param.getParameter<std::vector<int> >("pdgIdTP")
48  );
49 }
50 
51 
52 /* Destructor */
54 {
55  //do cleanup here
56 }
57 
58 
60  edm::Handle<edm::View<reco::Vertex> > & recoVertexes,
61  edm::Handle<TrackingVertexCollection> & trackingVertexes,
62  const edm::Event& event,
63  reco::RecoToSimCollection & associator
64 ) const
65 {
66  reco::VertexRecoToSimCollection outputCollection;
67 
68  std::map<TrackingVertexRef,std::pair<double, std::size_t> > matches;
69 
70  std::cout << "reco::VertexCollection size = " << recoVertexes->size();
71  std::cout << " ; TrackingVertexCollection size = " << trackingVertexes->size() << std::endl << std::endl;
72 
73  // Loop over RecoVertex
74  for (std::size_t recoIndex = 0; recoIndex < recoVertexes->size(); ++recoIndex)
75  {
76  reco::VertexBaseRef recoVertex(recoVertexes, recoIndex);
77 
78  matches.clear();
79 
80  double recoDaughterWeight = 0.;
81 
82  // Loop over daughter tracks of RecoVertex
83  for (
84  reco::Vertex::trackRef_iterator recoDaughter = recoVertex->tracks_begin();
85  recoDaughter != recoVertex->tracks_end();
86  ++recoDaughter
87  )
88  {
89  // Check the quality of the RecoDoughters
90  if ( !(*recoDaughter)->quality(trackQuality_) ) continue;
91 
92  // Check for association for the given RecoDaughter
93  if ( associator.numberOfAssociations(*recoDaughter) > 0 )
94  {
95  std::vector<std::pair<TrackingParticleRef,double> > associations = associator[*recoDaughter];
96 
97  // Loop over TrackingParticles associated with RecoDaughter
98  for (
99  std::vector<std::pair<TrackingParticleRef,double> >::const_iterator association = associations.begin();
100  association != associations.end();
101  ++association
102  )
103  {
104  // Get a reference to parent vertex of TrackingParticle associated to the RecoDaughter
105  TrackingVertexRef trackingVertex = association->first->parentVertex();
106  // Store matched RecoDaughter to the trackingVertex
107  matches[trackingVertex].first += recoVertex->trackWeight(*recoDaughter);
108  matches[trackingVertex].second++;
109  }
110  }
111 
112  recoDaughterWeight += recoVertex->trackWeight(*recoDaughter);
113  } //loop over daughter tracks of RecoVertex
114 
115  std::size_t assoIndex = 0;
116 
117  // Loop over map between TrackingVertexes and matched RecoDaugther
118  for (
119  std::map<TrackingVertexRef,std::pair<double, std::size_t> >::const_iterator match = matches.begin();
120  match != matches.end();
121  ++match
122  )
123  {
124  // Getting the TrackingVertex information
125  TrackingVertexRef trackingVertex = match->first;
126  double matchedDaughterWeight = match->second.first;
127  std::size_t matchedDaughterCounter = match->second.second;
128 
129  // Count for only reconstructible SimDaughters
130  std::size_t simDaughterCounter = 0;
131 
132  for (
133  TrackingVertex::tp_iterator simDaughter = trackingVertex->daughterTracks_begin();
134  simDaughter != trackingVertex->daughterTracks_end();
135  ++simDaughter
136  )
137  if ( selector_(**simDaughter) ) simDaughterCounter++;
138 
139  // Sanity condition in case that reconstructable condition is too tight
140  if ( simDaughterCounter < matchedDaughterCounter )
141  simDaughterCounter = matchedDaughterCounter;
142 
143  // Condition over S2RMatchedSimRatio
144  if ( (double)matchedDaughterCounter/simDaughterCounter < R2SMatchedSimRatio_ ) continue;
145 
146  double quality = (double)matchedDaughterWeight/recoDaughterWeight;
147 
148  // Condition over R2SMatchedRecoRatio
149  if (quality < R2SMatchedRecoRatio_) continue;
150 
151  outputCollection.insert(recoVertex, std::make_pair(trackingVertex, quality));
152 
153  std::cout << "R2S: INDEX " << assoIndex;
154  std::cout << " ; SimDaughterCounter = " << simDaughterCounter;
155  std::cout << " ; RecoDaughterWeight = " << recoDaughterWeight;
156  std::cout << " ; MatchedDaughterCounter = " << matchedDaughterCounter;
157  std::cout << " ; MatchedDaughterWeight = " << matchedDaughterWeight;
158  std::cout << " ; quality = " << quality << std::endl;
159 
160  assoIndex++;
161  }
162  } // Loop on RecoVertex
163 
164  std::cout << std::endl;
165  std::cout << "RecoToSim OUTPUT COLLECTION: outputCollection.size() = " << outputCollection.size() << std::endl << std::endl;
166 
167  return outputCollection;
168 }
169 
170 
171 
173  edm::Handle<edm::View<reco::Vertex> > & recoVertexes,
174  edm::Handle<TrackingVertexCollection> & trackingVertexes,
175  const edm::Event& event,
176  reco::SimToRecoCollection & associator
177 ) const
178 {
179  reco::VertexSimToRecoCollection outputCollection;
180 
181  // Loop over TrackingVertexes
182  std::map<std::size_t,std::pair<double, std::size_t> > matches;
183 
184  // Loop over TrackingVertexes
185  for (std::size_t simIndex = 0; simIndex < trackingVertexes->size(); ++simIndex)
186  {
187  TrackingVertexRef trackingVertex(trackingVertexes, simIndex);
188 
189  matches.clear();
190 
191  std::size_t simDaughterCounter = 0;
192 
193  // Loop over daughter tracks of TrackingVertex
194  for (
195  TrackingVertex::tp_iterator simDaughter = trackingVertex->daughterTracks_begin();
196  simDaughter != trackingVertex->daughterTracks_end();
197  ++simDaughter
198  )
199  {
200  // Select only reconstructible SimDaughters
201  if ( !selector_(**simDaughter) ) continue;
202 
203  // Check for association for the given RecoDaughter
204  if ( associator.numberOfAssociations(*simDaughter) > 0 )
205  {
206  std::vector<std::pair<reco::TrackBaseRef, double> > associations = associator[*simDaughter];
207 
208  // Loop over RecoTracks associated with TrackingParticle
209  for (
210  std::vector<std::pair<reco::TrackBaseRef,double> >::const_iterator association = associations.begin();
211  association != associations.end();
212  ++association
213  )
214  {
215  reco::TrackBaseRef recoTrack = association->first;
216 
217  for (std::size_t recoIndex = 0; recoIndex < recoVertexes->size(); ++recoIndex)
218  {
219  reco::VertexBaseRef recoVertex(recoVertexes, recoIndex);
220 
221  for (
222  reco::Vertex::trackRef_iterator recoDaughter = recoVertex->tracks_begin();
223  recoDaughter != recoVertex->tracks_end();
224  ++recoDaughter
225  )
226  {
227  // Store matched RecoDaughter to the RecoVertex
228  if (
229  recoDaughter->id() == recoTrack.id() &&
230  recoDaughter->key() == recoTrack.key()
231  )
232  {
233  matches[recoIndex].first += recoVertex->trackWeight(*recoDaughter);
234  matches[recoIndex].second++;
235  }
236  }
237  }
238  } // loop a recotracks
239  }
240 
241  simDaughterCounter++;
242  } // loop a simDaughter
243 
244  std::size_t assoIndex = 0;
245 
246  // Loop over map, set score, add to outputCollection
247  for (
248  std::map<std::size_t,std::pair<double,std::size_t> >::const_iterator match = matches.begin();
249  match != matches.end();
250  ++match
251  )
252  {
253  // Getting the TrackingVertex information
254  reco::VertexBaseRef recoVertex(recoVertexes, match->first);
255  double matchedDaughterWeight = match->second.first;
256  std::size_t matchedDaughterCounter = match->second.second;
257 
258  double recoDaughterWeight = 0.;
259 
260  // Weighted count of those tracks with a given quality of the RecoDoughters
261  for (
262  reco::Vertex::trackRef_iterator recoDaughter = recoVertex->tracks_begin();
263  recoDaughter != recoVertex->tracks_end();
264  ++recoDaughter
265  )
266  if ( (*recoDaughter)->quality(trackQuality_) )
267  recoDaughterWeight += recoVertex->trackWeight(*recoDaughter);
268 
269  // Sanity condition in case that track quality condition is too tight
270  if ( recoDaughterWeight < matchedDaughterWeight )
271  recoDaughterWeight = matchedDaughterWeight;
272 
273  // Condition over S2RMatchedRecoRatio
274  if ( matchedDaughterWeight/recoDaughterWeight < S2RMatchedRecoRatio_ ) continue;
275 
276  double quality = (double)matchedDaughterCounter/simDaughterCounter;
277 
278  // Condition over S2RMatchedSimRatio
279  if (quality < S2RMatchedSimRatio_) continue;
280 
281  outputCollection.insert(trackingVertex, std::make_pair(recoVertex, quality));
282 
283  std::cout << "R2S: INDEX " << assoIndex;
284  std::cout << " ; SimDaughterCounter = " << simDaughterCounter;
285  std::cout << " ; RecoDaughterWeight = " << recoDaughterWeight;
286  std::cout << " ; MatchedDaughterCounter = " << matchedDaughterCounter;
287  std::cout << " ; MatchedDaughterWeight = " << matchedDaughterWeight;
288  std::cout << " ; quality = " << quality << std::endl;
289 
290  assoIndex++;
291  }
292  } // loop over TrackingVertexes
293 
294  std::cout << "SimToReco OUTPUT COLLECTION: outputCollection.size() = " << outputCollection.size() << std::endl << std::endl;
295 
296  return outputCollection;
297 }
298 
T getParameter(std::string const &) const
trackRef_iterator tracks_end() const
last iterator over tracks
Definition: Vertex.cc:45
reco::TrackBase::TrackQuality trackQuality_
ProductID id() const
Definition: RefToBase.h:220
TrackingParticleSelector selector_
reco::VertexSimToRecoCollection associateSimToReco(edm::Handle< edm::View< reco::Vertex > > &, edm::Handle< TrackingVertexCollection > &, const edm::Event &, reco::SimToRecoCollection &) const
size_type numberOfAssociations(const key_type &k) const
number of associations to a key
reco::VertexRecoToSimCollection associateRecoToSim(edm::Handle< edm::View< reco::Vertex > > &, edm::Handle< TrackingVertexCollection > &, const edm::Event &, reco::RecoToSimCollection &) const
SingleObjectSelector< TrackingParticleCollection,::TrackingParticleSelector > TrackingParticleSelector
float trackWeight(const TrackBaseRef &r) const
returns the weight with which a Track has contributed to the vertex-fit.
How EventSelector::AcceptEvent() decides whether to accept an event for output otherwise it is excluding the probing of A single or multiple positive and the trigger will pass if any such matching triggers are PASS or EXCEPTION[A criterion thatmatches no triggers at all is detected and causes a throw.] A single negative with an expectation of appropriate bit checking in the decision and the trigger will pass if any such matching triggers are FAIL or EXCEPTION A wildcarded negative criterion that matches more than one trigger in the trigger but the state exists so we define the behavior If all triggers are the negative crieriion will lead to accepting the event(this again matches the behavior of"!*"before the partial wildcard feature was incorporated).The per-event"cost"of each negative criterion with multiple relevant triggers is about the same as!*was in the past
size_t key() const
Definition: RefToBase.h:228
static TrackQuality qualityByName(const std::string &name)
Definition: TrackBase.cc:46
size_type size() const
map size
void insert(const key_type &k, const data_type &v)
insert an association
VertexAssociatorByTracks(const edm::ParameterSet &)
std::vector< TrackBaseRef >::const_iterator trackRef_iterator
The iteratator for the vector&lt;TrackRef&gt;
Definition: Vertex.h:38
tuple cout
Definition: gather_cfg.py:41
std::pair< typename Association::data_type::first_type, double > match(Reference key, Association association, bool bestMatchByMaxValue)
Generic matching function.
Definition: Utils.h:6
trackRef_iterator tracks_begin() const
first iterator over tracks
Definition: Vertex.cc:40