Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-09-08 03:21:54

0001 #include "DQMOffline/RecoB/plugins/PrimaryVertexMonitor.h"
0002 #include "FWCore/ServiceRegistry/interface/Service.h"
0003 #include "FWCore/Framework/interface/MakerMacros.h"
0004 #include "FWCore/Utilities/interface/isFinite.h"
0005 
0006 #include "CommonTools/Statistics/interface/ChiSquaredProbability.h"
0007 
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 
0010 #include "TMath.h"
0011 
0012 using namespace reco;
0013 using namespace edm;
0014 
0015 PrimaryVertexMonitor::PrimaryVertexMonitor(const edm::ParameterSet& pSet)
0016     : conf_(pSet),
0017       TopFolderName_(pSet.getParameter<std::string>("TopFolderName")),
0018       AlignmentLabel_(pSet.getParameter<std::string>("AlignmentLabel")),
0019       ndof_(pSet.getParameter<int>("ndof")),
0020       useHPfoAlignmentPlots_(pSet.getParameter<bool>("useHPforAlignmentPlots")),
0021       errorPrinted_(false),
0022       nbvtx(nullptr),
0023       bsX(nullptr),
0024       bsY(nullptr),
0025       bsZ(nullptr),
0026       bsSigmaZ(nullptr),
0027       bsDxdz(nullptr),
0028       bsDydz(nullptr),
0029       bsBeamWidthX(nullptr),
0030       bsBeamWidthY(nullptr),
0031       bsType(nullptr),
0032       sumpt(nullptr),
0033       ntracks(nullptr),
0034       weight(nullptr),
0035       chi2ndf(nullptr),
0036       chi2prob(nullptr),
0037       dxy(nullptr),
0038       dxy2(nullptr),
0039       dz(nullptr),
0040       dxyErr(nullptr),
0041       dzErr(nullptr),
0042       phi_pt1(nullptr),
0043       eta_pt1(nullptr),
0044       phi_pt10(nullptr),
0045       eta_pt10(nullptr),
0046       dxyVsPhi_pt1(nullptr),
0047       dzVsPhi_pt1(nullptr),
0048       dxyVsEta_pt1(nullptr),
0049       dzVsEta_pt1(nullptr),
0050       dxyVsEtaVsPhi_pt1(nullptr),
0051       dzVsEtaVsPhi_pt1(nullptr),
0052       dxyVsPhi_pt10(nullptr),
0053       dzVsPhi_pt10(nullptr),
0054       dxyVsEta_pt10(nullptr),
0055       dzVsEta_pt10(nullptr),
0056       dxyVsEtaVsPhi_pt10(nullptr),
0057       dzVsEtaVsPhi_pt10(nullptr) {
0058   vertexInputTag_ = pSet.getParameter<InputTag>("vertexLabel");
0059   beamSpotInputTag_ = pSet.getParameter<InputTag>("beamSpotLabel");
0060   vertexToken_ = consumes<reco::VertexCollection>(vertexInputTag_);
0061   scoreToken_ = consumes<VertexScore>(vertexInputTag_);
0062   beamspotToken_ = consumes<reco::BeamSpot>(beamSpotInputTag_);
0063 }
0064 
0065 // -- BeginRun
0066 //---------------------------------------------------------------------------------//
0067 void PrimaryVertexMonitor::bookHistograms(DQMStore::IBooker& iBooker, edm::Run const&, edm::EventSetup const&) {
0068   std::string dqmLabel = "";
0069 
0070   //
0071   // Book all histograms.
0072   //
0073 
0074   //  get the store
0075   dqmLabel = TopFolderName_ + "/" + vertexInputTag_.label();
0076   iBooker.setCurrentFolder(dqmLabel);
0077 
0078   //   xPos = iBooker.book1D ("xPos","x Coordinate" ,100, -0.1, 0.1);
0079 
0080   nbvtx = iBooker.book1D("vtxNbr", "Reconstructed Vertices in Event", 80, -0.5, 79.5);
0081   nbgvtx = iBooker.book1D("goodvtxNbr", "Reconstructed Good Vertices in Event", 80, -0.5, 79.5);
0082 
0083   // to be configured each year...
0084   auto vposx = conf_.getParameter<double>("Xpos");
0085   auto vposy = conf_.getParameter<double>("Ypos");
0086 
0087   nbtksinvtx[0] = iBooker.book1D("otherVtxTrksNbr", "Reconstructed Tracks in Vertex (other Vtx)", 40, -0.5, 99.5);
0088   ntracksVsZ[0] = iBooker.bookProfile(
0089       "otherVtxTrksVsZ", "Reconstructed Tracks in Vertex (other Vtx) vs Z", 80, -20., 20., 50, 0, 100, "");
0090   ntracksVsZ[0]->setAxisTitle("z-bs", 1);
0091   ntracksVsZ[0]->setAxisTitle("#tracks", 2);
0092 
0093   score[0] = iBooker.book1D("otherVtxScore", "sqrt(score) (other Vtx)", 100, 0., 400.);
0094   trksWeight[0] = iBooker.book1D("otherVtxTrksWeight", "Total weight of Tracks in Vertex (other Vtx)", 40, 0, 100.);
0095   vtxchi2[0] = iBooker.book1D("otherVtxChi2", "#chi^{2} (other Vtx)", 100, 0., 200.);
0096   vtxndf[0] = iBooker.book1D("otherVtxNdf", "ndof (other Vtx)", 100, 0., 200.);
0097   vtxprob[0] = iBooker.book1D("otherVtxProb", "#chi^{2} probability (other Vtx)", 100, 0., 1.);
0098   nans[0] = iBooker.book1D("otherVtxNans", "Illegal values for x,y,z,xx,xy,xz,yy,yz,zz (other Vtx)", 9, 0.5, 9.5);
0099 
0100   nbtksinvtx[1] = iBooker.book1D("tagVtxTrksNbr", "Reconstructed Tracks in Vertex (tagged Vtx)", 100, -0.5, 99.5);
0101   ntracksVsZ[1] = iBooker.bookProfile(
0102       "tagVtxTrksVsZ", "Reconstructed Tracks in Vertex (tagged Vtx) vs Z", 80, -20., 20., 50, 0, 100, "");
0103   ntracksVsZ[1]->setAxisTitle("z-bs", 1);
0104   ntracksVsZ[1]->setAxisTitle("#tracks", 2);
0105 
0106   score[1] = iBooker.book1D("tagVtxScore", "sqrt(score) (tagged Vtx)", 100, 0., 400.);
0107   trksWeight[1] = iBooker.book1D("tagVtxTrksWeight", "Total weight of Tracks in Vertex (tagged Vtx)", 100, 0, 100.);
0108   vtxchi2[1] = iBooker.book1D("tagVtxChi2", "#chi^{2} (tagged Vtx)", 100, 0., 200.);
0109   vtxndf[1] = iBooker.book1D("tagVtxNdf", "ndof (tagged Vtx)", 100, 0., 200.);
0110   vtxprob[1] = iBooker.book1D("tagVtxProb", "#chi^{2} probability (tagged Vtx)", 100, 0., 1.);
0111   nans[1] = iBooker.book1D("tagVtxNans", "Illegal values for x,y,z,xx,xy,xz,yy,yz,zz (tagged Vtx)", 9, 0.5, 9.5);
0112 
0113   xrec[0] = iBooker.book1D("otherPosX", "Position x Coordinate (other Vtx)", 100, vposx - 0.1, vposx + 0.1);
0114   yrec[0] = iBooker.book1D("otherPosY", "Position y Coordinate (other Vtx)", 100, vposy - 0.1, vposy + 0.1);
0115   zrec[0] = iBooker.book1D("otherPosZ", "Position z Coordinate (other Vtx)", 100, -20., 20.);
0116   xDiff[0] = iBooker.book1D("otherDiffX", "X distance from BeamSpot (other Vtx)", 100, -500, 500);
0117   yDiff[0] = iBooker.book1D("otherDiffY", "Y distance from BeamSpot (other Vtx)", 100, -500, 500);
0118   xerr[0] = iBooker.book1D("otherErrX", "Uncertainty x Coordinate (other Vtx)", 100, 0., 100);
0119   yerr[0] = iBooker.book1D("otherErrY", "Uncertainty y Coordinate (other Vtx)", 100, 0., 100);
0120   zerr[0] = iBooker.book1D("otherErrZ", "Uncertainty z Coordinate (other Vtx)", 100, 0., 100);
0121   xerrVsTrks[0] = iBooker.book2D(
0122       "otherErrVsWeightX", "Uncertainty x Coordinate vs. track weight (other Vtx)", 100, 0, 100., 100, 0., 100);
0123   yerrVsTrks[0] = iBooker.book2D(
0124       "otherErrVsWeightY", "Uncertainty y Coordinate vs. track weight (other Vtx)", 100, 0, 100., 100, 0., 100);
0125   zerrVsTrks[0] = iBooker.book2D(
0126       "otherErrVsWeightZ", "Uncertainty z Coordinate vs. track weight (other Vtx)", 100, 0, 100., 100, 0., 100);
0127 
0128   xrec[1] = iBooker.book1D("tagPosX", "Position x Coordinate (tagged Vtx)", 100, vposx - 0.1, vposx + 0.1);
0129   yrec[1] = iBooker.book1D("tagPosY", "Position y Coordinate (tagged Vtx)", 100, vposy - 0.1, vposy + 0.1);
0130   zrec[1] = iBooker.book1D("tagPosZ", "Position z Coordinate (tagged Vtx)", 100, -20., 20.);
0131   xDiff[1] = iBooker.book1D("tagDiffX", "X distance from BeamSpot (tagged Vtx)", 100, -500, 500);
0132   yDiff[1] = iBooker.book1D("tagDiffY", "Y distance from BeamSpot (tagged Vtx)", 100, -500, 500);
0133   xerr[1] = iBooker.book1D("tagErrX", "Uncertainty x Coordinate (tagged Vtx)", 100, 0., 100);
0134   yerr[1] = iBooker.book1D("tagErrY", "Uncertainty y Coordinate (tagged Vtx)", 100, 0., 100);
0135   zerr[1] = iBooker.book1D("tagErrZ", "Uncertainty z Coordinate (tagged Vtx)", 100, 0., 100);
0136   xerrVsTrks[1] = iBooker.book2D(
0137       "tagErrVsWeightX", "Uncertainty x Coordinate vs. track weight (tagged Vtx)", 100, 0, 100., 100, 0., 100);
0138   yerrVsTrks[1] = iBooker.book2D(
0139       "tagErrVsWeightY", "Uncertainty y Coordinate vs. track weight (tagged Vtx)", 100, 0, 100., 100, 0., 100);
0140   zerrVsTrks[1] = iBooker.book2D(
0141       "tagErrVsWeightZ", "Uncertainty z Coordinate vs. track weight (tagged Vtx)", 100, 0, 100., 100, 0., 100);
0142 
0143   type[0] = iBooker.book1D("otherType", "Vertex type (other Vtx)", 3, -0.5, 2.5);
0144   type[1] = iBooker.book1D("tagType", "Vertex type (tagged Vtx)", 3, -0.5, 2.5);
0145   for (int i = 0; i < 2; ++i) {
0146     type[i]->setBinLabel(1, "Valid, real");
0147     type[i]->setBinLabel(2, "Valid, fake");
0148     type[i]->setBinLabel(3, "Invalid");
0149   }
0150 
0151   //  get the store
0152   dqmLabel = TopFolderName_ + "/" + beamSpotInputTag_.label();
0153   iBooker.setCurrentFolder(dqmLabel);
0154 
0155   bsX = iBooker.book1D("bsX", "BeamSpot x0", 100, vposx - 0.1, vposx + 0.1);
0156   bsY = iBooker.book1D("bsY", "BeamSpot y0", 100, vposy - 0.1, vposy + 0.1);
0157   bsZ = iBooker.book1D("bsZ", "BeamSpot z0", 100, -2., 2.);
0158   bsSigmaZ = iBooker.book1D("bsSigmaZ", "BeamSpot sigmaZ", 100, 0., 10.);
0159   bsDxdz = iBooker.book1D("bsDxdz", "BeamSpot dxdz", 100, -0.0003, 0.0003);
0160   bsDydz = iBooker.book1D("bsDydz", "BeamSpot dydz", 100, -0.0003, 0.0003);
0161   bsBeamWidthX = iBooker.book1D("bsBeamWidthX", "BeamSpot BeamWidthX", 100, 0., 100.);
0162   bsBeamWidthY = iBooker.book1D("bsBeamWidthY", "BeamSpot BeamWidthY", 100, 0., 100.);
0163   bsType = iBooker.book1D("bsType", "BeamSpot type", 4, -1.5, 2.5);
0164   bsType->setBinLabel(1, "Unknown");
0165   bsType->setBinLabel(2, "Fake");
0166   bsType->setBinLabel(3, "LHC");
0167   bsType->setBinLabel(4, "Tracker");
0168 
0169   //  get the store
0170   dqmLabel = TopFolderName_ + "/" + AlignmentLabel_;
0171   iBooker.setCurrentFolder(dqmLabel);
0172 
0173   int TKNoBin = conf_.getParameter<int>("TkSizeBin");
0174   double TKNoMin = conf_.getParameter<double>("TkSizeMin");
0175   double TKNoMax = conf_.getParameter<double>("TkSizeMax");
0176 
0177   int DxyBin = conf_.getParameter<int>("DxyBin");
0178   double DxyMin = conf_.getParameter<double>("DxyMin");
0179   double DxyMax = conf_.getParameter<double>("DxyMax");
0180 
0181   int DzBin = conf_.getParameter<int>("DzBin");
0182   double DzMin = conf_.getParameter<double>("DzMin");
0183   double DzMax = conf_.getParameter<double>("DzMax");
0184 
0185   int PhiBin = conf_.getParameter<int>("PhiBin");
0186   int PhiBin2D = conf_.getParameter<int>("PhiBin2D");
0187   double PhiMin = conf_.getParameter<double>("PhiMin");
0188   double PhiMax = conf_.getParameter<double>("PhiMax");
0189 
0190   int EtaBin = conf_.getParameter<int>("EtaBin");
0191   int EtaBin2D = conf_.getParameter<int>("EtaBin2D");
0192   double EtaMin = conf_.getParameter<double>("EtaMin");
0193   double EtaMax = conf_.getParameter<double>("EtaMax");
0194 
0195   ntracks = iBooker.book1D("ntracks", "number of PV tracks (p_{T} > 1 GeV)", TKNoBin, TKNoMin, TKNoMax);
0196   ntracks->setAxisTitle("Number of PV Tracks (p_{T} > 1 GeV) per Event", 1);
0197   ntracks->setAxisTitle("Number of Event", 2);
0198 
0199   weight = iBooker.book1D("weight", "weight of PV tracks (p_{T} > 1 GeV)", 100, 0., 1.);
0200   weight->setAxisTitle("weight of PV Tracks (p_{T} > 1 GeV) per Event", 1);
0201   weight->setAxisTitle("Number of Event", 2);
0202 
0203   sumpt = iBooker.book1D("sumpt", "#Sum p_{T} of PV tracks (p_{T} > 1 GeV)", 100, -0.5, 249.5);
0204   chi2ndf = iBooker.book1D("chi2ndf", "PV tracks (p_{T} > 1 GeV) #chi^{2}/ndof", 100, 0., 20.);
0205   chi2prob = iBooker.book1D("chi2prob", "PV tracks (p_{T} > 1 GeV) #chi^{2} probability", 100, 0., 1.);
0206 
0207   dxy = iBooker.book1D("dxy", "PV tracks (p_{T} > 1 GeV) d_{xy} (#mum)", DxyBin, DxyMin, DxyMax);
0208   dxy2 = iBooker.book1D("dxyzoom", "PV tracks (p_{T} > 1 GeV) d_{xy} (#mum)", DxyBin, DxyMin / 5., DxyMax / 5.);
0209   dxyErr = iBooker.book1D("dxyErr", "PV tracks (p_{T} > 1 GeV) d_{xy} error (#mum)", 100, 0., 2000.);
0210   dz = iBooker.book1D("dz", "PV tracks (p_{T} > 1 GeV) d_{z} (#mum)", DzBin, DzMin, DzMax);
0211   dzErr = iBooker.book1D("dzErr", "PV tracks (p_{T} > 1 GeV) d_{z} error(#mum)", 100, 0., 10000.);
0212 
0213   phi_pt1 = iBooker.book1D("phi_pt1", "PV tracks (p_{T} > 1 GeV) #phi; PV tracks #phi;#tracks", PhiBin, PhiMin, PhiMax);
0214   eta_pt1 = iBooker.book1D("eta_pt1", "PV tracks (p_{T} > 1 GeV) #eta; PV tracks #eta;#tracks", EtaBin, EtaMin, EtaMax);
0215   phi_pt10 =
0216       iBooker.book1D("phi_pt10", "PV tracks (p_{T} > 10 GeV) #phi; PV tracks #phi;#tracks", PhiBin, PhiMin, PhiMax);
0217   eta_pt10 =
0218       iBooker.book1D("eta_pt10", "PV tracks (p_{T} > 10 GeV) #phi; PV tracks #eta;#tracks", EtaBin, EtaMin, EtaMax);
0219 
0220   dxyVsPhi_pt1 = iBooker.bookProfile("dxyVsPhi_pt1",
0221                                      "PV tracks (p_{T} > 1 GeV) d_{xy} (#mum) VS track #phi",
0222                                      PhiBin,
0223                                      PhiMin,
0224                                      PhiMax,
0225                                      DxyBin,
0226                                      DxyMin,
0227                                      DxyMax,
0228                                      "");
0229   dxyVsPhi_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) #phi", 1);
0230   dxyVsPhi_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) d_{xy} (#mum)", 2);
0231 
0232   dzVsPhi_pt1 = iBooker.bookProfile("dzVsPhi_pt1",
0233                                     "PV tracks (p_{T} > 1 GeV) d_{z} (#mum) VS track #phi",
0234                                     PhiBin,
0235                                     PhiMin,
0236                                     PhiMax,
0237                                     DzBin,
0238                                     DzMin,
0239                                     DzMax,
0240                                     "");
0241   dzVsPhi_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) #phi", 1);
0242   dzVsPhi_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) d_{z} (#mum)", 2);
0243 
0244   dxyVsEta_pt1 = iBooker.bookProfile("dxyVsEta_pt1",
0245                                      "PV tracks (p_{T} > 1 GeV) d_{xy} (#mum) VS track #eta",
0246                                      EtaBin,
0247                                      EtaMin,
0248                                      EtaMax,
0249                                      DxyBin,
0250                                      DxyMin,
0251                                      DxyMax,
0252                                      "");
0253   dxyVsEta_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) #eta", 1);
0254   dxyVsEta_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) d_{xy} (#mum)", 2);
0255 
0256   dzVsEta_pt1 = iBooker.bookProfile("dzVsEta_pt1",
0257                                     "PV tracks (p_{T} > 1 GeV) d_{z} (#mum) VS track #eta",
0258                                     EtaBin,
0259                                     EtaMin,
0260                                     EtaMax,
0261                                     DzBin,
0262                                     DzMin,
0263                                     DzMax,
0264                                     "");
0265   dzVsEta_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) #eta", 1);
0266   dzVsEta_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) d_{z} (#mum)", 2);
0267 
0268   dxyVsEtaVsPhi_pt1 = iBooker.bookProfile2D("dxyVsEtaVsPhi_pt1",
0269                                             "PV tracks (p_{T} > 1 GeV) d_{xy} (#mum) VS track #eta VS track #phi",
0270                                             EtaBin2D,
0271                                             EtaMin,
0272                                             EtaMax,
0273                                             PhiBin2D,
0274                                             PhiMin,
0275                                             PhiMax,
0276                                             DxyBin,
0277                                             DxyMin,
0278                                             DxyMax,
0279                                             "");
0280   dxyVsEtaVsPhi_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) #eta", 1);
0281   dxyVsEtaVsPhi_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) #phi", 2);
0282   dxyVsEtaVsPhi_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) d_{xy} (#mum)", 3);
0283 
0284   dzVsEtaVsPhi_pt1 = iBooker.bookProfile2D("dzVsEtaVsPhi_pt1",
0285                                            "PV tracks (p_{T} > 1 GeV) d_{z} (#mum) VS track #eta VS track #phi",
0286                                            EtaBin2D,
0287                                            EtaMin,
0288                                            EtaMax,
0289                                            PhiBin2D,
0290                                            PhiMin,
0291                                            PhiMax,
0292                                            DzBin,
0293                                            DzMin,
0294                                            DzMax,
0295                                            "");
0296   dzVsEtaVsPhi_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) #eta", 1);
0297   dzVsEtaVsPhi_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) #phi", 2);
0298   dzVsEtaVsPhi_pt1->setAxisTitle("PV track (p_{T} > 1 GeV) d_{z} (#mum)", 3);
0299 
0300   dxyVsPhi_pt10 = iBooker.bookProfile("dxyVsPhi_pt10",
0301                                       "PV tracks (p_{T} > 10 GeV) d_{xy} (#mum) VS track #phi",
0302                                       PhiBin,
0303                                       PhiMin,
0304                                       PhiMax,
0305                                       DxyBin,
0306                                       DxyMin,
0307                                       DxyMax,
0308                                       "");
0309   dxyVsPhi_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) #phi", 1);
0310   dxyVsPhi_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) d_{xy} (#mum)", 2);
0311 
0312   dzVsPhi_pt10 = iBooker.bookProfile("dzVsPhi_pt10",
0313                                      "PV tracks (p_{T} > 10 GeV) d_{z} (#mum) VS track #phi",
0314                                      PhiBin,
0315                                      PhiMin,
0316                                      PhiMax,
0317                                      DzBin,
0318                                      DzMin,
0319                                      DzMax,
0320                                      "");
0321   dzVsPhi_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) #phi", 1);
0322   dzVsPhi_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) d_{z} (#mum)", 2);
0323 
0324   dxyVsEta_pt10 = iBooker.bookProfile("dxyVsEta_pt10",
0325                                       "PV tracks (p_{T} > 10 GeV) d_{xy} (#mum) VS track #eta",
0326                                       EtaBin,
0327                                       EtaMin,
0328                                       EtaMax,
0329                                       DxyBin,
0330                                       DxyMin,
0331                                       DxyMax,
0332                                       "");
0333   dxyVsEta_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) #eta", 1);
0334   dxyVsEta_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) d_{xy} (#mum)", 2);
0335 
0336   dzVsEta_pt10 = iBooker.bookProfile("dzVsEta_pt10",
0337                                      "PV tracks (p_{T} > 10 GeV) d_{z} (#mum) VS track #eta",
0338                                      EtaBin,
0339                                      EtaMin,
0340                                      EtaMax,
0341                                      DzBin,
0342                                      DzMin,
0343                                      DzMax,
0344                                      "");
0345   dzVsEta_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) #eta", 1);
0346   dzVsEta_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) d_{z} (#mum)", 2);
0347 
0348   dxyVsEtaVsPhi_pt10 = iBooker.bookProfile2D("dxyVsEtaVsPhi_pt10",
0349                                              "PV tracks (p_{T} > 10 GeV) d_{xy} (#mum) VS track #eta VS track #phi",
0350                                              EtaBin2D,
0351                                              EtaMin,
0352                                              EtaMax,
0353                                              PhiBin2D,
0354                                              PhiMin,
0355                                              PhiMax,
0356                                              DxyBin,
0357                                              DxyMin,
0358                                              DxyMax,
0359                                              "");
0360   dxyVsEtaVsPhi_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) #eta", 1);
0361   dxyVsEtaVsPhi_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) #phi", 2);
0362   dxyVsEtaVsPhi_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) d_{xy} (#mum)", 3);
0363 
0364   dzVsEtaVsPhi_pt10 = iBooker.bookProfile2D("dzVsEtaVsPhi_pt10",
0365                                             "PV tracks (p_{T} > 10 GeV) d_{z} (#mum) VS track #eta VS track #phi",
0366                                             EtaBin2D,
0367                                             EtaMin,
0368                                             EtaMax,
0369                                             PhiBin2D,
0370                                             PhiMin,
0371                                             PhiMax,
0372                                             DzBin,
0373                                             DzMin,
0374                                             DzMax,
0375                                             "");
0376   dzVsEtaVsPhi_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) #eta", 1);
0377   dzVsEtaVsPhi_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) #phi", 2);
0378   dzVsEtaVsPhi_pt10->setAxisTitle("PV track (p_{T} > 10 GeV) d_{z} (#mum)", 3);
0379 }
0380 
0381 void PrimaryVertexMonitor::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0382   Handle<reco::VertexCollection> recVtxs;
0383   iEvent.getByToken(vertexToken_, recVtxs);
0384 
0385   Handle<VertexScore> scores;
0386   iEvent.getByToken(scoreToken_, scores);
0387 
0388   edm::Handle<reco::BeamSpot> beamSpotHandle;
0389   iEvent.getByToken(beamspotToken_, beamSpotHandle);
0390 
0391   //
0392   // check for absent products and simply "return" in that case
0393   //
0394   if (recVtxs.isValid() == false || beamSpotHandle.isValid() == false) {
0395     edm::LogWarning("PrimaryVertexMonitor")
0396         << " Some products not available in the event: VertexCollection " << vertexInputTag_ << " " << recVtxs.isValid()
0397         << " BeamSpot " << beamSpotInputTag_ << " " << beamSpotHandle.isValid() << ". Skipping plots for this event";
0398     return;
0399   }
0400 
0401   // check upfront that refs to track are (likely) to be valid
0402   {
0403     bool ok = true;
0404     for (const auto& v : *recVtxs) {
0405       if (v.tracksSize() > 0) {
0406         const auto& ref = v.trackRefAt(0);
0407         if (ref.isNull() || !ref.isAvailable()) {
0408           if (!errorPrinted_)
0409             edm::LogWarning("PrimaryVertexMonitor")
0410                 << "Skipping vertex collection: " << vertexInputTag_
0411                 << " since likely the track collection the vertex has refs pointing to is missing (at least the first "
0412                    "TrackBaseRef is null or not available)";
0413           else
0414             errorPrinted_ = true;
0415           ok = false;
0416         }
0417       }
0418     }
0419     if (!ok)
0420       return;
0421   }
0422 
0423   BeamSpot beamSpot = *beamSpotHandle;
0424 
0425   nbvtx->Fill(recVtxs->size() * 1.);
0426   int ng = 0;
0427   for (auto const& vx : (*recVtxs))
0428     if (vx.isValid() && !vx.isFake() && vx.ndof() >= ndof_)
0429       ++ng;
0430   nbgvtx->Fill(ng * 1.);
0431 
0432   if (scores.isValid() && !(*scores).empty()) {
0433     auto pvScore = (*scores).get(0);
0434     score[1]->Fill(std::sqrt(pvScore));
0435     for (unsigned int i = 1; i < (*scores).size(); ++i)
0436       score[0]->Fill(std::sqrt((*scores).get(i)));
0437   }
0438 
0439   // fill PV tracks MEs (as now, for alignment)
0440   if (!recVtxs->empty()) {
0441     vertexPlots(recVtxs->front(), beamSpot, 1);
0442     pvTracksPlots(recVtxs->front());
0443 
0444     for (reco::VertexCollection::const_iterator v = recVtxs->begin() + 1; v != recVtxs->end(); ++v)
0445       vertexPlots(*v, beamSpot, 0);
0446   }
0447 
0448   // Beamline plots:
0449   bsX->Fill(beamSpot.x0());
0450   bsY->Fill(beamSpot.y0());
0451   bsZ->Fill(beamSpot.z0());
0452   bsSigmaZ->Fill(beamSpot.sigmaZ());
0453   bsDxdz->Fill(beamSpot.dxdz());
0454   bsDydz->Fill(beamSpot.dydz());
0455   bsBeamWidthX->Fill(beamSpot.BeamWidthX() * cmToUm);
0456   bsBeamWidthY->Fill(beamSpot.BeamWidthY() * cmToUm);
0457   bsType->Fill(beamSpot.type());
0458 }
0459 
0460 void PrimaryVertexMonitor::pvTracksPlots(const Vertex& v) {
0461   if (!v.isValid())
0462     return;
0463   if (v.isFake())
0464     return;
0465 
0466   if (v.tracksSize() == 0) {
0467     ntracks->Fill(0);
0468     return;
0469   }
0470 
0471   const math::XYZPoint myVertex(v.position().x(), v.position().y(), v.position().z());
0472 
0473   size_t nTracks = 0;
0474   float sumPT = 0.;
0475 
0476   for (reco::Vertex::trackRef_iterator t = v.tracks_begin(); t != v.tracks_end(); t++) {
0477     bool isHighPurity = (**t).quality(reco::TrackBase::highPurity);
0478     if (!isHighPurity && useHPfoAlignmentPlots_)
0479       continue;
0480 
0481     float pt = (**t).pt();
0482     if (pt < 1.)
0483       continue;
0484 
0485     nTracks++;
0486 
0487     float eta = (**t).eta();
0488     float phi = (**t).phi();
0489 
0490     float w = v.trackWeight(*t);
0491     float chi2NDF = (**t).normalizedChi2();
0492     float chi2Prob = TMath::Prob((**t).chi2(), (int)(**t).ndof());
0493     float Dxy = (**t).dxy(myVertex) * cmToUm;
0494     float Dz = (**t).dz(myVertex) * cmToUm;
0495     float DxyErr = (**t).dxyError() * cmToUm;
0496     float DzErr = (**t).dzError() * cmToUm;
0497 
0498     sumPT += pt * pt;
0499 
0500     // fill MEs
0501     weight->Fill(w);
0502     chi2ndf->Fill(chi2NDF);
0503     chi2prob->Fill(chi2Prob);
0504     dxy->Fill(Dxy);
0505     dxy2->Fill(Dxy);
0506     dz->Fill(Dz);
0507     dxyErr->Fill(DxyErr);
0508     dzErr->Fill(DzErr);
0509     phi_pt1->Fill(phi);
0510     eta_pt1->Fill(eta);
0511 
0512     dxyVsPhi_pt1->Fill(phi, Dxy);
0513     dzVsPhi_pt1->Fill(phi, Dz);
0514     dxyVsEta_pt1->Fill(eta, Dxy);
0515     dzVsEta_pt1->Fill(eta, Dz);
0516     dxyVsEtaVsPhi_pt1->Fill(eta, phi, Dxy);
0517     dzVsEtaVsPhi_pt1->Fill(eta, phi, Dz);
0518 
0519     if (pt < 10.)
0520       continue;
0521 
0522     phi_pt10->Fill(phi);
0523     eta_pt10->Fill(eta);
0524     dxyVsPhi_pt10->Fill(phi, Dxy);
0525     dzVsPhi_pt10->Fill(phi, Dz);
0526     dxyVsEta_pt10->Fill(eta, Dxy);
0527     dzVsEta_pt10->Fill(eta, Dz);
0528     dxyVsEtaVsPhi_pt10->Fill(eta, phi, Dxy);
0529     dzVsEtaVsPhi_pt10->Fill(eta, phi, Dz);
0530   }
0531   ntracks->Fill(float(nTracks));
0532   sumpt->Fill(sumPT);
0533 }
0534 
0535 void PrimaryVertexMonitor::vertexPlots(const Vertex& v, const BeamSpot& beamSpot, int i) {
0536   if (i < 0 || i > 1)
0537     return;
0538   if (!v.isValid())
0539     type[i]->Fill(2.);
0540   else if (v.isFake())
0541     type[i]->Fill(1.);
0542   else
0543     type[i]->Fill(0.);
0544 
0545   if (v.isValid() && !v.isFake()) {
0546     float weight = 0;
0547     for (reco::Vertex::trackRef_iterator t = v.tracks_begin(); t != v.tracks_end(); t++)
0548       weight += v.trackWeight(*t);
0549     trksWeight[i]->Fill(weight);
0550     nbtksinvtx[i]->Fill(v.tracksSize());
0551     ntracksVsZ[i]->Fill(v.position().z() - beamSpot.z0(), v.tracksSize());
0552 
0553     vtxchi2[i]->Fill(v.chi2());
0554     vtxndf[i]->Fill(v.ndof());
0555     vtxprob[i]->Fill(ChiSquaredProbability(v.chi2(), v.ndof()));
0556 
0557     xrec[i]->Fill(v.position().x());
0558     yrec[i]->Fill(v.position().y());
0559     zrec[i]->Fill(v.position().z());
0560 
0561     float xb = beamSpot.x0() + beamSpot.dxdz() * (v.position().z() - beamSpot.z0());
0562     float yb = beamSpot.y0() + beamSpot.dydz() * (v.position().z() - beamSpot.z0());
0563     xDiff[i]->Fill((v.position().x() - xb) * cmToUm);
0564     yDiff[i]->Fill((v.position().y() - yb) * cmToUm);
0565 
0566     xerr[i]->Fill(v.xError() * cmToUm);
0567     yerr[i]->Fill(v.yError() * cmToUm);
0568     zerr[i]->Fill(v.zError() * cmToUm);
0569     xerrVsTrks[i]->Fill(weight, v.xError() * cmToUm);
0570     yerrVsTrks[i]->Fill(weight, v.yError() * cmToUm);
0571     zerrVsTrks[i]->Fill(weight, v.zError() * cmToUm);
0572 
0573     nans[i]->Fill(1., edm::isNotFinite(v.position().x()) * 1.);
0574     nans[i]->Fill(2., edm::isNotFinite(v.position().y()) * 1.);
0575     nans[i]->Fill(3., edm::isNotFinite(v.position().z()) * 1.);
0576 
0577     int index = 3;
0578     for (int k = 0; k != 3; k++) {
0579       for (int j = k; j != 3; j++) {
0580         index++;
0581         nans[i]->Fill(index * 1., edm::isNotFinite(v.covariance(k, j)) * 1.);
0582         // in addition, diagonal element must be positive
0583         if (j == k && v.covariance(k, j) < 0) {
0584           nans[i]->Fill(index * 1., 1.);
0585         }
0586       }
0587     }
0588   }
0589 }
0590 
0591 //define this as a plug-in
0592 DEFINE_FWK_MODULE(PrimaryVertexMonitor);