CMS 3D CMS Logo
Primary4DVertexValidation.cc
Go to the documentation of this file.
1 #include <numeric>
2 #include <string>
3 #include <vector>
4 #include <map>
5 #include <algorithm>
6 
7 #include "FWCore/Framework/interface/Frameworkfwd.h"
8 #include "FWCore/Framework/interface/ESHandle.h"
9 #include "FWCore/Utilities/interface/isFinite.h"
10 
11 #include "DataFormats/Common/interface/ValidHandle.h"
12 #include "FWCore/Framework/interface/Event.h"
13 #include "FWCore/Framework/interface/EventSetup.h"
14 #include "FWCore/Framework/interface/MakerMacros.h"
15 #include "FWCore/ParameterSet/interface/ParameterSet.h"
16 #include "FWCore/MessageLogger/interface/MessageLogger.h"
17 #include "FWCore/Utilities/interface/isFinite.h"
18 
19 #include "DataFormats/Math/interface/LorentzVector.h"
20 #include "DataFormats/Math/interface/Point3D.h"
21 
22 // reco track and vertex
23 #include "DataFormats/TrackReco/interface/Track.h"
24 #include "DataFormats/TrackReco/interface/TrackFwd.h"
25 #include "DataFormats/VertexReco/interface/Vertex.h"
26 #include "DataFormats/VertexReco/interface/VertexFwd.h"
27 #include "RecoVertex/VertexPrimitives/interface/TransientVertex.h"
28 
29 // TrackingParticle
30 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
31 #include "SimDataFormats/Associations/interface/TrackToTrackingParticleAssociator.h"
32 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticleFwd.h"
33 #include "SimDataFormats/TrackingAnalysis/interface/TrackingVertexContainer.h"
34 
35 // Fastjet
36 #include <fastjet/internal/base.hh>
37 #include "fastjet/PseudoJet.hh"
38 #include "fastjet/JetDefinition.hh"
39 #include "fastjet/ClusterSequence.hh"
40 #include "fastjet/Selector.hh"
41 
42 // HepPDTRecord
43 #include "SimGeneral/HepPDTRecord/interface/ParticleDataTable.h"
44 
45 // pile-up
46 #include "SimDataFormats/PileupSummaryInfo/interface/PileupSummaryInfo.h"
47 
48 // associator
49 #include "SimTracker/VertexAssociation/interface/calculateVertexSharedTracks.h"
50 
51 // vertexing
52 #include "RecoVertex/PrimaryVertexProducer/interface/TrackFilterForPVFinding.h"
53 
54 // simulated vertex
55 #include "SimDataFormats/Associations/interface/VertexToTrackingVertexAssociator.h"
56 
57 // DQM
58 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
59 #include "DQMServices/Core/interface/DQMStore.h"
60 
61 #include "DataFormats/Math/interface/deltaR.h"
62 #include "DataFormats/Math/interface/GeantUnits.h"
63 #include "CLHEP/Units/PhysicalConstants.h"
64 
65 // class declaration
66 class Primary4DVertexValidation : public DQMEDAnalyzer {
67  typedef math::XYZTLorentzVector LorentzVector;
68 
69  // auxiliary class holding simulated vertices
70  struct simPrimaryVertex {
71  simPrimaryVertex(double x1, double y1, double z1, double t1)
72  : x(x1),
73  y(y1),
74  z(z1),
75  t(t1),
76  pt(0),
77  ptsq(0),
78  closest_vertex_distance_z(-1.),
79  nGenTrk(0),
80  num_matched_reco_tracks(0),
81  average_match_quality(0.0) {
82  ptot.setPx(0);
83  ptot.setPy(0);
84  ptot.setPz(0);
85  ptot.setE(0);
86  p4 = LorentzVector(0, 0, 0, 0);
87  r = sqrt(x * x + y * y);
88  };
89  double x, y, z, r, t;
90  HepMC::FourVector ptot;
91  LorentzVector p4;
92  double pt;
93  double ptsq;
94  double closest_vertex_distance_z;
95  int nGenTrk;
96  int num_matched_reco_tracks;
97  float average_match_quality;
98  EncodedEventId eventId;
99  TrackingVertexRef sim_vertex;
100  int OriginalIndex = -1;
101 
102  unsigned int nwosmatch = 0; // number of rec vertices dominated by this sim evt (by wos)
103  unsigned int nwntmatch = 0; // number of rec vertices dominated by this sim evt (by wnt)
104  std::vector<unsigned int> wos_dominated_recv; // list of dominated rec vertices (by wos, size==nwosmatch)
105 
106  std::map<unsigned int, double> wnt; // weighted number of tracks in recvtx (by index)
107  std::map<unsigned int, double> wos; // weight over sigma**2 in recvtx (by index)
108  double sumwos = 0; // sum of wos in any recvtx
109  double sumwnt = 0; // sum of weighted tracks
110  unsigned int rec = NOT_MATCHED; // best match (NO_MATCH if not matched)
111  unsigned int matchQuality = 0; // quality flag
112 
113  void addTrack(unsigned int irecv, double twos, double twt) {
114  sumwnt += twt;
115  if (wnt.find(irecv) == wnt.end()) {
116  wnt[irecv] = twt;
117  } else {
118  wnt[irecv] += twt;
119  }
120 
121  sumwos += twos;
122  if (wos.find(irecv) == wos.end()) {
123  wos[irecv] = twos;
124  } else {
125  wos[irecv] += twos;
126  }
127  }
128  };
129 
130  // auxiliary class holding reconstructed vertices
131  struct recoPrimaryVertex {
132  recoPrimaryVertex(double x1, double y1, double z1)
133  : x(x1),
134  y(y1),
135  z(z1),
136  pt(0),
137  ptsq(0),
138  closest_vertex_distance_z(-1.),
139  nRecoTrk(0),
140  num_matched_sim_tracks(0),
141  ndof(0.),
142  recVtx(nullptr) {
143  r = sqrt(x * x + y * y);
144  };
145  double x, y, z, r;
146  double pt;
147  double ptsq;
148  double closest_vertex_distance_z;
149  int nRecoTrk;
150  int num_matched_sim_tracks;
151  double ndof;
152  const reco::Vertex* recVtx;
153  reco::VertexBaseRef recVtxRef;
154  int OriginalIndex = -1;
155 
156  std::map<unsigned int, double> wos; // sim event -> wos
157  std::map<unsigned int, double> wnt; // sim event -> weighted number of truth matched tracks
158  unsigned int wosmatch = NOT_MATCHED; // index of the sim event providing the largest contribution to wos
159  unsigned int wntmatch = NOT_MATCHED; // index of the sim event providing the highest number of tracks
160  double sumwos = 0; // total sum of wos of all truth matched tracks
161  double sumwnt = 0; // total weighted number of truth matchted tracks
162  double maxwos = 0; // largest wos sum from one sim event (wosmatch)
163  double maxwnt = 0; // largest weighted number of tracks from one sim event (ntmatch)
164  unsigned int maxwosnt = 0; // number of tracks from the sim event with highest wos
165  unsigned int sim = NOT_MATCHED; // best match (NO_MATCH if not matched)
166  unsigned int matchQuality = 0; // quality flag
167 
168  bool is_real() { return (matchQuality > 0) && (matchQuality < 99); }
169 
170  bool is_fake() { return (matchQuality <= 0) || (matchQuality >= 99); }
171 
172  bool is_signal() { return (sim == 0); }
173 
174  int split_from() {
175  if (is_real())
176  return -1;
177  if ((maxwos > 0) && (maxwos > 0.3 * sumwos))
178  return wosmatch;
179  return -1;
180  }
181  bool other_fake() { return (is_fake() && (split_from() < 0)); }
182 
183  void addTrack(unsigned int iev, double twos, double wt) {
184  sumwnt += wt;
185  if (wnt.find(iev) == wnt.end()) {
186  wnt[iev] = wt;
187  } else {
188  wnt[iev] += wt;
189  }
190 
191  sumwos += twos;
192  if (wos.find(iev) == wos.end()) {
193  wos[iev] = twos;
194  } else {
195  wos[iev] += twos;
196  }
197  }
198  };
199 
200 public:
201  explicit Primary4DVertexValidation(const edm::ParameterSet&);
202  ~Primary4DVertexValidation() override;
203 
204  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
205  void analyze(const edm::Event&, const edm::EventSetup&) override;
206  void bookHistograms(DQMStore::IBooker& i, edm::Run const&, edm::EventSetup const&) override;
207 
208 private:
209  void matchReco2Sim(std::vector<recoPrimaryVertex>&,
210  std::vector<simPrimaryVertex>&,
211  const edm::ValueMap<float>&,
212  const edm::ValueMap<float>&,
213  const edm::Handle<reco::BeamSpot>&);
214  bool matchRecoTrack2SimSignal(const reco::TrackBaseRef&);
215  std::pair<const edm::Ref<std::vector<TrackingParticle>>*, int> getMatchedTP(const reco::TrackBaseRef&,
216  const TrackingVertexRef&);
217  double timeFromTrueMass(double, double, double, double);
218  bool select(const reco::Vertex&, int level = 0);
219  void observablesFromJets(const std::vector<reco::Track>&,
220  const std::vector<double>&,
221  const std::vector<int>&,
222  const std::string&,
223  unsigned int&,
224  double&,
225  double&,
226  double&,
227  double&);
228  void isParticle(const reco::TrackBaseRef&,
229  const edm::ValueMap<float>&,
230  const edm::ValueMap<float>&,
231  const edm::ValueMap<float>&,
232  const edm::ValueMap<float>&,
233  const edm::ValueMap<float>&,
234  unsigned int&,
235  bool&,
236  bool&,
237  bool&,
238  bool&);
239  void getWosWnt(const reco::Vertex&,
240  const reco::TrackBaseRef&,
241  const edm::ValueMap<float>&,
242  const edm::ValueMap<float>&,
243  const edm::Handle<reco::BeamSpot>&,
244  double&,
245  double&);
246  std::vector<Primary4DVertexValidation::simPrimaryVertex> getSimPVs(const edm::Handle<TrackingVertexCollection>&);
247  std::vector<Primary4DVertexValidation::recoPrimaryVertex> getRecoPVs(const edm::Handle<edm::View<reco::Vertex>>&);
248 
249  void printMatchedRecoTrackInfo(const reco::Vertex&,
250  const reco::TrackBaseRef&,
251  const TrackingParticleRef&,
252  const unsigned int&);
253  void printSimVtxRecoVtxInfo(const struct Primary4DVertexValidation::simPrimaryVertex&,
254  const struct Primary4DVertexValidation::recoPrimaryVertex&);
255  const bool trkTPSelLV(const TrackingParticle&);
256  const bool trkRecSel(const reco::TrackBase&);
257 
258  // ----------member data ---------------------------
259 
260  const std::string folder_;
261  static constexpr unsigned int NOT_MATCHED = 66666;
262  static constexpr double simUnit_ = 1e9; // sim time in s while reco time in ns
263  static constexpr double c_ = 2.99792458e1; // c in cm/ns
264  static constexpr double mvaL_ = 0.5; // MVA cuts for MVA categories
265  static constexpr double mvaH_ = 0.8;
266  static constexpr double selNdof_ = 4.;
267  static constexpr double maxRank_ = 8.;
268  static constexpr double maxTry_ = 10.;
269  static constexpr double zWosMatchMax_ = 1.;
270  static constexpr double etacutGEN_ = 4.; // |eta| < 4;
271  static constexpr double etacutREC_ = 3.; // |eta| < 3;
272  static constexpr double pTcut_ = 0.7; // PT > 0.7 GeV
273  static constexpr double deltaZcut_ = 0.1; // dz separation 1 mm
274  static constexpr double trackMaxBtlEta_ = 1.5;
275  static constexpr double trackMinEtlEta_ = 1.6;
276  static constexpr double trackMaxEtlEta_ = 3.;
277  static constexpr double tol_ = 1.e-4; // tolerance on reconstructed track time, [ns]
278  static constexpr double minThrSumWnt_ = 0.01; // min threshold for filling histograms with logarithmic scale
279  static constexpr double minThrSumWos_ = 0.1;
280  static constexpr double minThrSumPt_ = 0.01;
281  static constexpr double minThrSumPt2_ = 1.e-3;
282  static constexpr double minThrMetPt_ = 1.e-3;
283  static constexpr double minThrSumPz_ = 1.e-4;
284  static constexpr double rBTL_ = 110.0;
285  static constexpr double zETL_ = 290.0;
286 
287  static constexpr float c_cm_ns = geant_units::operators::convertMmToCm(CLHEP::c_light); // [mm/ns] -> [cm/ns]
288 
289  bool use_only_charged_tracks_;
290  bool optionalPlots_;
291  bool use3dNoTime_;
292 
293  const double minProbHeavy_;
294  const double trackweightTh_;
295  const double mvaTh_;
296  const reco::RecoToSimCollection* r2s_;
297  const reco::SimToRecoCollection* s2r_;
298 
299  edm::EDGetTokenT<reco::TrackCollection> RecTrackToken_;
300 
301  edm::EDGetTokenT<std::vector<PileupSummaryInfo>> vecPileupSummaryInfoToken_;
302 
303  edm::EDGetTokenT<TrackingParticleCollection> trackingParticleCollectionToken_;
304  edm::EDGetTokenT<TrackingVertexCollection> trackingVertexCollectionToken_;
305  edm::EDGetTokenT<reco::SimToRecoCollection> simToRecoAssociationToken_;
306  edm::EDGetTokenT<reco::RecoToSimCollection> recoToSimAssociationToken_;
307  edm::EDGetTokenT<reco::BeamSpot> RecBeamSpotToken_;
308  edm::EDGetTokenT<edm::View<reco::Vertex>> Rec4DVerToken_;
309 
310  edm::EDGetTokenT<edm::ValueMap<int>> trackAssocToken_;
311  edm::EDGetTokenT<edm::ValueMap<float>> pathLengthToken_;
312  edm::EDGetTokenT<edm::ValueMap<float>> momentumToken_;
313  edm::EDGetTokenT<edm::ValueMap<float>> timeToken_;
314 
315  edm::EDGetTokenT<edm::ValueMap<float>> t0PidToken_;
316  edm::EDGetTokenT<edm::ValueMap<float>> sigmat0PidToken_;
317  edm::EDGetTokenT<edm::ValueMap<float>> t0SafePidToken_;
318  edm::EDGetTokenT<edm::ValueMap<float>> sigmat0SafePidToken_;
319  edm::EDGetTokenT<edm::ValueMap<float>> trackMVAQualToken_;
320 
321  edm::EDGetTokenT<edm::ValueMap<float>> tmtdToken_;
322  edm::EDGetTokenT<edm::ValueMap<float>> tofPiToken_;
323  edm::EDGetTokenT<edm::ValueMap<float>> tofKToken_;
324  edm::EDGetTokenT<edm::ValueMap<float>> tofPToken_;
325  edm::EDGetTokenT<edm::ValueMap<float>> probPiToken_;
326  edm::EDGetTokenT<edm::ValueMap<float>> probKToken_;
327  edm::EDGetTokenT<edm::ValueMap<float>> probPToken_;
328  edm::ESGetToken<HepPDT::ParticleDataTable, edm::DefaultRecord> pdtToken_;
329 
330  // histogram declaration
331  MonitorElement* meTrackEffPtTot_;
332  MonitorElement* meTrackMatchedTPEffPtTot_;
333  MonitorElement* meTrackMatchedTPEffPtMtd_;
334  MonitorElement* meTrackEffEtaTot_;
335  MonitorElement* meTrackMatchedTPEffEtaTot_;
336  MonitorElement* meTrackMatchedTPEffEtaMtd_;
337  MonitorElement* meTrackMatchedTPResTot_;
338  MonitorElement* meTrackMatchedTPPullTot_;
339  MonitorElement* meTrackResTot_;
340  MonitorElement* meTrackPullTot_;
341  MonitorElement* meTrackRes_[3];
342  MonitorElement* meTrackPull_[3];
343  MonitorElement* meTrackResMass_[3];
344  MonitorElement* meTrackResMassTrue_[3];
345  MonitorElement* meTrackMatchedTPZposResTot_;
346  MonitorElement* meTrackZposResTot_;
347  MonitorElement* meTrackZposRes_[3];
348  MonitorElement* meTrack3DposRes_[3];
349  MonitorElement* meTimeRes_;
350  MonitorElement* meTimePull_;
351  MonitorElement* meTimeSignalRes_;
352  MonitorElement* meTimeSignalPull_;
353  MonitorElement* mePUvsRealV_;
354  MonitorElement* mePUvsFakeV_;
355  MonitorElement* mePUvsOtherFakeV_;
356  MonitorElement* mePUvsSplitV_;
357  MonitorElement* meMatchQual_;
358  MonitorElement* meDeltaZrealreal_;
359  MonitorElement* meDeltaZfakefake_;
360  MonitorElement* meDeltaZfakereal_;
361  MonitorElement* meDeltaTrealreal_;
362  MonitorElement* meDeltaTfakefake_;
363  MonitorElement* meDeltaTfakereal_;
364  MonitorElement* meRecoPosInSimCollection_;
365  MonitorElement* meRecoPosInRecoOrigCollection_;
366  MonitorElement* meSimPosInSimOrigCollection_;
367  MonitorElement* meRecoPVPosSignal_;
368  MonitorElement* meRecoPVPosSignalNotHighestPt_;
369  MonitorElement* meRecoVtxVsLineDensity_;
370  MonitorElement* meRecVerNumber_;
371  MonitorElement* meRecPVZ_;
372  MonitorElement* meRecPVT_;
373  MonitorElement* meSimVerNumber_;
374  MonitorElement* meSimPVZ_;
375  MonitorElement* meSimPVT_;
376  MonitorElement* meSimPVTvsZ_;
377 
378  MonitorElement* meVtxTrackMult_;
379  MonitorElement* meVtxTrackW_;
380  MonitorElement* meVtxTrackWnt_;
381  MonitorElement* meVtxTrackRecLVMult_;
382  MonitorElement* meVtxTrackRecLVW_;
383  MonitorElement* meVtxTrackRecLVWnt_;
384 
385  MonitorElement* mePUTrackMult_;
386  MonitorElement* mePUTrackRelMult_;
387  MonitorElement* meFakeTrackRelMult_;
388  MonitorElement* mePUTrackSumWnt_;
389  MonitorElement* mePUTrackRelSumWnt_;
390  MonitorElement* mePUTrackSumWos_;
391  MonitorElement* mePUTrackRelSumWos_;
392  MonitorElement* meSecTrackSumWos_;
393  MonitorElement* meSecTrackRelSumWos_;
394  MonitorElement* meFakeTrackRelSumWos_;
395  MonitorElement* mePUTrackWnt_;
396  MonitorElement* mePUTrackSumPt_;
397  MonitorElement* mePUTrackRelSumPt_;
398  MonitorElement* mePUTrackSumPt2_;
399  MonitorElement* mePUTrackRelSumPt2_;
400 
401  MonitorElement* mePUTrackRelMultvsMult_;
402  MonitorElement* meFakeTrackRelMultvsMult_;
403  MonitorElement* mePUTrackRelSumWntvsSumWnt_;
404  MonitorElement* mePUTrackRelSumWosvsSumWos_;
405  MonitorElement* meFakeTrackRelSumWosvsSumWos_;
406  MonitorElement* meSecTrackRelSumWosvsSumWos_;
407  MonitorElement* mePUTrackRelSumPtvsSumPt_;
408  MonitorElement* mePUTrackRelSumPt2vsSumPt2_;
409 
410  MonitorElement* mePUTrackRecLVMult_;
411  MonitorElement* mePUTrackRecLVRelMult_;
412  MonitorElement* meFakeTrackRecLVRelMult_;
413  MonitorElement* mePUTrackRecLVSumWnt_;
414  MonitorElement* mePUTrackRecLVRelSumWnt_;
415  MonitorElement* mePUTrackRecLVSumWos_;
416  MonitorElement* mePUTrackRecLVRelSumWos_;
417  MonitorElement* meSecTrackRecLVSumWos_;
418  MonitorElement* meSecTrackRecLVRelSumWos_;
419  MonitorElement* meFakeTrackRecLVRelSumWos_;
420  MonitorElement* mePUTrackRecLVWnt_;
421  MonitorElement* mePUTrackRecLVSumPt_;
422  MonitorElement* mePUTrackRecLVRelSumPt_;
423  MonitorElement* mePUTrackRecLVSumPt2_;
424  MonitorElement* mePUTrackRecLVRelSumPt2_;
425 
426  MonitorElement* mePUTrackRecLVRelMultvsMult_;
427  MonitorElement* meFakeTrackRecLVRelMultvsMult_;
428  MonitorElement* mePUTrackRecLVRelSumWntvsSumWnt_;
429  MonitorElement* mePUTrackRecLVRelSumWosvsSumWos_;
430  MonitorElement* meSecTrackRecLVRelSumWosvsSumWos_;
431  MonitorElement* meFakeTrackRecLVRelSumWosvsSumWos_;
432  MonitorElement* mePUTrackRecLVRelSumPtvsSumPt_;
433  MonitorElement* mePUTrackRecLVRelSumPt2vsSumPt2_;
434 
435  MonitorElement* meJetsPUMult_;
436  MonitorElement* meJetsPUHt_;
437  MonitorElement* meJetsPUSumPt2_;
438  MonitorElement* meJetsPUMetPt_;
439  MonitorElement* meJetsPUSumPz_;
440  MonitorElement* meJetsPURelMult_;
441  MonitorElement* meJetsPURelHt_;
442  MonitorElement* meJetsPURelSumPt2_;
443  MonitorElement* meJetsFakeRelSumPt2_;
444  MonitorElement* meJetsPURelMetPt_;
445  MonitorElement* meJetsPURelSumPz_;
446 
447  MonitorElement* meJetsPURelMultvsMult_;
448  MonitorElement* meJetsPURelHtvsHt_;
449  MonitorElement* meJetsPURelSumPt2vsSumPt2_;
450  MonitorElement* meJetsFakeRelSumPt2vsSumPt2_;
451  MonitorElement* meJetsPURelMetPtvsMetPt_;
452  MonitorElement* meJetsPURelSumPzvsSumPz_;
453 
454  MonitorElement* meJetsRecLVPUMult_;
455  MonitorElement* meJetsRecLVPUHt_;
456  MonitorElement* meJetsRecLVPUSumPt2_;
457  MonitorElement* meJetsRecLVPUMetPt_;
458  MonitorElement* meJetsRecLVPUSumPz_;
459  MonitorElement* meJetsRecLVPURelMult_;
460  MonitorElement* meJetsRecLVPURelHt_;
461  MonitorElement* meJetsRecLVPURelSumPt2_;
462  MonitorElement* meJetsRecLVFakeRelSumPt2_;
463  MonitorElement* meJetsRecLVPURelMetPt_;
464  MonitorElement* meJetsRecLVPURelSumPz_;
465 
466  MonitorElement* meJetsRecLVPURelMultvsMult_;
467  MonitorElement* meJetsRecLVPURelHtvsHt_;
468  MonitorElement* meJetsRecLVPURelSumPt2vsSumPt2_;
469  MonitorElement* meJetsRecLVFakeRelSumPt2vsSumPt2_;
470  MonitorElement* meJetsRecLVPURelMetPtvsMetPt_;
471  MonitorElement* meJetsRecLVPURelSumPzvsSumPz_;
472 
473  // some tests
474  MonitorElement* meTrackResLowPTot_;
475  MonitorElement* meTrackResHighPTot_;
476  MonitorElement* meTrackPullLowPTot_;
477  MonitorElement* meTrackPullHighPTot_;
478 
479  MonitorElement* meTrackResLowP_[3];
480  MonitorElement* meTrackResHighP_[3];
481  MonitorElement* meTrackPullLowP_[3];
482  MonitorElement* meTrackPullHighP_[3];
483 
484  MonitorElement* meTrackResMassProtons_[3];
485  MonitorElement* meTrackResMassTrueProtons_[3];
486  MonitorElement* meTrackResMassPions_[3];
487  MonitorElement* meTrackResMassTruePions_[3];
488 
489  MonitorElement* meBarrelPIDp_;
490  MonitorElement* meEndcapPIDp_;
491 
492  MonitorElement* meBarrelNoPIDtype_;
493  MonitorElement* meEndcapNoPIDtype_;
494 
495  MonitorElement* meBarrelTruePiNoPID_;
496  MonitorElement* meBarrelTrueKNoPID_;
497  MonitorElement* meBarrelTruePNoPID_;
498  MonitorElement* meEndcapTruePiNoPID_;
499  MonitorElement* meEndcapTrueKNoPID_;
500  MonitorElement* meEndcapTruePNoPID_;
501 
502  MonitorElement* meBarrelTruePiAsPi_;
503  MonitorElement* meBarrelTruePiAsK_;
504  MonitorElement* meBarrelTruePiAsP_;
505  MonitorElement* meEndcapTruePiAsPi_;
506  MonitorElement* meEndcapTruePiAsK_;
507  MonitorElement* meEndcapTruePiAsP_;
508 
509  MonitorElement* meBarrelTrueKAsPi_;
510  MonitorElement* meBarrelTrueKAsK_;
511  MonitorElement* meBarrelTrueKAsP_;
512  MonitorElement* meEndcapTrueKAsPi_;
513  MonitorElement* meEndcapTrueKAsK_;
514  MonitorElement* meEndcapTrueKAsP_;
515 
516  MonitorElement* meBarrelTruePAsPi_;
517  MonitorElement* meBarrelTruePAsK_;
518  MonitorElement* meBarrelTruePAsP_;
519  MonitorElement* meEndcapTruePAsPi_;
520  MonitorElement* meEndcapTruePAsK_;
521  MonitorElement* meEndcapTruePAsP_;
522 };
523 
524 // constructors and destructor
525 Primary4DVertexValidation::Primary4DVertexValidation(const edm::ParameterSet& iConfig)
526  : folder_(iConfig.getParameter<std::string>("folder")),
527  use_only_charged_tracks_(iConfig.getParameter<bool>("useOnlyChargedTracks")),
528  optionalPlots_(iConfig.getUntrackedParameter<bool>("optionalPlots")),
529  use3dNoTime_(iConfig.getParameter<bool>("use3dNoTime")),
530  minProbHeavy_(iConfig.getParameter<double>("minProbHeavy")),
531  trackweightTh_(iConfig.getParameter<double>("trackweightTh")),
532  mvaTh_(iConfig.getParameter<double>("mvaTh")),
533  pdtToken_(esConsumes<HepPDT::ParticleDataTable, edm::DefaultRecord>()) {
534  vecPileupSummaryInfoToken_ = consumes<std::vector<PileupSummaryInfo>>(edm::InputTag(std::string("addPileupInfo")));
535  trackingParticleCollectionToken_ =
536  consumes<TrackingParticleCollection>(iConfig.getParameter<edm::InputTag>("SimTag"));
537  trackingVertexCollectionToken_ = consumes<TrackingVertexCollection>(iConfig.getParameter<edm::InputTag>("SimTag"));
538  simToRecoAssociationToken_ =
539  consumes<reco::SimToRecoCollection>(iConfig.getParameter<edm::InputTag>("TPtoRecoTrackAssoc"));
540  recoToSimAssociationToken_ =
541  consumes<reco::RecoToSimCollection>(iConfig.getParameter<edm::InputTag>("TPtoRecoTrackAssoc"));
542  RecTrackToken_ = consumes<reco::TrackCollection>(iConfig.getParameter<edm::InputTag>("mtdTracks"));
543  RecBeamSpotToken_ = consumes<reco::BeamSpot>(iConfig.getParameter<edm::InputTag>("offlineBS"));
544  Rec4DVerToken_ = consumes<edm::View<reco::Vertex>>(iConfig.getParameter<edm::InputTag>("offline4DPV"));
545  trackAssocToken_ = consumes<edm::ValueMap<int>>(iConfig.getParameter<edm::InputTag>("trackAssocSrc"));
546  pathLengthToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("pathLengthSrc"));
547  momentumToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("momentumSrc"));
548  timeToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("timeSrc"));
549  t0PidToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("t0PID"));
550  sigmat0PidToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("sigmat0PID"));
551  t0SafePidToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("t0SafePID"));
552  sigmat0SafePidToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("sigmat0SafePID"));
553  trackMVAQualToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("trackMVAQual"));
554  tmtdToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("tmtd"));
555  tofPiToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("tofPi"));
556  tofKToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("tofK"));
557  tofPToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("tofP"));
558  probPiToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("probPi"));
559  probKToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("probK"));
560  probPToken_ = consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("probP"));
561 }
562 
563 Primary4DVertexValidation::~Primary4DVertexValidation() {}
564 
565 //
566 // member functions
567 //
568 void Primary4DVertexValidation::bookHistograms(DQMStore::IBooker& ibook,
569  edm::Run const& iRun,
570  edm::EventSetup const& iSetup) {
571  ibook.setCurrentFolder(folder_);
572  // --- histograms booking
573  meTrackEffPtTot_ = ibook.book1D("EffPtTot", "Pt of tracks associated to LV; track pt [GeV] ", 110, 0., 11.);
574  meTrackEffEtaTot_ = ibook.book1D("EffEtaTot", "Eta of tracks associated to LV; track eta ", 66, 0., 3.3);
575  meTrackMatchedTPEffPtTot_ =
576  ibook.book1D("MatchedTPEffPtTot", "Pt of tracks associated to LV matched to TP; track pt [GeV] ", 110, 0., 11.);
577  meTrackMatchedTPEffPtMtd_ = ibook.book1D(
578  "MatchedTPEffPtMtd", "Pt of tracks associated to LV matched to TP with time; track pt [GeV] ", 110, 0., 11.);
579  meTrackMatchedTPEffEtaTot_ =
580  ibook.book1D("MatchedTPEffEtaTot", "Eta of tracks associated to LV matched to TP; track eta ", 66, 0., 3.3);
581  meTrackMatchedTPEffEtaMtd_ = ibook.book1D(
582  "MatchedTPEffEtaMtd", "Eta of tracks associated to LV matched to TP with time; track eta ", 66, 0., 3.3);
583  meTrackMatchedTPResTot_ =
584  ibook.book1D("MatchedTPTrackRes",
585  "t_{rec} - t_{sim} for tracks associated to LV matched to TP; t_{rec} - t_{sim} [ns] ",
586  120,
587  -0.15,
588  0.15);
589  meTrackResTot_ = ibook.book1D("TrackRes", "t_{rec} - t_{sim} for tracks; t_{rec} - t_{sim} [ns] ", 120, -0.15, 0.15);
590  meTrackRes_[0] = ibook.book1D(
591  "TrackRes-LowMVA", "t_{rec} - t_{sim} for tracks with MVA < 0.5; t_{rec} - t_{sim} [ns] ", 100, -1., 1.);
592  meTrackRes_[1] = ibook.book1D(
593  "TrackRes-MediumMVA", "t_{rec} - t_{sim} for tracks with 0.5 < MVA < 0.8; t_{rec} - t_{sim} [ns] ", 100, -1., 1.);
594  meTrackRes_[2] = ibook.book1D(
595  "TrackRes-HighMVA", "t_{rec} - t_{sim} for tracks with MVA > 0.8; t_{rec} - t_{sim} [ns] ", 100, -1., 1.);
596  if (optionalPlots_) {
597  meTrackResMass_[0] = ibook.book1D(
598  "TrackResMass-LowMVA", "t_{rec} - t_{est} for tracks with MVA < 0.5; t_{rec} - t_{est} [ns] ", 100, -1., 1.);
599  meTrackResMass_[1] = ibook.book1D("TrackResMass-MediumMVA",
600  "t_{rec} - t_{est} for tracks with 0.5 < MVA < 0.8; t_{rec} - t_{est} [ns] ",
601  100,
602  -1.,
603  1.);
604  meTrackResMass_[2] = ibook.book1D(
605  "TrackResMass-HighMVA", "t_{rec} - t_{est} for tracks with MVA > 0.8; t_{rec} - t_{est} [ns] ", 100, -1., 1.);
606  meTrackResMassTrue_[0] = ibook.book1D(
607  "TrackResMassTrue-LowMVA", "t_{est} - t_{sim} for tracks with MVA < 0.5; t_{est} - t_{sim} [ns] ", 100, -1., 1.);
608  meTrackResMassTrue_[1] = ibook.book1D("TrackResMassTrue-MediumMVA",
609  "t_{est} - t_{sim} for tracks with 0.5 < MVA < 0.8; t_{est} - t_{sim} [ns] ",
610  100,
611  -1.,
612  1.);
613  meTrackResMassTrue_[2] = ibook.book1D("TrackResMassTrue-HighMVA",
614  "t_{est} - t_{sim} for tracks with MVA > 0.8; t_{est} - t_{sim} [ns] ",
615  100,
616  -1.,
617  1.);
618  }
619  meTrackMatchedTPPullTot_ = ibook.book1D(
620  "MatchedTPTrackPull", "Pull for tracks associated to LV matched to TP; (t_{rec}-t_{sim})/#sigma_{t}", 50, -5., 5.);
621  meTrackPullTot_ = ibook.book1D("TrackPull", "Pull for tracks; (t_{rec}-t_{sim})/#sigma_{t}", 100, -10., 10.);
622  meTrackPull_[0] =
623  ibook.book1D("TrackPull-LowMVA", "Pull for tracks with MVA < 0.5; (t_{rec}-t_{sim})/#sigma_{t}", 100, -10., 10.);
624  meTrackPull_[1] = ibook.book1D(
625  "TrackPull-MediumMVA", "Pull for tracks with 0.5 < MVA < 0.8; (t_{rec}-t_{sim})/#sigma_{t}", 100, -10., 10.);
626  meTrackPull_[2] =
627  ibook.book1D("TrackPull-HighMVA", "Pull for tracks with MVA > 0.8; (t_{rec}-t_{sim})/#sigma_{t}", 100, -10., 10.);
628  meTrackMatchedTPZposResTot_ =
629  ibook.book1D("MatchedTPTrackZposResTot",
630  "Z_{PCA} - Z_{sim} for tracks associated to LV matched to TP;Z_{PCA} - Z_{sim} [cm] ",
631  50,
632  -0.1,
633  0.1);
634  meTrackZposResTot_ =
635  ibook.book1D("TrackZposResTot", "Z_{PCA} - Z_{sim} for tracks;Z_{PCA} - Z_{sim} [cm] ", 50, -0.5, 0.5);
636  meTrackZposRes_[0] = ibook.book1D(
637  "TrackZposRes-LowMVA", "Z_{PCA} - Z_{sim} for tracks with MVA < 0.5;Z_{PCA} - Z_{sim} [cm] ", 50, -0.5, 0.5);
638  meTrackZposRes_[1] = ibook.book1D("TrackZposRes-MediumMVA",
639  "Z_{PCA} - Z_{sim} for tracks with 0.5 < MVA < 0.8 ;Z_{PCA} - Z_{sim} [cm] ",
640  50,
641  -0.5,
642  0.5);
643  meTrackZposRes_[2] = ibook.book1D(
644  "TrackZposRes-HighMVA", "Z_{PCA} - Z_{sim} for tracks with MVA > 0.8 ;Z_{PCA} - Z_{sim} [cm] ", 50, -0.5, 0.5);
645  meTrack3DposRes_[0] =
646  ibook.book1D("Track3DposRes-LowMVA",
647  "3dPos_{PCA} - 3dPos_{sim} for tracks with MVA < 0.5 ;3dPos_{PCA} - 3dPos_{sim} [cm] ",
648  50,
649  -0.5,
650  0.5);
651  meTrack3DposRes_[1] =
652  ibook.book1D("Track3DposRes-MediumMVA",
653  "3dPos_{PCA} - 3dPos_{sim} for tracks with 0.5 < MVA < 0.8 ;3dPos_{PCA} - 3dPos_{sim} [cm] ",
654  50,
655  -0.5,
656  0.5);
657  meTrack3DposRes_[2] =
658  ibook.book1D("Track3DposRes-HighMVA",
659  "3dPos_{PCA} - 3dPos_{sim} for tracks with MVA > 0.8;3dPos_{PCA} - 3dPos_{sim} [cm] ",
660  50,
661  -0.5,
662  0.5);
663  meTimeRes_ = ibook.book1D("TimeRes", "t_{rec} - t_{sim} ;t_{rec} - t_{sim} [ns] ", 100, -0.2, 0.2);
664  meTimePull_ = ibook.book1D("TimePull", "Pull; t_{rec} - t_{sim}/#sigma_{t rec}", 100, -10., 10.);
665  meTimeSignalRes_ =
666  ibook.book1D("TimeSignalRes", "t_{rec} - t_{sim} for signal ;t_{rec} - t_{sim} [ns] ", 50, -0.1, 0.1);
667  meTimeSignalPull_ =
668  ibook.book1D("TimeSignalPull", "Pull for signal; t_{rec} - t_{sim}/#sigma_{t rec}", 100, -10., 10.);
669  mePUvsRealV_ =
670  ibook.bookProfile("PUvsReal", "#PU vertices vs #real matched vertices;#PU;#real ", 100, 0, 300, 100, 0, 200);
671  mePUvsFakeV_ =
672  ibook.bookProfile("PUvsFake", "#PU vertices vs #fake matched vertices;#PU;#fake ", 100, 0, 300, 100, 0, 20);
673  mePUvsOtherFakeV_ = ibook.bookProfile(
674  "PUvsOtherFake", "#PU vertices vs #other fake matched vertices;#PU;#other fake ", 100, 0, 300, 100, 0, 20);
675  mePUvsSplitV_ =
676  ibook.bookProfile("PUvsSplit", "#PU vertices vs #split matched vertices;#PU;#split ", 100, 0, 300, 100, 0, 20);
677  meMatchQual_ = ibook.book1D("MatchQuality", "RECO-SIM vertex match quality; ", 8, 0, 8.);
678  meDeltaZrealreal_ = ibook.book1D("DeltaZrealreal", "#Delta Z real-real; |#Delta Z (r-r)| [cm]", 100, 0, 0.5);
679  meDeltaZfakefake_ = ibook.book1D("DeltaZfakefake", "#Delta Z fake-fake; |#Delta Z (f-f)| [cm]", 100, 0, 0.5);
680  meDeltaZfakereal_ = ibook.book1D("DeltaZfakereal", "#Delta Z fake-real; |#Delta Z (f-r)| [cm]", 100, 0, 0.5);
681  meDeltaTrealreal_ = ibook.book1D("DeltaTrealreal", "#Delta T real-real; |#Delta T (r-r)| [sigma]", 60, 0., 30.);
682  meDeltaTfakefake_ = ibook.book1D("DeltaTfakefake", "#Delta T fake-fake; |#Delta T (f-f)| [sigma]", 60, 0., 30.);
683  meDeltaTfakereal_ = ibook.book1D("DeltaTfakereal", "#Delta T fake-real; |#Delta T (f-r)| [sigma]", 60, 0., 30.);
684  if (optionalPlots_) {
685  meRecoPosInSimCollection_ = ibook.book1D(
686  "RecoPosInSimCollection", "Sim signal vertex index associated to Reco signal vertex; Sim PV index", 200, 0, 200);
687  meRecoPosInRecoOrigCollection_ =
688  ibook.book1D("RecoPosInRecoOrigCollection", "Reco signal index in OrigCollection; Reco index", 200, 0, 200);
689  meSimPosInSimOrigCollection_ =
690  ibook.book1D("SimPosInSimOrigCollection", "Sim signal index in OrigCollection; Sim index", 200, 0, 200);
691  }
692  meRecoPVPosSignal_ =
693  ibook.book1D("RecoPVPosSignal", "Position in reco collection of PV associated to sim signal", 20, 0, 20);
694  meRecoPVPosSignalNotHighestPt_ =
695  ibook.book1D("RecoPVPosSignalNotHighestPt",
696  "Position in reco collection of PV associated to sim signal not highest Pt",
697  20,
698  0,
699  20);
700  meRecVerNumber_ = ibook.book1D("RecVerNumber", "RECO Vertex Number: Number of vertices", 50, 0, 250);
701  meSimVerNumber_ = ibook.book1D("SimVerNumber", "SIM Vertex Number: Number of vertices", 50, 0, 250);
702  meRecPVZ_ = ibook.book1D("recPVZ", "#Rec vertices/10 mm", 30, -15., 15.);
703  meRecPVT_ = ibook.book1D("recPVT", "#Rec vertices/50 ps", 30, -0.75, 0.75);
704  meSimPVZ_ = ibook.book1D("simPVZ", "#Sim vertices/10 mm", 30, -15., 15.);
705  meSimPVT_ = ibook.book1D("simPVT", "#Sim vertices/50 ps", 30, -0.75, 0.75);
706  meSimPVTvsZ_ = ibook.bookProfile("simPVTvsZ", "PV Time vs Z", 30, -15., 15., 30, -0.75, 0.75);
707 
708  meVtxTrackMult_ = ibook.book1D("VtxTrackMult", "Log10(Vertex track multiplicity)", 80, 0.5, 2.5);
709  meVtxTrackW_ = ibook.book1D("VtxTrackW", "Vertex track weight (all)", 50, 0., 1.);
710  meVtxTrackWnt_ = ibook.book1D("VtxTrackWnt", "Vertex track Wnt", 50, 0., 1.);
711  meVtxTrackRecLVMult_ =
712  ibook.book1D("VtxTrackRecLVMult", "Log10(Vertex track multiplicity) for matched LV", 80, 0.5, 2.5);
713  meVtxTrackRecLVW_ = ibook.book1D("VtxTrackRecLVW", "Vertex track weight for matched LV (all)", 50, 0., 1.);
714  meVtxTrackRecLVWnt_ = ibook.book1D("VtxTrackRecLVWnt", "Vertex track Wnt for matched LV", 50, 0., 1.);
715 
716  mePUTrackRelMult_ = ibook.book1D(
717  "PUTrackRelMult", "Relative multiplicity of PU tracks for matched vertices; #PUTrks/#Trks", 50, 0., 1.);
718  meFakeTrackRelMult_ = ibook.book1D(
719  "FakeTrackRelMult", "Relative multiplicity of fake tracks for matched vertices; #fakeTrks/#Trks", 50, 0., 1.);
720  mePUTrackRelSumWnt_ =
721  ibook.book1D("PUTrackRelSumWnt",
722  "Relative Sum of wnt of PU tracks for matched vertices; PUSumW*min(Pt, 1.)/SumW*min(Pt, 1.)",
723  50,
724  0.,
725  1.);
726  mePUTrackRelSumWos_ = ibook.book1D(
727  "PUTrackRelSumWos", "Relative Sum of wos of PU tracks for matched vertices; PUSumWos/SumWos", 50, 0., 1.);
728  meSecTrackRelSumWos_ =
729  ibook.book1D("SecTrackRelSumWos",
730  "Relative Sum of wos of tracks from secondary vtx for matched vertices; SecSumWos/SumWos",
731  50,
732  0.,
733  1.);
734  meFakeTrackRelSumWos_ = ibook.book1D(
735  "FakeTrackRelSumWos", "Relative Sum of wos of fake tracks for matched vertices; FakeSumWos/SumWos", 50, 0., 1.);
736  mePUTrackRelSumPt_ = ibook.book1D(
737  "PUTrackRelSumPt", "Relative Sum of Pt of PU tracks for matched vertices; PUSumPt/SumPt", 50, 0., 1.);
738  mePUTrackRelSumPt2_ = ibook.book1D(
739  "PUTrackRelSumPt2", "Relative Sum of Pt2 for PU tracks for matched vertices; PUSumPt2/SumPt2", 50, 0., 1.);
740  mePUTrackRecLVRelMult_ = ibook.book1D(
741  "PUTrackRecLVRelMult", "Relative multiplicity of PU tracks for matched LV; #PUTrks/#Trks", 50, 0., 1.);
742  meFakeTrackRecLVRelMult_ = ibook.book1D(
743  "FakeTrackRecLVRelMult", "Relative multiplicity of fake tracks for matched LV; #FakeTrks/#Trks", 50, 0., 1.);
744  mePUTrackRecLVRelSumWnt_ =
745  ibook.book1D("PUTrackRecLVRelSumWnt",
746  "Relative Sum of Wnt of PU tracks for matched LV; PUSumW*min(Pt, 1.)/SumW*min(Pt, 1.)",
747  50,
748  0.,
749  1.);
750  mePUTrackRecLVRelSumWos_ = ibook.book1D(
751  "PUTrackRecLVRelSumWos", "Relative Sum of Wos of PU tracks for matched LV; PUSumWos/SumWos", 50, 0., 1.);
752  meSecTrackRecLVRelSumWos_ =
753  ibook.book1D("SecTrackRecLVRelSumWos",
754  "Relative Sum of wos of tracks from secondary vtx for matched LV; SecSumWos/SumWos",
755  50,
756  0.,
757  1.);
758  meFakeTrackRecLVRelSumWos_ = ibook.book1D(
759  "FakeTrackRecLVRelSumWos", "Relative Sum of wos of fake tracks for matched LV; FakeSumWos/SumWos", 50, 0., 1.);
760  mePUTrackRecLVRelSumPt_ =
761  ibook.book1D("PUTrackRecLVRelSumPt", "Relative Sum of Pt of PU tracks for matched LV; PUSumPt/SumPt", 50, 0., 1.);
762  mePUTrackRecLVRelSumPt2_ = ibook.book1D(
763  "PUTrackRecLVRelSumPt2", "Relative Sum of Pt2 of PU tracks for matched LV; PUSumPt2/SumPt2", 50, 0., 1.);
764 
765  if (optionalPlots_) {
766  mePUTrackMult_ = ibook.book1D("PUTrackMult", "Number of PU tracks for matched vertices; #PUTrks", 50, 0., 100.);
767  mePUTrackWnt_ = ibook.book1D("PUTrackWnt", "Wnt of PU tracks for matched vertices; PUTrkW*min(Pt, 1.)", 50, 0., 1.);
768  mePUTrackSumWnt_ = ibook.book1D(
769  "PUTrackSumWnt", "Sum of wnt of PU tracks for matched vertices; log10(PUSumW*min(Pt, 1.))", 50, -2., 3.);
770  mePUTrackSumWos_ =
771  ibook.book1D("PUTrackSumWos", "Sum of wos of PU tracks for matched vertices; log10(PUSumWos)", 50, -1., 7.);
772  meSecTrackSumWos_ = ibook.book1D(
773  "SecTrackSumWos", "Sum of wos of tracks from secondary vtx for matched vertices; log10(SecSumWos)", 50, -1., 7.);
774  mePUTrackSumPt_ =
775  ibook.book1D("PUTrackSumPt", "Sum of Pt of PU tracks for matched vertices; log10(PUSumPt)", 50, -2., 3.);
776  mePUTrackSumPt2_ =
777  ibook.book1D("PUTrackSumPt2", "Sum of Pt2 of PU tracks for matched vertices; log10(PUSumPt2)", 50, -3., 3.);
778 
779  mePUTrackRelMultvsMult_ =
780  ibook.bookProfile("PUTrackRelMultvsMult",
781  "Relative PU multiplicity vs Number of tracks for matched vertices; #Trks; #PUTrks/#Trks",
782  50,
783  0.,
784  200.,
785  0.,
786  1.,
787  "s");
788  meFakeTrackRelMultvsMult_ = ibook.bookProfile(
789  "FakeTrackRelMultvsMult",
790  "Relative multiplicity of fake tracks vs Number of tracks for matched vertices; #Trks; #FakeTrks/#Trks",
791  50,
792  0.,
793  200.,
794  0.,
795  1.,
796  "s");
797  mePUTrackRelSumWntvsSumWnt_ =
798  ibook.bookProfile("PUTrackRelSumWntvsSumWnt",
799  "Relative PU Sum of Wnt vs Sum of Wnt of tracks for matched vertices; log10(SumW*min(Pt, "
800  "1.)); PUSumW*min(Pt, 1.)/SumW*min(Pt, 1.)",
801  50,
802  0.,
803  2.5,
804  0.,
805  1.,
806  "s");
807  mePUTrackRelSumWosvsSumWos_ = ibook.bookProfile(
808  "PUTrackRelSumWosvsSumWos",
809  "Relative PU Sum of Wos vs Sum of Wos of tracks for matched vertices; log10(SumWos); PUSumWos/SumWos",
810  50,
811  2.5,
812  7.,
813  0.,
814  1.,
815  "s");
816  meSecTrackRelSumWosvsSumWos_ = ibook.bookProfile("SecTrackRelSumWosvsSumWos",
817  "Relative Sum of Wos of tracks from secondary vtx vs Sum of Wos "
818  "of tracks for matched vertices; log10(SumWos); SecSumWos/SumWos",
819  50,
820  2.,
821  7.,
822  0.,
823  1.,
824  "s");
825  meFakeTrackRelSumWosvsSumWos_ = ibook.bookProfile("FakeTrackRelSumWosvsSumWos",
826  "Relative Sum of Wos of fake tracks vs Sum of Wos of tracks for "
827  "matched vertices; log10(SumWos); FakeSumWos/SumWos",
828  50,
829  2.5,
830  7.5,
831  0.,
832  1.,
833  "s");
834  mePUTrackRelSumPtvsSumPt_ = ibook.bookProfile(
835  "PUTrackRelSumPtvsSumPt",
836  "Relative PU Sum of Pt vs Sum of Pt of tracks for matched vertices; log10(SumPt); PUSumPt/SumPt",
837  50,
838  0.,
839  3.,
840  0.,
841  1.,
842  "s");
843  mePUTrackRelSumPt2vsSumPt2_ = ibook.bookProfile(
844  "PUTrackRelSumPt2vsSumPt2",
845  "Relative PU Sum of Pt2 vs Sum of Pt2 of tracks for matched vertices; log10(SumPt2); PUSumPt2/SumPt2",
846  50,
847  0.,
848  4.,
849  0.,
850  1.,
851  "s");
852 
853  mePUTrackRecLVMult_ = ibook.book1D("PUTrackRecLVMult", "Number of PU tracks for matched LV; #PUTrks", 50, 0., 100.);
854  mePUTrackRecLVWnt_ =
855  ibook.book1D("PUTrackRecLVWnt", "Wnt of PU tracks for matched LV; PUTrkW*min(Pt, 1.)", 50, 0., 1.);
856  mePUTrackRecLVSumWnt_ = ibook.book1D(
857  "PUTrackRecLVSumWnt", "Sum of wnt of PU tracks for matched LV; log10(PUSumW*min(Pt, 1.))", 50, -2., 3.);
858  mePUTrackRecLVSumWos_ =
859  ibook.book1D("PUTrackRecLVSumWos", "Sum of wos of PU tracks for matched LV; log10(PUSumWos)", 50, -1., 7.);
860  meSecTrackRecLVSumWos_ = ibook.book1D(
861  "SecTrackRecLVSumWos", "Sum of wos of tracks from secondary vtx for matched LV; log10(SecSumWos)", 50, -1., 7.);
862  mePUTrackRecLVSumPt_ =
863  ibook.book1D("PUTrackRecLVSumPt", "Sum of Pt of PU tracks for matched LV; log10(PUSumPt)", 50, -2., 3.);
864  mePUTrackRecLVSumPt2_ =
865  ibook.book1D("PUTrackRecLVSumPt2", "Sum of Pt2 of PU tracks for matched LV; log10(PUSumPt2)", 50, -3., 3.);
866 
867  mePUTrackRecLVRelMultvsMult_ =
868  ibook.bookProfile("PUTrackRecLVRelMultvsMult",
869  "Relative PU multiplicity vs Number of tracks for matched LV; #Trks; #PUTrks/#Trks",
870  50,
871  0.,
872  200.,
873  0.,
874  1.,
875  "s");
876  meFakeTrackRecLVRelMultvsMult_ = ibook.bookProfile(
877  "FakeTrackRecLVRelMultvsMult",
878  "Relative multiplicity of fake tracks vs Number of tracks for matched LV; #Trks; #FakeTrks/#Trks",
879  50,
880  0.,
881  200.,
882  0.,
883  1.,
884  "s");
885  mePUTrackRecLVRelSumWntvsSumWnt_ =
886  ibook.bookProfile("PUTrackRecLVRelSumWntvsSumWnt",
887  "Relative PU Sum of Wnt vs Sum of Wnt of tracks for matched LV; log10(SumW*min(Pt, 1.)); "
888  "PUSumW*min(Pt, 1.)/SumW*min(Pt, 1.)",
889  50,
890  1.,
891  2.3,
892  0.,
893  1.,
894  "s");
895  mePUTrackRecLVRelSumWosvsSumWos_ = ibook.bookProfile(
896  "PUTrackRecLVRelSumWosvsSumWos",
897  "Relative PU Sum of Wos vs Sum of Wos of tracks for matched vertices; log10(SumWos); PUSumWos/SumWos",
898  50,
899  5.5,
900  7.,
901  0.,
902  1.,
903  "s");
904  meSecTrackRecLVRelSumWosvsSumWos_ =
905  ibook.bookProfile("SecTrackRecLVRelSumWosvsSumWos",
906  "Relative Sum of Wos of tracks from secondary vtx vs Sum of Wos of tracks for matched "
907  "vertices; log10(SumWos); SecSumWos/SumWos",
908  50,
909  5.5,
910  7.,
911  0.,
912  1.,
913  "s");
914  meFakeTrackRecLVRelSumWosvsSumWos_ = ibook.bookProfile("FakeTrackRecLVRelSumWosvsSumWos",
915  "Relative Sum of Wos of fake tracks vs Sum of Wos of tracks "
916  "for matched vertices; log10(SumWos); FakeSumWos/SumWos",
917  50,
918  5.5,
919  7.,
920  0.,
921  1.,
922  "s");
923  mePUTrackRecLVRelSumPtvsSumPt_ =
924  ibook.bookProfile("PUTrackRecLVRelSumPtvsSumPt",
925  "Relative PU Sum of Pt vs Sum of Pt of tracks for matched LV; log10(SumPt); PUSumPt/SumPt",
926  50,
927  1.4,
928  3.,
929  0.,
930  1.,
931  "s");
932  mePUTrackRecLVRelSumPt2vsSumPt2_ =
933  ibook.bookProfile("PUTrackRecLVRelSumPt2vsSumPt2",
934  "Relative PU Sum of Pt2 vs Sum of tracks for matched LV; log10(SumPt2); PUSumPt2/SumPt2",
935  50,
936  2.,
937  4.,
938  0.,
939  1.,
940  "s");
941  }
942 
943  meJetsPURelMult_ =
944  ibook.book1D("JetsPURelMult",
945  "Relative contribution of PU to jet multiplicity for matched vertices; (#Jets-#JetsNoPU)/#Jets",
946  50,
947  0.,
948  1.);
949  meJetsPURelHt_ =
950  ibook.book1D("JetsPURelHt",
951  "Relative contribution of PU to scalar sum of Et of jets for matched vertices; (Ht-HtNoPU)/HT",
952  50,
953  0.,
954  1.);
955  meJetsPURelSumPt2_ =
956  ibook.book1D("JetsPURelSumPt2",
957  "Relative contribution of PU to sum of Pt2 of jets for matched vertices; (SumPt2-SumPt2NoPU)/SumPt2",
958  50,
959  0.,
960  1.);
961  meJetsFakeRelSumPt2_ = ibook.book1D(
962  "JetsFakeRelSumPt2",
963  "Relative contribution of fake tracks to sum of Pt2 of jets for matched vertices; (SumPt2-SumPt2NoFake)/SumPt2",
964  50,
965  0.,
966  1.);
967  meJetsPURelMetPt_ =
968  ibook.book1D("JetsPURelMetPt",
969  "Relative contribution of PU to Missing Transverse Energy for matched vertices; (Met-MetNoPU)/Met",
970  50,
971  -1.,
972  1.);
973  meJetsPURelSumPz_ =
974  ibook.book1D("JetsPURelSumPz",
975  "Relative contribution of PU to sum of Pz of jets for matched vertices; (SumPz-SumPzNoPU)/SumPz",
976  50,
977  -1.,
978  1.);
979 
980  meJetsRecLVPURelMult_ =
981  ibook.book1D("JetsRecLVPURelMult",
982  "Relative contribution of PU to jet multiplicity for matched LV; (#Jets-#JetsNoPU)/#Jets",
983  50,
984  0.,
985  1.);
986  meJetsRecLVPURelHt_ =
987  ibook.book1D("JetsRecLVPURelHt",
988  "Relative contribution of PU to scalar sum of Et of jets for matched LV; (Ht-HtNoPU)/HT",
989  50,
990  0.,
991  1.);
992  meJetsRecLVPURelSumPt2_ =
993  ibook.book1D("JetsRecLVPURelSumPt2",
994  "Relative contribution of PU to sum of Pt2 of jets for matched LV; (SumPt2-SumPt2NoPU)/SumPt2",
995  50,
996  0.,
997  1.);
998  meJetsRecLVFakeRelSumPt2_ = ibook.book1D(
999  "JetsRecLVFakeRelSumPt2",
1000  "Relative contribution of fake tracks to sum of Pt2 of jets for matched LV; (SumPt2-SumPt2NoFake)/SumPt2",
1001  50,
1002  0.,
1003  1.);
1004  meJetsRecLVPURelMetPt_ =
1005  ibook.book1D("JetsRecLVPURelMetPt",
1006  "Relative contribution of PU to Missing Transverse Energy for matched LV; (Met-MetNoPU)/Met",
1007  50,
1008  -1.,
1009  1.);
1010  meJetsRecLVPURelSumPz_ =
1011  ibook.book1D("JetsRecLVPURelSumPz",
1012  "Relative contribution of PU to sum of Pz of jets for matched LV; (SumPz-SumPzNoPU)/SumPz",
1013  50,
1014  -1.,
1015  1.);
1016 
1017  if (optionalPlots_) {
1018  meJetsPUMult_ = ibook.book1D(
1019  "JetsPUMult", "Contribution of PU to jet multiplicity for matched vertices; #Jets-#JetsNoPU", 50, 0., 100.);
1020  meJetsPUHt_ = ibook.book1D("JetsPUHt",
1021  "Contribution of PU to scalar sum of Et of jets for matched vertices; log10(Ht-HtNoPU)",
1022  50,
1023  -2.,
1024  3.);
1025  meJetsPUSumPt2_ =
1026  ibook.book1D("JetsPUSumPt2",
1027  "Contribution of PU to sum of Pt2 of jets for matched vertices; log10(sumPt2-SumPt2NoPU)",
1028  50,
1029  -3.,
1030  3.);
1031  meJetsPUMetPt_ =
1032  ibook.book1D("JetsPUMetPt",
1033  "Contribution of PU to Missing Transverse Energy for matched vertices; log10(Met-MetNoPU)",
1034  50,
1035  -3.,
1036  2.);
1037  meJetsPUSumPz_ =
1038  ibook.book1D("JetsPUSumPz",
1039  "Contribution of PU to sum of Pz of jets for matched vertices; log10(abs(SumPz-SumPzNoPU))",
1040  50,
1041  -4.,
1042  3.);
1043 
1044  meJetsPURelMultvsMult_ = ibook.bookProfile("JetsPURelMultvsMult",
1045  "Relative contribution of PU to jet multiplicity vs number of jets for "
1046  "matched vertices; #Jets; (#Jets-#JetsNoPU)/#Jets",
1047  50,
1048  0.,
1049  120.,
1050  0.,
1051  1.,
1052  "s");
1053  meJetsPURelHtvsHt_ = ibook.bookProfile("JetsPURelHtvsHt",
1054  "Relative contribution of PU to scalar sum of Et of jets vs scalar sum of "
1055  "Et for matched vertices; log10(Ht); (Ht-HtNoPU)/HT",
1056  50,
1057  0.,
1058  3.,
1059  0.,
1060  1.,
1061  "s");
1062  meJetsPURelSumPt2vsSumPt2_ = ibook.bookProfile("JetsPURelSumPt2vsSumPt2",
1063  "Relative contribution of PU to sum of Pt2 of jets vs sum of Pt2 "
1064  "for matched vertices; log10(SumPt2); (SumPt2-SumPt2NoPU)/SumPt2",
1065  50,
1066  -1.,
1067  4.,
1068  0.,
1069  1.,
1070  "s");
1071  meJetsFakeRelSumPt2vsSumPt2_ =
1072  ibook.bookProfile("JetsFakeRelSumPt2vsSumPt2",
1073  "Relative contribution of fake tracks to sum of Pt2 of jets vs sum of Pt2 for matched "
1074  "vertices; log10(SumPt2); (SumPt2-SumPt2NoFake)/SumPt2",
1075  50,
1076  -1.,
1077  4.,
1078  0.,
1079  1.,
1080  "s");
1081  meJetsPURelMetPtvsMetPt_ = ibook.bookProfile("JetsPURelMetPtvsMetPt",
1082  "Relative contribution of PU to Missing Transverse Energy vs MET for "
1083  "matched vertices; log10(Met); (Met-MetNoPU)/Met",
1084  50,
1085  -1.,
1086  2.,
1087  -1.,
1088  1.,
1089  "s");
1090  meJetsPURelSumPzvsSumPz_ = ibook.bookProfile("JetsPURelSumPzvsSumPz",
1091  "Relative contribution of PU to sum of Pz of jets vs Sum of Pz for "
1092  "matched vertices; log10(abs SumPz); (SumPz-SumPzNoPU)/SumPz",
1093  50,
1094  -2.,
1095  3.,
1096  -1.,
1097  1.,
1098  "s");
1099 
1100  meJetsRecLVPUMult_ = ibook.book1D(
1101  "JetsRecLVPUMult", "Contribution of PU to jet multiplicity for matched LV; #Jets-#JetsNoPU", 50, 0., 100.);
1102  meJetsRecLVPUHt_ = ibook.book1D(
1103  "JetsRecLVPUHt", "Contribution of PU to scalar sum of Et of jets for matched LV; log10(Ht-HtNoPU)", 50, -2., 3.);
1104  meJetsRecLVPUSumPt2_ =
1105  ibook.book1D("JetsRecLVPUSumPt2",
1106  "Contribution of PU to sum of Pt2 of jets for matched LV; log10(sumPt2-SumPt2NoPU)",
1107  50,
1108  -3.,
1109  3.);
1110  meJetsRecLVPUMetPt_ =
1111  ibook.book1D("JetsRecLVPUMetPt",
1112  "Contribution of PU to Missing Transverse Energy for matched LV; log10(Met-MetNoPU)",
1113  50,
1114  -3.,
1115  2.);
1116  meJetsRecLVPUSumPz_ =
1117  ibook.book1D("JetsRecLVPUSumPz",
1118  "Contribution of PU to sum of Pz of jets for matched LV; log10(abs(SumPz-SumPzNoPU))",
1119  50,
1120  -4.,
1121  3.);
1122 
1123  meJetsRecLVPURelMultvsMult_ = ibook.bookProfile("JetsRecLVPURelMultvsMult",
1124  "Relative contribution of PU to jet multiplicity vs number of jets "
1125  "for matched vertices; #Jets; (#Jets-#JetsNoPU)/#Jets",
1126  50,
1127  0.,
1128  120.,
1129  0.,
1130  1.,
1131  "s");
1132  meJetsRecLVPURelHtvsHt_ = ibook.bookProfile("JetsRecLVPURelHtvsHt",
1133  "Relative contribution of PU to scalar sum of Et of jets vs scalar sum "
1134  "of Et for matched vertices; log10(Ht); (Ht-HtNoPU)/HT",
1135  50,
1136  1.5,
1137  3.,
1138  0.,
1139  1.,
1140  "s");
1141  meJetsRecLVPURelSumPt2vsSumPt2_ =
1142  ibook.bookProfile("JetsRecLVPURelSumPt2vsSumPt2",
1143  "Relative contribution of PU to sum of Pt2 of jets vs sum of Pt2 for matched vertices; "
1144  "log10(SumPt2); (SumPt2-SumPt2NoPU)/SumPt2",
1145  50,
1146  2.,
1147  5.,
1148  0.,
1149  1.,
1150  "s");
1151  meJetsRecLVFakeRelSumPt2vsSumPt2_ =
1152  ibook.bookProfile("JetsRecLVFakeRelSumPt2vsSumPt2",
1153  "Relative contribution of fake tracks to sum of Pt2 of jets vs sum of Pt2 for matched "
1154  "vertices; log10(SumPt2); (SumPt2-SumPt2NoFake)/SumPt2",
1155  50,
1156  2.,
1157  5.,
1158  0.,
1159  1.,
1160  "s");
1161  meJetsRecLVPURelMetPtvsMetPt_ = ibook.bookProfile("JetsRecLVPURelMetPtvsMetPt",
1162  "Relative contribution of PU to Missing Transverse Energy vs MET "
1163  "for matched vertices; log10(Met); (Met-MetNoPU)/Met",
1164  50,
1165  0.,
1166  2.5,
1167  -1.,
1168  1.,
1169  "s");
1170  meJetsRecLVPURelSumPzvsSumPz_ = ibook.bookProfile("JetsRecLVPURelSumPzvsSumPz",
1171  "Relative contribution of PU to sum of Pz of jets vs Sum of Pz "
1172  "for matched vertices; log10(abs SumPz); (SumPz-SumPzNoPU)/SumPz",
1173  50,
1174  0.5,
1175  3.5,
1176  -1.,
1177  1.,
1178  "s");
1179  }
1180 
1181  // some tests
1182  meTrackResLowPTot_ = ibook.book1D(
1183  "TrackResLowP", "t_{rec} - t_{sim} for tracks with p < 2 GeV; t_{rec} - t_{sim} [ns] ", 70, -0.15, 0.15);
1184  meTrackResLowP_[0] =
1185  ibook.book1D("TrackResLowP-LowMVA",
1186  "t_{rec} - t_{sim} for tracks with MVA < 0.5 and p < 2 GeV; t_{rec} - t_{sim} [ns] ",
1187  100,
1188  -1.,
1189  1.);
1190  meTrackResLowP_[1] =
1191  ibook.book1D("TrackResLowP-MediumMVA",
1192  "t_{rec} - t_{sim} for tracks with 0.5 < MVA < 0.8 and p < 2 GeV; t_{rec} - t_{sim} [ns] ",
1193  100,
1194  -1.,
1195  1.);
1196  meTrackResLowP_[2] =
1197  ibook.book1D("TrackResLowP-HighMVA",
1198  "t_{rec} - t_{sim} for tracks with MVA > 0.8 and p < 2 GeV; t_{rec} - t_{sim} [ns] ",
1199  100,
1200  -1.,
1201  1.);
1202  meTrackResHighPTot_ = ibook.book1D(
1203  "TrackResHighP", "t_{rec} - t_{sim} for tracks with p > 2 GeV; t_{rec} - t_{sim} [ns] ", 70, -0.15, 0.15);
1204  meTrackResHighP_[0] =
1205  ibook.book1D("TrackResHighP-LowMVA",
1206  "t_{rec} - t_{sim} for tracks with MVA < 0.5 and p > 2 GeV; t_{rec} - t_{sim} [ns] ",
1207  100,
1208  -1.,
1209  1.);
1210  meTrackResHighP_[1] =
1211  ibook.book1D("TrackResHighP-MediumMVA",
1212  "t_{rec} - t_{sim} for tracks with 0.5 < MVA < 0.8 and p > 2 GeV; t_{rec} - t_{sim} [ns] ",
1213  100,
1214  -1.,
1215  1.);
1216  meTrackResHighP_[2] =
1217  ibook.book1D("TrackResHighP-HighMVA",
1218  "t_{rec} - t_{sim} for tracks with MVA > 0.8 and p > 2 GeV; t_{rec} - t_{sim} [ns] ",
1219  100,
1220  -1.,
1221  1.);
1222  meTrackPullLowPTot_ =
1223  ibook.book1D("TrackPullLowP", "Pull for tracks with p < 2 GeV; (t_{rec}-t_{sim})/#sigma_{t}", 100, -10., 10.);
1224  meTrackPullLowP_[0] = ibook.book1D("TrackPullLowP-LowMVA",
1225  "Pull for tracks with MVA < 0.5 and p < 2 GeV; (t_{rec}-t_{sim})/#sigma_{t}",
1226  100,
1227  -10.,
1228  10.);
1229  meTrackPullLowP_[1] = ibook.book1D("TrackPullLowP-MediumMVA",
1230  "Pull for tracks with 0.5 < MVA < 0.8 and p < 2 GeV; (t_{rec}-t_{sim})/#sigma_{t}",
1231  100,
1232  -10.,
1233  10.);
1234  meTrackPullLowP_[2] = ibook.book1D("TrackPullLowP-HighMVA",
1235  "Pull for tracks with MVA > 0.8 and p < 2 GeV; (t_{rec}-t_{sim})/#sigma_{t}",
1236  100,
1237  -10.,
1238  10.);
1239  meTrackPullHighPTot_ =
1240  ibook.book1D("TrackPullHighP", "Pull for tracks with p > 2 GeV; (t_{rec}-t_{sim})/#sigma_{t}", 100, -10., 10.);
1241  meTrackPullHighP_[0] = ibook.book1D("TrackPullHighP-LowMVA",
1242  "Pull for tracks with MVA < 0.5 and p > 2 GeV; (t_{rec}-t_{sim})/#sigma_{t}",
1243  100,
1244  -10.,
1245  10.);
1246  meTrackPullHighP_[1] =
1247  ibook.book1D("TrackPullHighP-MediumMVA",
1248  "Pull for tracks with 0.5 < MVA < 0.8 and p > 2 GeV; (t_{rec}-t_{sim})/#sigma_{t}",
1249  100,
1250  -10.,
1251  10.);
1252  meTrackPullHighP_[2] = ibook.book1D("TrackPullHighP-HighMVA",
1253  "Pull for tracks with MVA > 0.8 and p > 2 GeV; (t_{rec}-t_{sim})/#sigma_{t}",
1254  100,
1255  -10.,
1256  10.);
1257  if (optionalPlots_) {
1258  meTrackResMassProtons_[0] =
1259  ibook.book1D("TrackResMass-Protons-LowMVA",
1260  "t_{rec} - t_{est} for proton tracks with MVA < 0.5; t_{rec} - t_{est} [ns] ",
1261  100,
1262  -1.,
1263  1.);
1264  meTrackResMassProtons_[1] =
1265  ibook.book1D("TrackResMass-Protons-MediumMVA",
1266  "t_{rec} - t_{est} for proton tracks with 0.5 < MVA < 0.8; t_{rec} - t_{est} [ns] ",
1267  100,
1268  -1.,
1269  1.);
1270  meTrackResMassProtons_[2] =
1271  ibook.book1D("TrackResMass-Protons-HighMVA",
1272  "t_{rec} - t_{est} for proton tracks with MVA > 0.8; t_{rec} - t_{est} [ns] ",
1273  100,
1274  -1.,
1275  1.);
1276  meTrackResMassTrueProtons_[0] =
1277  ibook.book1D("TrackResMassTrue-Protons-LowMVA",
1278  "t_{est} - t_{sim} for proton tracks with MVA < 0.5; t_{est} - t_{sim} [ns] ",
1279  100,
1280  -1.,
1281  1.);
1282  meTrackResMassTrueProtons_[1] =
1283  ibook.book1D("TrackResMassTrue-Protons-MediumMVA",
1284  "t_{est} - t_{sim} for proton tracks with 0.5 < MVA < 0.8; t_{est} - t_{sim} [ns] ",
1285  100,
1286  -1.,
1287  1.);
1288  meTrackResMassTrueProtons_[2] =
1289  ibook.book1D("TrackResMassTrue-Protons-HighMVA",
1290  "t_{est} - t_{sim} for proton tracks with MVA > 0.8; t_{est} - t_{sim} [ns] ",
1291  100,
1292  -1.,
1293  1.);
1294 
1295  meTrackResMassPions_[0] = ibook.book1D("TrackResMass-Pions-LowMVA",
1296  "t_{rec} - t_{est} for pion tracks with MVA < 0.5; t_{rec} - t_{est} [ns] ",
1297  100,
1298  -1.,
1299  1.);
1300  meTrackResMassPions_[1] =
1301  ibook.book1D("TrackResMass-Pions-MediumMVA",
1302  "t_{rec} - t_{est} for pion tracks with 0.5 < MVA < 0.8; t_{rec} - t_{est} [ns] ",
1303  100,
1304  -1.,
1305  1.);
1306  meTrackResMassPions_[2] = ibook.book1D("TrackResMass-Pions-HighMVA",
1307  "t_{rec} - t_{est} for pion tracks with MVA > 0.8; t_{rec} - t_{est} [ns] ",
1308  100,
1309  -1.,
1310  1.);
1311  meTrackResMassTruePions_[0] =
1312  ibook.book1D("TrackResMassTrue-Pions-LowMVA",
1313  "t_{est} - t_{sim} for pion tracks with MVA < 0.5; t_{est} - t_{sim} [ns] ",
1314  100,
1315  -1.,
1316  1.);
1317  meTrackResMassTruePions_[1] =
1318  ibook.book1D("TrackResMassTrue-Pions-MediumMVA",
1319  "t_{est} - t_{sim} for pion tracks with 0.5 < MVA < 0.8; t_{est} - t_{sim} [ns] ",
1320  100,
1321  -1.,
1322  1.);
1323  meTrackResMassTruePions_[2] =
1324  ibook.book1D("TrackResMassTrue-Pions-HighMVA",
1325  "t_{est} - t_{sim} for pion tracks with MVA > 0.8; t_{est} - t_{sim} [ns] ",
1326  100,
1327  -1.,
1328  1.);
1329  }
1330 
1331  meBarrelPIDp_ = ibook.book1D("BarrelPIDp", "PID track MTD momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1332  meEndcapPIDp_ = ibook.book1D("EndcapPIDp", "PID track with MTD momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1333 
1334  meBarrelNoPIDtype_ = ibook.book1D("BarrelNoPIDtype", "Barrel PID failure category", 4, 0., 4.);
1335  meEndcapNoPIDtype_ = ibook.book1D("EndcapNoPIDtype", "Endcap PID failure category", 4, 0., 4.);
1336 
1337  meBarrelTruePiNoPID_ =
1338  ibook.book1D("BarrelTruePiNoPID", "True pi NoPID momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1339  meBarrelTrueKNoPID_ =
1340  ibook.book1D("BarrelTrueKNoPID", "True k NoPID momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1341  meBarrelTruePNoPID_ =
1342  ibook.book1D("BarrelTruePNoPID", "True p NoPID momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1343  meEndcapTruePiNoPID_ =
1344  ibook.book1D("EndcapTruePiNoPID", "True pi NoPID momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1345  meEndcapTrueKNoPID_ =
1346  ibook.book1D("EndcapTrueKNoPID", "True k NoPID momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1347  meEndcapTruePNoPID_ =
1348  ibook.book1D("EndcapTruePNoPID", "True p NoPID momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1349 
1350  meBarrelTruePiAsPi_ =
1351  ibook.book1D("BarrelTruePiAsPi", "True pi as pi momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1352  meBarrelTruePiAsK_ =
1353  ibook.book1D("BarrelTruePiAsK", "True pi as k momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1354  meBarrelTruePiAsP_ =
1355  ibook.book1D("BarrelTruePiAsP", "True pi as p momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1356  meEndcapTruePiAsPi_ =
1357  ibook.book1D("EndcapTruePiAsPi", "True pi as pi momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1358  meEndcapTruePiAsK_ =
1359  ibook.book1D("EndcapTruePiAsK", "True pi as k momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1360  meEndcapTruePiAsP_ =
1361  ibook.book1D("EndcapTruePiAsP", "True pi as p momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1362 
1363  meBarrelTrueKAsPi_ =
1364  ibook.book1D("BarrelTrueKAsPi", "True k as pi momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1365  meBarrelTrueKAsK_ = ibook.book1D("BarrelTrueKAsK", "True k as k momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1366  meBarrelTrueKAsP_ = ibook.book1D("BarrelTrueKAsP", "True k as p momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1367  meEndcapTrueKAsPi_ =
1368  ibook.book1D("EndcapTrueKAsPi", "True k as pi momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1369  meEndcapTrueKAsK_ = ibook.book1D("EndcapTrueKAsK", "True k as k momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1370  meEndcapTrueKAsP_ = ibook.book1D("EndcapTrueKAsP", "True k as p momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1371 
1372  meBarrelTruePAsPi_ =
1373  ibook.book1D("BarrelTruePAsPi", "True p as pi momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1374  meBarrelTruePAsK_ = ibook.book1D("BarrelTruePAsK", "True p as k momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1375  meBarrelTruePAsP_ = ibook.book1D("BarrelTruePAsP", "True p as p momentum spectrum, |eta| < 1.5;p [GeV]", 25, 0., 10.);
1376  meEndcapTruePAsPi_ =
1377  ibook.book1D("EndcapTruePAsPi", "True p as pi momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1378  meEndcapTruePAsK_ = ibook.book1D("EndcapTruePAsK", "True p as k momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1379  meEndcapTruePAsP_ = ibook.book1D("EndcapTruePAsP", "True p as p momentum spectrum, |eta| > 1.6;p [GeV]", 25, 0., 10.);
1380 }
1381 
1382 bool Primary4DVertexValidation::matchRecoTrack2SimSignal(const reco::TrackBaseRef& recoTrack) {
1383  auto found = r2s_->find(recoTrack);
1384 
1385  // reco track not matched to any TP
1386  if (found == r2s_->end())
1387  return false;
1388 
1390  for (const auto& tp : found->val) {
1391  if (tp.first->eventId().bunchCrossing() == 0 && tp.first->eventId().event() == 0)
1392  return true;
1393  }
1394 
1395  // reco track not matched to any TP from signal vertex
1396  return false;
1397 }
1398 
1399 std::pair<const edm::Ref<std::vector<TrackingParticle>>*, int> Primary4DVertexValidation::getMatchedTP(
1400  const reco::TrackBaseRef& recoTrack, const TrackingVertexRef& vsim) {
1401  auto found = r2s_->find(recoTrack);
1402 
1403  // reco track not matched to any TP (fake tracks)
1404  if (found == r2s_->end())
1405  return std::make_pair(nullptr, -1);
1406 
1407  // matched TP equal to any TP of a given sim vertex
1408  for (const auto& tp : found->val) {
1409  if (std::find_if(vsim->daughterTracks_begin(), vsim->daughterTracks_end(), [&](const TrackingParticleRef& vtp) {
1410  return tp.first == vtp;
1411  }) != vsim->daughterTracks_end())
1412  return std::make_pair(&tp.first, 0);
1413  // matched TP not associated to any daughter track of a given sim vertex but having the same eventID (track from secondary vtx)
1414  else if (tp.first->eventId().bunchCrossing() == vsim->eventId().bunchCrossing() &&
1415  tp.first->eventId().event() == vsim->eventId().event()) {
1416  return std::make_pair(&tp.first, 1);
1417  }
1418  // matched TP not associated to any sim vertex of a given simulated event (PU track)
1419  else {
1420  return std::make_pair(&tp.first, 2);
1421  }
1422  }
1423 
1424  // reco track not matched to any TP from vertex
1425  return std::make_pair(nullptr, -1);
1426 }
1427 
1428 double Primary4DVertexValidation::timeFromTrueMass(double mass, double pathlength, double momentum, double time) {
1429  if (time > 0 && pathlength > 0 && mass > 0) {
1430  double gammasq = 1. + momentum * momentum / (mass * mass);
1431  double v = c_ * std::sqrt(1. - 1. / gammasq); // cm / ns
1432  double t_est = time - (pathlength / v);
1433 
1434  return t_est;
1435  } else {
1436  return -1;
1437  }
1438 }
1439 
1440 bool Primary4DVertexValidation::select(const reco::Vertex& v, int level) {
1441  /* level
1442  0 !isFake && ndof>4 (default)
1443  1 !isFake && ndof>4 && prob > 0.01
1444  2 !isFake && ndof>4 && prob > 0.01 && ptmax2 > 0.4
1445  */
1446  if (v.isFake())
1447  return false;
1448  if ((level == 0) && (v.ndof() > selNdof_))
1449  return true;
1450  /*if ((level == 1) && (v.ndof() > selNdof_) && (vertex_pxy(v) > 0.01))
1451  return true;
1452  if ((level == 2) && (v.ndof() > selNdof_) && (vertex_pxy(v) > 0.01) && (vertex_ptmax2(v) > 0.4))
1453  return true;
1454  if ((level == 3) && (v.ndof() > selNdof_) && (vertex_ptmax2(v) < 0.4))
1455  return true;*/
1456  return false;
1457 }
1458 
1459 void Primary4DVertexValidation::observablesFromJets(const std::vector<reco::Track>& reco_Tracks,
1460  const std::vector<double>& mass_Tracks,
1461  const std::vector<int>& category_Tracks,
1462  const std::string& skip_Tracks,
1463  unsigned int& n_Jets,
1464  double& sum_EtJets,
1465  double& sum_Pt2Jets,
1466  double& met_Pt,
1467  double& sum_PzJets) {
1468  double sum_PtJets = 0;
1469  n_Jets = 0;
1470  sum_EtJets = 0;
1471  sum_Pt2Jets = 0;
1472  met_Pt = 0;
1473  sum_PzJets = 0;
1474  auto met = LorentzVector(0, 0, 0, 0);
1475  std::vector<fastjet::PseudoJet> fjInputs_;
1476  fjInputs_.clear();
1477  size_t countScale0 = 0;
1478  for (size_t i = 0; i < reco_Tracks.size(); i++) {
1479  const auto recotr = reco_Tracks[i];
1480  const auto mass = mass_Tracks[i];
1481  float scale = 1.;
1482  if (recotr.charge() == 0) {
1483  continue;
1484  }
1485  // skip PU tracks in jet definition if skip_PU is required
1486  if (skip_Tracks == "skip_PU" && category_Tracks[i] == 2) {
1487  continue;
1488  }
1489  // skip fake tracks in jet definition if skip_Fake is required
1490  if (skip_Tracks == "skip_Fake" && category_Tracks[i] == -1) {
1491  continue;
1492  }
1493  if (recotr.pt() != 0) {
1494  scale = (recotr.pt() - recotr.ptError()) / recotr.pt();
1495  }
1496  if (edm::isNotFinite(scale)) {
1497  edm::LogWarning("Primary4DVertexValidation") << "Scaling is NAN ignoring this recotrack" << std::endl;
1498  scale = 0;
1499  }
1500  if (scale < 0) {
1501  scale = 0;
1502  countScale0++;
1503  }
1504  if (scale != 0) {
1505  fjInputs_.push_back(fastjet::PseudoJet(recotr.px() * scale,
1506  recotr.py() * scale,
1507  recotr.pz() * scale,
1508  std::sqrt(recotr.p() * recotr.p() + mass * mass) * scale));
1509  }
1510  }
1511  fastjet::ClusterSequence sequence(fjInputs_, fastjet::JetDefinition(fastjet::antikt_algorithm, 0.4));
1512  auto jets = fastjet::sorted_by_pt(sequence.inclusive_jets(0));
1513  for (const auto& pj : jets) {
1514  auto p4 = LorentzVector(pj.px(), pj.py(), pj.pz(), pj.e());
1515  sum_EtJets += std::sqrt(p4.e() * p4.e() - p4.P() * p4.P() + p4.pt() * p4.pt());
1516  sum_PtJets += p4.pt();
1517  sum_Pt2Jets += (p4.pt() * p4.pt() * 0.8 * 0.8);
1518  met += p4;
1519  sum_PzJets += p4.pz();
1520  n_Jets++;
1521  }
1522  met_Pt = met.pt();
1523  double metAbove = met_Pt - 2 * std::sqrt(sum_PtJets);
1524  if (metAbove > 0) {
1525  sum_Pt2Jets += (metAbove * metAbove);
1526  }
1527  if (countScale0 == reco_Tracks.size()) {
1528  sum_Pt2Jets = countScale0 * 0.01; //leave some epsilon value to sort vertices with unknown pt
1529  }
1530 }
1531 
1532 void Primary4DVertexValidation::isParticle(const reco::TrackBaseRef& recoTrack,
1533  const edm::ValueMap<float>& sigmat0,
1534  const edm::ValueMap<float>& sigmat0Safe,
1535  const edm::ValueMap<float>& probPi,
1536  const edm::ValueMap<float>& probK,
1537  const edm::ValueMap<float>& probP,
1538  unsigned int& no_PIDtype,
1539  bool& no_PID,
1540  bool& is_Pi,
1541  bool& is_K,
1542  bool& is_P) {
1543  no_PIDtype = 0;
1544  no_PID = false;
1545  is_Pi = false;
1546  is_K = false;
1547  is_P = false;
1548  if (probPi[recoTrack] == -1) {
1549  no_PIDtype = 1;
1550  } else if (edm::isNotFinite(probPi[recoTrack])) {
1551  no_PIDtype = 2;
1552  } else if (probPi[recoTrack] == 1 && probK[recoTrack] == 0 && probP[recoTrack] == 0 &&
1553  sigmat0[recoTrack] < sigmat0Safe[recoTrack]) {
1554  no_PIDtype = 3;
1555  }
1556  no_PID = no_PIDtype > 0;
1557  is_Pi = !no_PID && 1. - probPi[recoTrack] < minProbHeavy_;
1558  is_K = !no_PID && !is_Pi && probK[recoTrack] > probP[recoTrack];
1559  is_P = !no_PID && !is_Pi && !is_K;
1560 }
1561 
1562 void Primary4DVertexValidation::getWosWnt(const reco::Vertex& recoVtx,
1563  const reco::TrackBaseRef& recoTrk,
1564  const edm::ValueMap<float>& sigmat0,
1565  const edm::ValueMap<float>& mtdQualMVA,
1566  const edm::Handle<reco::BeamSpot>& BS,
1567  double& wos,
1568  double& wnt) {
1569  double dz2_beam = pow((*BS).BeamWidthX() * cos(recoTrk->phi()) / tan(recoTrk->theta()), 2) +
1570  pow((*BS).BeamWidthY() * sin(recoTrk->phi()) / tan(recoTrk->theta()), 2);
1571  double dz2 =
1572  pow(recoTrk->dzError(), 2) + dz2_beam + pow(0.0020, 2); // added 20 um, some tracks have crazy small resolutions
1573  wos = recoVtx.trackWeight(recoTrk) / dz2;
1574  wnt = recoVtx.trackWeight(recoTrk) *
1575  std::min(recoTrk->pt(), 1.0); // pt-weighted number of tracks (downweights pt < 1 GeV tracks)
1576 
1577  // If tracks have time information, give more weight to tracks with good resolution
1578  if (sigmat0[recoTrk] > 0 && mtdQualMVA[recoTrk] > mvaTh_) {
1579  double sigmaZ = (*BS).sigmaZ();
1580  double sigmaT = sigmaZ / c_; // c in cm/ns
1581  wos = wos / erf(sigmat0[recoTrk] / sigmaT);
1582  }
1583 }
1584 
1585 /* Extract information form TrackingParticles/TrackingVertex and fill
1586  * the helper class simPrimaryVertex with proper generation-level
1587  * information */
1588 std::vector<Primary4DVertexValidation::simPrimaryVertex> Primary4DVertexValidation::getSimPVs(
1589  const edm::Handle<TrackingVertexCollection>& tVC) {
1590  std::vector<Primary4DVertexValidation::simPrimaryVertex> simpv;
1591  int current_event = -1;
1592  int s = -1;
1593  for (TrackingVertexCollection::const_iterator v = tVC->begin(); v != tVC->end(); ++v) {
1594  // We keep only the first vertex from all the events at BX=0.
1595  if (v->eventId().bunchCrossing() != 0)
1596  continue;
1597  if (v->eventId().event() != current_event) {
1598  current_event = v->eventId().event();
1599  } else {
1600  continue;
1601  }
1602  s++;
1603  if (std::abs(v->position().z()) > 1000)
1604  continue; // skip junk vertices
1605 
1606  // could be a new vertex, check all primaries found so far to avoid multiple entries
1607  simPrimaryVertex sv(v->position().x(), v->position().y(), v->position().z(), v->position().t());
1608  sv.eventId = v->eventId();
1609  sv.sim_vertex = TrackingVertexRef(tVC, std::distance(tVC->begin(), v));
1610  sv.OriginalIndex = s;
1611 
1612  for (TrackingParticleRefVector::iterator iTrack = v->daughterTracks_begin(); iTrack != v->daughterTracks_end();
1613  ++iTrack) {
1614  assert((**iTrack).eventId().bunchCrossing() == 0);
1615  }
1616  simPrimaryVertex* vp = nullptr; // will become non-NULL if a vertex is found and then point to it
1617  for (std::vector<simPrimaryVertex>::iterator v0 = simpv.begin(); v0 != simpv.end(); v0++) {
1618  if ((sv.eventId == v0->eventId) && (std::abs(sv.x - v0->x) < 1e-5) && (std::abs(sv.y - v0->y) < 1e-5) &&
1619  (std::abs(sv.z - v0->z) < 1e-5)) {
1620  vp = &(*v0);
1621  break;
1622  }
1623  }
1624  if (!vp) {
1625  // this is a new vertex, add it to the list of sim-vertices
1626  simpv.push_back(sv);
1627  vp = &simpv.back();
1628  }
1629 
1630  // Loop over daughter track(s) as Tracking Particles
1631  for (TrackingVertex::tp_iterator iTP = v->daughterTracks_begin(); iTP != v->daughterTracks_end(); ++iTP) {
1632  auto momentum = (*(*iTP)).momentum();
1633  const reco::Track* matched_best_reco_track = nullptr;
1634  double match_quality = -1;
1635  if (use_only_charged_tracks_ && (**iTP).charge() == 0)
1636  continue;
1637  if (s2r_->find(*iTP) != s2r_->end()) {
1638  matched_best_reco_track = (*s2r_)[*iTP][0].first.get();
1639  match_quality = (*s2r_)[*iTP][0].second;
1640  }
1641 
1642  vp->ptot.setPx(vp->ptot.x() + momentum.x());
1643  vp->ptot.setPy(vp->ptot.y() + momentum.y());
1644  vp->ptot.setPz(vp->ptot.z() + momentum.z());
1645  vp->ptot.setE(vp->ptot.e() + (**iTP).energy());
1646  vp->pt += (**iTP).pt();
1647  vp->ptsq += ((**iTP).pt() * (**iTP).pt());
1648  vp->nGenTrk++;
1649 
1650  if (matched_best_reco_track) {
1651  vp->num_matched_reco_tracks++;
1652  vp->average_match_quality += match_quality;
1653  }
1654  } // End of for loop on daughters sim-particles
1655  if (vp->num_matched_reco_tracks) {
1656  vp->average_match_quality /= static_cast<float>(vp->num_matched_reco_tracks);
1657  }
1658  LogTrace("Primary4DVertexValidation")
1659  << "average number of associated tracks: " << vp->num_matched_reco_tracks / static_cast<float>(vp->nGenTrk)
1660  << " with average quality: " << vp->average_match_quality;
1661  } // End of for loop on tracking vertices
1662 
1663  // In case of no simulated vertices, break here
1664  if (simpv.empty())
1665  return simpv;
1666 
1667  // Now compute the closest distance in z between all simulated vertex
1668  // first initialize
1669  auto prev_z = simpv.back().z;
1670  for (simPrimaryVertex& vsim : simpv) {
1671  vsim.closest_vertex_distance_z = std::abs(vsim.z - prev_z);
1672  prev_z = vsim.z;
1673  }
1674  // then calculate
1675  for (std::vector<simPrimaryVertex>::iterator vsim = simpv.begin(); vsim != simpv.end(); vsim++) {
1676  std::vector<simPrimaryVertex>::iterator vsim2 = vsim;
1677  vsim2++;
1678  for (; vsim2 != simpv.end(); vsim2++) {
1679  double distance = std::abs(vsim->z - vsim2->z);
1680  // need both to be complete
1681  vsim->closest_vertex_distance_z = std::min(vsim->closest_vertex_distance_z, distance);
1682  vsim2->closest_vertex_distance_z = std::min(vsim2->closest_vertex_distance_z, distance);
1683  }
1684  }
1685  return simpv;
1686 }
1687 
1688 /* Extract information form reco Vertex and fill the helper class
1689  * recoPrimaryVertex with proper reco-level information */
1690 std::vector<Primary4DVertexValidation::recoPrimaryVertex> Primary4DVertexValidation::getRecoPVs(
1691  const edm::Handle<edm::View<reco::Vertex>>& tVC) {
1692  std::vector<Primary4DVertexValidation::recoPrimaryVertex> recopv;
1693  int r = -1;
1694  for (auto v = tVC->begin(); v != tVC->end(); ++v) {
1695  r++;
1696  // Skip junk vertices
1697  if (std::abs(v->z()) > 1000)
1698  continue;
1699  if (v->isFake() || !v->isValid())
1700  continue;
1701 
1702  recoPrimaryVertex sv(v->position().x(), v->position().y(), v->position().z());
1703  sv.recVtx = &(*v);
1704  sv.recVtxRef = reco::VertexBaseRef(tVC, std::distance(tVC->begin(), v));
1705 
1706  sv.OriginalIndex = r;
1707  sv.ndof = v->ndof();
1708  // this is a new vertex, add it to the list of reco-vertices
1709  recopv.push_back(sv);
1710  Primary4DVertexValidation::recoPrimaryVertex* vp = &recopv.back();
1711 
1712  // Loop over daughter track(s)
1713  for (auto iTrack = v->tracks_begin(); iTrack != v->tracks_end(); ++iTrack) {
1714  auto momentum = (*(*iTrack)).innerMomentum();
1715  if (momentum.mag2() == 0)
1716  momentum = (*(*iTrack)).momentum();
1717  vp->pt += std::sqrt(momentum.perp2());
1718  vp->ptsq += (momentum.perp2());
1719  vp->nRecoTrk++;
1720 
1721  auto matched = r2s_->find(*iTrack);
1722  if (matched != r2s_->end()) {
1723  vp->num_matched_sim_tracks++;
1724  }
1725 
1726  } // End of for loop on daughters reconstructed tracks
1727  } // End of for loop on tracking vertices
1728 
1729  // In case of no reco vertices, break here
1730  if (recopv.empty())
1731  return recopv;
1732 
1733  // Now compute the closest distance in z between all reconstructed vertex
1734  // first initialize
1735  auto prev_z = recopv.back().z;
1736  for (recoPrimaryVertex& vreco : recopv) {
1737  vreco.closest_vertex_distance_z = std::abs(vreco.z - prev_z);
1738  prev_z = vreco.z;
1739  }
1740  for (std::vector<recoPrimaryVertex>::iterator vreco = recopv.begin(); vreco != recopv.end(); vreco++) {
1741  std::vector<recoPrimaryVertex>::iterator vreco2 = vreco;
1742  vreco2++;
1743  for (; vreco2 != recopv.end(); vreco2++) {
1744  double distance = std::abs(vreco->z - vreco2->z);
1745  // need both to be complete
1746  vreco->closest_vertex_distance_z = std::min(vreco->closest_vertex_distance_z, distance);
1747  vreco2->closest_vertex_distance_z = std::min(vreco2->closest_vertex_distance_z, distance);
1748  }
1749  }
1750  return recopv;
1751 }
1752 
1753 // ------------ method called to produce the data ------------
1754 void Primary4DVertexValidation::matchReco2Sim(std::vector<recoPrimaryVertex>& recopv,
1755  std::vector<simPrimaryVertex>& simpv,
1756  const edm::ValueMap<float>& sigmat0,
1757  const edm::ValueMap<float>& MVA,
1758  const edm::Handle<reco::BeamSpot>& BS) {
1759  // Initialization (clear wos and wnt)
1760  for (auto vv : simpv) {
1761  vv.wnt.clear();
1762  vv.wos.clear();
1763  }
1764  for (auto rv : recopv) {
1765  rv.wnt.clear();
1766  rv.wos.clear();
1767  }
1768 
1769  // Filling infos for matching rec and sim vertices
1770  for (unsigned int iv = 0; iv < recopv.size(); iv++) {
1771  const reco::Vertex* vertex = recopv.at(iv).recVtx;
1772  LogTrace("Primary4DVertexValidation") << "iv (rec): " << iv;
1773 
1774  for (unsigned int iev = 0; iev < simpv.size(); iev++) {
1775  double wnt = 0;
1776  double wos = 0;
1777  double evwnt = 0; // weighted number of tracks from sim event iev in the current recvtx
1778  double evwos = 0; // weight over sigma**2 of sim event iev in the current recvtx
1779  unsigned int evnt = 0; // number of tracks from sim event iev in the current recvtx
1780 
1781  for (auto iTrack = vertex->tracks_begin(); iTrack != vertex->tracks_end(); ++iTrack) {
1782  if (vertex->trackWeight(*iTrack) < trackweightTh_) {
1783  continue;
1784  }
1785 
1786  auto tp_info = getMatchedTP(*iTrack, simpv.at(iev).sim_vertex).first;
1787  int matchCategory = getMatchedTP(*iTrack, simpv.at(iev).sim_vertex).second;
1788  // matched TP equal to any TP of a given sim vertex
1789  if (tp_info != nullptr && matchCategory == 0) {
1790  getWosWnt(*vertex, *iTrack, MVA, sigmat0, BS, wos, wnt);
1791  simpv.at(iev).addTrack(iv, wos, wnt);
1792  recopv.at(iv).addTrack(iev, wos, wnt);
1793  evwos += wos;
1794  evwnt += wnt;
1795  evnt++;
1796  }
1797  } // RecoTracks loop
1798 
1799  // require 2 tracks for a wos-match
1800  if ((evwos > 0) && (evwos > recopv.at(iv).maxwos) && (evnt > 1)) {
1801  recopv.at(iv).wosmatch = iev;
1802  recopv.at(iv).maxwos = evwos;
1803  recopv.at(iv).maxwosnt = evnt;
1804  LogTrace("Primary4DVertexValidation") << "dominating sim event (iev): " << iev << " evwos = " << evwos;
1805  }
1806 
1807  // weighted track counting match, require at least one track
1808  if ((evnt > 0) && (evwnt > recopv.at(iv).maxwnt)) {
1809  recopv.at(iv).wntmatch = iev;
1810  recopv.at(iv).maxwnt = evwnt;
1811  }
1812  } // TrackingVertex loop
1813  if (recopv.at(iv).maxwos > 0) {
1814  simpv.at(recopv.at(iv).wosmatch).wos_dominated_recv.push_back(iv);
1815  simpv.at(recopv.at(iv).wosmatch).nwosmatch++; // count the rec vertices dominated by a sim vertex using wos
1816  assert(iv < recopv.size());
1817  }
1818  LogTrace("Primary4DVertexValidation") << "largest contribution to wos: wosmatch (iev) = " << recopv.at(iv).wosmatch
1819  << " maxwos = " << recopv.at(iv).maxwos;
1820  if (recopv.at(iv).maxwnt > 0) {
1821  simpv.at(recopv.at(iv).wntmatch).nwntmatch++; // count the rec vertices dominated by a sim vertex using wnt
1822  }
1823  } // RecoPrimaryVertex
1824 
1825  // reset
1826  for (auto& vrec : recopv) {
1827  vrec.sim = NOT_MATCHED;
1828  vrec.matchQuality = 0;
1829  }
1830  unsigned int iev = 0;
1831  for (auto& vv : simpv) {
1832  LogTrace("Primary4DVertexValidation") << "iev (sim): " << iev;
1833  LogTrace("Primary4DVertexValidation") << "wos_dominated_recv.size: " << vv.wos_dominated_recv.size();
1834  for (unsigned int i = 0; i < vv.wos_dominated_recv.size(); i++) {
1835  auto recov = vv.wos_dominated_recv.at(i);
1836  LogTrace("Primary4DVertexValidation")
1837  << "index of reco vertex: " << recov << " that has a wos: " << vv.wos.at(recov) << " at position " << i;
1838  }
1839  vv.rec = NOT_MATCHED;
1840  vv.matchQuality = 0;
1841  iev++;
1842  }
1843  // after filling infos, goes for the sim-reco match
1844  // this tries a one-to-one match, taking simPV with highest wos if there are > 1 simPV candidates
1845  for (unsigned int rank = 1; rank < maxRank_; rank++) {
1846  for (unsigned int iev = 0; iev < simpv.size(); iev++) { //loop on SimPV
1847  if (simpv.at(iev).rec != NOT_MATCHED) {
1848  continue; // only sim vertices not already matched
1849  }
1850  if (simpv.at(iev).nwosmatch == 0) {
1851  continue; // the sim vertex does not dominate any reco vertex
1852  }
1853  if (simpv.at(iev).nwosmatch > rank) {
1854  continue; // start with sim vertices dominating one rec vertex (rank 1), then go with the ones dominating more
1855  }
1856  unsigned int iv = NOT_MATCHED; // select a rec vertex index
1857  for (unsigned int k = 0; k < simpv.at(iev).wos_dominated_recv.size(); k++) {
1858  unsigned int rec = simpv.at(iev).wos_dominated_recv.at(k); //candidate rec vertex index
1859  auto vrec = recopv.at(rec);
1860  if (vrec.sim != NOT_MATCHED) {
1861  continue; // already matched
1862  }
1863  if (std::abs(simpv.at(iev).z - vrec.z) > zWosMatchMax_) {
1864  continue; // insanely far away
1865  }
1866  if ((iv == NOT_MATCHED) || simpv.at(iev).wos.at(rec) > simpv.at(iev).wos.at(iv)) {
1867  iv = rec;
1868  }
1869  }
1870  // if a viable candidate was found, make the link
1871  if (iv !=
1872  NOT_MATCHED) { // if the rec vertex has already been associated is possible that iv remains NOT_MATCHED at this point
1873  recopv.at(iv).sim = iev;
1874  simpv.at(iev).rec = iv;
1875  recopv.at(iv).matchQuality = rank;
1876  simpv.at(iev).matchQuality = rank;
1877  }
1878  } // iev
1879  } // rank
1880 
1881  // Reco vertices that are not necessarily dominated by a sim vertex, or whose dominating sim-vertex
1882  // has been matched already to another overlapping reco vertex, can still be associated to a specific
1883  // sim vertex (without being classified as dominating).
1884  // In terms of fitting 1d-distributions, this corresponds to a small peak close to a much larger nearby peak
1885  unsigned int ntry = 0;
1886  while (ntry++ < maxTry_) {
1887  unsigned nmatch = 0;
1888  for (unsigned int iev = 0; iev < simpv.size(); iev++) {
1889  if ((simpv.at(iev).rec != NOT_MATCHED) || (simpv.at(iev).wos.empty())) {
1890  continue;
1891  }
1892  // find a rec vertex for the NOT_MATCHED sim vertex
1893  unsigned int rec = NOT_MATCHED;
1894  for (auto rv : simpv.at(iev).wos) {
1895  if ((rec == NOT_MATCHED) || (rv.second > simpv.at(iev).wos.at(rec))) {
1896  rec = rv.first;
1897  }
1898  }
1899 
1900  if (rec == NOT_MATCHED) { // try with wnt match
1901  for (auto rv : simpv.at(iev).wnt) {
1902  if ((rec == NOT_MATCHED) || (rv.second > simpv.at(iev).wnt.at(rec))) {
1903  rec = rv.first;
1904  }
1905  }
1906  }
1907 
1908  if (rec == NOT_MATCHED) {
1909  continue; // should not happen
1910  }
1911  if (recopv.at(rec).sim != NOT_MATCHED) {
1912  continue; // already gone
1913  }
1914 
1915  // check if the rec vertex can be matched
1916  unsigned int rec2sim = NOT_MATCHED;
1917  for (auto sv : recopv.at(rec).wos) {
1918  if (simpv.at(sv.first).rec != NOT_MATCHED) {
1919  continue; // already used
1920  }
1921  if ((rec2sim == NOT_MATCHED) || (sv.second > recopv.at(rec).wos.at(rec2sim))) {
1922  rec2sim = sv.first;
1923  }
1924  }
1925  if (iev == rec2sim) {
1926  // do the match and assign lowest quality (i.e. max rank)
1927  recopv.at(rec).sim = iev;
1928  recopv.at(rec).matchQuality = maxRank_;
1929  simpv.at(iev).rec = rec;
1930  simpv.at(iev).matchQuality = maxRank_;
1931  nmatch++;
1932  }
1933  } // sim loop
1934  if (nmatch == 0) {
1935  break;
1936  }
1937  } // ntry
1938 
1939 // Debugging
1940 #ifdef EDM_ML_DEBUG
1941  unsigned int nmatch_tot = 0, n_dzgtsz = 0;
1942  unsigned int n_rank1 = 0, n_rank2 = 0, n_rank3 = 0, n_rank8 = 0;
1943 
1944  for (unsigned int iev = 0; iev < simpv.size(); iev++) { //loop on SimPV
1945  if (simpv.at(iev).rec != NOT_MATCHED) {
1946  unsigned int rec = simpv.at(iev).rec;
1947  unsigned int wosmatch = recopv.at(rec).wosmatch;
1948  LogTrace("Primary4DVertexValidation")
1949  << "Final match: iev (sim) = " << std::setw(4) << iev << " sim.rec = " << std::setw(4) << rec
1950  << " rec.wosmatch = " << std::setw(5) << wosmatch << " dZ/sigmaZ = " << std::setw(6) << std::setprecision(2)
1951  << std::abs((recopv.at(rec).z - simpv.at(iev).z) / recopv.at(rec).recVtx->zError())
1952  << " match qual = " << std::setw(1) << recopv.at(rec).matchQuality;
1953  nmatch_tot++;
1954  if (std::abs((recopv.at(rec).z - simpv.at(iev).z) / recopv.at(rec).recVtx->zError()) > 3.) {
1955  n_dzgtsz++;
1956  }
1957  if (recopv.at(rec).matchQuality == 1) {
1958  n_rank1++;
1959  }
1960  if (recopv.at(rec).matchQuality == 2) {
1961  n_rank2++;
1962  }
1963  if (recopv.at(rec).matchQuality == 3) {
1964  n_rank3++;
1965  }
1966  if (recopv.at(rec).matchQuality == 8) {
1967  n_rank8++;
1968  }
1969  }
1970  }
1971  LogTrace("Primary4DVertexValidation") << "n_sim = " << simpv.size() << " n_rec = " << recopv.size()
1972  << " nmatch_tot = " << nmatch_tot << " n(dZ>sigmaZ) = " << n_dzgtsz
1973  << " n_rank1 = " << n_rank1 << " n_rank2 = " << n_rank2
1974  << " n_rank3 = " << n_rank3 << " n_rank8 = " << n_rank8;
1975 #endif
1976 }
1977 
1978 void Primary4DVertexValidation::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
1979  using edm::Handle;
1980  using edm::View;
1981  using std::cout;
1982  using std::endl;
1983  using std::vector;
1984  using namespace reco;
1985 
1986  std::vector<float> pileUpInfo_z;
1987 
1988  // get the pileup information
1989  edm::Handle<std::vector<PileupSummaryInfo>> puinfoH;
1990  if (iEvent.getByToken(vecPileupSummaryInfoToken_, puinfoH)) {
1991  for (auto const& pu_info : *puinfoH.product()) {
1992  if (pu_info.getBunchCrossing() == 0) {
1993  pileUpInfo_z = pu_info.getPU_zpositions();
1994  break;
1995  }
1996  }
1997  }
1998 
1999  edm::Handle<TrackingParticleCollection> TPCollectionH;
2000  iEvent.getByToken(trackingParticleCollectionToken_, TPCollectionH);
2001  if (!TPCollectionH.isValid())
2002  edm::LogWarning("Primary4DVertexValidation") << "TPCollectionH is not valid";
2003 
2004  edm::Handle<TrackingVertexCollection> TVCollectionH;
2005  iEvent.getByToken(trackingVertexCollectionToken_, TVCollectionH);
2006  if (!TVCollectionH.isValid())
2007  edm::LogWarning("Primary4DVertexValidation") << "TVCollectionH is not valid";
2008 
2009  edm::Handle<reco::SimToRecoCollection> simToRecoH;
2010  iEvent.getByToken(simToRecoAssociationToken_, simToRecoH);
2011  if (simToRecoH.isValid())
2012  s2r_ = simToRecoH.product();
2013  else
2014  edm::LogWarning("Primary4DVertexValidation") << "simToRecoH is not valid";
2015 
2016  edm::Handle<reco::RecoToSimCollection> recoToSimH;
2017  iEvent.getByToken(recoToSimAssociationToken_, recoToSimH);
2018  if (recoToSimH.isValid())
2019  r2s_ = recoToSimH.product();
2020  else
2021  edm::LogWarning("Primary4DVertexValidation") << "recoToSimH is not valid";
2022 
2023  edm::Handle<reco::BeamSpot> BeamSpotH;
2024  iEvent.getByToken(RecBeamSpotToken_, BeamSpotH);
2025  if (!BeamSpotH.isValid())
2026  edm::LogWarning("Primary4DVertexValidation") << "BeamSpotH is not valid";
2027 
2028  std::vector<simPrimaryVertex> simpv; // a list of simulated primary MC vertices
2029  simpv = getSimPVs(TVCollectionH);
2030  // this bool check if first vertex in that with highest pT
2031  bool signal_is_highest_pt =
2032  std::max_element(simpv.begin(), simpv.end(), [](const simPrimaryVertex& lhs, const simPrimaryVertex& rhs) {
2033  return lhs.ptsq < rhs.ptsq;
2034  }) == simpv.begin();
2035 
2036  std::vector<recoPrimaryVertex> recopv; // a list of reconstructed primary MC vertices
2037  edm::Handle<edm::View<reco::Vertex>> recVtxs;
2038  iEvent.getByToken(Rec4DVerToken_, recVtxs);
2039  if (!recVtxs.isValid())
2040  edm::LogWarning("Primary4DVertexValidation") << "recVtxs is not valid";
2041  recopv = getRecoPVs(recVtxs);
2042 
2043  const auto& trackAssoc = iEvent.get(trackAssocToken_);
2044  const auto& pathLength = iEvent.get(pathLengthToken_);
2045  const auto& momentum = iEvent.get(momentumToken_);
2046  const auto& time = iEvent.get(timeToken_);
2047  const auto& t0Pid = iEvent.get(t0PidToken_);
2048  const auto& sigmat0 = iEvent.get(sigmat0PidToken_);
2049  const auto& t0Safe = iEvent.get(t0SafePidToken_);
2050  const auto& sigmat0Safe = iEvent.get(sigmat0SafePidToken_);
2051  const auto& mtdQualMVA = iEvent.get(trackMVAQualToken_);
2052  const auto& tMtd = iEvent.get(tmtdToken_);
2053  const auto& tofPi = iEvent.get(tofPiToken_);
2054  const auto& tofK = iEvent.get(tofKToken_);
2055  const auto& tofP = iEvent.get(tofPToken_);
2056  const auto& probPi = iEvent.get(probPiToken_);
2057  const auto& probK = iEvent.get(probKToken_);
2058  const auto& probP = iEvent.get(probPToken_);
2059  const auto& fPDGTable = iSetup.getHandle(pdtToken_);
2060 
2061  // I have simPV and recoPV collections
2062  matchReco2Sim(recopv, simpv, sigmat0Safe, mtdQualMVA, BeamSpotH);
2063 
2064  // Loop on tracks
2065  for (unsigned int iv = 0; iv < recopv.size(); iv++) {
2066  if (recopv.at(iv).ndof > selNdof_) {
2067  const reco::Vertex* vertex = recopv.at(iv).recVtx;
2068 
2069  for (unsigned int iev = 0; iev < simpv.size(); iev++) {
2070  auto vsim = simpv.at(iev).sim_vertex;
2071 
2072  bool selectedVtxMatching = recopv.at(iv).sim == iev && simpv.at(iev).rec == iv;
2073  bool selectedLV = simpv.at(iev).eventId.bunchCrossing() == 0 && simpv.at(iev).eventId.event() == 0 &&
2074  recopv.at(iv).OriginalIndex == 0;
2075  bool selectedLVMatching = selectedVtxMatching && selectedLV; // bool for reco vtx leading match
2076  if (selectedLVMatching && !recopv.at(iv).is_signal()) {
2077  edm::LogWarning("Primary4DVertexValidation")
2078  << "Reco vtx leading match inconsistent: BX/ID " << simpv.at(iev).eventId.bunchCrossing() << " "
2079  << simpv.at(iev).eventId.event();
2080  }
2081 #ifdef EDM_ML_DEBUG
2082  if (selectedLVMatching) {
2083  printSimVtxRecoVtxInfo(simpv.at(iev), recopv.at(iv));
2084  }
2085 #endif
2086  double vzsim = simpv.at(iev).z;
2087  double vtsim = simpv.at(iev).t * simUnit_;
2088 
2089  double wnt = 0, wos = 0;
2090  double PUsumWnt = 0, PUsumWos = 0, SecsumWos = 0, FakesumWos = 0, PUsumPt = 0, PUsumPt2 = 0;
2091  double sumWnt = 0, sumWos = 0, sumPt = 0, sumPt2 = 0;
2092  unsigned int nt = 0, PUnt = 0, Fakent = 0;
2093 
2094  std::vector<double> massVector;
2095  std::vector<reco::Track> recotracks;
2096  std::vector<int> categoryVector;
2097  double sumEtJets = 0, sumPt2Jets = 0, metPt = 0, sumPzJets = 0;
2098  double sumEtJetsnoPU = 0, sumPt2JetsnoPU = 0, metPtnoPU = 0, sumPzJetsnoPU = 0;
2099  double sumEtJetsnoFake = 0, sumPt2JetsnoFake = 0, metPtnoFake = 0, sumPzJetsnoFake = 0;
2100  unsigned int nJets = 0, nJetsnoPU = 0, nJetsnoFake = 0;
2101  for (auto iTrack = vertex->tracks_begin(); iTrack != vertex->tracks_end(); ++iTrack) {
2102  if (trackAssoc[*iTrack] == -1) {
2103  edm::LogWarning("mtdTracks") << "Extended track not associated";
2104  continue;
2105  }
2106 
2107  // monitor all track weights associated to a vertex before selection on it
2108  if (selectedVtxMatching) {
2109  meVtxTrackW_->Fill(vertex->trackWeight(*iTrack));
2110  if (selectedLV) {
2111  meVtxTrackRecLVW_->Fill(vertex->trackWeight(*iTrack));
2112  }
2113  }
2114 
2115  if (vertex->trackWeight(*iTrack) < trackweightTh_)
2116  continue;
2117  bool noCrack = std::abs((*iTrack)->eta()) < trackMaxBtlEta_ || std::abs((*iTrack)->eta()) > trackMinEtlEta_;
2118 
2119  bool selectRecoTrk = trkRecSel(**iTrack);
2120  if (selectedLVMatching && selectRecoTrk) {
2121  if (noCrack) {
2122  meTrackEffPtTot_->Fill((*iTrack)->pt());
2123  }
2124  meTrackEffEtaTot_->Fill(std::abs((*iTrack)->eta()));
2125  }
2126 
2127  auto tp_info = getMatchedTP(*iTrack, vsim).first;
2128  int matchCategory = getMatchedTP(*iTrack, vsim).second;
2129 
2130  // PU, fake and secondary tracks
2131  if (selectedVtxMatching) {
2132  unsigned int no_PIDtype = 0;
2133  bool no_PID, is_Pi, is_K, is_P;
2134  int PartID = 211; // pion
2135  isParticle(*iTrack, sigmat0, sigmat0Safe, probPi, probK, probP, no_PIDtype, no_PID, is_Pi, is_K, is_P);
2136  if (!use3dNoTime_) {
2137  if (no_PID || is_Pi) {
2138  PartID = 211;
2139  } else if (is_K) {
2140  PartID = 321;
2141  } else if (is_P) {
2142  PartID = 2212;
2143  }
2144  }
2145  const HepPDT::ParticleData* PData = fPDGTable->particle(HepPDT::ParticleID(abs(PartID)));
2146  double mass = PData->mass().value();
2147  massVector.push_back(mass);
2148  recotracks.push_back(**iTrack);
2149  getWosWnt(*vertex, *iTrack, sigmat0Safe, mtdQualMVA, BeamSpotH, wos, wnt);
2150  meVtxTrackWnt_->Fill(wnt);
2151  if (selectedLV) {
2152  meVtxTrackRecLVWnt_->Fill(wnt);
2153  }
2154  // reco track matched to any TP
2155  if (tp_info != nullptr) {
2156 #ifdef EDM_ML_DEBUG
2157  if (selectedLV) {
2158  printMatchedRecoTrackInfo(*vertex, *iTrack, *tp_info, matchCategory);
2159  }
2160 #endif
2161  // matched TP not associated to any daughter track of a given sim vertex but having the same eventID (track from secondary vtx)
2162  if (matchCategory == 1) {
2163  categoryVector.push_back(matchCategory);
2164  SecsumWos += wos;
2165  }
2166  // matched TP not associated to any sim vertex of a given simulated event (PU track)
2167  if (matchCategory == 2) {
2168  if (optionalPlots_) {
2169  mePUTrackWnt_->Fill(wnt);
2170  if (selectedLV) {
2171  mePUTrackRecLVWnt_->Fill(wnt);
2172  }
2173  }
2174  PUsumWnt += wnt;
2175  PUsumWos += wos;
2176  PUsumPt += (*iTrack)->pt();
2177  PUsumPt2 += ((*iTrack)->pt() * (*iTrack)->pt());
2178  PUnt++;
2179  categoryVector.push_back(2);
2180  }
2181  }
2182  // reco track not matched to any TP (fake tracks)
2183  else {
2184  categoryVector.push_back(matchCategory);
2185  FakesumWos += wos;
2186  Fakent++;
2187  }
2188  nt++;
2189  sumWnt += wnt;
2190  sumWos += wos;
2191  sumPt += (*iTrack)->pt();
2192  sumPt2 += ((*iTrack)->pt() * (*iTrack)->pt());
2193  }
2194 
2195  // matched TP equal to any TP of a given sim vertex
2196  if (tp_info != nullptr && matchCategory == 0) {
2197  categoryVector.push_back(matchCategory);
2198  double mass = (*tp_info)->mass();
2199  double tsim = (*tp_info)->parentVertex()->position().t() * simUnit_;
2200  double tEst = timeFromTrueMass(mass, pathLength[*iTrack], momentum[*iTrack], time[*iTrack]);
2201 
2202  double xsim = (*tp_info)->parentVertex()->position().x();
2203  double ysim = (*tp_info)->parentVertex()->position().y();
2204  double zsim = (*tp_info)->parentVertex()->position().z();
2205  double xPCA = (*iTrack)->vx();
2206  double yPCA = (*iTrack)->vy();
2207  double zPCA = (*iTrack)->vz();
2208 
2209  double dZ = zPCA - zsim;
2210  double d3D = std::sqrt((xPCA - xsim) * (xPCA - xsim) + (yPCA - ysim) * (yPCA - ysim) + dZ * dZ);
2211  // orient d3D according to the projection of RECO - SIM onto simulated momentum
2212  if ((xPCA - xsim) * ((*tp_info)->px()) + (yPCA - ysim) * ((*tp_info)->py()) + dZ * ((*tp_info)->pz()) <
2213  0.) {
2214  d3D = -d3D;
2215  }
2216 
2217  // select TPs associated to the signal event
2218  bool selectTP = trkTPSelLV(**tp_info);
2219 
2220  if (selectedLVMatching && selectRecoTrk && selectTP) {
2221  meTrackMatchedTPZposResTot_->Fill((*iTrack)->vz() - vzsim);
2222  if (noCrack) {
2223  meTrackMatchedTPEffPtTot_->Fill((*iTrack)->pt());
2224  }
2225  meTrackMatchedTPEffEtaTot_->Fill(std::abs((*iTrack)->eta()));
2226  }
2227 
2228  if (sigmat0Safe[*iTrack] == -1)
2229  continue;
2230 
2231  if (selectedLVMatching && selectRecoTrk && selectTP) {
2232  meTrackMatchedTPResTot_->Fill(t0Safe[*iTrack] - vtsim);
2233  meTrackMatchedTPPullTot_->Fill((t0Safe[*iTrack] - vtsim) / sigmat0Safe[*iTrack]);
2234  if (noCrack) {
2235  meTrackMatchedTPEffPtMtd_->Fill((*iTrack)->pt());
2236  }
2237  meTrackMatchedTPEffEtaMtd_->Fill(std::abs((*iTrack)->eta()));
2238 
2239  unsigned int noPIDtype = 0;
2240  bool noPID = false, isPi = false, isK = false, isP = false;
2241  isParticle(*iTrack, sigmat0, sigmat0Safe, probPi, probK, probP, noPIDtype, noPID, isPi, isK, isP);
2242 
2243  if ((isPi && std::abs(tMtd[*iTrack] - tofPi[*iTrack] - t0Pid[*iTrack]) > tol_) ||
2244  (isK && std::abs(tMtd[*iTrack] - tofK[*iTrack] - t0Pid[*iTrack]) > tol_) ||
2245  (isP && std::abs(tMtd[*iTrack] - tofP[*iTrack] - t0Pid[*iTrack]) > tol_)) {
2246  edm::LogWarning("Primary4DVertexValidation")
2247  << "No match between mass hyp. and time: " << std::abs((*tp_info)->pdgId()) << " mass hyp pi/k/p "
2248  << isPi << " " << isK << " " << isP << " t0/t0safe " << t0Pid[*iTrack] << " " << t0Safe[*iTrack]
2249  << " tMtd - tof pi/K/p " << tMtd[*iTrack] - tofPi[*iTrack] << " " << tMtd[*iTrack] - tofK[*iTrack]
2250  << " " << tMtd[*iTrack] - tofP[*iTrack] << " Prob pi/K/p " << probPi[*iTrack] << " "
2251  << probK[*iTrack] << " " << probP[*iTrack];
2252  }
2253 
2254  if (std::abs((*iTrack)->eta()) < trackMaxBtlEta_) {
2255  meBarrelPIDp_->Fill((*iTrack)->p());
2256  meBarrelNoPIDtype_->Fill(noPIDtype + 0.5);
2257  if (std::abs((*tp_info)->pdgId()) == 211) {
2258  if (noPID) {
2259  meBarrelTruePiNoPID_->Fill((*iTrack)->p());
2260  } else if (isPi) {
2261  meBarrelTruePiAsPi_->Fill((*iTrack)->p());
2262  } else if (isK) {
2263  meBarrelTruePiAsK_->Fill((*iTrack)->p());
2264  } else if (isP) {
2265  meBarrelTruePiAsP_->Fill((*iTrack)->p());
2266  } else {
2267  edm::LogWarning("Primary4DVertexValidation")
2268  << "No PID class: " << std::abs((*tp_info)->pdgId()) << " t0/t0safe " << t0Pid[*iTrack] << " "
2269  << t0Safe[*iTrack] << " Prob pi/K/p " << probPi[*iTrack] << " " << probK[*iTrack] << " "
2270  << probP[*iTrack];
2271  }
2272  } else if (std::abs((*tp_info)->pdgId()) == 321) {
2273  if (noPID) {
2274  meBarrelTrueKNoPID_->Fill((*iTrack)->p());
2275  } else if (isPi) {
2276  meBarrelTrueKAsPi_->Fill((*iTrack)->p());
2277  } else if (isK) {
2278  meBarrelTrueKAsK_->Fill((*iTrack)->p());
2279  } else if (isP) {
2280  meBarrelTrueKAsP_->Fill((*iTrack)->p());
2281  } else {
2282  edm::LogWarning("Primary4DVertexValidation")
2283  << "No PID class: " << std::abs((*tp_info)->pdgId()) << " t0/t0safe " << t0Pid[*iTrack] << " "
2284  << t0Safe[*iTrack] << " Prob pi/K/p " << probPi[*iTrack] << " " << probK[*iTrack] << " "
2285  << probP[*iTrack];
2286  }
2287  } else if (std::abs((*tp_info)->pdgId()) == 2212) {
2288  if (noPID) {
2289  meBarrelTruePNoPID_->Fill((*iTrack)->p());
2290  } else if (isPi) {
2291  meBarrelTruePAsPi_->Fill((*iTrack)->p());
2292  } else if (isK) {
2293  meBarrelTruePAsK_->Fill((*iTrack)->p());
2294  } else if (isP) {
2295  meBarrelTruePAsP_->Fill((*iTrack)->p());
2296  } else {
2297  edm::LogWarning("Primary4DVertexValidation")
2298  << "No PID class: " << std::abs((*tp_info)->pdgId()) << " t0/t0safe " << t0Pid[*iTrack] << " "
2299  << t0Safe[*iTrack] << " Prob pi/K/p " << probPi[*iTrack] << " " << probK[*iTrack] << " "
2300  << probP[*iTrack];
2301  }
2302  }
2303  } else if (std::abs((*iTrack)->eta()) > trackMinEtlEta_ && std::abs((*iTrack)->eta()) < trackMaxEtlEta_) {
2304  meEndcapPIDp_->Fill((*iTrack)->p());
2305  meEndcapNoPIDtype_->Fill(noPIDtype + 0.5);
2306  if (std::abs((*tp_info)->pdgId()) == 211) {
2307  if (noPID) {
2308  meEndcapTruePiNoPID_->Fill((*iTrack)->p());
2309  } else if (isPi) {
2310  meEndcapTruePiAsPi_->Fill((*iTrack)->p());
2311  } else if (isK) {
2312  meEndcapTruePiAsK_->Fill((*iTrack)->p());
2313  } else if (isP) {
2314  meEndcapTruePiAsP_->Fill((*iTrack)->p());
2315  } else {
2316  edm::LogWarning("Primary4DVertexValidation")
2317  << "No PID class: " << std::abs((*tp_info)->pdgId()) << " t0/t0safe " << t0Pid[*iTrack] << " "
2318  << t0Safe[*iTrack] << " Prob pi/K/p " << probPi[*iTrack] << " " << probK[*iTrack] << " "
2319  << probP[*iTrack];
2320  }
2321  } else if (std::abs((*tp_info)->pdgId()) == 321) {
2322  if (noPID) {
2323  meEndcapTrueKNoPID_->Fill((*iTrack)->p());
2324  } else if (isPi) {
2325  meEndcapTrueKAsPi_->Fill((*iTrack)->p());
2326  } else if (isK) {
2327  meEndcapTrueKAsK_->Fill((*iTrack)->p());
2328  } else if (isP) {
2329  meEndcapTrueKAsP_->Fill((*iTrack)->p());
2330  } else {
2331  edm::LogWarning("Primary4DVertexValidation")
2332  << "No PID class: " << std::abs((*tp_info)->pdgId()) << " t0/t0safe " << t0Pid[*iTrack] << " "
2333  << t0Safe[*iTrack] << " Prob pi/K/p " << probPi[*iTrack] << " " << probK[*iTrack] << " "
2334  << probP[*iTrack];
2335  }
2336  } else if (std::abs((*tp_info)->pdgId()) == 2212) {
2337  if (noPID) {
2338  meEndcapTruePNoPID_->Fill((*iTrack)->p());
2339  } else if (isPi) {
2340  meEndcapTruePAsPi_->Fill((*iTrack)->p());
2341  } else if (isK) {
2342  meEndcapTruePAsK_->Fill((*iTrack)->p());
2343  } else if (isP) {
2344  meEndcapTruePAsP_->Fill((*iTrack)->p());
2345  } else {
2346  edm::LogWarning("Primary4DVertexValidation")
2347  << "No PID class: " << std::abs((*tp_info)->pdgId()) << " t0/t0safe " << t0Pid[*iTrack] << " "
2348  << t0Safe[*iTrack] << " Prob pi/K/p " << probPi[*iTrack] << " " << probK[*iTrack] << " "
2349  << probP[*iTrack];
2350  }
2351  }
2352  }
2353  }
2354  meTrackResTot_->Fill(t0Safe[*iTrack] - tsim);
2355  meTrackPullTot_->Fill((t0Safe[*iTrack] - tsim) / sigmat0Safe[*iTrack]);
2356  meTrackZposResTot_->Fill(dZ);
2357  if ((*iTrack)->p() <= 2) {
2358  meTrackResLowPTot_->Fill(t0Safe[*iTrack] - tsim);
2359  meTrackPullLowPTot_->Fill((t0Safe[*iTrack] - tsim) / sigmat0Safe[*iTrack]);
2360  } else {
2361  meTrackResHighPTot_->Fill(t0Safe[*iTrack] - tsim);
2362  meTrackPullHighPTot_->Fill((t0Safe[*iTrack] - tsim) / sigmat0Safe[*iTrack]);
2363  }
2364 
2365  if (mtdQualMVA[(*iTrack)] < mvaL_) {
2366  meTrackZposRes_[0]->Fill(dZ);
2367  meTrack3DposRes_[0]->Fill(d3D);
2368  meTrackRes_[0]->Fill(t0Safe[*iTrack] - tsim);
2369  meTrackPull_[0]->Fill((t0Safe[*iTrack] - tsim) / sigmat0Safe[*iTrack]);
2370 
2371  if (optionalPlots_) {
2372  meTrackResMass_[0]->Fill(t0Safe[*iTrack] - tEst);
2373  meTrackResMassTrue_[0]->Fill(tEst - tsim);
2374  }
2375 
2376  if ((*iTrack)->p() <= 2) {
2377  meTrackResLowP_[0]->Fill(t0Safe[*iTrack] - tsim);
2378  meTrackPullLowP_[0]->Fill((t0Safe[*iTrack] - tsim) / sigmat0Safe[*iTrack]);
2379  } else if ((*iTrack)->p() > 2) {
2380  meTrackResHighP_[0]->Fill(t0Safe[*iTrack] - tsim);
2381  meTrackPullHighP_[0]->Fill((t0Safe[*iTrack] - tsim) / sigmat0Safe[*iTrack]);
2382  }
2383 
2384  if (optionalPlots_) {
2385  if (std::abs((*tp_info)->pdgId()) == 2212) {
2386  meTrackResMassProtons_[0]->Fill(t0Safe[*iTrack] - tEst);
2387  meTrackResMassTrueProtons_[0]->Fill(tEst - tsim);
2388  } else if (std::abs((*tp_info)->pdgId()) == 211) {
2389  meTrackResMassPions_[0]->Fill(t0Safe[*iTrack] - tEst);
2390  meTrackResMassTruePions_[0]->Fill(tEst - tsim);
2391  }
2392  }
2393 
2394  } else if (mtdQualMVA[(*iTrack)] > mvaL_ && mtdQualMVA[(*iTrack)] < mvaH_) {
2395  meTrackZposRes_[1]->Fill(dZ);
2396  meTrack3DposRes_[1]->Fill(d3D);
2397  meTrackRes_[1]->Fill(t0Safe[*iTrack] - tsim);
2398  meTrackPull_[1]->Fill((t0Safe[*iTrack] - tsim) / sigmat0Safe[*iTrack]);
2399 
2400  if (optionalPlots_) {
2401  meTrackResMass_[1]->Fill(t0Safe[*iTrack] - tEst);
2402  meTrackResMassTrue_[1]->Fill(tEst - tsim);
2403  }
2404 
2405  if ((*iTrack)->p() <= 2) {
2406  meTrackResLowP_[1]->Fill(t0Safe[*iTrack] - tsim);
2407  meTrackPullLowP_[1]->Fill((t0Safe[*iTrack] - tsim) / sigmat0Safe[*iTrack]);
2408  } else if ((*iTrack)->p() > 2) {
2409  meTrackResHighP_[1]->Fill(t0Safe[*iTrack] - tsim);
2410  meTrackPullHighP_[1]->Fill((t0Safe[*iTrack] - tsim) / sigmat0Safe[*iTrack]);
2411  }
2412 
2413  if (optionalPlots_) {
2414  if (std::abs((*tp_info)->pdgId()) == 2212) {
2415  meTrackResMassProtons_[1]->Fill(t0Safe[*iTrack] - tEst);
2416  meTrackResMassTrueProtons_[1]->Fill(tEst - tsim);
2417  } else if (std::abs((*tp_info)->pdgId()) == 211) {
2418  meTrackResMassPions_[1]->Fill(t0Safe[*iTrack] - tEst);
2419  meTrackResMassTruePions_[1]->Fill(tEst - tsim);
2420  }
2421  }
2422 
2423  } else if (mtdQualMVA[(*iTrack)] > mvaH_) {
2424  meTrackZposRes_[2]->Fill(dZ);
2425  meTrack3DposRes_[2]->Fill(d3D);
2426  meTrackRes_[2]->Fill(t0Safe[*iTrack] - tsim);
2427  meTrackPull_[2]->Fill((t0Safe[*iTrack] - tsim) / sigmat0Safe[*iTrack]);
2428 
2429  if (optionalPlots_) {
2430  meTrackResMass_[2]->Fill(t0Safe[*iTrack] - tEst);
2431  meTrackResMassTrue_[2]->Fill(tEst - tsim);
2432  }
2433 
2434  if ((*iTrack)->p() <= 2) {
2435  meTrackResLowP_[2]->Fill(t0Safe[*iTrack] - tsim);
2436  meTrackPullLowP_[2]->Fill((t0Safe[*iTrack] - tsim) / sigmat0Safe[*iTrack]);
2437  } else if ((*iTrack)->p() > 2) {
2438  meTrackResHighP_[2]->Fill(t0Safe[*iTrack] - tsim);
2439  meTrackPullHighP_[2]->Fill((t0Safe[*iTrack] - tsim) / sigmat0Safe[*iTrack]);
2440  }
2441 
2442  if (optionalPlots_) {
2443  if (std::abs((*tp_info)->pdgId()) == 2212) {
2444  meTrackResMassProtons_[2]->Fill(t0Safe[*iTrack] - tEst);
2445  meTrackResMassTrueProtons_[2]->Fill(tEst - tsim);
2446  } else if (std::abs((*tp_info)->pdgId()) == 211) {
2447  meTrackResMassPions_[2]->Fill(t0Safe[*iTrack] - tEst);
2448  meTrackResMassTruePions_[2]->Fill(tEst - tsim);
2449  }
2450  }
2451  }
2452  } // if tp_info != nullptr && MatchCategory == 0
2453  } // loop on reco tracks
2454  if (selectedVtxMatching) {
2455  meVtxTrackMult_->Fill(log10(nt));
2456  mePUTrackRelMult_->Fill(static_cast<double>(PUnt) / nt);
2457  meFakeTrackRelMult_->Fill(static_cast<double>(Fakent) / nt);
2458  mePUTrackRelSumWnt_->Fill(PUsumWnt / sumWnt);
2459  mePUTrackRelSumWos_->Fill(PUsumWos / sumWos);
2460  meSecTrackRelSumWos_->Fill(SecsumWos / sumWos);
2461  meFakeTrackRelSumWos_->Fill(FakesumWos / sumWos);
2462  mePUTrackRelSumPt_->Fill(PUsumPt / sumPt);
2463  mePUTrackRelSumPt2_->Fill(PUsumPt2 / sumPt2);
2464 
2465  observablesFromJets(
2466  recotracks, massVector, categoryVector, "use_allTracks", nJets, sumEtJets, sumPt2Jets, metPt, sumPzJets);
2467  observablesFromJets(recotracks,
2468  massVector,
2469  categoryVector,
2470  "skip_Fake",
2471  nJetsnoFake,
2472  sumEtJetsnoFake,
2473  sumPt2JetsnoFake,
2474  metPtnoFake,
2475  sumPzJetsnoFake);
2476  observablesFromJets(recotracks,
2477  massVector,
2478  categoryVector,
2479  "skip_PU",
2480  nJetsnoPU,
2481  sumEtJetsnoPU,
2482  sumPt2JetsnoPU,
2483  metPtnoPU,
2484  sumPzJetsnoPU);
2485 
2486  meJetsPURelMult_->Fill(static_cast<double>(nJets - nJetsnoPU) / nJets);
2487  meJetsPURelHt_->Fill((sumEtJets - sumEtJetsnoPU) / sumEtJets);
2488  meJetsPURelSumPt2_->Fill((sumPt2Jets - sumPt2JetsnoPU) / sumPt2Jets);
2489  meJetsFakeRelSumPt2_->Fill((sumPt2Jets - sumPt2JetsnoFake) / sumPt2Jets);
2490  meJetsPURelMetPt_->Fill((metPt - metPtnoPU) / metPt);
2491  meJetsPURelSumPz_->Fill((sumPzJets - sumPzJetsnoPU) / sumPzJets);
2492 
2493  if (optionalPlots_) {
2494  mePUTrackMult_->Fill(PUnt);
2495  mePUTrackSumWnt_->Fill(log10(std::max(minThrSumWnt_, PUsumWnt)));
2496  mePUTrackSumWos_->Fill(log10(std::max(minThrSumWos_, PUsumWos)));
2497  meSecTrackSumWos_->Fill(log10(std::max(minThrSumWos_, SecsumWos)));
2498  mePUTrackSumPt_->Fill(log10(std::max(minThrSumPt_, PUsumPt)));
2499  mePUTrackSumPt2_->Fill(log10(std::max(minThrSumPt2_, PUsumPt2)));
2500 
2501  mePUTrackRelMultvsMult_->Fill(nt, static_cast<double>(PUnt) / nt);
2502  meFakeTrackRelMultvsMult_->Fill(nt, static_cast<double>(Fakent) / nt);
2503  mePUTrackRelSumWntvsSumWnt_->Fill(log10(std::max(minThrSumWnt_, sumWnt)), PUsumWnt / sumWnt);
2504  mePUTrackRelSumWosvsSumWos_->Fill(log10(std::max(minThrSumWos_, sumWos)), PUsumWos / sumWos);
2505  meSecTrackRelSumWosvsSumWos_->Fill(log10(std::max(minThrSumWos_, sumWos)), SecsumWos / sumWos);
2506  meFakeTrackRelSumWosvsSumWos_->Fill(log10(std::max(minThrSumWos_, sumWos)), FakesumWos / sumWos);
2507  mePUTrackRelSumPtvsSumPt_->Fill(log10(std::max(minThrSumPt_, sumPt)), PUsumPt / sumPt);
2508  mePUTrackRelSumPt2vsSumPt2_->Fill(log10(std::max(minThrSumPt2_, sumPt2)), PUsumPt2 / sumPt2);
2509 
2510  meJetsPUMult_->Fill(nJets - nJetsnoPU);
2511  meJetsPUHt_->Fill(log10(std::max(minThrSumPt_, sumEtJets - sumEtJetsnoPU)));
2512  meJetsPUSumPt2_->Fill(log10(std::max(minThrSumPt2_, sumPt2Jets - sumPt2JetsnoPU)));
2513  meJetsPUMetPt_->Fill(log10(std::max(minThrMetPt_, metPt - metPtnoPU)));
2514  meJetsPUSumPz_->Fill(log10(std::max(minThrSumPz_, std::abs(sumPzJets - sumPzJetsnoPU))));
2515 
2516  meJetsPURelMultvsMult_->Fill(nJets, static_cast<double>(nJets - nJetsnoPU) / nJets);
2517  meJetsPURelHtvsHt_->Fill(log10(std::max(minThrSumPt_, sumEtJets)), (sumEtJets - sumEtJetsnoPU) / sumEtJets);
2518  meJetsPURelSumPt2vsSumPt2_->Fill(log10(std::max(minThrSumPt2_, sumPt2Jets)),
2519  (sumPt2Jets - sumPt2JetsnoPU) / sumPt2Jets);
2520  meJetsFakeRelSumPt2vsSumPt2_->Fill(log10(std::max(minThrSumPt2_, sumPt2Jets)),
2521  (sumPt2Jets - sumPt2JetsnoFake) / sumPt2Jets);
2522  meJetsPURelMetPtvsMetPt_->Fill(log10(std::max(minThrMetPt_, metPt)), (metPt - metPtnoPU) / metPt);
2523  meJetsPURelSumPzvsSumPz_->Fill(log10(std::max(minThrSumPz_, std::abs(sumPzJets))),
2524  (sumPzJets - sumPzJetsnoPU) / sumPzJets);
2525  }
2526  if (selectedLV) {
2527  meVtxTrackRecLVMult_->Fill(log10(nt));
2528  mePUTrackRecLVRelMult_->Fill(static_cast<double>(PUnt) / nt);
2529  meFakeTrackRecLVRelMult_->Fill(static_cast<double>(Fakent) / nt);
2530  mePUTrackRecLVRelSumWnt_->Fill(PUsumWnt / sumWnt);
2531  mePUTrackRecLVRelSumWos_->Fill(PUsumWos / sumWos);
2532  meSecTrackRecLVRelSumWos_->Fill(SecsumWos / sumWos);
2533  meFakeTrackRecLVRelSumWos_->Fill(FakesumWos / sumWos);
2534  mePUTrackRecLVRelSumPt_->Fill(PUsumPt / sumPt);
2535  mePUTrackRecLVRelSumPt2_->Fill(PUsumPt2 / sumPt2);
2536 
2537  meJetsRecLVPURelMult_->Fill(static_cast<double>(nJets - nJetsnoPU) / nJets);
2538  meJetsRecLVPURelHt_->Fill((sumEtJets - sumEtJetsnoPU) / sumEtJets);
2539  meJetsRecLVPURelSumPt2_->Fill((sumPt2Jets - sumPt2JetsnoPU) / sumPt2Jets);
2540  meJetsRecLVFakeRelSumPt2_->Fill((sumPt2Jets - sumPt2JetsnoFake) / sumPt2Jets);
2541  meJetsRecLVPURelMetPt_->Fill((metPt - metPtnoPU) / metPt);
2542  meJetsRecLVPURelSumPz_->Fill((sumPzJets - sumPzJetsnoPU) / sumPzJets);
2543 
2544  LogTrace("Primary4DVertexValidation")
2545  << "#PUTrks = " << PUnt << " #Trks = " << nt << " PURelMult = " << std::setprecision(3)
2546  << static_cast<double>(PUnt) / nt;
2547  LogTrace("Primary4DVertexValidation")
2548  << "PUsumWnt = " << std::setprecision(3) << PUsumWnt << " sumWnt = " << std::setprecision(3) << sumWnt
2549  << " PURelsumWnt = " << std::setprecision(3) << PUsumWnt / sumWnt;
2550  LogTrace("Primary4DVertexValidation")
2551  << "PUsumWos = " << std::setprecision(3) << PUsumWos << " sumWos = " << std::setprecision(3) << sumWos
2552  << " PURelsumWos = " << std::setprecision(3) << PUsumWos / sumWos;
2553  LogTrace("Primary4DVertexValidation")
2554  << "PuSumPt = " << std::setprecision(3) << PUsumPt << " SumPt = " << std::setprecision(4) << sumPt
2555  << " PURelSumPt = " << std::setprecision(3) << PUsumPt / sumPt;
2556  LogTrace("Primary4DVertexValidation")
2557  << "PuSumPt2 = " << std::setprecision(3) << PUsumPt2 << " SumPt2 = " << std::setprecision(4) << sumPt2
2558  << " PURelSumPt2 = " << std::setprecision(3) << PUsumPt2 / sumPt2;
2559  if (optionalPlots_) {
2560  mePUTrackRecLVMult_->Fill(PUnt);
2561  mePUTrackRecLVSumWnt_->Fill(log10(std::max(minThrSumWnt_, PUsumWnt)));
2562  mePUTrackRecLVSumWos_->Fill(log10(std::max(minThrSumWos_, PUsumWos)));
2563  meSecTrackRecLVSumWos_->Fill(log10(std::max(minThrSumWos_, PUsumWos)));
2564  mePUTrackRecLVSumPt_->Fill(log10(std::max(minThrSumPt_, PUsumPt)));
2565  mePUTrackRecLVSumPt2_->Fill(log10(std::max(minThrSumPt2_, PUsumPt2)));
2566 
2567  mePUTrackRecLVRelMultvsMult_->Fill(nt, static_cast<double>(PUnt) / nt);
2568  meFakeTrackRecLVRelMultvsMult_->Fill(nt, static_cast<double>(Fakent) / nt);
2569  mePUTrackRecLVRelSumWntvsSumWnt_->Fill(log10(std::max(minThrSumWnt_, sumWnt)), PUsumWnt / sumWnt);
2570  mePUTrackRecLVRelSumWosvsSumWos_->Fill(log10(std::max(minThrSumWos_, sumWos)), PUsumWos / sumWos);
2571  meSecTrackRecLVRelSumWosvsSumWos_->Fill(log10(std::max(minThrSumWos_, sumWos)), SecsumWos / sumWos);
2572  meFakeTrackRecLVRelSumWosvsSumWos_->Fill(log10(std::max(minThrSumWos_, sumWos)), FakesumWos / sumWos);
2573  mePUTrackRecLVRelSumPtvsSumPt_->Fill(log10(std::max(minThrSumPt_, sumPt)), PUsumPt / sumPt);
2574  mePUTrackRecLVRelSumPt2vsSumPt2_->Fill(log10(std::max(minThrSumPt2_, sumPt2)), PUsumPt2 / sumPt2);
2575 
2576  meJetsRecLVPUMult_->Fill(nJets - nJetsnoPU);
2577  meJetsRecLVPUHt_->Fill(log10(std::max(minThrSumPt_, sumEtJets - sumEtJetsnoPU)));
2578  meJetsRecLVPUSumPt2_->Fill(log10(std::max(minThrSumPt2_, sumPt2Jets - sumPt2JetsnoPU)));
2579  meJetsRecLVPUMetPt_->Fill(log10(std::max(minThrMetPt_, metPt - metPtnoPU)));
2580  meJetsRecLVPUSumPz_->Fill(log10(std::max(minThrSumPz_, std::abs(sumPzJets - sumPzJetsnoPU))));
2581 
2582  meJetsRecLVPURelMultvsMult_->Fill(nJets, static_cast<double>(nJets - nJetsnoPU) / nJets);
2583  meJetsRecLVPURelHtvsHt_->Fill(log10(std::max(minThrSumPt_, sumEtJets)),
2584  (sumEtJets - sumEtJetsnoPU) / sumEtJets);
2585  meJetsRecLVPURelSumPt2vsSumPt2_->Fill(log10(std::max(minThrSumPt2_, sumPt2Jets)),
2586  (sumPt2Jets - sumPt2JetsnoPU) / sumPt2Jets);
2587  meJetsRecLVFakeRelSumPt2vsSumPt2_->Fill(log10(std::max(minThrSumPt2_, sumPt2Jets)),
2588  (sumPt2Jets - sumPt2JetsnoFake) / sumPt2Jets);
2589  meJetsRecLVPURelMetPtvsMetPt_->Fill(log10(std::max(minThrMetPt_, metPt)), (metPt - metPtnoPU) / metPt);
2590  meJetsRecLVPURelSumPzvsSumPz_->Fill(log10(std::max(minThrSumPz_, std::abs(sumPzJets))),
2591  (sumPzJets - sumPzJetsnoPU) / sumPzJets);
2592  }
2593  }
2594  }
2595  } // loop on simpv
2596  } // ndof
2597  } // loop on recopv
2598 
2599  int real = 0;
2600  int fake = 0;
2601  int other_fake = 0;
2602  int split = 0;
2603 
2604  meRecVerNumber_->Fill(recopv.size());
2605  for (unsigned int ir = 0; ir < recopv.size(); ir++) {
2606  if (recopv.at(ir).ndof > selNdof_) {
2607  meRecPVZ_->Fill(recopv.at(ir).z);
2608  if (recopv.at(ir).recVtx->tError() > 0.) {
2609  meRecPVT_->Fill(recopv.at(ir).recVtx->t());
2610  }
2611  LogTrace("Primary4DVertexValidation") << "************* IR: " << ir;
2612  LogTrace("Primary4DVertexValidation") << "is_real: " << recopv.at(ir).is_real();
2613  LogTrace("Primary4DVertexValidation") << "is_fake: " << recopv.at(ir).is_fake();
2614  LogTrace("Primary4DVertexValidation") << "is_signal: " << recopv.at(ir).is_signal();
2615  LogTrace("Primary4DVertexValidation") << "split_from: " << recopv.at(ir).split_from();
2616  LogTrace("Primary4DVertexValidation") << "other fake: " << recopv.at(ir).other_fake();
2617  if (recopv.at(ir).is_real()) {
2618  real++;
2619  }
2620  if (recopv.at(ir).is_fake()) {
2621  fake++;
2622  }
2623  if (recopv.at(ir).other_fake()) {
2624  other_fake++;
2625  }
2626  if (recopv.at(ir).split_from() != -1) {
2627  split++;
2628  }
2629  } // ndof
2630  }
2631 
2632  LogTrace("Primary4DVertexValidation") << "is_real: " << real;
2633  LogTrace("Primary4DVertexValidation") << "is_fake: " << fake;
2634  LogTrace("Primary4DVertexValidation") << "split_from: " << split;
2635  LogTrace("Primary4DVertexValidation") << "other fake: " << other_fake;
2636  mePUvsRealV_->Fill(simpv.size(), real);
2637  mePUvsFakeV_->Fill(simpv.size(), fake);
2638  mePUvsOtherFakeV_->Fill(simpv.size(), other_fake);
2639  mePUvsSplitV_->Fill(simpv.size(), split);
2640 
2641  // fill vertices histograms here in a new loop
2642  meSimVerNumber_->Fill(simpv.size());
2643  for (unsigned int is = 0; is < simpv.size(); is++) {
2644  meSimPVZ_->Fill(simpv.at(is).z);
2645  meSimPVT_->Fill(simpv.at(is).t * simUnit_);
2646  meSimPVTvsZ_->Fill(simpv.at(is).z, simpv.at(is).t * simUnit_);
2647  if (is == 0 && optionalPlots_) {
2648  meSimPosInSimOrigCollection_->Fill(simpv.at(is).OriginalIndex);
2649  }
2650 
2651  if (simpv.at(is).rec == NOT_MATCHED) {
2652  LogTrace("Primary4DVertexValidation") << "sim vertex: " << is << " is not matched with any reco";
2653  continue;
2654  }
2655 
2656  for (unsigned int ir = 0; ir < recopv.size(); ir++) {
2657  if (recopv.at(ir).ndof > selNdof_) {
2658  if (recopv.at(ir).sim == is && simpv.at(is).rec == ir) {
2659  meTimeRes_->Fill(recopv.at(ir).recVtx->t() - simpv.at(is).t * simUnit_);
2660  meTimePull_->Fill((recopv.at(ir).recVtx->t() - simpv.at(is).t * simUnit_) / recopv.at(ir).recVtx->tError());
2661  meMatchQual_->Fill(recopv.at(ir).matchQuality - 0.5);
2662  if (ir == 0) { // signal vertex plots
2663  meTimeSignalRes_->Fill(recopv.at(ir).recVtx->t() - simpv.at(is).t * simUnit_);
2664  meTimeSignalPull_->Fill((recopv.at(ir).recVtx->t() - simpv.at(is).t * simUnit_) /
2665  recopv.at(ir).recVtx->tError());
2666  if (optionalPlots_) {
2667  meRecoPosInSimCollection_->Fill(recopv.at(ir).sim);
2668  meRecoPosInRecoOrigCollection_->Fill(recopv.at(ir).OriginalIndex);
2669  }
2670  }
2671  if (simpv.at(is).eventId.bunchCrossing() == 0 && simpv.at(is).eventId.event() == 0) {
2672  if (!recopv.at(ir).is_signal()) {
2673  edm::LogWarning("Primary4DVertexValidation")
2674  << "Reco vtx leading match inconsistent: BX/ID " << simpv.at(is).eventId.bunchCrossing() << " "
2675  << simpv.at(is).eventId.event();
2676  }
2677  meRecoPVPosSignal_->Fill(
2678  recopv.at(ir).OriginalIndex); // position in reco vtx correction associated to sim signal
2679  if (!signal_is_highest_pt) {
2680  meRecoPVPosSignalNotHighestPt_->Fill(
2681  recopv.at(ir).OriginalIndex); // position in reco vtx correction associated to sim signal
2682  }
2683  }
2684  LogTrace("Primary4DVertexValidation") << "*** Matching RECO: " << ir << "with SIM: " << is << " ***";
2685  LogTrace("Primary4DVertexValidation") << "Match Quality is " << recopv.at(ir).matchQuality;
2686  LogTrace("Primary4DVertexValidation") << "****";
2687  }
2688  } // ndof
2689  }
2690  }
2691 
2692  // dz histos
2693  for (unsigned int iv = 0; iv < recVtxs->size() - 1; iv++) {
2694  if (recVtxs->at(iv).ndof() > selNdof_) {
2695  double mindistance_realreal = 1e10;
2696 
2697  for (unsigned int jv = iv; jv < recVtxs->size(); jv++) {
2698  if ((!(jv == iv)) && select(recVtxs->at(jv))) {
2699  double dz = recVtxs->at(iv).z() - recVtxs->at(jv).z();
2700  double dtsigma = std::sqrt(recVtxs->at(iv).covariance(3, 3) + recVtxs->at(jv).covariance(3, 3));
2701  double dt = (std::abs(dz) <= deltaZcut_ && dtsigma > 0.)
2702  ? (recVtxs->at(iv).t() - recVtxs->at(jv).t()) / dtsigma
2703  : -9999.;
2704  if (recopv.at(iv).is_real() && recopv.at(jv).is_real()) {
2705  meDeltaZrealreal_->Fill(std::abs(dz));
2706  if (dt != -9999.) {
2707  meDeltaTrealreal_->Fill(std::abs(dt));
2708  }
2709  if (std::abs(dz) < std::abs(mindistance_realreal)) {
2710  mindistance_realreal = dz;
2711  }
2712  } else if (recopv.at(iv).is_fake() && recopv.at(jv).is_fake()) {
2713  meDeltaZfakefake_->Fill(std::abs(dz));
2714  if (dt != -9999.) {
2715  meDeltaTfakefake_->Fill(std::abs(dt));
2716  }
2717  }
2718  }
2719  }
2720 
2721  double mindistance_fakereal = 1e10;
2722  double mindistance_realfake = 1e10;
2723  for (unsigned int jv = 0; jv < recVtxs->size(); jv++) {
2724  if ((!(jv == iv)) && select(recVtxs->at(jv))) {
2725  double dz = recVtxs->at(iv).z() - recVtxs->at(jv).z();
2726  double dtsigma = std::sqrt(recVtxs->at(iv).covariance(3, 3) + recVtxs->at(jv).covariance(3, 3));
2727  double dt = (std::abs(dz) <= deltaZcut_ && dtsigma > 0.)
2728  ? (recVtxs->at(iv).t() - recVtxs->at(jv).t()) / dtsigma
2729  : -9999.;
2730  if (recopv.at(iv).is_fake() && recopv.at(jv).is_real()) {
2731  meDeltaZfakereal_->Fill(std::abs(dz));
2732  if (dt != -9999.) {
2733  meDeltaTfakereal_->Fill(std::abs(dt));
2734  }
2735  if (std::abs(dz) < std::abs(mindistance_fakereal)) {
2736  mindistance_fakereal = dz;
2737  }
2738  }
2739 
2740  if (recopv.at(iv).is_real() && recopv.at(jv).is_fake() && (std::abs(dz) < std::abs(mindistance_realfake))) {
2741  mindistance_realfake = dz;
2742  }
2743  }
2744  }
2745  } // ndof
2746  }
2747 
2748 } // end of analyze
2749 
2750 void Primary4DVertexValidation::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
2751  edm::ParameterSetDescription desc;
2752 
2753  desc.add<std::string>("folder", "MTD/Vertices");
2754  desc.add<edm::InputTag>("TPtoRecoTrackAssoc", edm::InputTag("trackingParticleRecoTrackAsssociation"));
2755  desc.add<edm::InputTag>("mtdTracks", edm::InputTag("trackExtenderWithMTD"));
2756  desc.add<edm::InputTag>("SimTag", edm::InputTag("mix", "MergedTrackTruth"));
2757  desc.add<edm::InputTag>("offlineBS", edm::InputTag("offlineBeamSpot"));
2758  desc.add<edm::InputTag>("offline4DPV", edm::InputTag("offlinePrimaryVertices4D"));
2759  desc.add<edm::InputTag>("trackAssocSrc", edm::InputTag("trackExtenderWithMTD:generalTrackassoc"))
2760  ->setComment("Association between General and MTD Extended tracks");
2761  desc.add<edm::InputTag>("pathLengthSrc", edm::InputTag("trackExtenderWithMTD:generalTrackPathLength"));
2762  desc.add<edm::InputTag>("momentumSrc", edm::InputTag("trackExtenderWithMTD:generalTrackp"));
2763  desc.add<edm::InputTag>("tmtd", edm::InputTag("trackExtenderWithMTD:generalTracktmtd"));
2764  desc.add<edm::InputTag>("timeSrc", edm::InputTag("trackExtenderWithMTD:generalTracktmtd"));
2765  desc.add<edm::InputTag>("sigmaSrc", edm::InputTag("trackExtenderWithMTD:generalTracksigmatmtd"));
2766  desc.add<edm::InputTag>("t0PID", edm::InputTag("tofPID:t0"));
2767  desc.add<edm::InputTag>("sigmat0PID", edm::InputTag("tofPID:sigmat0"));
2768  desc.add<edm::InputTag>("t0SafePID", edm::InputTag("tofPID:t0safe"));
2769  desc.add<edm::InputTag>("sigmat0SafePID", edm::InputTag("tofPID:sigmat0safe"));
2770  desc.add<edm::InputTag>("trackMVAQual", edm::InputTag("mtdTrackQualityMVA:mtdQualMVA"));
2771  desc.add<edm::InputTag>("tofPi", edm::InputTag("trackExtenderWithMTD:generalTrackTofPi"));
2772  desc.add<edm::InputTag>("tofK", edm::InputTag("trackExtenderWithMTD:generalTrackTofK"));
2773  desc.add<edm::InputTag>("tofP", edm::InputTag("trackExtenderWithMTD:generalTrackTofP"));
2774  desc.add<edm::InputTag>("probPi", edm::InputTag("tofPID:probPi"));
2775  desc.add<edm::InputTag>("probK", edm::InputTag("tofPID:probK"));
2776  desc.add<edm::InputTag>("probP", edm::InputTag("tofPID:probP"));
2777  desc.add<bool>("useOnlyChargedTracks", true);
2778  desc.addUntracked<bool>("optionalPlots", false);
2779  desc.add<bool>("use3dNoTime", false);
2780  desc.add<double>("trackweightTh", 0.5);
2781  desc.add<double>("mvaTh", 0.8);
2782  desc.add<double>("minProbHeavy", 0.75);
2783  descriptions.add("vertices4DValid", desc);
2784 }
2785 
2786 void Primary4DVertexValidation::printMatchedRecoTrackInfo(const reco::Vertex& vtx,
2787  const reco::TrackBaseRef& trk,
2788  const TrackingParticleRef& tp,
2789  const unsigned int& categ) {
2790  std::string strTrk;
2791  switch (categ) {
2792  case 0:
2793  strTrk = "Reco_Track:";
2794  break;
2795  case 1:
2796  strTrk = "SecRecoTrk:";
2797  break;
2798  case 2:
2799  strTrk = "PU_RecoTrk:";
2800  break;
2801  }
2802  LogTrace("Primary4DVertexValidation") << strTrk << " w =" << std::setw(6) << std::setprecision(2)
2803  << vtx.trackWeight(trk) << " pt =" << std::setw(6) << std::setprecision(2)
2804  << trk->pt() << " eta =" << std::setw(6) << std::setprecision(2) << trk->eta()
2805  << " MatchedTP: Pt =" << std::setw(6) << std::setprecision(2) << tp->pt()
2806  << " eta =" << std::setw(6) << std::setprecision(2) << tp->eta()
2807  << " Parent vtx: z =" << std::setw(8) << std::setprecision(4)
2808  << tp->parentVertex()->position().z() << " t =" << std::setw(8)
2809  << std::setprecision(4) << tp->parentVertex()->position().t() * simUnit_
2810  << " BX =" << tp->parentVertex()->eventId().bunchCrossing()
2811  << " ev =" << tp->parentVertex()->eventId().event() << std::endl;
2812 }
2813 
2814 void Primary4DVertexValidation::printSimVtxRecoVtxInfo(
2815  const struct Primary4DVertexValidation::simPrimaryVertex& simpVtx,
2816  const struct Primary4DVertexValidation::recoPrimaryVertex& recopVtx) {
2817  LogTrace("Primary4DVertexValidation") << "Sim vtx (x,y,z,t) = (" << std::setprecision(4) << simpVtx.x << ","
2818  << std::setprecision(4) << simpVtx.y << "," << std::setprecision(4) << simpVtx.z
2819  << "," << std::setprecision(4) << simpVtx.t * simUnit_ << ")"
2820  << " Simvtx.rec = " << simpVtx.rec;
2821  LogTrace("Primary4DVertexValidation") << "Sim vtx: pt = " << std::setprecision(4) << simpVtx.pt
2822  << " ptsq = " << std::setprecision(6) << simpVtx.ptsq
2823  << " nGenTrk = " << simpVtx.nGenTrk
2824  << " nmatch recotrks = " << simpVtx.num_matched_reco_tracks;
2825  LogTrace("Primary4DVertexValidation") << "Reco vtx (x,y,z) = (" << std::setprecision(4) << recopVtx.x << ","
2826  << std::setprecision(4) << recopVtx.y << "," << std::setprecision(4)
2827  << recopVtx.z << ")"
2828  << " Recovtx.sim = " << recopVtx.sim;
2829  LogTrace("Primary4DVertexValidation") << "Reco vtx: pt = " << std::setprecision(4) << recopVtx.pt
2830  << " ptsq = " << std::setprecision(6) << recopVtx.ptsq
2831  << " nrecotrks = " << recopVtx.nRecoTrk
2832  << " nmatch simtrks = " << recopVtx.num_matched_sim_tracks;
2833  LogTrace("Primary4DVertexValidation") << "wnt " << recopVtx.sumwnt << " wos = " << recopVtx.sumwos;
2834  for (auto iTP = simpVtx.sim_vertex->daughterTracks_begin(); iTP != simpVtx.sim_vertex->daughterTracks_end(); ++iTP) {
2835  if (use_only_charged_tracks_ && (**iTP).charge() == 0) {
2836  continue;
2837  }
2838  LogTrace("Primary4DVertexValidation")
2839  << "Daughter track of sim vertex: pt =" << std::setw(6) << std::setprecision(2) << (*iTP)->pt()
2840  << " eta =" << std::setw(6) << std::setprecision(2) << (*iTP)->eta();
2841  }
2842 }
2843 
2844 const bool Primary4DVertexValidation::trkTPSelLV(const TrackingParticle& tp) {
2845  bool match = false;
2846  if (tp.status() != 1) {
2847  return match;
2848  }
2849 
2850  auto x_pv = tp.parentVertex()->position().x();
2851  auto y_pv = tp.parentVertex()->position().y();
2852  auto z_pv = tp.parentVertex()->position().z();
2853 
2854  auto r_pv = std::sqrt(x_pv * x_pv + y_pv * y_pv);
2855 
2856  match =
2857  tp.charge() != 0 && std::abs(tp.eta()) < etacutGEN_ && tp.pt() > pTcut_ && r_pv < rBTL_ && std::abs(z_pv) < zETL_;
2858  return match;
2859 }
2860 
2861 const bool Primary4DVertexValidation::trkRecSel(const reco::TrackBase& trk) {
2862  bool match = false;
2863  match = std::abs(trk.eta()) <= etacutREC_ && trk.pt() > pTcut_;
2864  return match;
2865 }
2866 
2867 DEFINE_FWK_MODULE(Primary4DVertexValidation);
Primary4DVertexValidation::trackMaxBtlEta_
static constexpr double trackMaxBtlEta_
Definition: Primary4DVertexValidation.cc:274
dt
float dt
Definition: AMPTWrapper.h:136
Primary4DVertexValidation::Rec4DVerToken_
edm::EDGetTokenT< edm::View< reco::Vertex > > Rec4DVerToken_
Definition: Primary4DVertexValidation.cc:308
ProducerED_cfi.InputTag
InputTag
Definition: ProducerED_cfi.py:5
deDxTools::esConsumes
ESGetTokenH3DDVariant esConsumes(std::string const &Record, edm::ConsumesCollector &)
Definition: DeDxTools.cc:283
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
Primary4DVertexValidation::trackingVertexCollectionToken_
edm::EDGetTokenT< TrackingVertexCollection > trackingVertexCollectionToken_
Definition: Primary4DVertexValidation.cc:304
bphysicsOniaDQM_cfi.vertex
vertex
Definition: bphysicsOniaDQM_cfi.py:7
Primary4DVertexValidation::trackMaxEtlEta_
static constexpr double trackMaxEtlEta_
Definition: Primary4DVertexValidation.cc:276
genfragment_ptgun_cfg.PartID
PartID
Definition: genfragment_ptgun_cfg.py:6
mps_fire.i
i
Definition: mps_fire.py:429
Primary4DVertexValidation::etacutREC_
static constexpr double etacutREC_
Definition: Primary4DVertexValidation.cc:271
folder_
std::string folder_
Definition: DQMEDAnalyzer.cc:60
muonTagProbeFilters_cff.matched
matched
Definition: muonTagProbeFilters_cff.py:62
Primary4DVertexValidation::mePUTrackRecLVRelSumWntvsSumWnt_
MonitorElement * mePUTrackRecLVRelSumWntvsSumWnt_
Definition: Primary4DVertexValidation.cc:428
Primary4DVertexValidation::trackingParticleCollectionToken_
edm::EDGetTokenT< TrackingParticleCollection > trackingParticleCollectionToken_
Definition: Primary4DVertexValidation.cc:303
dqm::implementation::NavigatorBase::setCurrentFolder
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:36
Primary4DVertexValidation::sigmat0SafePidToken_
edm::EDGetTokenT< edm::ValueMap< float > > sigmat0SafePidToken_
Definition: Primary4DVertexValidation.cc:318
ParticleDataTable
HepPDT::ParticleDataTable ParticleDataTable
Definition: ParticleDataTable.h:8
Primary4DVertexValidation::getSimPVs
std::vector< Primary4DVertexValidation::simPrimaryVertex > getSimPVs(const edm::Handle< TrackingVertexCollection > &)
Definition: Primary4DVertexValidation.cc:1588
Primary4DVertexValidation::RecBeamSpotToken_
edm::EDGetTokenT< reco::BeamSpot > RecBeamSpotToken_
Definition: Primary4DVertexValidation.cc:307
edm::isNotFinite
constexpr bool isNotFinite(T x)
Definition: isFinite.h:9
Primary4DVertexValidation::trkTPSelLV
const bool trkTPSelLV(const TrackingParticle &)
Definition: Primary4DVertexValidation.cc:2844
WZElectronSkims53X_cff.max
max
Definition: WZElectronSkims53X_cff.py:197
Primary4DVertexValidation::trkRecSel
const bool trkRecSel(const reco::TrackBase &)
Definition: Primary4DVertexValidation.cc:2861
Primary4DVertexValidation::select
bool select(const reco::Vertex &, int level=0)
Definition: Primary4DVertexValidation.cc:1440
Primary4DVertexValidation::matchReco2Sim
void matchReco2Sim(std::vector< recoPrimaryVertex > &, std::vector< simPrimaryVertex > &, const edm::ValueMap< float > &, const edm::ValueMap< float > &, const edm::Handle< reco::BeamSpot > &)
Definition: Primary4DVertexValidation.cc:1754
Primary4DVertexValidation::observablesFromJets
void observablesFromJets(const std::vector< reco::Track > &, const std::vector< double > &, const std::vector< int > &, const std::string &, unsigned int &, double &, double &, double &, double &)
Definition: Primary4DVertexValidation.cc:1459
Primary4DVertexValidation::probKToken_
edm::EDGetTokenT< edm::ValueMap< float > > probKToken_
Definition: Primary4DVertexValidation.cc:326
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
reco::VertexBaseRef
edm::RefToBase< reco::Vertex > VertexBaseRef
persistent reference to a Vertex, using views
Definition: VertexFwd.h:21
edm::Handle::product
T const * product() const
Definition: Handle.h:70
L1BJetProducer_cff.vtx
vtx
Definition: L1BJetProducer_cff.py:10
Primary4DVertexValidation::matchRecoTrack2SimSignal
bool matchRecoTrack2SimSignal(const reco::TrackBaseRef &)
Definition: Primary4DVertexValidation.cc:1382
caHitNtupletGeneratorKernels::iev
TkSoAView< TrackerTraits > HitToTuple< TrackerTraits > const *__restrict__ int32_t int32_t int iev
Definition: CAHitNtupletGeneratorKernelsImpl.h:860
ALPAKA_ACCELERATOR_NAMESPACE::vertexFinder::iv
int32_t *__restrict__ iv
Definition: fitVertices.h:39
Primary4DVertexValidation::etacutGEN_
static constexpr double etacutGEN_
Definition: Primary4DVertexValidation.cc:270
Primary4DVertexValidation::sigmat0PidToken_
edm::EDGetTokenT< edm::ValueMap< float > > sigmat0PidToken_
Definition: Primary4DVertexValidation.cc:316
l1ctLayer1_cff.dZ
dZ
Definition: l1ctLayer1_cff.py:69
PDWG_EXODelayedJetMET_cff.jets
jets
Definition: PDWG_EXODelayedJetMET_cff.py:14
Primary4DVertexValidation::minThrSumWos_
static constexpr double minThrSumWos_
Definition: Primary4DVertexValidation.cc:279
Primary4DVertexValidation::meFakeTrackRecLVRelSumWosvsSumWos_
MonitorElement * meFakeTrackRecLVRelSumWosvsSumWos_
Definition: Primary4DVertexValidation.cc:431
MillePedeFileConverter_cfg.e
e
Definition: MillePedeFileConverter_cfg.py:37
pfDeepBoostedJetPreprocessParams_cfi.sv
sv
Definition: pfDeepBoostedJetPreprocessParams_cfi.py:352
Primary4DVertexValidation::meJetsRecLVPURelSumPt2vsSumPt2_
MonitorElement * meJetsRecLVPURelSumPt2vsSumPt2_
Definition: Primary4DVertexValidation.cc:468
ALPAKA_ACCELERATOR_NAMESPACE::brokenline::constexpr
if constexpr(n > 3)
Definition: BrokenLine.h:164
nJets
static constexpr int nJets
Definition: Phase2L1CaloJetEmulator.h:15
cms::cuda::assert
assert(be >=bs)
Primary4DVertexValidation::analyze
void analyze(const edm::Event &, const edm::EventSetup &) override
Definition: Primary4DVertexValidation.cc:1978
Primary4DVertexValidation::timeToken_
edm::EDGetTokenT< edm::ValueMap< float > > timeToken_
Definition: Primary4DVertexValidation.cc:313
Primary4DVertexValidation::getRecoPVs
std::vector< Primary4DVertexValidation::recoPrimaryVertex > getRecoPVs(const edm::Handle< edm::View< reco::Vertex >> &)
Definition: Primary4DVertexValidation.cc:1690
hcalRecHitTable_cff.time
time
Definition: hcalRecHitTable_cff.py:14
Primary4DVertexValidation::trackMinEtlEta_
static constexpr double trackMinEtlEta_
Definition: Primary4DVertexValidation.cc:275
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:256
edm::AssociationMap::find
const_iterator find(const key_type &k) const
find element with specified reference key
Definition: AssociationMap.h:175
LHEGenericFilter_cfi.ParticleID
ParticleID
Definition: LHEGenericFilter_cfi.py:6
LogTrace
#define LogTrace(id)
Definition: MessageLogger.h:242
Primary4DVertexValidation::meTrackResMassProtons_
MonitorElement * meTrackResMassProtons_[3]
Definition: Primary4DVertexValidation.cc:484
Primary4DVertexValidation::meTrackResMassTrueProtons_
MonitorElement * meTrackResMassTrueProtons_[3]
Definition: Primary4DVertexValidation.cc:485
edm::AssociationMap::end
const_iterator end() const
last iterator over the map (read only)
Definition: AssociationMap.h:173
dqm::impl::MonitorElement::Fill
void Fill(long long x)
Definition: MonitorElement.h:292
edm::second
U second(std::pair< T, U > const &p)
Definition: ParameterSet.cc:222
math::XYZTLorentzVector
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
reco::TrackBase::pt
double pt() const
track transverse momentum
Definition: TrackBase.h:637
Primary4DVertexValidation::timeFromTrueMass
double timeFromTrueMass(double, double, double, double)
Definition: Primary4DVertexValidation.cc:1428
iEvent
int iEvent
Definition: GenABIO.cc:224
gammaJetAnalysis_CHSJECs_cff.level
level
Definition: gammaJetAnalysis_CHSJECs_cff.py:6
Primary4DVertexValidation::LorentzVector
math::XYZTLorentzVector LorentzVector
Definition: Primary4DVertexValidation.cc:67
Primary4DVertexValidation::meSecTrackRecLVRelSumWosvsSumWos_
MonitorElement * meSecTrackRecLVRelSumWosvsSumWos_
Definition: Primary4DVertexValidation.cc:430
Primary4DVertexValidation::simPrimaryVertex::addTrack
void addTrack(unsigned int irecv, double twos, double twt)
Definition: Primary4DVertexValidation.cc:113
dqm::implementation::IBooker::bookProfile
MonitorElement * bookProfile(TString const &name, TString const &title, int nchX, double lowX, double highX, int, double lowY, double highY, char const *option="s", FUNC onbooking=NOOP())
Definition: DQMStore.h:408
reco::TrackBase::dzError
double dzError() const
error on dz
Definition: TrackBase.h:778
V0Monitor_cff.v0
v0
Definition: V0Monitor_cff.py:7
Primary4DVertexValidation::fillDescriptions
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
Definition: Primary4DVertexValidation.cc:2750
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:23
Primary4DVertexValidation::Primary4DVertexValidation
Primary4DVertexValidation(const edm::ParameterSet &)
Definition: Primary4DVertexValidation.cc:525
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
reco::Vertex::trackWeight
float trackWeight(const TREF &r) const
returns the weight with which a Track has contributed to the vertex-fit.
Definition: Vertex.h:97
Primary4DVertexValidation::probPiToken_
edm::EDGetTokenT< edm::ValueMap< float > > probPiToken_
Definition: Primary4DVertexValidation.cc:325
Primary4DVertexValidation::minThrMetPt_
static constexpr double minThrMetPt_
Definition: Primary4DVertexValidation.cc:282
funct::tan
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
Primary4DVertexValidation::mePUTrackRecLVRelSumPt2vsSumPt2_
MonitorElement * mePUTrackRecLVRelSumPt2vsSumPt2_
Definition: Primary4DVertexValidation.cc:433
Primary4DVertexValidation::probPToken_
edm::EDGetTokenT< edm::ValueMap< float > > probPToken_
Definition: Primary4DVertexValidation.cc:327
Primary4DVertexValidation::meJetsRecLVFakeRelSumPt2vsSumPt2_
MonitorElement * meJetsRecLVFakeRelSumPt2vsSumPt2_
Definition: Primary4DVertexValidation.cc:469
reco::TrackBase::phi
double phi() const
azimuthal angle of momentum vector
Definition: TrackBase.h:649
DEFINE_FWK_MODULE
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
Primary4DVertexValidation::momentumToken_
edm::EDGetTokenT< edm::ValueMap< float > > momentumToken_
Definition: Primary4DVertexValidation.cc:312
Primary4DVertexValidation::meTrackResMassTruePions_
MonitorElement * meTrackResMassTruePions_[3]
Definition: Primary4DVertexValidation.cc:487
Primary4DVertexValidation::tofKToken_
edm::EDGetTokenT< edm::ValueMap< float > > tofKToken_
Definition: Primary4DVertexValidation.cc:323
Primary4DVertexValidation::RecTrackToken_
edm::EDGetTokenT< reco::TrackCollection > RecTrackToken_
Definition: Primary4DVertexValidation.cc:299
HLT_2024v14_cff.distance
distance
Definition: HLT_2024v14_cff.py:5740
ParticleData
HepPDT::ParticleData ParticleData
Definition: ParticleDataTable.h:9
edm::EventSetup::getHandle
ESHandle< T > getHandle(const ESGetToken< T, R > &iToken) const
Definition: EventSetup.h:130
Primary4DVertexValidation::r2s_
const reco::RecoToSimCollection * r2s_
Definition: Primary4DVertexValidation.cc:296
Primary4DVertexValidation::vecPileupSummaryInfoToken_
edm::EDGetTokenT< std::vector< PileupSummaryInfo > > vecPileupSummaryInfoToken_
Definition: Primary4DVertexValidation.cc:301
nt
int nt
Definition: AMPTWrapper.h:42
Primary4DVertexValidation::simPrimaryVertex::simPrimaryVertex
simPrimaryVertex(double x1, double y1, double z1, double t1)
Definition: Primary4DVertexValidation.cc:71
Primary4DVertexValidation::deltaZcut_
static constexpr double deltaZcut_
Definition: Primary4DVertexValidation.cc:273
EgHLTOffHistBins_cfi.mass
mass
Definition: EgHLTOffHistBins_cfi.py:34
Primary4DVertexValidation::t0PidToken_
edm::EDGetTokenT< edm::ValueMap< float > > t0PidToken_
Definition: Primary4DVertexValidation.cc:315
TrackingVertexRef
edm::Ref< TrackingVertexCollection > TrackingVertexRef
Definition: TrackingVertexContainer.h:9
Primary4DVertexValidation::NOT_MATCHED
static constexpr unsigned int NOT_MATCHED
Definition: Primary4DVertexValidation.cc:261
reco::TrackBase::eta
double eta() const
pseudorapidity of momentum vector
Definition: TrackBase.h:652
Primary4DVertexValidation::simToRecoAssociationToken_
edm::EDGetTokenT< reco::SimToRecoCollection > simToRecoAssociationToken_
Definition: Primary4DVertexValidation.cc:305
Primary4DVertexValidation::trackAssocToken_
edm::EDGetTokenT< edm::ValueMap< int > > trackAssocToken_
Definition: Primary4DVertexValidation.cc:310
angle_units::operators::convertMmToCm
constexpr NumType convertMmToCm(NumType millimeters)
Definition: angle_units.h:44
Primary4DVertexValidation::minThrSumPz_
static constexpr double minThrSumPz_
Definition: Primary4DVertexValidation.cc:283
Primary4DVertexValidation::isParticle
void isParticle(const reco::TrackBaseRef &, const edm::ValueMap< float > &, const edm::ValueMap< float > &, const edm::ValueMap< float > &, const edm::ValueMap< float > &, const edm::ValueMap< float > &, unsigned int &, bool &, bool &, bool &, bool &)
Definition: Primary4DVertexValidation.cc:1532
edm::ConfigurationDescriptions::add
void add(std::string const &label, ParameterSetDescription const &psetDescription)
Definition: ConfigurationDescriptions.cc:57
nano_mu_local_reco_cff.bool
bool
Definition: nano_mu_local_reco_cff.py:14
Primary4DVertexValidation::minThrSumPt2_
static constexpr double minThrSumPt2_
Definition: Primary4DVertexValidation.cc:281
testProducerWithPsetDescEmpty_cfi.y1
y1
Definition: testProducerWithPsetDescEmpty_cfi.py:29
Primary4DVertexValidation::recoPrimaryVertex::addTrack
void addTrack(unsigned int iev, double twos, double wt)
Definition: Primary4DVertexValidation.cc:183
edm::HandleBase::isValid
bool isValid() const
Definition: HandleBase.h:70
Primary4DVertexValidation::printSimVtxRecoVtxInfo
void printSimVtxRecoVtxInfo(const struct Primary4DVertexValidation::simPrimaryVertex &, const struct Primary4DVertexValidation::recoPrimaryVertex &)
Definition: Primary4DVertexValidation.cc:2814
BTaggingMonitor_cfi.met
met
Definition: BTaggingMonitor_cfi.py:74
reco::TrackBase::theta
double theta() const
polar angle
Definition: TrackBase.h:602
Primary4DVertexValidation::getMatchedTP
std::pair< const edm::Ref< std::vector< TrackingParticle > > *, int > getMatchedTP(const reco::TrackBaseRef &, const TrackingVertexRef &)
Definition: Primary4DVertexValidation.cc:1399
reco
fixed size matrix
Definition: AlignmentAlgorithmBase.h:46
edm
HLT enums.
Definition: AlignableModifier.h:19
Primary4DVertexValidation::s2r_
const reco::SimToRecoCollection * s2r_
Definition: Primary4DVertexValidation.cc:297
ppsModifySingularModes_cfi.z1
z1
Definition: ppsModifySingularModes_cfi.py:8
Primary4DVertexValidation::tofPiToken_
edm::EDGetTokenT< edm::ValueMap< float > > tofPiToken_
Definition: Primary4DVertexValidation.cc:322
TrackingParticle
Monte Carlo truth information used for tracking validation.
Definition: TrackingParticle.h:29
Primary4DVertexValidation::bookHistograms
void bookHistograms(DQMStore::IBooker &i, edm::Run const &, edm::EventSetup const &) override
Definition: Primary4DVertexValidation.cc:568
Primary4DVertexValidation::recoToSimAssociationToken_
edm::EDGetTokenT< reco::RecoToSimCollection > recoToSimAssociationToken_
Definition: Primary4DVertexValidation.cc:306
match
std::pair< typename Association::data_type::first_type, double > match(Reference key, Association association, bool bestMatchByMaxValue)
Generic matching function.
Definition: Utils.h:10
gather_cfg.cout
cout
Definition: gather_cfg.py:144
Primary4DVertexValidation::pdtToken_
edm::ESGetToken< HepPDT::ParticleDataTable, edm::DefaultRecord > pdtToken_
Definition: Primary4DVertexValidation.cc:328
edm::LogWarning
Log< level::Warning, false > LogWarning
Definition: MessageLogger.h:128
dqm::implementation::IBooker::book1D
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:98
Primary4DVertexValidation::t0SafePidToken_
edm::EDGetTokenT< edm::ValueMap< float > > t0SafePidToken_
Definition: Primary4DVertexValidation.cc:317
Primary4DVertexValidation::tmtdToken_
edm::EDGetTokenT< edm::ValueMap< float > > tmtdToken_
Definition: Primary4DVertexValidation.cc:321
testProducerWithPsetDescEmpty_cfi.x1
x1
Definition: testProducerWithPsetDescEmpty_cfi.py:33
Primary4DVertexValidation::zWosMatchMax_
static constexpr double zWosMatchMax_
Definition: Primary4DVertexValidation.cc:269
Primary4DVertexValidation::mePUTrackRecLVRelSumWosvsSumWos_
MonitorElement * mePUTrackRecLVRelSumWosvsSumWos_
Definition: Primary4DVertexValidation.cc:429
Primary4DVertexValidation::minThrSumPt_
static constexpr double minThrSumPt_
Definition: Primary4DVertexValidation.cc:280
Primary4DVertexValidation::tofPToken_
edm::EDGetTokenT< edm::ValueMap< float > > tofPToken_
Definition: Primary4DVertexValidation.cc:324
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29
TtFullHadEvtBuilder_cfi.sumPt
sumPt
Definition: TtFullHadEvtBuilder_cfi.py:38
Primary4DVertexValidation::minThrSumWnt_
static constexpr double minThrSumWnt_
Definition: Primary4DVertexValidation.cc:278
Primary4DVertexValidation::getWosWnt
void getWosWnt(const reco::Vertex &, const reco::TrackBaseRef &, const edm::ValueMap< float > &, const edm::ValueMap< float > &, const edm::Handle< reco::BeamSpot > &, double &, double &)
Definition: Primary4DVertexValidation.cc:1562
edm::Run
Definition: Run.h:45
Primary4DVertexValidation::trackMVAQualToken_
edm::EDGetTokenT< edm::ValueMap< float > > trackMVAQualToken_
Definition: Primary4DVertexValidation.cc:319
Primary4DVertexValidation::printMatchedRecoTrackInfo
void printMatchedRecoTrackInfo(const reco::Vertex &, const reco::TrackBaseRef &, const TrackingParticleRef &, const unsigned int &)
Definition: Primary4DVertexValidation.cc:2786
Primary4DVertexValidation::pathLengthToken_
edm::EDGetTokenT< edm::ValueMap< float > > pathLengthToken_
Definition: Primary4DVertexValidation.cc:311
Generated for CMSSW Reference Manual by doxygen 1.8.14