Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-09-07 04:36:33

0001 #include <iostream>
0002 #include <sstream>
0003 #include <map>
0004 #include <string>
0005 #include <vector>
0006 
0007 #include "FWCore/Framework/interface/Frameworkfwd.h"
0008 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0009 #include "FWCore/Framework/interface/Event.h"
0010 #include "FWCore/Framework/interface/EventSetup.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/IdealGeometryRecord.h"
0017 #include "Geometry/HGCalGeometry/interface/HGCalGeometry.h"
0018 #include "DataFormats/ForwardDetId/interface/ForwardSubdetector.h"
0019 #include "DataFormats/ForwardDetId/interface/HGCalDetId.h"
0020 #include "DataFormats/ForwardDetId/interface/HGCScintillatorDetId.h"
0021 #include "DataFormats/ForwardDetId/interface/HGCSiliconDetId.h"
0022 
0023 class HGCalValidScintTest : public edm::one::EDAnalyzer<> {
0024 public:
0025   explicit HGCalValidScintTest(const edm::ParameterSet&);
0026   ~HGCalValidScintTest() override = default;
0027 
0028   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0029   void analyze(edm::Event const& iEvent, edm::EventSetup const&) override;
0030 
0031 private:
0032   const edm::ESGetToken<HGCalGeometry, IdealGeometryRecord> geomToken_;
0033   struct layerInfo {
0034     int ringMin, ringMax;
0035     double rMin, rMax;
0036     layerInfo(int minR = 100, double rMn = 0, int maxR = 0, double rMx = 0)
0037         : ringMin(minR), ringMax(maxR), rMin(rMn), rMax(rMx) {}
0038   };
0039 };
0040 
0041 HGCalValidScintTest::HGCalValidScintTest(const edm::ParameterSet& iC)
0042     : geomToken_(esConsumes<HGCalGeometry, IdealGeometryRecord>(edm::ESInputTag{"", "HGCalHEScintillatorSensitive"})) {}
0043 
0044 void HGCalValidScintTest::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0045   edm::ParameterSetDescription desc;
0046   descriptions.add("hgcalValidScintTest", desc);
0047 }
0048 
0049 void HGCalValidScintTest::analyze(const edm::Event&, const edm::EventSetup& iSetup) {
0050   const auto& geom = &iSetup.getData(geomToken_);
0051   DetId::Detector det = DetId::HGCalHSc;
0052   edm::LogVerbatim("HGCalGeom") << "Perform test for HGCalHEScintillatorSensitive Detector " << det << " Mode "
0053                                 << geom->topology().dddConstants().geomMode();
0054 
0055   int firstLayer = geom->topology().dddConstants().firstLayer();
0056   int lastLayer = geom->topology().dddConstants().lastLayer(true);
0057   const std::vector<DetId>& ids = geom->getValidDetIds();
0058   edm::LogVerbatim("HGCalGeom") << "doTest: " << ids.size() << " valid ids for " << geom->cellElement();
0059 
0060   int zside(1);
0061   std::map<int, layerInfo> layerMap;
0062   for (int layer = firstLayer; layer <= lastLayer; ++layer) {
0063     std::vector<std::pair<int, int> > done;
0064     for (auto const& id : ids) {
0065       HGCScintillatorDetId hid(id);
0066       if ((hid.zside() != zside) || (hid.layer() != layer))
0067         continue;
0068       std::pair<int, int> ring = std::make_pair(hid.ring(), 0);
0069       if (std::find(done.begin(), done.end(), ring) != done.end())
0070         continue;
0071       done.emplace_back(ring);
0072       edm::LogVerbatim("HGCalGeom") << "Corners for " << hid;
0073 
0074       const auto cor = geom->getNewCorners(id);
0075       std::ostringstream st1;
0076       for (unsigned int k = 0; k < cor.size(); ++k)
0077         st1 << " [" << k << "] " << std::setprecision(4) << cor[k];
0078       edm::LogVerbatim("HGCalGeom") << st1.str();
0079 
0080       double r = cor[0].perp();
0081       auto itr = layerMap.find(layer);
0082       if (itr == layerMap.end()) {
0083         layerInfo info;
0084         info.ringMin = info.ringMax = ring.first;
0085         info.rMin = info.rMax = r;
0086         layerMap[layer] = info;
0087       } else {
0088         layerInfo info = itr->second;
0089         if (info.ringMin > ring.first) {
0090           info.ringMin = ring.first;
0091           info.rMin = r;
0092         } else if (info.ringMax < ring.first) {
0093           info.ringMax = ring.first;
0094           info.rMax = r;
0095         }
0096         layerMap[layer] = info;
0097       }
0098     }
0099   }
0100   edm::LogVerbatim("HGCalGeom") << "\n\nSummary of " << layerMap.size() << " Scintillator layers";
0101   for (auto itr = layerMap.begin(); itr != layerMap.end(); ++itr)
0102     edm::LogVerbatim("HGCalGeom") << "Layer " << itr->first << " lowest Ring " << (itr->second).ringMin << ":"
0103                                   << (itr->second).rMin << " largest Ring " << (itr->second).ringMax << ":"
0104                                   << (itr->second).rMax;
0105 }
0106 
0107 //define this as a plug-in
0108 DEFINE_FWK_MODULE(HGCalValidScintTest);