Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-10-16 05:03:57

0001 #include <iostream>
0002 #include <sstream>
0003 #include <string>
0004 #include <vector>
0005 
0006 #include "FWCore/Framework/interface/Frameworkfwd.h"
0007 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/Framework/interface/EventSetup.h"
0010 #include "FWCore/Framework/interface/ESTransientHandle.h"
0011 #include "FWCore/Framework/interface/MakerMacros.h"
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0015 
0016 #include "Geometry/Records/interface/HcalRecNumberingRecord.h"
0017 #include "Geometry/CaloTopology/interface/HcalTopology.h"
0018 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0019 #include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
0020 
0021 class HcalTopologyTester : public edm::one::EDAnalyzer<edm::one::WatchRuns> {
0022 public:
0023   explicit HcalTopologyTester(const edm::ParameterSet&);
0024 
0025   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0026 
0027 private:
0028   void analyze(edm::Event const&, edm::EventSetup const&) override;
0029   void beginJob() override {}
0030   void beginRun(edm::Run const&, edm::EventSetup const&) override {}
0031   void endRun(edm::Run const&, edm::EventSetup const&) override {}
0032   void doTest(const HcalTopology& topology);
0033 
0034   // ----------member data ---------------------------
0035   const edm::ESGetToken<HcalTopology, HcalRecNumberingRecord> tokTopo_;
0036 };
0037 
0038 HcalTopologyTester::HcalTopologyTester(const edm::ParameterSet&)
0039     : tokTopo_{esConsumes<HcalTopology, HcalRecNumberingRecord>(edm::ESInputTag{})} {}
0040 
0041 void HcalTopologyTester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0042   edm::ParameterSetDescription desc;
0043   desc.setUnknown();
0044   descriptions.add("hcalTopologyTester", desc);
0045 }
0046 
0047 void HcalTopologyTester::analyze(edm::Event const&, edm::EventSetup const& iSetup) { doTest(iSetup.getData(tokTopo_)); }
0048 
0049 void HcalTopologyTester::doTest(const HcalTopology& topology) {
0050   // Total number of valid cells
0051   edm::LogVerbatim("HCalGeom") << "Total number of cells in HB:" << topology.getHBSize()
0052                                << " HE: " << topology.getHESize() << " HF: " << topology.getHFSize()
0053                                << " HO: " << topology.getHOSize() << " HT: " << topology.getHTSize()
0054                                << " Calib: " << topology.getCALIBSize() << " Overall: " << topology.ncells();
0055   std::vector<std::string> dets = {"HB", "HE", "HO", "HF"};
0056   for (int det = 1; det <= 4; ++det)
0057     edm::LogVerbatim("HCalGeom") << "Valid cells for " << dets[det - 1] << " = " << topology.ncells(det);
0058 
0059   // First test on movements along eta/phi directions
0060   edm::LogVerbatim("HCalGeom") << "\nTest on movements along eta/phi directions"
0061                                << "\n==========================================";
0062   for (int idet = 0; idet < 4; idet++) {
0063     HcalSubdetector subdet = HcalBarrel;
0064     if (idet == 1)
0065       subdet = HcalOuter;
0066     else if (idet == 2)
0067       subdet = HcalEndcap;
0068     else if (idet == 3)
0069       subdet = HcalForward;
0070     for (int depth = 1; depth < 4; ++depth) {
0071       for (int ieta = -41; ieta <= 41; ieta++) {
0072         for (int iphi = 1; iphi <= 72; iphi++) {
0073           const HcalDetId id(subdet, ieta, iphi, depth);
0074           if (topology.valid(id)) {
0075             std::vector<DetId> idE = topology.east(id);
0076             std::vector<DetId> idW = topology.west(id);
0077             std::vector<DetId> idN = topology.north(id);
0078             std::vector<DetId> idS = topology.south(id);
0079             std::vector<DetId> idU = topology.up(id);
0080             edm::LogVerbatim("HCalGeom") << "Neighbours for : Tower " << id;
0081             std::ostringstream st1;
0082             st1 << "          " << idE.size() << " sets along East:";
0083             for (auto& i : idE)
0084               st1 << " " << (HcalDetId)(i());
0085             edm::LogVerbatim("HCalGeom") << st1.str();
0086             std::ostringstream st2;
0087             st2 << "          " << idW.size() << " sets along West:";
0088             for (auto& i : idW)
0089               st2 << " " << (HcalDetId)(i());
0090             edm::LogVerbatim("HCalGeom") << st2.str();
0091             std::ostringstream st3;
0092             st3 << "          " << idN.size() << " sets along North:";
0093             for (auto& i : idN)
0094               st3 << " " << (HcalDetId)(i());
0095             edm::LogVerbatim("HCalGeom") << st3.str();
0096             std::ostringstream st4;
0097             st4 << "          " << idS.size() << " sets along South:";
0098             for (auto& i : idS)
0099               st4 << " " << (HcalDetId)(i());
0100             edm::LogVerbatim("HCalGeom") << st4.str();
0101             std::ostringstream st5;
0102             st5 << "          " << idU.size() << " sets up in depth:";
0103             for (auto& i : idU)
0104               st5 << " " << (HcalDetId)(i());
0105             edm::LogVerbatim("HCalGeom") << st5.str();
0106           }
0107         }
0108       }
0109     }
0110   }
0111 
0112   // Check on Dense Index
0113   edm::LogVerbatim("HCalGeom") << "\nCheck on Dense Index"
0114                                << "\n=====================";
0115   int maxDepthHB = topology.maxDepthHB();
0116   int maxDepthHE = topology.maxDepthHE();
0117   for (int det = 1; det <= HcalForward; det++) {
0118     for (int eta = -HcalDetId::kHcalEtaMask2; eta <= (int)(HcalDetId::kHcalEtaMask2); eta++) {
0119       for (unsigned int phi = 0; phi <= HcalDetId::kHcalPhiMask2; phi++) {
0120         for (int depth = 1; depth < maxDepthHB + maxDepthHE; depth++) {
0121           HcalDetId cell((HcalSubdetector)det, eta, phi, depth);
0122           if (topology.valid(cell)) {
0123             unsigned int dense = topology.detId2denseId(DetId(cell));
0124             DetId id = topology.denseId2detId(dense);
0125             std::string cherr = (cell != HcalDetId(id)) ? " **** ERROR *****" : "";
0126             edm::LogVerbatim("HCalGeom") << cell << " Dense " << std::hex << dense << std::dec << " o/p "
0127                                          << HcalDetId(id) << cherr;
0128           }
0129         }
0130       }
0131     }
0132   }
0133 
0134   // Check list of depths
0135   edm::LogVerbatim("HCalGeom") << "\nCheck list of Depths"
0136                                << "\n====================";
0137   for (int eta = topology.lastHERing() - 2; eta <= topology.lastHERing(); ++eta) {
0138     for (unsigned int phi = 0; phi <= HcalDetId::kHcalPhiMask2; phi++) {
0139       for (int depth = 1; depth <= maxDepthHE; depth++) {
0140         HcalDetId cell(HcalEndcap, eta, phi, depth);
0141         if (topology.valid(cell)) {
0142           std::vector<int> depths = topology.mergedDepthList29(cell);
0143           std::ostringstream st1;
0144           st1 << cell << " is with merge depth flag " << topology.mergedDepth29(cell) << " having " << depths.size()
0145               << " merged depths:";
0146           for (unsigned int k = 0; k < depths.size(); ++k)
0147             st1 << " [" << k << "]:" << depths[k];
0148           edm::LogVerbatim("HCalGeom") << st1.str();
0149         }
0150       }
0151     }
0152   }
0153 }
0154 
0155 //define this as a plug-in
0156 DEFINE_FWK_MODULE(HcalTopologyTester);