Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-10-01 01:02:31

0001 // system include files
0002 #include <cmath>
0003 #include <iostream>
0004 #include <fstream>
0005 #include <map>
0006 #include <mutex>
0007 #include <string>
0008 #include <vector>
0009 
0010 #include "DataFormats/ForwardDetId/interface/HFNoseDetId.h"
0011 #include "DataFormats/ForwardDetId/interface/HGCSiliconDetId.h"
0012 #include "DataFormats/ForwardDetId/interface/HGCScintillatorDetId.h"
0013 #include "DataFormats/HGCDigi/interface/HGCDigiCollections.h"
0014 #include "DataFormats/HGCRecHit/interface/HGCRecHitCollections.h"
0015 
0016 #include "FWCore/Framework/interface/Frameworkfwd.h"
0017 #include "FWCore/Framework/interface/stream/EDAnalyzer.h"
0018 #include "FWCore/Framework/interface/Event.h"
0019 #include "FWCore/Framework/interface/EventSetup.h"
0020 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0021 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0022 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0023 #include "FWCore/Utilities/interface/InputTag.h"
0024 #include "FWCore/Utilities/interface/ESGetToken.h"
0025 
0026 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0027 #include "Geometry/HGCalCommonData/interface/HGCalDDDConstants.h"
0028 #include "Geometry/HGCalCommonData/interface/HGCalWaferIndex.h"
0029 
0030 #include "SimDataFormats/CaloHit/interface/PCaloHitContainer.h"
0031 #include "SimDataFormats/CaloHit/interface/PCaloHit.h"
0032 
0033 namespace HGCalValidSimhitCheck {
0034   struct Counters {
0035     Counters() {
0036       badTypes_.clear();
0037       occupancy_.clear();
0038       goodChannels_.clear();
0039     }
0040     CMS_THREAD_GUARD(mtx_) mutable std::map<int, int> badTypes_, occupancy_;
0041     CMS_THREAD_GUARD(mtx_) mutable std::vector<int> goodChannels_;
0042     mutable std::mutex mtx_;
0043   };
0044 }  // namespace HGCalValidSimhitCheck
0045 
0046 class HGCalWaferHitCheck : public edm::stream::EDAnalyzer<edm::GlobalCache<HGCalValidSimhitCheck::Counters> > {
0047 public:
0048   explicit HGCalWaferHitCheck(const edm::ParameterSet&, const HGCalValidSimhitCheck::Counters* count);
0049 
0050   static std::unique_ptr<HGCalValidSimhitCheck::Counters> initializeGlobalCache(edm::ParameterSet const& iConfig) {
0051     return std::make_unique<HGCalValidSimhitCheck::Counters>();
0052   }
0053 
0054   void analyze(edm::Event const&, edm::EventSetup const&) override;
0055   void endStream() override;
0056   static void globalEndJob(const HGCalValidSimhitCheck::Counters* counters);
0057   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0058 
0059 protected:
0060   void beginRun(edm::Run const&, edm::EventSetup const&) override;
0061   void endRun(edm::Run const&, edm::EventSetup const&) override {}
0062 
0063 private:
0064   template <class T>
0065   void analyzeHits(const std::string&, const T&);
0066 
0067   // ----------member data ---------------------------
0068   enum inputType { Sim = 1, Digi = 2, Reco = 3 };
0069   const std::string nameDetector_, caloHitSource_;
0070   const edm::InputTag source_;
0071   const int inpType_;
0072   const int verbosity_;
0073   const bool ifNose_;
0074   const edm::ESGetToken<HGCalDDDConstants, IdealGeometryRecord> geomToken_;
0075   const edm::EDGetTokenT<HGCalDigiCollection> digiSource_;
0076   const edm::EDGetTokenT<HGCRecHitCollection> recHitSource_;
0077   const edm::EDGetTokenT<edm::PCaloHitContainer> tok_hit_;
0078   const HGCalDDDConstants* hgcons_;
0079   std::map<int, int> badTypes_, occupancy_;
0080   std::vector<int> goodChannels_;
0081   static const int waferMax = 12;
0082   static const int layerMax = 28;
0083 };
0084 
0085 HGCalWaferHitCheck::HGCalWaferHitCheck(const edm::ParameterSet& iConfig, const HGCalValidSimhitCheck::Counters*)
0086     : nameDetector_(iConfig.getParameter<std::string>("detectorName")),
0087       caloHitSource_(iConfig.getParameter<std::string>("caloHitSource")),
0088       source_(iConfig.getParameter<edm::InputTag>("source")),
0089       inpType_(iConfig.getParameter<int>("inputType")),
0090       verbosity_(iConfig.getUntrackedParameter<int>("verbosity", 0)),
0091       ifNose_(iConfig.getUntrackedParameter<bool>("ifNose", false)),
0092       geomToken_(esConsumes<HGCalDDDConstants, IdealGeometryRecord, edm::Transition::BeginRun>(
0093           edm::ESInputTag{"", nameDetector_})),
0094       digiSource_(consumes<HGCalDigiCollection>(source_)),
0095       recHitSource_(consumes<HGCRecHitCollection>(source_)),
0096       tok_hit_(consumes<edm::PCaloHitContainer>(edm::InputTag("g4SimHits", caloHitSource_))) {
0097   edm::LogVerbatim("HGCalValidation") << "Validation for input type " << inpType_ << " for " << nameDetector_;
0098 }
0099 
0100 void HGCalWaferHitCheck::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0101   edm::ParameterSetDescription desc;
0102   desc.add<std::string>("detectorName", "HGCalEESensitive");
0103   desc.add<std::string>("caloHitSource", "HGCHitsEE");
0104   desc.add<edm::InputTag>("source", edm::InputTag("simHGCalUnsuppressedDigis", "EE"));
0105   desc.add<int>("inputType", 1);
0106   desc.addUntracked<int>("verbosity", 0);
0107   desc.addUntracked<bool>("ifNose", false);
0108   descriptions.add("hgcalWaferHitCheckEE", desc);
0109 }
0110 
0111 void HGCalWaferHitCheck::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0112   //Get collection depending on inputType
0113   if (inpType_ <= Sim) {
0114     const edm::Handle<edm::PCaloHitContainer>& theCaloHitContainer = iEvent.getHandle(tok_hit_);
0115     if (theCaloHitContainer.isValid()) {
0116       if (verbosity_ > 0)
0117         edm::LogVerbatim("HGCalValidation") << nameDetector_ << " with " << theCaloHitContainer->size() << " SimHits";
0118       analyzeHits(nameDetector_, *(theCaloHitContainer.product()));
0119     } else if (verbosity_ > 0) {
0120       edm::LogVerbatim("HGCalValidation") << "PCaloHitContainer does not exist for " << nameDetector_;
0121     }
0122   } else if (inpType_ == Digi) {
0123     const edm::Handle<HGCalDigiCollection>& theDigiContainer = iEvent.getHandle(digiSource_);
0124     if (theDigiContainer.isValid()) {
0125       if (verbosity_ > 0)
0126         edm::LogVerbatim("HGCalValidation") << nameDetector_ << " with " << theDigiContainer->size() << " Digis";
0127       analyzeHits(nameDetector_, *(theDigiContainer.product()));
0128     } else if (verbosity_ > 0) {
0129       edm::LogVerbatim("HGCalValidation") << "DigiContainer does not exist for " << nameDetector_;
0130     }
0131   } else {
0132     const edm::Handle<HGCRecHitCollection>& theRecHitContainer = iEvent.getHandle(recHitSource_);
0133     if (theRecHitContainer.isValid()) {
0134       if (verbosity_ > 0)
0135         edm::LogVerbatim("HGCalValidation") << nameDetector_ << " with " << theRecHitContainer->size() << " hits";
0136       analyzeHits(nameDetector_, *(theRecHitContainer.product()));
0137     } else if (verbosity_ > 0) {
0138       edm::LogVerbatim("HGCalValidation") << "RecHitContainer does not exist for " << nameDetector_;
0139     }
0140   }
0141 }
0142 
0143 template <class T>
0144 void HGCalWaferHitCheck::analyzeHits(std::string const& name, T const& hits) {
0145   for (auto const& hit : hits) {
0146     uint32_t id = hit.id();
0147     int zside, type, layer, waferU, waferV;
0148     if (ifNose_) {
0149       HFNoseDetId detId = HFNoseDetId(id);
0150       waferU = detId.waferU();
0151       waferV = detId.waferV();
0152       type = detId.type();
0153       layer = detId.layer();
0154       zside = detId.zside();
0155     } else {
0156       HGCSiliconDetId detId = HGCSiliconDetId(id);
0157       waferU = detId.waferU();
0158       waferV = detId.waferV();
0159       type = detId.type();
0160       layer = detId.layer();
0161       zside = detId.zside();
0162     }
0163     int index = HGCalWaferIndex::waferIndex(layer, waferU, waferV, false);
0164     int typef = hgcons_->waferType(layer, waferU, waferV, true);
0165     if (zside < 0)
0166       index *= -1;
0167     auto itr = occupancy_.find(index);
0168     if (itr == occupancy_.end())
0169       occupancy_[index] = 1;
0170     else
0171       ++occupancy_[index];
0172     if (type != typef) {
0173       auto ktr = badTypes_.find(index);
0174       if (ktr == badTypes_.end())
0175         badTypes_[index] = 1;
0176       else
0177         ++badTypes_[index];
0178       if (verbosity_ == 1)
0179         edm::LogVerbatim("HGCalValidation")
0180             << "Detector " << name << " zside = " << zside << " layer = " << layer << " type = " << type << ":" << typef
0181             << " wafer = " << waferU << ":" << waferV << " index " << index;
0182     }
0183     if (verbosity_ > 1)
0184       edm::LogVerbatim("HGCalValidation")
0185           << "Detector " << name << " zside = " << zside << " layer = " << layer << " type = " << type << ":" << typef
0186           << " wafer = " << waferU << ":" << waferV;
0187   }
0188 }
0189 
0190 // ------------ method called when starting to processes a run  ------------
0191 void HGCalWaferHitCheck::beginRun(const edm::Run&, const edm::EventSetup& iSetup) {
0192   hgcons_ = &iSetup.getData(geomToken_);
0193   if (verbosity_ > 0)
0194     edm::LogVerbatim("HGCalValidation") << nameDetector_ << " defined with " << hgcons_->layers(false)
0195                                         << " Layers with " << (hgcons_->firstLayer() - 1) << " in front";
0196 
0197   goodChannels_.clear();
0198   for (int iz = 0; iz < 2; ++iz) {
0199     int zside = iz * 2 - 1;
0200     for (int layer = 1; layer <= layerMax; ++layer) {
0201       for (int waferU = -waferMax; waferU <= waferMax; ++waferU) {
0202         for (int waferV = -waferMax; waferV <= waferMax; ++waferV) {
0203           int index = zside * HGCalWaferIndex::waferIndex(layer, waferU, waferV, false);
0204           if (hgcons_->isValidHex8(layer, waferU, waferV, true))
0205             goodChannels_.emplace_back(index);
0206         }
0207       }
0208     }
0209   }
0210 }
0211 
0212 void HGCalWaferHitCheck::endStream() {
0213   std::scoped_lock lock(globalCache()->mtx_);
0214   for (auto [id, count] : occupancy_) {
0215     if (globalCache()->occupancy_.find(id) == globalCache()->occupancy_.end())
0216       globalCache()->occupancy_[id] = count;
0217     else
0218       globalCache()->occupancy_[id] += count;
0219   }
0220   for (auto [id, count] : badTypes_) {
0221     if (globalCache()->badTypes_.find(id) == globalCache()->badTypes_.end())
0222       globalCache()->badTypes_[id] = count;
0223     else
0224       globalCache()->badTypes_[id] += count;
0225   }
0226 
0227   globalCache()->goodChannels_ = goodChannels_;
0228   globalCache()->mtx_.unlock();
0229 }
0230 
0231 void HGCalWaferHitCheck::globalEndJob(const HGCalValidSimhitCheck::Counters* count) {
0232   int allbad(0), nocc(0);
0233   for (auto const& index : count->goodChannels_) {
0234     int zside = (index < 0) ? -1 : 1;
0235     int layer = HGCalWaferIndex::waferLayer(std::abs(index));
0236     int waferU = HGCalWaferIndex::waferU(std::abs(index));
0237     int waferV = HGCalWaferIndex::waferV(std::abs(index));
0238     int occ = (count->occupancy_.find(index) == count->occupancy_.end()) ? 0 : count->occupancy_[index];
0239     int bad = (count->badTypes_.find(index) == count->badTypes_.end()) ? 0 : count->badTypes_[index];
0240     if (occ == 0)
0241       ++nocc;
0242     if (bad > 0)
0243       ++allbad;
0244     if (occ == 0 || bad > 0) {
0245       edm::LogVerbatim("HGCalValidation") << "ZS:Layer:u:v:index " << zside << ":" << layer << ":" << waferU << ":"
0246                                           << waferV << ":" << index << " Occ " << occ << " bad " << bad;
0247     }
0248   }
0249   edm::LogVerbatim("HGCalValidation") << "\n\n"
0250                                       << allbad << " channels with bad types among " << count->goodChannels_.size()
0251                                       << " channels and " << nocc << " channels with zero occupancy\n\n";
0252 }
0253 
0254 #include "FWCore/Framework/interface/MakerMacros.h"
0255 //define this as a plug-in
0256 DEFINE_FWK_MODULE(HGCalWaferHitCheck);