Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 09:43:07

0001 // -*- C++ -*-
0002 //
0003 // Package:    SiPixelMonitorRecHits
0004 // Class:      SiPixelRecHitSource
0005 //
0006 /**\class
0007 
0008  Description: Pixel DQM source for RecHits
0009 
0010  Implementation:
0011      Originally based on the code for Digis, adapted
0012         to read RecHits and create relevant histograms
0013 */
0014 //
0015 // Original Author:  Vincenzo Chiochia
0016 //         Created:
0017 //
0018 //
0019 // Adapted by:  Keith Rose
0020 //      For use in SiPixelMonitorClient for RecHits
0021 // Updated by: Lukas Wehrli
0022 // for pixel offline DQM
0023 
0024 #include "DQM/SiPixelMonitorRecHit/interface/SiPixelRecHitSource.h"
0025 // Framework
0026 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0027 #include "FWCore/ServiceRegistry/interface/Service.h"
0028 #include "FWCore/Framework/interface/MakerMacros.h"
0029 // DQM Framework
0030 #include "DQM/SiPixelCommon/interface/SiPixelFolderOrganizer.h"
0031 #include "DQMServices/Core/interface/DQMStore.h"
0032 // Geometry
0033 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0034 // DataFormats
0035 #include "DataFormats/DetId/interface/DetId.h"
0036 #include "DataFormats/TrackerCommon/interface/PixelBarrelName.h"
0037 #include "DataFormats/SiPixelDetId/interface/PixelBarrelNameUpgrade.h"
0038 #include "DataFormats/TrackerCommon/interface/PixelEndcapName.h"
0039 #include "DataFormats/SiPixelDetId/interface/PixelEndcapNameUpgrade.h"
0040 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0041 
0042 //
0043 #include <cstdlib>
0044 #include <iostream>
0045 #include <string>
0046 using namespace std;
0047 using namespace edm;
0048 
0049 SiPixelRecHitSource::SiPixelRecHitSource(const edm::ParameterSet &iConfig)
0050     : conf_(iConfig),
0051       src_(consumes<SiPixelRecHitCollection>(conf_.getParameter<edm::InputTag>("src"))),
0052       trackerTopoTokenBeginRun_(esConsumes<TrackerTopology, TrackerTopologyRcd, edm::Transition::BeginRun>()),
0053       trackerGeomTokenBeginRun_(esConsumes<TrackerGeometry, TrackerDigiGeometryRecord, edm::Transition::BeginRun>()),
0054       saveFile(conf_.getUntrackedParameter<bool>("saveFile", false)),
0055       isPIB(conf_.getUntrackedParameter<bool>("isPIB", false)),
0056       slowDown(conf_.getUntrackedParameter<bool>("slowDown", false)),
0057       modOn(conf_.getUntrackedParameter<bool>("modOn", true)),
0058       twoDimOn(conf_.getUntrackedParameter<bool>("twoDimOn", true)),
0059       reducedSet(conf_.getUntrackedParameter<bool>("reducedSet", false)),
0060       ladOn(conf_.getUntrackedParameter<bool>("ladOn", false)),
0061       layOn(conf_.getUntrackedParameter<bool>("layOn", false)),
0062       phiOn(conf_.getUntrackedParameter<bool>("phiOn", false)),
0063       ringOn(conf_.getUntrackedParameter<bool>("ringOn", false)),
0064       bladeOn(conf_.getUntrackedParameter<bool>("bladeOn", false)),
0065       diskOn(conf_.getUntrackedParameter<bool>("diskOn", false)),
0066       isUpgrade(conf_.getUntrackedParameter<bool>("isUpgrade", false)) {
0067   firstRun = true;
0068   LogInfo("PixelDQM") << "SiPixelRecHitSource::SiPixelRecHitSource: Got DQM BackEnd interface" << endl;
0069   topFolderName_ = conf_.getParameter<std::string>("TopFolderName");
0070 }
0071 
0072 SiPixelRecHitSource::~SiPixelRecHitSource() {
0073   // do anything here that needs to be done at desctruction time
0074   // (e.g. close files, deallocate resources etc.)
0075   LogInfo("PixelDQM") << "SiPixelRecHitSource::~SiPixelRecHitSource: Destructor" << endl;
0076   std::map<uint32_t, SiPixelRecHitModule *>::iterator struct_iter;
0077   for (struct_iter = thePixelStructure.begin(); struct_iter != thePixelStructure.end(); struct_iter++) {
0078     delete struct_iter->second;
0079     struct_iter->second = nullptr;
0080   }
0081 }
0082 
0083 void SiPixelRecHitSource::dqmBeginRun(const edm::Run &r, const edm::EventSetup &iSetup) {
0084   LogInfo("PixelDQM") << " SiPixelRecHitSource::beginJob - Initialisation ... " << std::endl;
0085   LogInfo("PixelDQM") << "Mod/Lad/Lay/Phi " << modOn << "/" << ladOn << "/" << layOn << "/" << phiOn << std::endl;
0086   LogInfo("PixelDQM") << "Blade/Disk/Ring" << bladeOn << "/" << diskOn << "/" << ringOn << std::endl;
0087   LogInfo("PixelDQM") << "2DIM IS " << twoDimOn << "\n";
0088 
0089   if (firstRun) {
0090     eventNo = 0;
0091     // Build map
0092     buildStructure(iSetup);
0093     // Book Monitoring Elements
0094     firstRun = false;
0095   }
0096 }
0097 
0098 void SiPixelRecHitSource::bookHistograms(DQMStore::IBooker &iBooker, edm::Run const &, const edm::EventSetup &iSetup) {
0099   bookMEs(iBooker, iSetup);
0100 }
0101 
0102 //------------------------------------------------------------------
0103 // Method called for every event
0104 //------------------------------------------------------------------
0105 void SiPixelRecHitSource::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0106   eventNo++;
0107   // cout << eventNo << endl;
0108   // get input data
0109   edm::Handle<SiPixelRecHitCollection> recHitColl;
0110   iEvent.getByToken(src_, recHitColl);
0111 
0112   std::map<uint32_t, SiPixelRecHitModule *>::iterator struct_iter;
0113   for (struct_iter = thePixelStructure.begin(); struct_iter != thePixelStructure.end(); struct_iter++) {
0114     uint32_t TheID = (*struct_iter).first;
0115 
0116     SiPixelRecHitCollection::const_iterator match = recHitColl->find(TheID);
0117 
0118     // if( pixelrechitRangeIteratorBegin == pixelrechitRangeIteratorEnd) {cout
0119     // << "oops" << endl;}
0120     float rechit_x = 0;
0121     float rechit_y = 0;
0122     int rechit_count = 0;
0123 
0124     if (match != recHitColl->end()) {
0125       SiPixelRecHitCollection::DetSet pixelrechitRange = *match;
0126       SiPixelRecHitCollection::DetSet::const_iterator pixelrechitRangeIteratorBegin = pixelrechitRange.begin();
0127       SiPixelRecHitCollection::DetSet::const_iterator pixelrechitRangeIteratorEnd = pixelrechitRange.end();
0128       SiPixelRecHitCollection::DetSet::const_iterator pixeliter = pixelrechitRangeIteratorBegin;
0129 
0130       for (; pixeliter != pixelrechitRangeIteratorEnd; pixeliter++) {
0131         rechit_count++;
0132         // cout << TheID << endl;
0133         SiPixelRecHit::ClusterRef const &clust = pixeliter->cluster();
0134         int sizeX = (*clust).sizeX();
0135         // cout << sizeX << endl;
0136         int sizeY = (*clust).sizeY();
0137         // cout << sizeY << endl;
0138         LocalPoint lp = pixeliter->localPosition();
0139         rechit_x = lp.x();
0140         rechit_y = lp.y();
0141 
0142         LocalError lerr = pixeliter->localPositionError();
0143         float lerr_x = sqrt(lerr.xx());
0144         float lerr_y = sqrt(lerr.yy());
0145         // std::cout << "errors " << lerr_x << " " << lerr_y << std::endl;
0146         (*struct_iter)
0147             .second->fill(rechit_x,
0148                           rechit_y,
0149                           sizeX,
0150                           sizeY,
0151                           lerr_x,
0152                           lerr_y,
0153                           modOn,
0154                           ladOn,
0155                           layOn,
0156                           phiOn,
0157                           bladeOn,
0158                           diskOn,
0159                           ringOn,
0160                           twoDimOn,
0161                           reducedSet);
0162       }
0163     }
0164     if (rechit_count > 0)
0165       (*struct_iter).second->nfill(rechit_count, modOn, ladOn, layOn, phiOn, bladeOn, diskOn, ringOn);
0166   }
0167 
0168   // slow down...
0169   if (slowDown)
0170     usleep(10000);
0171 }
0172 
0173 //------------------------------------------------------------------
0174 // Build data structure
0175 //------------------------------------------------------------------
0176 void SiPixelRecHitSource::buildStructure(const edm::EventSetup &iSetup) {
0177   LogInfo("PixelDQM") << " SiPixelRecHitSource::buildStructure";
0178 
0179   edm::ESHandle<TrackerGeometry> pDD = iSetup.getHandle(trackerGeomTokenBeginRun_);
0180 
0181   edm::ESHandle<TrackerTopology> tTopoHandle = iSetup.getHandle(trackerTopoTokenBeginRun_);
0182   const TrackerTopology *pTT = tTopoHandle.product();
0183 
0184   LogVerbatim("PixelDQM") << " *** Geometry node for TrackerGeom is  " << &(*pDD) << std::endl;
0185   LogVerbatim("PixelDQM") << " *** I have " << pDD->dets().size() << " detectors" << std::endl;
0186   LogVerbatim("PixelDQM") << " *** I have " << pDD->detTypes().size() << " types" << std::endl;
0187 
0188   for (TrackerGeometry::DetContainer::const_iterator it = pDD->dets().begin(); it != pDD->dets().end(); it++) {
0189     if (dynamic_cast<PixelGeomDetUnit const *>((*it)) != nullptr) {
0190       DetId detId = (*it)->geographicalId();
0191 
0192       if ((detId.subdetId() == static_cast<int>(PixelSubdetector::PixelBarrel)) ||
0193           (detId.subdetId() == static_cast<int>(PixelSubdetector::PixelEndcap))) {
0194         uint32_t id = detId();
0195 
0196         if (detId.subdetId() == static_cast<int>(PixelSubdetector::PixelBarrel)) {
0197           if (isPIB)
0198             continue;
0199           LogDebug("PixelDQM") << " ---> Adding Barrel Module " << detId.rawId() << endl;
0200           SiPixelRecHitModule *theModule = new SiPixelRecHitModule(id);
0201           thePixelStructure.insert(pair<uint32_t, SiPixelRecHitModule *>(id, theModule));
0202 
0203         } else if ((detId.subdetId() == static_cast<int>(PixelSubdetector::PixelEndcap))) {
0204           LogDebug("PixelDQM") << " ---> Adding Endcap Module " << detId.rawId() << endl;
0205 
0206           PixelEndcapName::HalfCylinder side = PixelEndcapName(DetId(id), pTT, isUpgrade).halfCylinder();
0207           int disk = PixelEndcapName(DetId(id), pTT, isUpgrade).diskName();
0208           int blade = PixelEndcapName(DetId(id), pTT, isUpgrade).bladeName();
0209           int panel = PixelEndcapName(DetId(id), pTT, isUpgrade).pannelName();
0210           int module = PixelEndcapName(DetId(id), pTT, isUpgrade).plaquetteName();
0211 
0212           char sside[80];
0213           sprintf(sside, "HalfCylinder_%i", side);
0214           char sdisk[80];
0215           sprintf(sdisk, "Disk_%i", disk);
0216           char sblade[80];
0217           sprintf(sblade, "Blade_%02i", blade);
0218           char spanel[80];
0219           sprintf(spanel, "Panel_%i", panel);
0220           char smodule[80];
0221           sprintf(smodule, "Module_%i", module);
0222           std::string side_str = sside;
0223           std::string disk_str = sdisk;
0224           bool mask = side_str.find("HalfCylinder_1") != string::npos ||
0225                       side_str.find("HalfCylinder_2") != string::npos ||
0226                       side_str.find("HalfCylinder_4") != string::npos || disk_str.find("Disk_2") != string::npos;
0227           if (isPIB && mask)
0228             continue;
0229 
0230           SiPixelRecHitModule *theModule = new SiPixelRecHitModule(id);
0231           thePixelStructure.insert(pair<uint32_t, SiPixelRecHitModule *>(id, theModule));
0232         }
0233       }
0234     }
0235   }  // FOR_LOOP
0236 
0237   LogInfo("PixelDQM") << " *** Pixel Structure Size " << thePixelStructure.size() << endl;
0238 }
0239 //------------------------------------------------------------------
0240 // Book MEs
0241 //------------------------------------------------------------------
0242 void SiPixelRecHitSource::bookMEs(DQMStore::IBooker &iBooker, const edm::EventSetup &iSetup) {
0243   std::map<uint32_t, SiPixelRecHitModule *>::iterator struct_iter;
0244 
0245   SiPixelFolderOrganizer theSiPixelFolder(false);
0246 
0247   edm::ESHandle<TrackerTopology> tTopoHandle = iSetup.getHandle(trackerTopoTokenBeginRun_);
0248   const TrackerTopology *pTT = tTopoHandle.product();
0249 
0250   for (struct_iter = thePixelStructure.begin(); struct_iter != thePixelStructure.end(); struct_iter++) {
0251     /// Create folder tree and book histograms
0252     if (modOn) {
0253       if (theSiPixelFolder.setModuleFolder(iBooker, (*struct_iter).first, 0, isUpgrade)) {
0254         (*struct_iter).second->book(conf_, iBooker, pTT, 0, twoDimOn, reducedSet, isUpgrade);
0255       } else {
0256         if (!isPIB)
0257           throw cms::Exception("LogicError") << "[SiPixelDigiSource::bookMEs] Creation of DQM folder failed";
0258       }
0259     }
0260     if (ladOn) {
0261       if (theSiPixelFolder.setModuleFolder(iBooker, (*struct_iter).first, 1, isUpgrade)) {
0262         (*struct_iter).second->book(conf_, iBooker, pTT, 1, twoDimOn, reducedSet, isUpgrade);
0263       } else {
0264         LogDebug("PixelDQM") << "PROBLEM WITH LADDER-FOLDER\n";
0265       }
0266     }
0267     if (layOn) {
0268       if (theSiPixelFolder.setModuleFolder(iBooker, (*struct_iter).first, 2, isUpgrade)) {
0269         (*struct_iter).second->book(conf_, iBooker, pTT, 2, twoDimOn, reducedSet, isUpgrade);
0270       } else {
0271         LogDebug("PixelDQM") << "PROBLEM WITH LAYER-FOLDER\n";
0272       }
0273     }
0274     if (phiOn) {
0275       if (theSiPixelFolder.setModuleFolder(iBooker, (*struct_iter).first, 3, isUpgrade)) {
0276         (*struct_iter).second->book(conf_, iBooker, pTT, 3, twoDimOn, reducedSet, isUpgrade);
0277       } else {
0278         LogDebug("PixelDQM") << "PROBLEM WITH PHI-FOLDER\n";
0279       }
0280     }
0281     if (bladeOn) {
0282       if (theSiPixelFolder.setModuleFolder(iBooker, (*struct_iter).first, 4, isUpgrade)) {
0283         (*struct_iter).second->book(conf_, iBooker, pTT, 4, twoDimOn, reducedSet, isUpgrade);
0284       } else {
0285         LogDebug("PixelDQM") << "PROBLEM WITH BLADE-FOLDER\n";
0286       }
0287     }
0288     if (diskOn) {
0289       if (theSiPixelFolder.setModuleFolder(iBooker, (*struct_iter).first, 5, isUpgrade)) {
0290         (*struct_iter).second->book(conf_, iBooker, pTT, 5, twoDimOn, reducedSet, isUpgrade);
0291       } else {
0292         LogDebug("PixelDQM") << "PROBLEM WITH DISK-FOLDER\n";
0293       }
0294     }
0295     if (ringOn) {
0296       if (theSiPixelFolder.setModuleFolder(iBooker, (*struct_iter).first, 6, isUpgrade)) {
0297         (*struct_iter).second->book(conf_, iBooker, pTT, 6, twoDimOn, reducedSet, isUpgrade);
0298       } else {
0299         LogDebug("PixelDQM") << "PROBLEM WITH RING-FOLDER\n";
0300       }
0301     }
0302   }
0303 }
0304 
0305 // define this as a plug-in
0306 DEFINE_FWK_MODULE(SiPixelRecHitSource);