Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-04-26 01:55:02

0001 // -*- C++ -*-
0002 //
0003 // Package:    Alignment/GeneralPurposeVertexAnalyzer
0004 // Class:      GeneralPurposeVertexAnalyzer
0005 //
0006 /**\class GeneralPurposeVertexAnalyzer GeneralPurposeVertexAnalyzer.cc Alignment/GeneralPurposeVertexAnalyzer/plugins/GeneralPurposeVertexAnalyzer.cc
0007    Description: monitor vertex properties for alignment purposes, largely copied from DQMOffline/RecoB/plugins/PrimaryVertexMonitor.cc
0008 
0009 */
0010 //
0011 // Original Author:  Marco Musich
0012 //         Created:  Thu, 13 Apr 2023 14:16:43 GMT
0013 //
0014 //
0015 
0016 // ROOT includes files
0017 #include "TMath.h"
0018 #include "TFile.h"
0019 #include "TH1I.h"
0020 #include "TH1D.h"
0021 #include "TH2D.h"
0022 #include "TProfile.h"
0023 #include "TProfile2D.h"
0024 
0025 // system include files
0026 #include <memory>
0027 #include <fmt/format.h>
0028 #include <fmt/printf.h>
0029 
0030 // user include files
0031 #include "CommonTools/Statistics/interface/ChiSquaredProbability.h"
0032 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0033 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0034 #include "DataFormats/Common/interface/Association.h"
0035 #include "DataFormats/Common/interface/ValueMap.h"
0036 #include "DataFormats/TrackReco/interface/Track.h"
0037 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0038 #include "DataFormats/VertexReco/interface/Vertex.h"
0039 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0040 #include "FWCore/ServiceRegistry/interface/Service.h"
0041 #include "FWCore/Framework/interface/Event.h"
0042 #include "FWCore/Framework/interface/Frameworkfwd.h"
0043 #include "FWCore/Framework/interface/MakerMacros.h"
0044 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0045 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0046 #include "FWCore/Utilities/interface/InputTag.h"
0047 #include "FWCore/Utilities/interface/isFinite.h"
0048 
0049 //
0050 // class declaration
0051 //
0052 
0053 using reco::TrackCollection;
0054 
0055 class GeneralPurposeVertexAnalyzer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0056 public:
0057   explicit GeneralPurposeVertexAnalyzer(const edm::ParameterSet &);
0058   ~GeneralPurposeVertexAnalyzer() override = default;
0059 
0060   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0061 
0062 private:
0063   void pvTracksPlots(const reco::Vertex &v);
0064   void vertexPlots(const reco::Vertex &v, const reco::BeamSpot &beamSpot, int i);
0065   template <typename T, typename... Args>
0066   T *book(const Args &...args) const;
0067   void beginJob() override;
0068   void analyze(const edm::Event &, const edm::EventSetup &) override;
0069 
0070   // ----------member data ---------------------------
0071   edm::Service<TFileService> fs_;
0072 
0073   const int ndof_;
0074   bool errorPrinted_;
0075 
0076   const edm::InputTag vertexInputTag_, beamSpotInputTag_;
0077   const edm::EDGetTokenT<reco::VertexCollection> vertexToken_;
0078   using VertexScore = edm::ValueMap<float>;
0079   const edm::EDGetTokenT<VertexScore> scoreToken_;
0080   const edm::EDGetTokenT<reco::BeamSpot> beamspotToken_;
0081 
0082   static constexpr int cmToUm = 10000;
0083 
0084   const double vposx_;
0085   const double vposy_;
0086   const int tkNoBin_;
0087   const double tkNoMin_;
0088   const double tkNoMax_;
0089 
0090   const int dxyBin_;
0091   const double dxyMin_;
0092   const double dxyMax_;
0093 
0094   const int dzBin_;
0095   const double dzMin_;
0096   const double dzMax_;
0097 
0098   const int phiBin_;
0099   const int phiBin2D_;
0100   const double phiMin_;
0101   const double phiMax_;
0102 
0103   const int etaBin_;
0104   const int etaBin2D_;
0105   const double etaMin_;
0106   const double etaMax_;
0107 
0108   // the histos
0109   TH1I *nbvtx, *nbgvtx;
0110   TH1D *nbtksinvtx[2], *trksWeight[2], *score[2];
0111   TH1D *tt[2];
0112   TH1D *xrec[2], *yrec[2], *zrec[2], *xDiff[2], *yDiff[2], *xerr[2], *yerr[2], *zerr[2];
0113   TH2D *xerrVsTrks[2], *yerrVsTrks[2], *zerrVsTrks[2];
0114   TH1D *ntracksVsZ[2];
0115   TH1D *vtxchi2[2], *vtxndf[2], *vtxprob[2], *nans[2];
0116   TH1D *type[2];
0117   TH1D *bsX, *bsY, *bsZ, *bsSigmaZ, *bsDxdz, *bsDydz, *bsBeamWidthX, *bsBeamWidthY, *bsType;
0118 
0119   TH1D *sumpt, *ntracks, *weight, *chi2ndf, *chi2prob;
0120   TH1D *dxy, *dxy2, *dz, *dxyErr, *dzErr;
0121   TH1D *phi_pt1, *eta_pt1;
0122   TH1D *phi_pt10, *eta_pt10;
0123   TProfile *dxyVsPhi_pt1, *dzVsPhi_pt1;
0124   TProfile *dxyVsEta_pt1, *dzVsEta_pt1;
0125   TProfile2D *dxyVsEtaVsPhi_pt1, *dzVsEtaVsPhi_pt1;
0126   TProfile *dxyVsPhi_pt10, *dzVsPhi_pt10;
0127   TProfile *dxyVsEta_pt10, *dzVsEta_pt10;
0128   TProfile2D *dxyVsEtaVsPhi_pt10, *dzVsEtaVsPhi_pt10;
0129 };
0130 
0131 //
0132 // constructors and destructor
0133 //
0134 GeneralPurposeVertexAnalyzer::GeneralPurposeVertexAnalyzer(const edm::ParameterSet &iConfig)
0135     : ndof_(iConfig.getParameter<int>("ndof")),
0136       errorPrinted_(false),
0137       vertexInputTag_(iConfig.getParameter<edm::InputTag>("vertexLabel")),
0138       beamSpotInputTag_(iConfig.getParameter<edm::InputTag>("beamSpotLabel")),
0139       vertexToken_(consumes<reco::VertexCollection>(vertexInputTag_)),
0140       scoreToken_(consumes<VertexScore>(vertexInputTag_)),
0141       beamspotToken_(consumes<reco::BeamSpot>(beamSpotInputTag_)),
0142       // to be configured for each year...
0143       vposx_(iConfig.getParameter<double>("Xpos")),
0144       vposy_(iConfig.getParameter<double>("Ypos")),
0145       tkNoBin_(iConfig.getParameter<int>("TkSizeBin")),
0146       tkNoMin_(iConfig.getParameter<double>("TkSizeMin")),
0147       tkNoMax_(iConfig.getParameter<double>("TkSizeMax")),
0148       dxyBin_(iConfig.getParameter<int>("DxyBin")),
0149       dxyMin_(iConfig.getParameter<double>("DxyMin")),
0150       dxyMax_(iConfig.getParameter<double>("DxyMax")),
0151       dzBin_(iConfig.getParameter<int>("DzBin")),
0152       dzMin_(iConfig.getParameter<double>("DzMin")),
0153       dzMax_(iConfig.getParameter<double>("DzMax")),
0154       phiBin_(iConfig.getParameter<int>("PhiBin")),
0155       phiBin2D_(iConfig.getParameter<int>("PhiBin2D")),
0156       phiMin_(iConfig.getParameter<double>("PhiMin")),
0157       phiMax_(iConfig.getParameter<double>("PhiMax")),
0158       etaBin_(iConfig.getParameter<int>("EtaBin")),
0159       etaBin2D_(iConfig.getParameter<int>("EtaBin2D")),
0160       etaMin_(iConfig.getParameter<double>("EtaMin")),
0161       etaMax_(iConfig.getParameter<double>("EtaMax")),
0162       // histograms
0163       nbvtx(nullptr),
0164       bsX(nullptr),
0165       bsY(nullptr),
0166       bsZ(nullptr),
0167       bsSigmaZ(nullptr),
0168       bsDxdz(nullptr),
0169       bsDydz(nullptr),
0170       bsBeamWidthX(nullptr),
0171       bsBeamWidthY(nullptr),
0172       bsType(nullptr),
0173       sumpt(nullptr),
0174       ntracks(nullptr),
0175       weight(nullptr),
0176       chi2ndf(nullptr),
0177       chi2prob(nullptr),
0178       dxy(nullptr),
0179       dxy2(nullptr),
0180       dz(nullptr),
0181       dxyErr(nullptr),
0182       dzErr(nullptr),
0183       phi_pt1(nullptr),
0184       eta_pt1(nullptr),
0185       phi_pt10(nullptr),
0186       eta_pt10(nullptr),
0187       dxyVsPhi_pt1(nullptr),
0188       dzVsPhi_pt1(nullptr),
0189       dxyVsEta_pt1(nullptr),
0190       dzVsEta_pt1(nullptr),
0191       dxyVsEtaVsPhi_pt1(nullptr),
0192       dzVsEtaVsPhi_pt1(nullptr),
0193       dxyVsPhi_pt10(nullptr),
0194       dzVsPhi_pt10(nullptr),
0195       dxyVsEta_pt10(nullptr),
0196       dzVsEta_pt10(nullptr),
0197       dxyVsEtaVsPhi_pt10(nullptr),
0198       dzVsEtaVsPhi_pt10(nullptr) {
0199   usesResource(TFileService::kSharedResource);
0200 }
0201 
0202 //
0203 // member functions
0204 //
0205 
0206 // ------------ method called for each event  ------------
0207 void GeneralPurposeVertexAnalyzer::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0208   using namespace edm;
0209 
0210   const auto &recVtxs = iEvent.getHandle(vertexToken_);
0211   const auto &scores = iEvent.getHandle(scoreToken_);
0212   const auto &beamSpotHandle = iEvent.getHandle(beamspotToken_);
0213   reco::BeamSpot beamSpot = *beamSpotHandle;
0214 
0215   // check for absent products and simply "return" in that case
0216   if (recVtxs.isValid() == false || beamSpotHandle.isValid() == false) {
0217     edm::LogWarning("GeneralPurposeVertexAnalyzer")
0218         << " Some products not available in the event: VertexCollection " << vertexInputTag_ << " " << recVtxs.isValid()
0219         << " BeamSpot " << beamSpotInputTag_ << " " << beamSpotHandle.isValid() << ". Skipping plots for this event";
0220     return;
0221   }
0222 
0223   // check upfront that refs to track are (likely) to be valid
0224   bool ok{true};
0225   for (const auto &v : *recVtxs) {
0226     if (v.tracksSize() > 0) {
0227       const auto &ref = v.trackRefAt(0);
0228       if (ref.isNull() || !ref.isAvailable()) {
0229         if (!errorPrinted_) {
0230           edm::LogWarning("GeneralPurposeVertexAnalyzer")
0231               << "Skipping vertex collection: " << vertexInputTag_
0232               << " since likely the track collection the vertex has refs pointing to is missing (at least the first "
0233                  "TrackBaseRef is null or not available)";
0234         } else {
0235           errorPrinted_ = true;
0236         }
0237         ok = false;
0238       }
0239     }
0240   }
0241 
0242   if (!ok) {
0243     return;
0244   }
0245 
0246   nbvtx->Fill(recVtxs->size());
0247   int ng = 0;
0248   for (auto const &vx : (*recVtxs)) {
0249     if (vx.isValid() && !vx.isFake() && vx.ndof() >= ndof_) {
0250       ++ng;
0251     }
0252   }
0253   nbgvtx->Fill(ng);
0254 
0255   if (scores.isValid() && !(*scores).empty()) {
0256     auto pvScore = (*scores).get(0);
0257     score[1]->Fill(std::sqrt(pvScore));
0258     for (unsigned int i = 1; i < (*scores).size(); ++i) {
0259       score[0]->Fill(std::sqrt((*scores).get(i)));
0260     }
0261   }
0262 
0263   // fill PV tracks MEs (as now, for alignment)
0264   if (!recVtxs->empty()) {
0265     vertexPlots(recVtxs->front(), beamSpot, 1);
0266     pvTracksPlots(recVtxs->front());
0267 
0268     for (reco::VertexCollection::const_iterator v = recVtxs->begin() + 1; v != recVtxs->end(); ++v) {
0269       vertexPlots(*v, beamSpot, 0);
0270     }
0271   }
0272 
0273   // Beamline plots:
0274   bsX->Fill(beamSpot.x0());
0275   bsY->Fill(beamSpot.y0());
0276   bsZ->Fill(beamSpot.z0());
0277   bsSigmaZ->Fill(beamSpot.sigmaZ());
0278   bsDxdz->Fill(beamSpot.dxdz());
0279   bsDydz->Fill(beamSpot.dydz());
0280   bsBeamWidthX->Fill(beamSpot.BeamWidthX() * cmToUm);
0281   bsBeamWidthY->Fill(beamSpot.BeamWidthY() * cmToUm);
0282   bsType->Fill(beamSpot.type());
0283 }
0284 
0285 void GeneralPurposeVertexAnalyzer::pvTracksPlots(const reco::Vertex &v) {
0286   if (!v.isValid())
0287     return;
0288   if (v.isFake())
0289     return;
0290 
0291   if (v.tracksSize() == 0) {
0292     ntracks->Fill(0);
0293     return;
0294   }
0295 
0296   const math::XYZPoint myVertex(v.position().x(), v.position().y(), v.position().z());
0297 
0298   float sumPT = 0.;
0299   for (const auto &t : v.tracks()) {
0300     const bool isHighPurity = t->quality(reco::TrackBase::highPurity);
0301     if (!isHighPurity) {
0302       continue;
0303     }
0304 
0305     const float pt = t->pt();
0306     if (pt < 1.f) {
0307       continue;
0308     }
0309 
0310     const float pt2 = pt * pt;
0311     const float eta = t->eta();
0312     const float phi = t->phi();
0313     const float w = v.trackWeight(t);
0314     const float chi2NDF = t->normalizedChi2();
0315     const float chi2Prob = TMath::Prob(t->chi2(), static_cast<int>(t->ndof()));
0316     const float Dxy = t->dxy(myVertex) * cmToUm;
0317     const float Dz = t->dz(myVertex) * cmToUm;
0318     const float DxyErr = t->dxyError() * cmToUm;
0319     const float DzErr = t->dzError() * cmToUm;
0320 
0321     sumPT += pt2;
0322 
0323     weight->Fill(w);
0324     chi2ndf->Fill(chi2NDF);
0325     chi2prob->Fill(chi2Prob);
0326     dxy->Fill(Dxy);
0327     dxy2->Fill(Dxy);
0328     dz->Fill(Dz);
0329     dxyErr->Fill(DxyErr);
0330     dzErr->Fill(DzErr);
0331     phi_pt1->Fill(phi);
0332     eta_pt1->Fill(eta);
0333     dxyVsPhi_pt1->Fill(phi, Dxy);
0334     dzVsPhi_pt1->Fill(phi, Dz);
0335     dxyVsEta_pt1->Fill(eta, Dxy);
0336     dzVsEta_pt1->Fill(eta, Dz);
0337     dxyVsEtaVsPhi_pt1->Fill(eta, phi, Dxy);
0338     dzVsEtaVsPhi_pt1->Fill(eta, phi, Dz);
0339 
0340     if (pt >= 10.f) {
0341       phi_pt10->Fill(phi);
0342       eta_pt10->Fill(eta);
0343       dxyVsPhi_pt10->Fill(phi, Dxy);
0344       dzVsPhi_pt10->Fill(phi, Dz);
0345       dxyVsEta_pt10->Fill(eta, Dxy);
0346       dzVsEta_pt10->Fill(eta, Dz);
0347       dxyVsEtaVsPhi_pt10->Fill(eta, phi, Dxy);
0348       dzVsEtaVsPhi_pt10->Fill(eta, phi, Dz);
0349     }
0350   }
0351 
0352   ntracks->Fill(static_cast<float>(v.tracks().size()));
0353   sumpt->Fill(sumPT);
0354 }
0355 
0356 void GeneralPurposeVertexAnalyzer::vertexPlots(const reco::Vertex &v, const reco::BeamSpot &beamSpot, int i) {
0357   if (i < 0 || i > 1)
0358     return;
0359   if (!v.isValid())
0360     type[i]->Fill(2.);
0361   else if (v.isFake())
0362     type[i]->Fill(1.);
0363   else
0364     type[i]->Fill(0.);
0365 
0366   if (v.isValid() && !v.isFake()) {
0367     float weight = 0;
0368     for (reco::Vertex::trackRef_iterator t = v.tracks_begin(); t != v.tracks_end(); t++)
0369       weight += v.trackWeight(*t);
0370     trksWeight[i]->Fill(weight);
0371     nbtksinvtx[i]->Fill(v.tracksSize());
0372     ntracksVsZ[i]->Fill(v.position().z() - beamSpot.z0(), v.tracksSize());
0373 
0374     vtxchi2[i]->Fill(v.chi2());
0375     vtxndf[i]->Fill(v.ndof());
0376     vtxprob[i]->Fill(ChiSquaredProbability(v.chi2(), v.ndof()));
0377 
0378     xrec[i]->Fill(v.position().x());
0379     yrec[i]->Fill(v.position().y());
0380     zrec[i]->Fill(v.position().z());
0381 
0382     float xb = beamSpot.x0() + beamSpot.dxdz() * (v.position().z() - beamSpot.z0());
0383     float yb = beamSpot.y0() + beamSpot.dydz() * (v.position().z() - beamSpot.z0());
0384     xDiff[i]->Fill((v.position().x() - xb) * cmToUm);
0385     yDiff[i]->Fill((v.position().y() - yb) * cmToUm);
0386 
0387     xerr[i]->Fill(v.xError() * cmToUm);
0388     yerr[i]->Fill(v.yError() * cmToUm);
0389     zerr[i]->Fill(v.zError() * cmToUm);
0390     xerrVsTrks[i]->Fill(weight, v.xError() * cmToUm);
0391     yerrVsTrks[i]->Fill(weight, v.yError() * cmToUm);
0392     zerrVsTrks[i]->Fill(weight, v.zError() * cmToUm);
0393 
0394     nans[i]->Fill(1., edm::isNotFinite(v.position().x()) * 1.);
0395     nans[i]->Fill(2., edm::isNotFinite(v.position().y()) * 1.);
0396     nans[i]->Fill(3., edm::isNotFinite(v.position().z()) * 1.);
0397 
0398     int index = 3;
0399     for (int k = 0; k != 3; k++) {
0400       for (int j = k; j != 3; j++) {
0401         index++;
0402         nans[i]->Fill(index * 1., edm::isNotFinite(v.covariance(k, j)) * 1.);
0403         // in addition, diagonal element must be positive
0404         if (j == k && v.covariance(k, j) < 0) {
0405           nans[i]->Fill(index * 1., 1.);
0406         }
0407       }
0408     }
0409   }
0410 }
0411 
0412 template <typename T, typename... Args>
0413 T *GeneralPurposeVertexAnalyzer::book(const Args &...args) const {
0414   T *t = fs_->make<T>(args...);
0415   return t;
0416 }
0417 
0418 // ------------ method called once each job just before starting event loop  ------------
0419 void GeneralPurposeVertexAnalyzer::beginJob() {
0420   nbvtx = book<TH1I>("vtxNbr", "Reconstructed Vertices in Event", 80, -0.5, 79.5);
0421   nbgvtx = book<TH1I>("goodvtxNbr", "Reconstructed Good Vertices in Event", 80, -0.5, 79.5);
0422 
0423   nbtksinvtx[0] = book<TH1D>("otherVtxTrksNbr", "Reconstructed Tracks in Vertex (other Vtx)", 40, -0.5, 99.5);
0424   ntracksVsZ[0] =
0425       book<TProfile>("otherVtxTrksVsZ", "Reconstructed Tracks in Vertex (other Vtx) vs Z", 80, -20., 20., 0., 100., "");
0426   ntracksVsZ[0]->SetXTitle("z-bs");
0427   ntracksVsZ[0]->SetYTitle("#tracks");
0428 
0429   score[0] = book<TH1D>("otherVtxScore", "sqrt(score) (other Vtx)", 100, 0., 400.);
0430   trksWeight[0] = book<TH1D>("otherVtxTrksWeight", "Total weight of Tracks in Vertex (other Vtx)", 40, 0, 100.);
0431   vtxchi2[0] = book<TH1D>("otherVtxChi2", "#chi^{2} (other Vtx)", 100, 0., 200.);
0432   vtxndf[0] = book<TH1D>("otherVtxNdf", "ndof (other Vtx)", 100, 0., 200.);
0433   vtxprob[0] = book<TH1D>("otherVtxProb", "#chi^{2} probability (other Vtx)", 100, 0., 1.);
0434   nans[0] = book<TH1D>("otherVtxNans", "Illegal values for x,y,z,xx,xy,xz,yy,yz,zz (other Vtx)", 9, 0.5, 9.5);
0435 
0436   nbtksinvtx[1] = book<TH1D>("tagVtxTrksNbr", "Reconstructed Tracks in Vertex (tagged Vtx)", 100, -0.5, 99.5);
0437   ntracksVsZ[1] =
0438       book<TProfile>("tagVtxTrksVsZ", "Reconstructed Tracks in Vertex (tagged Vtx) vs Z", 80, -20., 20., 0., 100., "");
0439   ntracksVsZ[1]->SetXTitle("z-bs");
0440   ntracksVsZ[1]->SetYTitle("#tracks");
0441 
0442   score[1] = book<TH1D>("tagVtxScore", "sqrt(score) (tagged Vtx)", 100, 0., 400.);
0443   trksWeight[1] = book<TH1D>("tagVtxTrksWeight", "Total weight of Tracks in Vertex (tagged Vtx)", 100, 0, 100.);
0444   vtxchi2[1] = book<TH1D>("tagVtxChi2", "#chi^{2} (tagged Vtx)", 100, 0., 200.);
0445   vtxndf[1] = book<TH1D>("tagVtxNdf", "ndof (tagged Vtx)", 100, 0., 200.);
0446   vtxprob[1] = book<TH1D>("tagVtxProb", "#chi^{2} probability (tagged Vtx)", 100, 0., 1.);
0447   nans[1] = book<TH1D>("tagVtxNans", "Illegal values for x,y,z,xx,xy,xz,yy,yz,zz (tagged Vtx)", 9, 0.5, 9.5);
0448 
0449   xrec[0] = book<TH1D>("otherPosX", "Position x Coordinate (other Vtx)", 100, vposx_ - 0.1, vposx_ + 0.1);
0450   yrec[0] = book<TH1D>("otherPosY", "Position y Coordinate (other Vtx)", 100, vposy_ - 0.1, vposy_ + 0.1);
0451   zrec[0] = book<TH1D>("otherPosZ", "Position z Coordinate (other Vtx)", 100, -20., 20.);
0452   xDiff[0] = book<TH1D>("otherDiffX", "X distance from BeamSpot (other Vtx)", 100, -500, 500);
0453   yDiff[0] = book<TH1D>("otherDiffY", "Y distance from BeamSpot (other Vtx)", 100, -500, 500);
0454   xerr[0] = book<TH1D>("otherErrX", "Uncertainty x Coordinate (other Vtx)", 100, 0., 100);
0455   yerr[0] = book<TH1D>("otherErrY", "Uncertainty y Coordinate (other Vtx)", 100, 0., 100);
0456   zerr[0] = book<TH1D>("otherErrZ", "Uncertainty z Coordinate (other Vtx)", 100, 0., 100);
0457   xerrVsTrks[0] = book<TH2D>(
0458       "otherErrVsWeightX", "Uncertainty x Coordinate vs. track weight (other Vtx)", 100, 0, 100., 100, 0., 100);
0459   yerrVsTrks[0] = book<TH2D>(
0460       "otherErrVsWeightY", "Uncertainty y Coordinate vs. track weight (other Vtx)", 100, 0, 100., 100, 0., 100);
0461   zerrVsTrks[0] = book<TH2D>(
0462       "otherErrVsWeightZ", "Uncertainty z Coordinate vs. track weight (other Vtx)", 100, 0, 100., 100, 0., 100);
0463 
0464   xrec[1] = book<TH1D>("tagPosX", "Position x Coordinate (tagged Vtx)", 100, vposx_ - 0.1, vposx_ + 0.1);
0465   yrec[1] = book<TH1D>("tagPosY", "Position y Coordinate (tagged Vtx)", 100, vposy_ - 0.1, vposy_ + 0.1);
0466   zrec[1] = book<TH1D>("tagPosZ", "Position z Coordinate (tagged Vtx)", 100, -20., 20.);
0467   xDiff[1] = book<TH1D>("tagDiffX", "X distance from BeamSpot (tagged Vtx)", 100, -500, 500);
0468   yDiff[1] = book<TH1D>("tagDiffY", "Y distance from BeamSpot (tagged Vtx)", 100, -500, 500);
0469   xerr[1] = book<TH1D>("tagErrX", "Uncertainty x Coordinate (tagged Vtx)", 100, 0., 100);
0470   yerr[1] = book<TH1D>("tagErrY", "Uncertainty y Coordinate (tagged Vtx)", 100, 0., 100);
0471   zerr[1] = book<TH1D>("tagErrZ", "Uncertainty z Coordinate (tagged Vtx)", 100, 0., 100);
0472   xerrVsTrks[1] = book<TH2D>(
0473       "tagErrVsWeightX", "Uncertainty x Coordinate vs. track weight (tagged Vtx)", 100, 0, 100., 100, 0., 100);
0474   yerrVsTrks[1] = book<TH2D>(
0475       "tagErrVsWeightY", "Uncertainty y Coordinate vs. track weight (tagged Vtx)", 100, 0, 100., 100, 0., 100);
0476   zerrVsTrks[1] = book<TH2D>(
0477       "tagErrVsWeightZ", "Uncertainty z Coordinate vs. track weight (tagged Vtx)", 100, 0, 100., 100, 0., 100);
0478 
0479   type[0] = book<TH1D>("otherType", "Vertex type (other Vtx)", 3, -0.5, 2.5);
0480   type[1] = book<TH1D>("tagType", "Vertex type (tagged Vtx)", 3, -0.5, 2.5);
0481 
0482   for (int i = 0; i < 2; ++i) {
0483     type[i]->GetXaxis()->SetBinLabel(1, "Valid, real");
0484     type[i]->GetXaxis()->SetBinLabel(2, "Valid, fake");
0485     type[i]->GetXaxis()->SetBinLabel(3, "Invalid");
0486   }
0487 
0488   bsX = book<TH1D>("bsX", "BeamSpot x0", 100, vposx_ - 0.1, vposx_ + 0.1);
0489   bsY = book<TH1D>("bsY", "BeamSpot y0", 100, vposy_ - 0.1, vposy_ + 0.1);
0490   bsZ = book<TH1D>("bsZ", "BeamSpot z0", 100, -2., 2.);
0491   bsSigmaZ = book<TH1D>("bsSigmaZ", "BeamSpot sigmaZ", 100, 0., 10.);
0492   bsDxdz = book<TH1D>("bsDxdz", "BeamSpot dxdz", 100, -0.0003, 0.0003);
0493   bsDydz = book<TH1D>("bsDydz", "BeamSpot dydz", 100, -0.0003, 0.0003);
0494   bsBeamWidthX = book<TH1D>("bsBeamWidthX", "BeamSpot BeamWidthX", 100, 0., 100.);
0495   bsBeamWidthY = book<TH1D>("bsBeamWidthY", "BeamSpot BeamWidthY", 100, 0., 100.);
0496   bsType = book<TH1D>("bsType", "BeamSpot type", 4, -1.5, 2.5);
0497   bsType->GetXaxis()->SetBinLabel(1, "Unknown");
0498   bsType->GetXaxis()->SetBinLabel(2, "Fake");
0499   bsType->GetXaxis()->SetBinLabel(3, "LHC");
0500   bsType->GetXaxis()->SetBinLabel(4, "Tracker");
0501 
0502   // repeated strings in titles
0503   std::string s_1 = "PV Tracks (p_{T} > 1 GeV)";
0504   std::string s_10 = "PV Tracks (p_{T} > 10 GeV)";
0505 
0506   ntracks = book<TH1D>("ntracks", fmt::sprintf("number of %s", s_1).c_str(), tkNoBin_, tkNoMin_, tkNoMax_);
0507   ntracks->SetXTitle(fmt::sprintf("Number of %s per Event", s_1).c_str());
0508   ntracks->SetYTitle("Number of Events");
0509 
0510   weight = book<TH1D>("weight", fmt::sprintf("weight of %s", s_1).c_str(), 100, 0., 1.);
0511   weight->SetXTitle(fmt::sprintf("weight of %s per Event", s_1).c_str());
0512   weight->SetYTitle("Number of Events");
0513 
0514   sumpt = book<TH1D>("sumpt", fmt::sprintf("#Sum p_{T} of %s", s_1).c_str(), 100, -0.5, 249.5);
0515   chi2ndf = book<TH1D>("chi2ndf", fmt::sprintf("%s #chi^{2}/ndof", s_1).c_str(), 100, 0., 20.);
0516   chi2prob = book<TH1D>("chi2prob", fmt::sprintf("%s #chi^{2} probability", s_1).c_str(), 100, 0., 1.);
0517   dxy = book<TH1D>("dxy", fmt::sprintf("%s d_{xy} (#mum)", s_1).c_str(), dxyBin_, dxyMin_, dxyMax_);
0518   dxy2 = book<TH1D>("dxyzoom", fmt::sprintf("%s d_{xy} (#mum)", s_1).c_str(), dxyBin_, dxyMin_ / 5., dxyMax_ / 5.);
0519   dxyErr = book<TH1D>("dxyErr", fmt::sprintf("%s d_{xy} error (#mum)", s_1).c_str(), 100, 0., 2000.);
0520   dz = book<TH1D>("dz", fmt::sprintf("%s d_{z} (#mum)", s_1).c_str(), dzBin_, dzMin_, dzMax_);
0521   dzErr = book<TH1D>("dzErr", fmt::sprintf("%s d_{z} error(#mum)", s_1).c_str(), 100, 0., 10000.);
0522 
0523   phi_pt1 =
0524       book<TH1D>("phi_pt1", fmt::sprintf("%s #phi; PV tracks #phi;#tracks", s_1).c_str(), phiBin_, phiMin_, phiMax_);
0525   eta_pt1 =
0526       book<TH1D>("eta_pt1", fmt::sprintf("%s #eta; PV tracks #eta;#tracks", s_1).c_str(), etaBin_, etaMin_, etaMax_);
0527   phi_pt10 =
0528       book<TH1D>("phi_pt10", fmt::sprintf("%s #phi; PV tracks #phi;#tracks", s_10).c_str(), phiBin_, phiMin_, phiMax_);
0529   eta_pt10 =
0530       book<TH1D>("eta_pt10", fmt::sprintf("%s #phi; PV tracks #eta;#tracks", s_10).c_str(), etaBin_, etaMin_, etaMax_);
0531 
0532   dxyVsPhi_pt1 = book<TProfile>("dxyVsPhi_pt1",
0533                                 fmt::sprintf("%s d_{xy} (#mum) VS track #phi", s_1).c_str(),
0534                                 phiBin_,
0535                                 phiMin_,
0536                                 phiMax_,
0537                                 dxyMin_,
0538                                 dxyMax_);
0539   dxyVsPhi_pt1->SetXTitle("PV track (p_{T} > 1 GeV) #phi");
0540   dxyVsPhi_pt1->SetYTitle("PV track (p_{T} > 1 GeV) d_{xy} (#mum)");
0541 
0542   dzVsPhi_pt1 = book<TProfile>("dzVsPhi_pt1",
0543                                fmt::sprintf("%s d_{z} (#mum) VS track #phi", s_1).c_str(),
0544                                phiBin_,
0545                                phiMin_,
0546                                phiMax_,
0547                                dzMin_,
0548                                dzMax_);
0549   dzVsPhi_pt1->SetXTitle("PV track (p_{T} > 1 GeV) #phi");
0550   dzVsPhi_pt1->SetYTitle("PV track (p_{T} > 1 GeV) d_{z} (#mum)");
0551 
0552   dxyVsEta_pt1 = book<TProfile>("dxyVsEta_pt1",
0553                                 fmt::sprintf("%s d_{xy} (#mum) VS track #eta", s_1).c_str(),
0554                                 etaBin_,
0555                                 etaMin_,
0556                                 etaMax_,
0557                                 dxyMin_,
0558                                 dxyMax_);
0559   dxyVsEta_pt1->SetXTitle("PV track (p_{T} > 1 GeV) #eta");
0560   dxyVsEta_pt1->SetYTitle("PV track (p_{T} > 1 GeV) d_{xy} (#mum)");
0561 
0562   dzVsEta_pt1 = book<TProfile>("dzVsEta_pt1",
0563                                fmt::sprintf("%s d_{z} (#mum) VS track #eta", s_1).c_str(),
0564                                etaBin_,
0565                                etaMin_,
0566                                etaMax_,
0567                                dzMin_,
0568                                dzMax_);
0569   dzVsEta_pt1->SetXTitle("PV track (p_{T} > 1 GeV) #eta");
0570   dzVsEta_pt1->SetYTitle("PV track (p_{T} > 1 GeV) d_{z} (#mum)");
0571 
0572   dxyVsEtaVsPhi_pt1 = book<TProfile2D>("dxyVsEtaVsPhi_pt1",
0573                                        fmt::sprintf("%s d_{xy} (#mum) VS track #eta VS track #phi", s_1).c_str(),
0574                                        etaBin2D_,
0575                                        etaMin_,
0576                                        etaMax_,
0577                                        phiBin2D_,
0578                                        phiMin_,
0579                                        phiMax_,
0580                                        dxyMin_,
0581                                        dxyMax_);
0582   dxyVsEtaVsPhi_pt1->SetXTitle("PV track (p_{T} > 1 GeV) #eta");
0583   dxyVsEtaVsPhi_pt1->SetYTitle("PV track (p_{T} > 1 GeV) #phi");
0584   dxyVsEtaVsPhi_pt1->SetZTitle("PV track (p_{T} > 1 GeV) d_{xy} (#mum)");
0585 
0586   dzVsEtaVsPhi_pt1 = book<TProfile2D>("dzVsEtaVsPhi_pt1",
0587                                       fmt::sprintf("%s d_{z} (#mum) VS track #eta VS track #phi", s_1).c_str(),
0588                                       etaBin2D_,
0589                                       etaMin_,
0590                                       etaMax_,
0591                                       phiBin2D_,
0592                                       phiMin_,
0593                                       phiMax_,
0594                                       dzMin_,
0595                                       dzMax_);
0596   dzVsEtaVsPhi_pt1->SetXTitle("PV track (p_{T} > 1 GeV) #eta");
0597   dzVsEtaVsPhi_pt1->SetYTitle("PV track (p_{T} > 1 GeV) #phi");
0598   dzVsEtaVsPhi_pt1->SetZTitle("PV track (p_{T} > 1 GeV) d_{z} (#mum)");
0599 
0600   dxyVsPhi_pt10 = book<TProfile>("dxyVsPhi_pt10",
0601                                  fmt::sprintf("%s d_{xy} (#mum) VS track #phi", s_10).c_str(),
0602                                  phiBin_,
0603                                  phiMin_,
0604                                  phiMax_,
0605                                  dxyMin_,
0606                                  dxyMax_);
0607   dxyVsPhi_pt10->SetXTitle("PV track (p_{T} > 10 GeV) #phi");
0608   dxyVsPhi_pt10->SetYTitle("PV track (p_{T} > 10 GeV) d_{xy} (#mum)");
0609 
0610   dzVsPhi_pt10 = book<TProfile>("dzVsPhi_pt10",
0611                                 fmt::sprintf("%s d_{z} (#mum) VS track #phi", s_10).c_str(),
0612                                 phiBin_,
0613                                 phiMin_,
0614                                 phiMax_,
0615                                 dzMin_,
0616                                 dzMax_);
0617   dzVsPhi_pt10->SetXTitle("PV track (p_{T} > 10 GeV) #phi");
0618   dzVsPhi_pt10->SetYTitle("PV track (p_{T} > 10 GeV) d_{z} (#mum)");
0619 
0620   dxyVsEta_pt10 = book<TProfile>("dxyVsEta_pt10",
0621                                  fmt::sprintf("%s d_{xy} (#mum) VS track #eta", s_10).c_str(),
0622                                  etaBin_,
0623                                  etaMin_,
0624                                  etaMax_,
0625                                  dxyMin_,
0626                                  dxyMax_);
0627   dxyVsEta_pt10->SetXTitle("PV track (p_{T} > 10 GeV) #eta");
0628   dxyVsEta_pt10->SetYTitle("PV track (p_{T} > 10 GeV) d_{xy} (#mum)");
0629 
0630   dzVsEta_pt10 = book<TProfile>("dzVsEta_pt10",
0631                                 fmt::sprintf("%s d_{z} (#mum) VS track #eta", s_10).c_str(),
0632                                 etaBin_,
0633                                 etaMin_,
0634                                 etaMax_,
0635                                 dzMin_,
0636                                 dzMax_);
0637   dzVsEta_pt10->SetXTitle("PV track (p_{T} > 10 GeV) #eta");
0638   dzVsEta_pt10->SetYTitle("PV track (p_{T} > 10 GeV) d_{z} (#mum)");
0639 
0640   dxyVsEtaVsPhi_pt10 = book<TProfile2D>("dxyVsEtaVsPhi_pt10",
0641                                         fmt::sprintf("%s d_{xy} (#mum) VS track #eta VS track #phi", s_10).c_str(),
0642                                         etaBin2D_,
0643                                         etaMin_,
0644                                         etaMax_,
0645                                         phiBin2D_,
0646                                         phiMin_,
0647                                         phiMax_,
0648                                         dxyMin_,
0649                                         dxyMax_);
0650   dxyVsEtaVsPhi_pt10->SetXTitle("PV track (p_{T} > 10 GeV) #eta");
0651   dxyVsEtaVsPhi_pt10->SetYTitle("PV track (p_{T} > 10 GeV) #phi");
0652   dxyVsEtaVsPhi_pt10->SetZTitle("PV track (p_{T} > 10 GeV) d_{xy} (#mum)");
0653 
0654   dzVsEtaVsPhi_pt10 = book<TProfile2D>("dzVsEtaVsPhi_pt10",
0655                                        fmt::sprintf("%s d_{z} (#mum) VS track #eta VS track #phi", s_10).c_str(),
0656                                        etaBin2D_,
0657                                        etaMin_,
0658                                        etaMax_,
0659                                        phiBin2D_,
0660                                        phiMin_,
0661                                        phiMax_,
0662                                        dzMin_,
0663                                        dzMax_);
0664   dzVsEtaVsPhi_pt10->SetXTitle("PV track (p_{T} > 10 GeV) #eta");
0665   dzVsEtaVsPhi_pt10->SetYTitle("PV track (p_{T} > 10 GeV) #phi");
0666   dzVsEtaVsPhi_pt10->SetZTitle("PV track (p_{T} > 10 GeV) d_{z} (#mum)");
0667 }
0668 
0669 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0670 void GeneralPurposeVertexAnalyzer::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0671   edm::ParameterSetDescription desc;
0672   desc.add<int>("ndof", 4);
0673   desc.add<edm::InputTag>("vertexLabel", edm::InputTag("offlinePrimaryVertices"));
0674   desc.add<edm::InputTag>("beamSpotLabel", edm::InputTag("offlineBeamSpot"));
0675   desc.add<double>("Xpos", 0.1);
0676   desc.add<double>("Ypos", 0.0);
0677   desc.add<int>("TkSizeBin", 100);
0678   desc.add<double>("TkSizeMin", 499.5);
0679   desc.add<double>("TkSizeMax", -0.5);
0680   desc.add<int>("DxyBin", 100);
0681   desc.add<double>("DxyMin", 5000.);
0682   desc.add<double>("DxyMax", -5000.);
0683   desc.add<int>("DzBin", 100);
0684   desc.add<double>("DzMin", -2000.0);
0685   desc.add<double>("DzMax", 2000.0);
0686   desc.add<int>("PhiBin", 32);
0687   desc.add<int>("PhiBin2D", 12);
0688   desc.add<double>("PhiMin", -M_PI);
0689   desc.add<double>("PhiMax", M_PI);
0690   desc.add<int>("EtaBin", 26);
0691   desc.add<int>("EtaBin2D", 8);
0692   desc.add<double>("EtaMin", -2.7);
0693   desc.add<double>("EtaMax", 2.7);
0694   descriptions.addWithDefaultLabel(desc);
0695 }
0696 
0697 //define this as a plug-in
0698 DEFINE_FWK_MODULE(GeneralPurposeVertexAnalyzer);