Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-07-07 22:33:15

0001 #include <iostream>
0002 #include <string>
0003 #include <vector>
0004 
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 #include "FWCore/Framework/interface/Frameworkfwd.h"
0007 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0008 
0009 #include "FWCore/Framework/interface/Event.h"
0010 #include "FWCore/Framework/interface/EventSetup.h"
0011 #include "FWCore/Framework/interface/ESTransientHandle.h"
0012 #include "FWCore/Framework/interface/MakerMacros.h"
0013 
0014 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0015 #include "Geometry/HGCalGeometry/interface/HGCalGeometry.h"
0016 #include "DataFormats/ForwardDetId/interface/ForwardSubdetector.h"
0017 #include "DataFormats/ForwardDetId/interface/HFNoseDetId.h"
0018 #include "CoralBase/Exception.h"
0019 
0020 class HFNoseGeometryTester : public edm::one::EDAnalyzer<> {
0021 public:
0022   explicit HFNoseGeometryTester(const edm::ParameterSet&);
0023   ~HFNoseGeometryTester() override {}
0024 
0025   void beginJob() override {}
0026   void analyze(edm::Event const& iEvent, edm::EventSetup const&) override;
0027   void endJob() override {}
0028 
0029 private:
0030   void doTestWafer(const HGCalGeometry* geom);
0031 
0032   const std::string name_;
0033   const edm::ESGetToken<HGCalGeometry, IdealGeometryRecord> tokGeom_;
0034 };
0035 
0036 HFNoseGeometryTester::HFNoseGeometryTester(const edm::ParameterSet& iC)
0037     : name_(iC.getParameter<std::string>("Detector")),
0038       tokGeom_(esConsumes<HGCalGeometry, IdealGeometryRecord>(edm::ESInputTag{"", name_})) {}
0039 
0040 void HFNoseGeometryTester::analyze(const edm::Event&, const edm::EventSetup& iSetup) {
0041   const HGCalGeometry* geom = &iSetup.getData(tokGeom_);
0042   if (geom->topology().isHFNose()) {
0043     doTestWafer(geom);
0044   } else {
0045     std::cout << name_ << " is not a valid name for HFNose Detecor" << std::endl;
0046   }
0047 }
0048 
0049 void HFNoseGeometryTester::doTestWafer(const HGCalGeometry* geom) {
0050   const std::vector<DetId>& ids = geom->getValidDetIds();
0051   std::cout << "doTestWafer:: " << ids.size() << " valid ids for " << geom->cellElement() << std::endl;
0052   int layers[] = {1, 4, 7};
0053   int zsides[] = {1, -1};
0054   int cells[] = {1, 4, 7};
0055   int wafers[] = {7, 5, 3, -3, -5, -7};
0056   for (int zside : zsides) {
0057     for (int layer : layers) {
0058       for (int waferU : wafers) {
0059         for (int waferV : wafers) {
0060           int type = geom->topology().dddConstants().getTypeHex(layer, waferU, waferV);
0061           std::cout << "zside " << zside << " layer " << layer << " wafer " << waferU << ":" << waferV << " type "
0062                     << type << std::endl;
0063           for (int cellU : cells) {
0064             for (int cellV : cells) {
0065               std::cout << " cell " << cellU << ":" << cellV << std::endl;
0066               DetId id1 = (DetId)(HFNoseDetId(zside, type, layer, waferU, waferV, cellU, cellV));
0067               std::cout << HFNoseDetId(id1) << std::endl;
0068               if (geom->topology().valid(id1)) {
0069                 auto icell1 = geom->getGeometry(id1);
0070                 GlobalPoint global1 = geom->getPosition(id1);
0071                 DetId idc1 = geom->getClosestCell(global1);
0072                 std::cout << "DetId (" << zside << ":" << type << ":" << layer << ":" << waferU << ":" << waferV << ":"
0073                           << cellU << ":" << cellV << ") Geom " << icell1 << " position (" << global1.x() << ", "
0074                           << global1.y() << ", " << global1.z() << ") ids " << std::hex << id1.rawId() << ":"
0075                           << idc1.rawId() << std::dec << ":" << HFNoseDetId(id1) << ":" << HFNoseDetId(idc1)
0076                           << " parameter[3] = " << icell1->param()[2] << ":" << icell1->param()[2];
0077                 if (id1.rawId() != idc1.rawId())
0078                   std::cout << "***** ERROR *****" << std::endl;
0079                 else
0080                   std::cout << std::endl;
0081                 std::vector<GlobalPoint> corners = geom->getCorners(idc1);
0082                 std::cout << corners.size() << " corners";
0083                 for (auto const& cor : corners)
0084                   std::cout << " [" << cor.x() << "," << cor.y() << "," << cor.z() << "]";
0085                 std::cout << std::endl;
0086               }
0087             }
0088           }
0089         }
0090       }
0091     }
0092   }
0093 }
0094 
0095 //define this as a plug-in
0096 DEFINE_FWK_MODULE(HFNoseGeometryTester);