File indexing completed on 2023-03-17 10:56:11
0001 #include "DQM/SiPixelCommon/interface/SiPixelHistogramId.h"
0002 #include "DQM/SiPixelMonitorRecHit/interface/SiPixelRecHitModule.h"
0003 #include "DQMServices/Core/interface/DQMStore.h"
0004
0005 #include "FWCore/ServiceRegistry/interface/Service.h"
0006
0007 #include <cstdlib>
0008 #include <iostream>
0009 #include <memory>
0010 #include <string>
0011 #include <vector>
0012
0013
0014 #include "DataFormats/DetId/interface/DetId.h"
0015 #include "DataFormats/TrackerCommon/interface/PixelBarrelName.h"
0016 #include "DataFormats/SiPixelDetId/interface/PixelBarrelNameUpgrade.h"
0017 #include "DataFormats/TrackerCommon/interface/PixelEndcapName.h"
0018 #include "DataFormats/SiPixelDetId/interface/PixelEndcapNameUpgrade.h"
0019 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0020
0021
0022
0023 SiPixelRecHitModule::SiPixelRecHitModule() : id_(0) {}
0024
0025 SiPixelRecHitModule::SiPixelRecHitModule(const uint32_t &id) : id_(id) {}
0026
0027
0028
0029
0030 SiPixelRecHitModule::~SiPixelRecHitModule() {}
0031
0032
0033
0034 void SiPixelRecHitModule::book(const edm::ParameterSet &iConfig,
0035 DQMStore::IBooker &iBooker,
0036 const TrackerTopology *pTT,
0037 int type,
0038 bool twoD,
0039 bool reducedSet,
0040 bool isUpgrade) {
0041 bool barrel = DetId(id_).subdetId() == static_cast<int>(PixelSubdetector::PixelBarrel);
0042 bool endcap = DetId(id_).subdetId() == static_cast<int>(PixelSubdetector::PixelEndcap);
0043 bool isHalfModule = false;
0044
0045 if (barrel) {
0046
0047 isHalfModule = PixelBarrelName(DetId(id_), pTT, isUpgrade).isHalfModule();
0048
0049
0050
0051 }
0052
0053 std::string hid;
0054
0055 edm::InputTag src = iConfig.getParameter<edm::InputTag>("src");
0056
0057
0058 if (type == 0) {
0059 SiPixelHistogramId *theHistogramId = new SiPixelHistogramId(src.label());
0060 if (!reducedSet) {
0061 if (twoD) {
0062
0063 hid = theHistogramId->setHistoId("xypos", id_);
0064 meXYPos_ = iBooker.book2D(hid, "XY Position", 100, -1., 1, 100, -4, 4);
0065 meXYPos_->setAxisTitle("X Position", 1);
0066 meXYPos_->setAxisTitle("Y Position", 2);
0067 } else {
0068
0069 hid = theHistogramId->setHistoId("xypos", id_);
0070 meXYPos_px_ = iBooker.book1D(hid + "_px", "X Position", 100, -1., 1);
0071 meXYPos_px_->setAxisTitle("X Position", 1);
0072 meXYPos_py_ = iBooker.book1D(hid + "_py", "Y Position", 100, -4, 4);
0073 meXYPos_py_->setAxisTitle("Y Position", 1);
0074 }
0075 }
0076 hid = theHistogramId->setHistoId("ClustX", id_);
0077 meClustX_ = iBooker.book1D(hid, "RecHit X size", 10, 0., 10.);
0078 meClustX_->setAxisTitle("RecHit size X dimension", 1);
0079 hid = theHistogramId->setHistoId("ClustY", id_);
0080 meClustY_ = iBooker.book1D(hid, "RecHit Y size", 15, 0., 15.);
0081 meClustY_->setAxisTitle("RecHit size Y dimension", 1);
0082
0083 hid = theHistogramId->setHistoId("ErrorX", id_);
0084 meErrorX_ = iBooker.book1D(hid, "RecHit error X", 100, 0., 0.02);
0085 meErrorX_->setAxisTitle("RecHit error X", 1);
0086 hid = theHistogramId->setHistoId("ErrorY", id_);
0087 meErrorY_ = iBooker.book1D(hid, "RecHit error Y", 100, 0., 0.02);
0088 meErrorY_->setAxisTitle("RecHit error Y", 1);
0089
0090
0091
0092
0093
0094 delete theHistogramId;
0095 }
0096
0097 if (type == 1 && barrel) {
0098 uint32_t DBladder;
0099
0100
0101 DBladder = PixelBarrelName(DetId(id_), pTT, isUpgrade).ladderName();
0102 char sladder[80];
0103 sprintf(sladder, "Ladder_%02i", DBladder);
0104 hid = src.label() + "_" + sladder;
0105 if (isHalfModule)
0106 hid += "H";
0107 else
0108 hid += "F";
0109 if (!reducedSet) {
0110 if (twoD) {
0111 meXYPosLad_ = iBooker.book2D("xypos_" + hid, "XY Position", 100, -1., 1, 100, -4, 4);
0112 meXYPosLad_->setAxisTitle("X Position", 1);
0113 meXYPosLad_->setAxisTitle("Y Position", 2);
0114 } else {
0115
0116 meXYPosLad_px_ = iBooker.book1D("xypos_" + hid + "_px", "X Position", 100, -1., 1);
0117 meXYPosLad_px_->setAxisTitle("X Position", 1);
0118 meXYPosLad_py_ = iBooker.book1D("xypos_" + hid + "_py", "Y Position", 100, -4, 4);
0119 meXYPosLad_py_->setAxisTitle("Y Position", 1);
0120 }
0121 }
0122 meClustXLad_ = iBooker.book1D("ClustX_" + hid, "RecHit X size", 10, 0., 10.);
0123 meClustXLad_->setAxisTitle("RecHit size X dimension", 1);
0124 meClustYLad_ = iBooker.book1D("ClustY_" + hid, "RecHit Y size", 15, 0., 15.);
0125 meClustYLad_->setAxisTitle("RecHit size Y dimension", 1);
0126 meErrorXLad_ = iBooker.book1D("ErrorX_" + hid, "RecHit error X", 100, 0., 0.02);
0127 meErrorXLad_->setAxisTitle("RecHit error X", 1);
0128 meErrorYLad_ = iBooker.book1D("ErrorY_" + hid, "RecHit error Y", 100, 0., 0.02);
0129 meErrorYLad_->setAxisTitle("RecHit error Y", 1);
0130 menRecHitsLad_ = iBooker.book1D("nRecHits_" + hid, "# of rechits in this module", 8, 0, 8);
0131 menRecHitsLad_->setAxisTitle("number of rechits", 1);
0132 }
0133
0134 if (type == 2 && barrel) {
0135 uint32_t DBlayer;
0136
0137
0138 DBlayer = PixelBarrelName(DetId(id_), pTT, isUpgrade).layerName();
0139 char slayer[80];
0140 sprintf(slayer, "Layer_%i", DBlayer);
0141 hid = src.label() + "_" + slayer;
0142
0143 if (!reducedSet) {
0144 if (twoD) {
0145 meXYPosLay_ = iBooker.book2D("xypos_" + hid, "XY Position", 100, -1., 1, 100, -4, 4);
0146 meXYPosLay_->setAxisTitle("X Position", 1);
0147 meXYPosLay_->setAxisTitle("Y Position", 2);
0148 } else {
0149
0150 meXYPosLay_px_ = iBooker.book1D("xypos_" + hid + "_px", "X Position", 100, -1., 1);
0151 meXYPosLay_px_->setAxisTitle("X Position", 1);
0152 meXYPosLay_py_ = iBooker.book1D("xypos_" + hid + "_py", "Y Position", 100, -4, 4);
0153 meXYPosLay_py_->setAxisTitle("Y Position", 1);
0154 }
0155 }
0156
0157 meClustXLay_ = iBooker.book1D("ClustX_" + hid, "RecHit X size", 10, 0., 10.);
0158 meClustXLay_->setAxisTitle("RecHit size X dimension", 1);
0159 meClustYLay_ = iBooker.book1D("ClustY_" + hid, "RecHit Y size", 15, 0., 15.);
0160 meClustYLay_->setAxisTitle("RecHit size Y dimension", 1);
0161 meErrorXLay_ = iBooker.book1D("ErrorX_" + hid, "RecHit error X", 100, 0., 0.02);
0162 meErrorXLay_->setAxisTitle("RecHit error X", 1);
0163 meErrorYLay_ = iBooker.book1D("ErrorY_" + hid, "RecHit error Y", 100, 0., 0.02);
0164 meErrorYLay_->setAxisTitle("RecHit error Y", 1);
0165 menRecHitsLay_ = iBooker.book1D("nRecHits_" + hid, "# of rechits in this module", 8, 0, 8);
0166 menRecHitsLay_->setAxisTitle("number of rechits", 1);
0167 }
0168
0169 if (type == 3 && barrel) {
0170 uint32_t DBmodule;
0171
0172
0173 DBmodule = PixelBarrelName(DetId(id_), pTT, isUpgrade).moduleName();
0174 char smodule[80];
0175 sprintf(smodule, "Ring_%i", DBmodule);
0176 hid = src.label() + "_" + smodule;
0177
0178 if (!reducedSet) {
0179 if (twoD) {
0180 meXYPosPhi_ = iBooker.book2D("xypos_" + hid, "XY Position", 100, -1., 1, 100, -4, 4);
0181 meXYPosPhi_->setAxisTitle("X Position", 1);
0182 meXYPosPhi_->setAxisTitle("Y Position", 2);
0183 } else {
0184
0185 meXYPosPhi_px_ = iBooker.book1D("xypos_" + hid + "_px", "X Position", 100, -1., 1);
0186 meXYPosPhi_px_->setAxisTitle("X Position", 1);
0187 meXYPosPhi_py_ = iBooker.book1D("xypos_" + hid + "_py", "Y Position", 100, -4, 4);
0188 meXYPosPhi_py_->setAxisTitle("Y Position", 1);
0189 }
0190 }
0191 meClustXPhi_ = iBooker.book1D("ClustX_" + hid, "RecHit X size", 10, 0., 10.);
0192 meClustXPhi_->setAxisTitle("RecHit size X dimension", 1);
0193 meClustYPhi_ = iBooker.book1D("ClustY_" + hid, "RecHit Y size", 15, 0., 15.);
0194 meClustYPhi_->setAxisTitle("RecHit size Y dimension", 1);
0195 meErrorXPhi_ = iBooker.book1D("ErrorX_" + hid, "RecHit error X", 100, 0., 0.02);
0196 meErrorXPhi_->setAxisTitle("RecHit error X", 1);
0197 meErrorYPhi_ = iBooker.book1D("ErrorY_" + hid, "RecHit error Y", 100, 0., 0.02);
0198 meErrorYPhi_->setAxisTitle("RecHit error Y", 1);
0199 menRecHitsPhi_ = iBooker.book1D("nRecHits_" + hid, "# of rechits in this module", 8, 0, 8);
0200 menRecHitsPhi_->setAxisTitle("number of rechits", 1);
0201 }
0202
0203 if (type == 4 && endcap) {
0204 uint32_t blade;
0205
0206
0207 blade = PixelEndcapName(DetId(id_), pTT, isUpgrade).bladeName();
0208
0209 char sblade[80];
0210 sprintf(sblade, "Blade_%02i", blade);
0211 hid = src.label() + "_" + sblade;
0212
0213
0214
0215
0216 meClustXBlade_ = iBooker.book1D("ClustX_" + hid, "RecHit X size", 10, 0., 10.);
0217 meClustXBlade_->setAxisTitle("RecHit size X dimension", 1);
0218 meClustYBlade_ = iBooker.book1D("ClustY_" + hid, "RecHit Y size", 15, 0., 15.);
0219 meClustYBlade_->setAxisTitle("RecHit size Y dimension", 1);
0220 meErrorXBlade_ = iBooker.book1D("ErrorX_" + hid, "RecHit error X", 100, 0., 0.02);
0221 meErrorXBlade_->setAxisTitle("RecHit error X", 1);
0222 meErrorYBlade_ = iBooker.book1D("ErrorY_" + hid, "RecHit error Y", 100, 0., 0.02);
0223 meErrorYBlade_->setAxisTitle("RecHit error Y", 1);
0224 menRecHitsBlade_ = iBooker.book1D("nRecHits_" + hid, "# of rechits in this module", 8, 0, 8);
0225 menRecHitsBlade_->setAxisTitle("number of rechits", 1);
0226 }
0227 if (type == 5 && endcap) {
0228 uint32_t disk;
0229
0230
0231 disk = PixelEndcapName(DetId(id_), pTT, isUpgrade).diskName();
0232
0233 char sdisk[80];
0234 sprintf(sdisk, "Disk_%i", disk);
0235 hid = src.label() + "_" + sdisk;
0236
0237
0238
0239
0240 meClustXDisk_ = iBooker.book1D("ClustX_" + hid, "RecHit X size", 10, 0., 10.);
0241 meClustXDisk_->setAxisTitle("RecHit size X dimension", 1);
0242 meClustYDisk_ = iBooker.book1D("ClustY_" + hid, "RecHit Y size", 15, 0., 15.);
0243 meClustYDisk_->setAxisTitle("RecHit size Y dimension", 1);
0244 meErrorXDisk_ = iBooker.book1D("ErrorX_" + hid, "RecHit error X", 100, 0., 0.02);
0245 meErrorXDisk_->setAxisTitle("RecHit error X", 1);
0246 meErrorYDisk_ = iBooker.book1D("ErrorY_" + hid, "RecHit error Y", 100, 0., 0.02);
0247 meErrorYDisk_->setAxisTitle("RecHit error Y", 1);
0248 menRecHitsDisk_ = iBooker.book1D("nRecHits_" + hid, "# of rechits in this module", 8, 0, 8);
0249 menRecHitsDisk_->setAxisTitle("number of rechits", 1);
0250 }
0251
0252 if (type == 6 && endcap) {
0253 uint32_t panel;
0254 uint32_t module;
0255
0256
0257
0258
0259
0260
0261
0262 panel = PixelEndcapName(DetId(id_), pTT, isUpgrade).pannelName();
0263 module = PixelEndcapName(DetId(id_), pTT, isUpgrade).plaquetteName();
0264
0265 char slab[80];
0266 sprintf(slab, "Panel_%i_Ring_%i", panel, module);
0267 hid = src.label() + "_" + slab;
0268
0269 if (!reducedSet) {
0270 if (twoD) {
0271 meXYPosRing_ = iBooker.book2D("xypos_" + hid, "XY Position", 100, -1., 1, 100, -4, 4);
0272 meXYPosRing_->setAxisTitle("X Position", 1);
0273 meXYPosRing_->setAxisTitle("Y Position", 2);
0274 } else {
0275
0276 meXYPosRing_px_ = iBooker.book1D("xypos_" + hid + "_px", "X Position", 100, -1., 1);
0277 meXYPosRing_px_->setAxisTitle("X Position", 1);
0278 meXYPosRing_py_ = iBooker.book1D("xypos_" + hid + "_py", "Y Position", 100, -4, 4);
0279 meXYPosRing_py_->setAxisTitle("Y Position", 1);
0280 }
0281 }
0282 meClustXRing_ = iBooker.book1D("ClustX_" + hid, "RecHit X size", 10, 0., 10.);
0283 meClustXRing_->setAxisTitle("RecHit size X dimension", 1);
0284 meClustYRing_ = iBooker.book1D("ClustY_" + hid, "RecHit Y size", 15, 0., 15.);
0285 meClustYRing_->setAxisTitle("RecHit size Y dimension", 1);
0286 meErrorXRing_ = iBooker.book1D("ErrorX_" + hid, "RecHit error X", 100, 0., 0.02);
0287 meErrorXRing_->setAxisTitle("RecHit error X", 1);
0288 meErrorYRing_ = iBooker.book1D("ErrorY_" + hid, "RecHit error Y", 100, 0., 0.02);
0289 meErrorYRing_->setAxisTitle("RecHit error Y", 1);
0290 menRecHitsRing_ = iBooker.book1D("nRecHits_" + hid, "# of rechits in this module", 8, 0, 8);
0291 menRecHitsRing_->setAxisTitle("number of rechits", 1);
0292 }
0293 }
0294
0295
0296
0297 void SiPixelRecHitModule::fill(const float &rechit_x,
0298 const float &rechit_y,
0299 const int &sizeX,
0300 const int &sizeY,
0301 const float &lerr_x,
0302 const float &lerr_y,
0303 bool modon,
0304 bool ladon,
0305 bool layon,
0306 bool phion,
0307 bool bladeon,
0308 bool diskon,
0309 bool ringon,
0310 bool twoD,
0311 bool reducedSet) {
0312 bool barrel = DetId(id_).subdetId() == static_cast<int>(PixelSubdetector::PixelBarrel);
0313 bool endcap = DetId(id_).subdetId() == static_cast<int>(PixelSubdetector::PixelEndcap);
0314
0315
0316
0317 if (modon) {
0318 meClustX_->Fill(sizeX);
0319 meClustY_->Fill(sizeY);
0320 meErrorX_->Fill(lerr_x);
0321 meErrorY_->Fill(lerr_y);
0322 if (!reducedSet) {
0323 if (twoD)
0324 meXYPos_->Fill(rechit_x, rechit_y);
0325 else {
0326 meXYPos_px_->Fill(rechit_x);
0327 meXYPos_py_->Fill(rechit_y);
0328 }
0329 }
0330 }
0331
0332
0333 if (ladon && barrel) {
0334 meClustXLad_->Fill(sizeX);
0335 meClustYLad_->Fill(sizeY);
0336 meErrorXLad_->Fill(lerr_x);
0337 meErrorYLad_->Fill(lerr_y);
0338 if (!reducedSet) {
0339 if (twoD)
0340 meXYPosLad_->Fill(rechit_x, rechit_y);
0341 else {
0342 meXYPosLad_px_->Fill(rechit_x);
0343 meXYPosLad_py_->Fill(rechit_y);
0344 }
0345 }
0346 }
0347
0348 if (layon && barrel) {
0349 meClustXLay_->Fill(sizeX);
0350 meClustYLay_->Fill(sizeY);
0351 meErrorXLay_->Fill(lerr_x);
0352 meErrorYLay_->Fill(lerr_y);
0353 if (!reducedSet) {
0354 if (twoD)
0355 meXYPosLay_->Fill(rechit_x, rechit_y);
0356 else {
0357 meXYPosLay_px_->Fill(rechit_x);
0358 meXYPosLay_py_->Fill(rechit_y);
0359 }
0360 }
0361 }
0362
0363 if (phion && barrel) {
0364 meClustXPhi_->Fill(sizeX);
0365 meClustYPhi_->Fill(sizeY);
0366 meErrorXPhi_->Fill(lerr_x);
0367 meErrorYPhi_->Fill(lerr_y);
0368 if (!reducedSet) {
0369 if (twoD)
0370 meXYPosPhi_->Fill(rechit_x, rechit_y);
0371 else {
0372 meXYPosPhi_px_->Fill(rechit_x);
0373 meXYPosPhi_py_->Fill(rechit_y);
0374 }
0375 }
0376 }
0377
0378 if (bladeon && endcap) {
0379
0380 meClustXBlade_->Fill(sizeX);
0381 meClustYBlade_->Fill(sizeY);
0382 meErrorXBlade_->Fill(lerr_x);
0383 meErrorYBlade_->Fill(lerr_y);
0384 }
0385
0386 if (diskon && endcap) {
0387
0388 meClustXDisk_->Fill(sizeX);
0389 meClustYDisk_->Fill(sizeY);
0390 meErrorXDisk_->Fill(lerr_x);
0391 meErrorYDisk_->Fill(lerr_y);
0392 }
0393
0394 if (ringon && endcap) {
0395 meClustXRing_->Fill(sizeX);
0396 meClustYRing_->Fill(sizeY);
0397 meErrorXRing_->Fill(lerr_x);
0398 meErrorYRing_->Fill(lerr_y);
0399 if (!reducedSet) {
0400 if (twoD)
0401 meXYPosRing_->Fill(rechit_x, rechit_y);
0402 else {
0403 meXYPosRing_px_->Fill(rechit_x);
0404 meXYPosRing_py_->Fill(rechit_y);
0405 }
0406 }
0407 }
0408 }
0409
0410 void SiPixelRecHitModule::nfill(
0411 const int &nrec, bool modon, bool ladon, bool layon, bool phion, bool bladeon, bool diskon, bool ringon) {
0412 bool barrel = DetId(id_).subdetId() == static_cast<int>(PixelSubdetector::PixelBarrel);
0413 bool endcap = DetId(id_).subdetId() == static_cast<int>(PixelSubdetector::PixelEndcap);
0414
0415
0416
0417 if (ladon && barrel)
0418 menRecHitsLad_->Fill(nrec);
0419 if (layon && barrel)
0420 menRecHitsLay_->Fill(nrec);
0421 if (phion && barrel)
0422 menRecHitsPhi_->Fill(nrec);
0423
0424 if (bladeon && endcap)
0425 menRecHitsBlade_->Fill(nrec);
0426 if (diskon && endcap)
0427 menRecHitsDisk_->Fill(nrec);
0428 if (ringon && endcap)
0429 menRecHitsRing_->Fill(nrec);
0430 }