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