CMS 3D CMS Logo

recoBSVTagInfoValidationAnalyzer.cc
Go to the documentation of this file.
2 #include "TH1D.h"
3 #include "TH1F.h"
4 #include "TH2D.h"
5 #include "TMath.h"
6 #include <Math/Functions.h>
7 #include <Math/SMatrix.h>
8 #include <Math/SVector.h>
9 #include <cmath>
10 #include <map>
11 #include <string>
12 
20 
24 #include "Math/VectorUtil.h"
26 #include "TROOT.h"
27 #include <Math/GenVector/PxPyPzE4D.h>
28 #include <Math/GenVector/PxPyPzM4D.h>
29 #include <TVector3.h>
30 //
31 // class decleration
32 //
33 using namespace reco;
34 using namespace std;
35 using namespace edm;
36 
38 public:
40 
41 private:
42  void analyze(const edm::Event &, const edm::EventSetup &) override;
43  void endJob() override;
44  // Member data
45 
47 
51 
54 
55  Int_t n_event;
57  Int_t rs_total_nsv;
58  Int_t rs_total_nbv;
60  Int_t rs_total_ncv;
61  Int_t rs_total_nlv;
62  Int_t total_nfake;
63 
65  Int_t sr_total_nsv;
66  Int_t sr_total_nbv;
68  Int_t sr_total_ncv;
69  Int_t sr_total_nlv;
70  Int_t total_nmiss;
71 
72  // Bookeeping of all the histograms per category
73  void bookRecoToSim(std::string const &);
74  void bookSimToReco(std::string const &);
75 
76  // Fill all histogram per category
77  void fillRecoToSim(std::string const &, reco::Vertex const &, TrackingVertexRef const &);
78  void fillSimToReco(std::string const &, reco::VertexBaseRef const &, TrackingVertexRef const &);
79 
80  // Histogram handlers
81  std::map<std::string, MonitorElement *> HistIndex_;
82 
83  // consumes
86 };
87 
89  : classifier_(config, consumesCollector()) {
90  // Initialize counters
91  n_event = 0;
92  rs_total_nall = 0;
93  rs_total_nsv = 0;
94  rs_total_nbv = 0;
95  rs_total_nbsv = 0;
96  rs_total_ncv = 0;
97  rs_total_nlv = 0;
98  total_nfake = 0;
99 
100  sr_total_nall = 0;
101  sr_total_nsv = 0;
102  sr_total_nbv = 0;
103  sr_total_nbsv = 0;
104  sr_total_ncv = 0;
105  sr_total_nlv = 0;
106  total_nmiss = 0;
107 
108  // get the store
110  dqmLabel = "SVValidation/";
111  dqmStore_->setCurrentFolder(dqmLabel);
112 
113  // Get the track collection
114  svInfoToken = consumes<reco::SecondaryVertexTagInfoCollection>(config.getParameter<InputTag>("svTagInfoProducer"));
115  // Name of the traking pariticle collection
116  tvToken = consumes<TrackingVertexCollection>(config.getParameter<InputTag>("trackingTruth"));
117  // Number of track categories
119 
120  // Define histogram for counting categories
121  HistIndex_["VertexClassifier"] = dqmStore_->book1D("VertexClassifier",
122  "Frequency for the different track categories",
124  -0.5,
126 
127  //--- RecoToSim
128  HistIndex_["rs_All_MatchQuality"] = dqmStore_->book1D("rs_All_MatchQuality", "Quality of Match", 51, -0.01, 1.01);
129  HistIndex_["rs_All_FlightDistance2d"] =
130  dqmStore_->book1D("rs_All_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5);
131  HistIndex_["rs_SecondaryVertex_FlightDistance2d"] =
132  dqmStore_->book1D("rs_SecondaryVertex_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5);
133  HistIndex_["rs_BSV_FlightDistance2d"] =
134  dqmStore_->book1D("rs_BSV_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5);
135  HistIndex_["rs_BWeakDecay_FlightDistance2d"] =
136  dqmStore_->book1D("rs_BWeakDecay_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5);
137  HistIndex_["rs_CWeakDecay_FlightDistance2d"] =
138  dqmStore_->book1D("rs_CWeakDecay_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5);
139  HistIndex_["rs_Light_FlightDistance2d"] =
140  dqmStore_->book1D("rs_Light_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5);
141 
142  HistIndex_["rs_All_nRecVtx"] = dqmStore_->book1D("rs_All_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
143  HistIndex_["rs_SecondaryVertex_nRecVtx"] =
144  dqmStore_->book1D("rs_SecondaryVertex_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
145  HistIndex_["rs_BSV_nRecVtx"] = dqmStore_->book1D("rs_BSV_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
146  HistIndex_["rs_BWeakDecay_nRecVtx"] =
147  dqmStore_->book1D("rs_BWeakDecay_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
148  HistIndex_["rs_CWeakDecay_nRecVtx"] =
149  dqmStore_->book1D("rs_CWeakDecay_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
150  HistIndex_["rs_Light_nRecVtx"] =
151  dqmStore_->book1D("rs_Light_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
152 
153  //--- SimToReco
154  HistIndex_["sr_All_MatchQuality"] = dqmStore_->book1D("sr_All_MatchQuality", "Quality of Match", 51, -0.01, 1.01);
155  HistIndex_["sr_All_nRecVtx"] = dqmStore_->book1D("sr_All_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
156  HistIndex_["sr_SecondaryVertex_nRecVtx"] =
157  dqmStore_->book1D("sr_SecondaryVertex_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
158  HistIndex_["sr_BSV_nRecVtx"] = dqmStore_->book1D("sr_BSV_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
159  HistIndex_["sr_BWeakDecay_nRecVtx"] =
160  dqmStore_->book1D("sr_BWeakDecay_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
161  HistIndex_["sr_CWeakDecay_nRecVtx"] =
162  dqmStore_->book1D("sr_CWeakDecay_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
163  HistIndex_["sr_Light_nRecVtx"] =
164  dqmStore_->book1D("sr_Light_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
165 
166  // Set the proper categories names
167  for (Int_t i = 0; i < numberVertexClassifier_; ++i)
168  HistIndex_["VertexClassifier"]->getTH1F()->GetXaxis()->SetBinLabel(i + 1, VertexCategories::Names[i]);
169 
170  // book histograms
171  bookRecoToSim("rs_All");
172  bookRecoToSim("rs_SecondaryVertex");
173  bookRecoToSim("rs_BSV");
174  bookRecoToSim("rs_BWeakDecay");
175  bookRecoToSim("rs_CWeakDecay");
176  bookRecoToSim("rs_Light");
177 
178  bookSimToReco("sr_All");
179  bookSimToReco("sr_SecondaryVertex");
180  bookSimToReco("sr_BSV");
181  bookSimToReco("sr_BWeakDecay");
182  bookSimToReco("sr_CWeakDecay");
183  bookSimToReco("sr_Light");
184 }
185 
187  ++n_event;
188 
189  std::cout << "*** Analyzing " << event.id() << " n_event = " << n_event << std::endl << std::endl;
190 
191  // Set the classifier for a new event
192  classifier_.newEvent(event, setup);
193 
194  // Vertex collection
196  event.getByToken(svInfoToken, svTagInfoCollection);
197  // Get a constant reference to the track history associated to the classifier
198  VertexHistory const &tracer = classifier_.history();
199 
200  cout << "* Event " << n_event << " ; svTagInfoCollection->size() = " << svTagInfoCollection->size() << endl;
201 
202  int rs_nall = 0;
203  int rs_nsv = 0;
204  int rs_nbv = 0;
205  int rs_nbsv = 0;
206  int rs_ncv = 0;
207  int rs_nlv = 0;
208  int nfake = 0;
209 
210  int sr_nall = 0;
211  int sr_nsv = 0;
212  int sr_nbv = 0;
213  int sr_nbsv = 0;
214  int sr_ncv = 0;
215  int sr_nlv = 0;
216  int nmiss = 0;
217 
218  // Loop over the svTagInfo collection.
219  for (std::size_t index = 0; index < svTagInfoCollection->size(); ++index) {
220  reco::SecondaryVertexTagInfoRef svTagInfo(svTagInfoCollection, index);
221 
222  // Loop over the vertexes in svTagInfo
223  for (std::size_t vindex = 0; vindex < svTagInfo->nVertices(); ++vindex) {
224  // Classify the vertices
225  classifier_.evaluate(svTagInfo, vindex);
226 
227  // quality of the match
228  double rs_quality = tracer.quality();
229 
230  // Fill the histogram with the categories
231  for (Int_t i = 0; i != numberVertexClassifier_; ++i) {
233  HistIndex_["VertexClassifier"]->Fill(i);
234  }
235  }
237  HistIndex_["rs_All_MatchQuality"]->Fill(rs_quality);
238  fillRecoToSim("rs_All", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
239  HistIndex_["rs_All_FlightDistance2d"]->Fill(svTagInfo->flightDistance(vindex, true).value());
240  rs_nall++;
241 
243  fillRecoToSim("rs_SecondaryVertex", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
244  HistIndex_["rs_SecondaryVertex_FlightDistance2d"]->Fill(svTagInfo->flightDistance(vindex, true).value());
245  rs_nsv++;
246  }
247 
249  fillRecoToSim("rs_BWeakDecay", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
250  HistIndex_["rs_BWeakDecay_FlightDistance2d"]->Fill(svTagInfo->flightDistance(vindex, true).value());
251  rs_nbv++;
252 
254  fillRecoToSim("rs_BSV", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
255  HistIndex_["rs_BSV_FlightDistance2d"]->Fill(svTagInfo->flightDistance(vindex, true).value());
256  rs_nbsv++;
257  }
258  } // BWeakDecay
259 
261  fillRecoToSim("rs_CWeakDecay", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
262  HistIndex_["rs_CWeakDecay_FlightDistance2d"]->Fill(svTagInfo->flightDistance(vindex, true).value());
263  rs_ncv++;
264 
265  } else {
266  fillRecoToSim("rs_Light", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
267  HistIndex_["rs_Light_FlightDistance2d"]->Fill(svTagInfo->flightDistance(vindex, true).value());
268  rs_nlv++;
269  }
270  } // end if classifier
271 
272  else {
273  cout << " VertexCategories::Fake!!" << endl;
274  nfake++;
275  }
276 
277  } // end loop over vertices in svTagInfo
278 
279  } // loop over svTagInfo
280 
281  HistIndex_["rs_All_nRecVtx"]->Fill(rs_nall);
282  HistIndex_["rs_SecondaryVertex_nRecVtx"]->Fill(rs_nsv);
283  HistIndex_["rs_BWeakDecay_nRecVtx"]->Fill(rs_nbv);
284  HistIndex_["rs_BSV_nRecVtx"]->Fill(rs_nbsv);
285  HistIndex_["rs_CWeakDecay_nRecVtx"]->Fill(rs_ncv);
286  HistIndex_["rs_Light_nRecVtx"]->Fill(rs_nlv);
287  cout << endl;
288 
289  //----------------------------------------------------------------
290  // SIM TO RECO!
291 
292  // Vertex collection
294  event.getByToken(tvToken, TVCollection);
295  // Loop over the TV collection.
296  for (std::size_t index = 0; index < TVCollection->size(); ++index) {
297  TrackingVertexRef trackingVertex(TVCollection, index);
298 
299  classifier_.evaluate(trackingVertex);
300 
301  double sr_quality = tracer.quality();
302 
304  HistIndex_["sr_All_MatchQuality"]->Fill(sr_quality);
305  fillSimToReco("sr_All", tracer.recoVertex(), trackingVertex);
306  sr_nall++;
307 
309  fillSimToReco("sr_SecondaryVertex", tracer.recoVertex(), trackingVertex);
310  sr_nsv++;
311  }
312 
314  fillSimToReco("sr_BWeakDecay", tracer.recoVertex(), trackingVertex);
315  sr_nbv++;
316 
318  fillSimToReco("sr_BSV", tracer.recoVertex(), trackingVertex);
319  sr_nbsv++;
320  }
321 
322  } // BWeakDecay
323 
325  fillSimToReco("sr_CWeakDecay", tracer.recoVertex(), trackingVertex);
326  sr_ncv++;
327  }
328 
329  else {
330  fillSimToReco("sr_Light", tracer.recoVertex(), trackingVertex);
331  sr_nlv++;
332  }
333 
334  } // Reconstructed
335  else {
336  // cout << "##### Not reconstructed!" << endl;
337  nmiss++;
338  }
339 
340  } // TVCollection.size()
341 
342  HistIndex_["sr_All_nRecVtx"]->Fill(sr_nall);
343  HistIndex_["sr_SecondaryVertex_nRecVtx"]->Fill(sr_nsv);
344  HistIndex_["sr_BWeakDecay_nRecVtx"]->Fill(sr_nbv);
345  HistIndex_["sr_BSV_nRecVtx"]->Fill(sr_nbsv);
346  HistIndex_["sr_CWeakDecay_nRecVtx"]->Fill(sr_ncv);
347  HistIndex_["sr_Light_nRecVtx"]->Fill(rs_nlv);
348 
349  rs_total_nall += rs_nall;
350  rs_total_nsv += rs_nsv;
351  rs_total_nbv += rs_nbv;
352  rs_total_nbsv += rs_nbsv;
353  rs_total_ncv += rs_ncv;
354  rs_total_nlv += rs_nlv;
355  total_nfake += nfake;
356 
357  sr_total_nall += sr_nall;
358  sr_total_nsv += sr_nsv;
359  sr_total_nbv += sr_nbv;
360  sr_total_nbsv += sr_nbsv;
361  sr_total_ncv += sr_ncv;
362  sr_total_nlv += sr_nlv;
363  total_nmiss += nmiss;
364 }
365 
367  // Book pull histograms
368 
369  std::string name = prefix + "_Pullx";
370  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -10., 10.);
371  name = prefix + "_Pully";
372  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -10., 10.);
373  name = prefix + "_Pullz";
374  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -10., 10.);
375 
376  name = prefix + "_Resx";
377  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.05, 0.05);
378  name = prefix + "_Resy";
379  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.05, 0.05);
380  name = prefix + "_Resz";
381  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.05, 0.05);
382 
383  name = prefix + "_Chi2Norm";
384  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, 0, 10.);
385  name = prefix + "_Chi2Prob";
386  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, 0., 1.);
387 
388  name = prefix + "_nRecTrks";
389  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 501, -0.5, 500.5);
390 
391  name = prefix + "_AverageTrackWeight";
392  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.1, 1.1);
393 
394  name = prefix + "_Mass";
395  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 65, 0., 6.5);
396 
397  name = prefix + "_RecPt";
398  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 2000, 0., 1000.);
399 
400  name = prefix + "_RecEta";
401  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
402 
403  name = prefix + "_RecCharge";
404  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 21, -0.5, 20.5);
405 
406  name = prefix + "_RecTrackPt";
407  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 2000, 0., 1000.);
408 
409  name = prefix + "_RecTrackEta";
410  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
411 
412  name = prefix + "_nSimTrks";
413  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 501, -0.5, 500.5);
414 
415  name = prefix + "_SimPt";
416  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 2000, 0., 1000.);
417 
418  name = prefix + "_SimEta";
419  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
420 
421  name = prefix + "_SimCharge";
422  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 21, -0.5, 20.5);
423 
424  name = prefix + "_SimTrackPt";
425  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 500, 0., 500.);
426 
427  name = prefix + "_SimTrackEta";
428  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
429 }
430 
432  // Book pull histograms
433 
434  std::string name = prefix + "_Pullx";
435  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -10., 10.);
436  name = prefix + "_Pully";
437  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -10., 10.);
438  name = prefix + "_Pullz";
439  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -10., 10.);
440 
441  name = prefix + "_Resx";
442  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.05, 0.05);
443  name = prefix + "_Resy";
444  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.05, 0.05);
445  name = prefix + "_Resz";
446  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.05, 0.05);
447 
448  name = prefix + "_Chi2Norm";
449  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, 0, 10.);
450  name = prefix + "_Chi2Prob";
451  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, 0., 1.);
452 
453  name = prefix + "_nRecTrks";
454  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 501, -0.5, 500.5);
455 
456  name = prefix + "_AverageTrackWeight";
457  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.1, 1.1);
458 
459  name = prefix + "_Mass";
460  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 65, 0., 6.5);
461 
462  name = prefix + "_RecPt";
463  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 2000, 0., 1000.);
464 
465  name = prefix + "_RecEta";
466  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
467 
468  name = prefix + "_RecCharge";
469  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 21, -0.5, 20.5);
470 
471  name = prefix + "_RecTrackPt";
472  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 2000, 0., 1000.);
473 
474  name = prefix + "_RecTrackEta";
475  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
476 
477  name = prefix + "_nSimTrks";
478  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 501, -0.5, 500.5);
479 
480  name = prefix + "_SimPt";
481  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 2000, 0., 1000.);
482 
483  name = prefix + "_SimEta";
484  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
485 
486  name = prefix + "_SimCharge";
487  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 21, -0.5, 20.5);
488 
489  name = prefix + "_SimTrackPt";
490  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 500, 0., 500.);
491 
492  name = prefix + "_SimTrackEta";
493  HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
494 }
495 
497  reco::Vertex const &vertex,
498  TrackingVertexRef const &simVertex) {
499  double pullx = (vertex.x() - simVertex->position().x()) / vertex.xError();
500  double pully = (vertex.y() - simVertex->position().y()) / vertex.yError();
501  double pullz = (vertex.z() - simVertex->position().z()) / vertex.zError();
502 
503  double resx = vertex.x() - simVertex->position().x();
504  double resy = vertex.y() - simVertex->position().y();
505  double resz = vertex.z() - simVertex->position().z();
506 
507  double chi2norm = vertex.normalizedChi2();
508  double chi2prob = ChiSquaredProbability(vertex.chi2(), vertex.ndof());
509 
510  double sum_weight = 0.;
511  double weight = 0.;
512  double tracksize = vertex.tracksSize();
513  math::XYZVector momentum;
515  int charge = 0;
516  double thePiMass = 0.13957;
517  for (reco::Vertex::trackRef_iterator recDaughter = vertex.tracks_begin(); recDaughter != vertex.tracks_end();
518  ++recDaughter) {
519  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzM4D<double>> vec;
520 
521  vec.SetPx((**recDaughter).px());
522  vec.SetPy((**recDaughter).py());
523  vec.SetPz((**recDaughter).pz());
524  vec.SetM(thePiMass);
525 
526  sum += vec;
527 
528  weight = vertex.trackWeight(*recDaughter);
529  sum_weight += weight;
530 
532  p = (**recDaughter).momentum();
533  momentum += p;
534 
535  charge += (*recDaughter)->charge();
536 
537  HistIndex_[prefix + "_RecTrackPt"]->Fill((*recDaughter)->pt());
538  HistIndex_[prefix + "_RecTrackEta"]->Fill((*recDaughter)->eta());
539  } // end loop to recDaughters
540  // cout << " average sum of weights = " <<
541  // sum_weight/tracksize << endl;
542 
543  double mass = sum.M();
544  double pt = sqrt(momentum.Perp2());
545  double eta = momentum.Eta();
546 
547  math::XYZVector simmomentum;
548  int simcharge = 0;
549  for (TrackingVertex::tp_iterator simDaughter = simVertex->daughterTracks_begin();
550  simDaughter != simVertex->daughterTracks_end();
551  ++simDaughter) {
553  p = (**simDaughter).momentum();
554  simmomentum += p;
555 
556  simcharge += (*simDaughter)->charge();
557 
558  HistIndex_[prefix + "_SimTrackPt"]->Fill((*simDaughter)->pt());
559  HistIndex_[prefix + "_SimTrackEta"]->Fill((*simDaughter)->eta());
560  }
561 
562  double simpt = sqrt(simmomentum.Perp2());
563  double simeta = simmomentum.Eta();
564 
565  // cout << "[fillRecoToSim] vertex.tracksSize() = " << vertex.tracksSize() <<
566  // " ; simVertex->nDaughterTracks() = " << simVertex->nDaughterTracks() <<
567  // endl;
568 
569  HistIndex_[prefix + "_nRecTrks"]->Fill(vertex.tracksSize());
570  HistIndex_[prefix + "_nSimTrks"]->Fill(simVertex->nDaughterTracks());
571  HistIndex_[prefix + "_Pullx"]->Fill(pullx);
572  HistIndex_[prefix + "_Pully"]->Fill(pully);
573  HistIndex_[prefix + "_Pullz"]->Fill(pullz);
574  HistIndex_[prefix + "_Resx"]->Fill(resx);
575  HistIndex_[prefix + "_Resy"]->Fill(resy);
576  HistIndex_[prefix + "_Resz"]->Fill(resz);
577  HistIndex_[prefix + "_AverageTrackWeight"]->Fill(sum_weight / tracksize);
578  HistIndex_[prefix + "_Chi2Norm"]->Fill(chi2norm);
579  HistIndex_[prefix + "_Chi2Prob"]->Fill(chi2prob);
580  HistIndex_[prefix + "_RecPt"]->Fill(pt);
581  HistIndex_[prefix + "_RecEta"]->Fill(eta);
582  HistIndex_[prefix + "_RecCharge"]->Fill(charge);
583  HistIndex_[prefix + "_Mass"]->Fill(mass);
584  HistIndex_[prefix + "_SimPt"]->Fill(simpt);
585  HistIndex_[prefix + "_SimEta"]->Fill(simeta);
586  HistIndex_[prefix + "_SimCharge"]->Fill(simcharge);
587 }
588 
590  reco::VertexBaseRef const &vertex,
591  TrackingVertexRef const &simVertex) {
592  double pullx = (vertex->x() - simVertex->position().x()) / vertex->xError();
593  double pully = (vertex->y() - simVertex->position().y()) / vertex->yError();
594  double pullz = (vertex->z() - simVertex->position().z()) / vertex->zError();
595 
596  double resx = vertex->x() - simVertex->position().x();
597  double resy = vertex->y() - simVertex->position().y();
598  double resz = vertex->z() - simVertex->position().z();
599 
600  double chi2norm = vertex->normalizedChi2();
601  double chi2prob = ChiSquaredProbability(vertex->chi2(), vertex->ndof());
602 
603  double sum_weight = 0.;
604  double weight = 0.;
605  double tracksize = vertex->tracksSize();
606  math::XYZVector momentum;
608  int charge = 0;
609  double thePiMass = 0.13957;
610  for (reco::Vertex::trackRef_iterator recDaughter = vertex->tracks_begin(); recDaughter != vertex->tracks_end();
611  ++recDaughter) {
612  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzM4D<double>> vec;
613 
614  vec.SetPx((**recDaughter).px());
615  vec.SetPy((**recDaughter).py());
616  vec.SetPz((**recDaughter).pz());
617  vec.SetM(thePiMass);
618 
619  sum += vec;
620 
621  weight = vertex->trackWeight(*recDaughter);
622  sum_weight += weight;
623 
625  p = (**recDaughter).momentum();
626  momentum += p;
627 
628  charge += (*recDaughter)->charge();
629 
630  HistIndex_[prefix + "_RecTrackPt"]->Fill((*recDaughter)->pt());
631  HistIndex_[prefix + "_RecTrackEta"]->Fill((*recDaughter)->eta());
632  }
633  // cout << " average sum of weights = " <<
634  // sum_weight/tracksize << endl;
635 
636  double mass = sum.M();
637  double pt = sqrt(momentum.Perp2());
638  double eta = momentum.Eta();
639 
640  math::XYZVector simmomentum;
641  int simcharge = 0;
642  for (TrackingVertex::tp_iterator simDaughter = simVertex->daughterTracks_begin();
643  simDaughter != simVertex->daughterTracks_end();
644  ++simDaughter) {
646  p = (**simDaughter).momentum();
647  simmomentum += p;
648 
649  simcharge += (*simDaughter)->charge();
650 
651  HistIndex_[prefix + "_SimTrackPt"]->Fill((*simDaughter)->pt());
652  HistIndex_[prefix + "_SimTrackEta"]->Fill((*simDaughter)->eta());
653  }
654 
655  double simpt = sqrt(simmomentum.Perp2());
656  double simeta = simmomentum.Eta();
657 
658  // cout << "[fillSimToReco] vertex->tracksSize() = " << vertex->tracksSize()
659  // << " ; simVertex->nDaughterTracks() = " << simVertex->nDaughterTracks() <<
660  // endl;
661 
662  HistIndex_[prefix + "_nRecTrks"]->Fill(vertex->tracksSize());
663  HistIndex_[prefix + "_nSimTrks"]->Fill(simVertex->nDaughterTracks());
664  HistIndex_[prefix + "_Pullx"]->Fill(pullx);
665  HistIndex_[prefix + "_Pully"]->Fill(pully);
666  HistIndex_[prefix + "_Pullz"]->Fill(pullz);
667  HistIndex_[prefix + "_Resx"]->Fill(resx);
668  HistIndex_[prefix + "_Resy"]->Fill(resy);
669  HistIndex_[prefix + "_Resz"]->Fill(resz);
670  HistIndex_[prefix + "_AverageTrackWeight"]->Fill(sum_weight / tracksize);
671  HistIndex_[prefix + "_Chi2Norm"]->Fill(chi2norm);
672  HistIndex_[prefix + "_Chi2Prob"]->Fill(chi2prob);
673  HistIndex_[prefix + "_RecPt"]->Fill(pt);
674  HistIndex_[prefix + "_RecEta"]->Fill(eta);
675  HistIndex_[prefix + "_RecCharge"]->Fill(charge);
676  HistIndex_[prefix + "_Mass"]->Fill(mass);
677  HistIndex_[prefix + "_SimPt"]->Fill(simpt);
678  HistIndex_[prefix + "_SimEta"]->Fill(simeta);
679  HistIndex_[prefix + "_SimCharge"]->Fill(simcharge);
680 }
681 
683  std::cout << std::endl;
684  std::cout << " ====== Total Number of analyzed events: " << n_event << " ====== " << std::endl;
685  std::cout << " ====== Total Number of R2S All: " << rs_total_nall << " ====== " << std::endl;
686  std::cout << " ====== Total Number of R2S SecondaryVertex: " << rs_total_nsv << " ====== " << std::endl;
687  std::cout << " ====== Total Number of R2S BWeakDecay: " << rs_total_nbv << " ====== " << std::endl;
688  std::cout << " ====== Total Number of R2S BWeakDecay::SecondaryVertex: " << rs_total_nbsv << " ====== " << std::endl;
689  std::cout << " ====== Total Number of R2S CWeakDecay: " << rs_total_ncv << " ====== " << std::endl;
690  std::cout << " ====== Total Number of R2S Light: " << rs_total_nlv << " ====== " << std::endl;
691  std::cout << std::endl;
692  std::cout << " ====== Total Number of S2R All: " << sr_total_nall << " ====== " << std::endl;
693  std::cout << " ====== Total Number of S2R SecondaryVertex: " << sr_total_nsv << " ====== " << std::endl;
694  std::cout << " ====== Total Number of S2R BWeakDecay: " << sr_total_nbv << " ====== " << std::endl;
695  std::cout << " ====== Total Number of S2R BWeakDecay::SecondaryVertex: " << sr_total_nbsv << " ====== " << std::endl;
696  std::cout << " ====== Total Number of S2R CWeakDecay: " << sr_total_ncv << " ====== " << std::endl;
697  std::cout << " ====== Total Number of S2R Light: " << sr_total_nlv << " ====== " << std::endl;
698  std::cout << std::endl;
699  std::cout << " ====== Total Number of Fake Vertices: " << total_nfake << " ====== " << std::endl;
700 }
701 
T getParameter(std::string const &) const
This class traces the simulated and generated history of a given track.
Definition: VertexHistory.h:18
VertexHistory const & history() const
Returns a reference to the vertex history used in the classification.
trackRef_iterator tracks_end() const
last iterator over tracks
Definition: Vertex.cc:81
double zError() const
error on z
Definition: Vertex.h:123
double y() const
y coordinate
Definition: Vertex.h:113
def setup(process, global_tag, zero_tesla=False)
Definition: GeneralSetup.py:2
Definition: weight.py:1
Definition: config.py:1
const reco::VertexBaseRef & recoVertex() const
Return a reference to the reconstructed track.
Definition: VertexHistory.h:56
example_stream void analyze(const edm::Event &, const edm::EventSetup &) override
void newEvent(edm::Event const &event, edm::EventSetup const &config) override
Pre-process event information (for accessing reconstraction information).
void fillSimToReco(std::string const &, reco::VertexBaseRef const &, TrackingVertexRef const &)
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
Category
Categories available to vertexes.
T sqrt(T t)
Definition: SSEVec.h:18
recoBSVTagInfoValidationAnalyzer(const edm::ParameterSet &)
VertexClassifierByProxy< reco::SecondaryVertexTagInfoCollection > classifier_
edm::EDGetTokenT< TrackingVertexCollection > tvToken
double chi2() const
chi-squares
Definition: Vertex.h:98
double z() const
z coordinate
Definition: Vertex.h:115
std::map< std::string, MonitorElement * > HistIndex_
void analyze(const edm::Event &, const edm::EventSetup &) override
float ChiSquaredProbability(double chiSquared, double nrDOF)
float trackWeight(const TREF &r) const
returns the weight with which a Track has contributed to the vertex-fit.
Definition: Vertex.h:81
bool is(Category category) const
Returns track flag for a given category.
double ndof() const
Definition: Vertex.h:105
const TrackingVertexRef & simVertex() const
Return the initial tracking vertex from the history.
Definition: HistoryBase.h:70
double x() const
x coordinate
Definition: Vertex.h:111
double xError() const
error on x
Definition: Vertex.h:119
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:30
VertexClassifierByProxy< Collection > const & evaluate(TrackingVertexRef const &vertex)
Classify the TrackingVertex in categories.
fixed size matrix
HLT enums.
std::vector< TrackBaseRef >::const_iterator trackRef_iterator
The iteratator for the vector<TrackRef>
Definition: Vertex.h:37
double normalizedChi2() const
chi-squared divided by n.d.o.f.
Definition: Vertex.h:107
double quality() const
Return the quality of the match.
Definition: VertexHistory.h:59
static const char *const Names[]
Name of the different categories.
trackRef_iterator tracks_begin() const
first iterator over tracks
Definition: Vertex.cc:76
void fillRecoToSim(std::string const &, reco::Vertex const &, TrackingVertexRef const &)
edm::EDGetTokenT< reco::SecondaryVertexTagInfoCollection > svInfoToken
double yError() const
error on y
Definition: Vertex.h:121
size_t tracksSize() const
number of tracks
Definition: Vertex.cc:71
Definition: event.py:1