GEMGeometry/test/GEMGeometryAnalyzer.cc

0001
0002 /** Derived from DTGeometryAnalyzer by Nicola Amapane
0003  *
0004  *  \author M. Maggi - INFN Bari
0005  */
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/ParameterSet/interface/ParameterSet.h"
0012
0013 #include "Geometry/GEMGeometry/interface/GEMGeometry.h"
0014 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0015 #include "Geometry/GEMGeometry/interface/GEMEtaPartitionSpecs.h"
0016 #include "Geometry/CommonTopologies/interface/StripTopology.h"
0017
0018 #include "DataFormats/Math/interface/deltaPhi.h"
0019
0020 #include <memory>
0021 #include <fstream>
0022 #include <string>
0023 #include <cmath>
0024 #include <vector>
0025 #include <iomanip>
0026 #include <set>
0027
0028 class GEMGeometryAnalyzer : public edm::one::EDAnalyzer<> {
0029 public:
0030   GEMGeometryAnalyzer(const edm::ParameterSet& pset);
0031
0032   ~GEMGeometryAnalyzer() override;
0033
0034   void beginJob() override {}
0035   void analyze(edm::Event const& iEvent, edm::EventSetup const&) override;
0036   void endJob() override {}
0037
0038 private:
0039   const std::string& myName() { return myName_; }
0040
0041   const int dashedLineWidth_;
0042   const std::string dashedLine_;
0043   const std::string myName_;
0044   const edm::ESGetToken<GEMGeometry, MuonGeometryRecord> tokGeom_;
0045   std::ofstream ofos;
0046 };
0047
0048 using namespace std;
0049 GEMGeometryAnalyzer::GEMGeometryAnalyzer(const edm::ParameterSet& /*iConfig*/)
0050     : dashedLineWidth_(104),
0051       dashedLine_(std::string(dashedLineWidth_, '-')),
0052       myName_("GEMGeometryAnalyzer"),
0053       tokGeom_{esConsumes<GEMGeometry, MuonGeometryRecord>(edm::ESInputTag{})} {
0054   ofos.open("GEMtestOutput.out");
0055   ofos << "======================== Opening output file" << std::endl;
0056 }
0057
0058 namespace {
0059   bool compareSupChm(const GEMSuperChamber* schm1, const GEMSuperChamber* schm2) {
0060     return (schm1->id().v12Form() < schm2->id().v12Form());
0061   }
0062 }  // namespace
0063
0064 GEMGeometryAnalyzer::~GEMGeometryAnalyzer() {
0065   ofos.close();
0066   ofos << "======================== Closing output file" << std::endl;
0067 }
0068
0069 void GEMGeometryAnalyzer::analyze(const edm::Event& /*iEvent*/, const edm::EventSetup& iSetup) {
0070   const auto& pDD = iSetup.getData(tokGeom_);
0071
0072   ofos << myName() << ": Analyzer..." << std::endl;
0073   ofos << "start " << dashedLine_ << std::endl;
0074
0075   ofos << " Geometry node for GEMGeom is  " << &pDD << endl;
0076   ofos << " detTypes       \t" << pDD.detTypes().size() << endl;
0077   ofos << " GeomDetUnit       \t" << pDD.detUnits().size() << endl;
0078   ofos << " GeomDet           \t" << pDD.dets().size() << endl;
0079   ofos << " GeomDetUnit DetIds\t" << pDD.detUnitIds().size() << endl;
0080   ofos << " eta partitions \t" << pDD.etaPartitions().size() << endl;
0081   ofos << " chambers       \t" << pDD.chambers().size() << endl;
0082   ofos << " super chambers  \t" << pDD.superChambers().size() << endl;
0083   ofos << " rings  \t\t" << pDD.rings().size() << endl;
0084   ofos << " stations  \t\t" << pDD.stations().size() << endl;
0085   ofos << " regions  \t\t" << pDD.regions().size() << endl;
0086
0087   // checking uniqueness of roll detIds
0088   bool flagNonUniqueRollID = false;
0089   bool flagNonUniqueRollRawID = false;
0090   int nstrips = 0;
0091   int npads = 0;
0092   for (auto roll1 : pDD.etaPartitions()) {
0093     nstrips += roll1->nstrips();
0094     npads += roll1->npads();
0095     for (auto roll2 : pDD.etaPartitions()) {
0096       if (roll1 != roll2) {
0097         if (roll1->id() == roll2->id())
0098           flagNonUniqueRollID = true;
0099         if (roll1->id().rawId() == roll2->id().rawId())
0100           flagNonUniqueRollRawID = true;
0101       }
0102     }
0103   }
0104   // checking the number of strips and pads
0105   ofos << " total number of strips\t" << nstrips << endl;
0106   ofos << " total number of pads  \t" << npads << endl;
0107   if (flagNonUniqueRollID or flagNonUniqueRollRawID)
0108     ofos << " -- WARNING: non unique roll Ids!!!" << endl;
0109
0110   // checking uniqueness of chamber detIds
0111   bool flagNonUniqueChID = false;
0112   bool flagNonUniqueChRawID = false;
0113   for (auto ch1 : pDD.chambers()) {
0114     for (auto ch2 : pDD.chambers()) {
0115       if (ch1 != ch2) {
0116         if (ch1->id() == ch2->id())
0117           flagNonUniqueChID = true;
0118         if (ch1->id().rawId() == ch2->id().rawId())
0119           flagNonUniqueChRawID = true;
0120       }
0121     }
0122   }
0123   if (flagNonUniqueChID or flagNonUniqueChRawID)
0124     ofos << " -- WARNING: non unique chamber Ids!!!" << endl;
0125
0126   ofos << myName() << ": Begin iteration over geometry..." << endl;
0127   ofos << "iter " << dashedLine_ << endl;
0128
0129   //----------------------- Global GEMGeometry TEST -------------------------------------------------------
0130   ofos << myName() << "Begin GEMGeometry structure TEST" << endl;
0131
0132   for (auto region : pDD.regions()) {
0133     ofos << "  GEMRegion " << region->region() << " has " << region->nStations() << " stations." << endl;
0134     for (auto station : region->stations()) {
0135       ofos << "    GEMStation " << station->getName() << " has " << station->nRings() << " rings." << endl;
0136       for (auto ring : station->rings()) {
0137         ofos << "      GEMRing " << ring->region() << " " << ring->station() << " " << ring->ring() << " has "
0138              << ring->nSuperChambers() << " super chambers." << endl;
0139         int i = 1;
0140         auto supChmSort = ring->superChambers();
0141         std::sort(supChmSort.begin(), supChmSort.end(), compareSupChm);
0142         for (auto sch : supChmSort) {
0143           GEMDetId schId(sch->id());
0144           ofos << "        GEMSuperChamber " << i << ", GEMDetId = " << schId.rawId() << ", " << schId << " has "
0145                << sch->nChambers() << " chambers." << endl;
0146           // checking the dimensions of each partition & chamber
0147           int j = 1;
0148           for (auto ch : sch->chambers()) {
0149             GEMDetId chId(ch->id());
0150             int nRolls(ch->nEtaPartitions());
0151             ofos << "          GEMChamber " << j << ", GEMDetId = " << chId.rawId() << ", " << chId << " has " << nRolls
0152                  << " eta partitions." << endl;
0153
0154             int k = 1;
0155             auto& rolls(ch->etaPartitions());
0156
0157             /*
0158          * possible checklist for an eta partition:
0159          *   base_bottom, base_top, height, strips, pads
0160          *   cx, cy, cz, ceta, cphi
0161          *   tx, ty, tz, teta, tphi
0162          *   bx, by, bz, beta, bphi
0163          *   pitch center, pitch bottom, pitch top
0164          *   deta, dphi
0165          *   gap thickness
0166          *   sum of all dx + gap = chamber height
0167          */
0168
0169             for (auto roll : rolls) {
0170               GEMDetId rId(roll->id());
0171               ofos << "            GEMEtaPartition " << k << ", GEMDetId = " << rId.rawId() << ", " << rId << endl;
0172
0173               const BoundPlane& bSurface(roll->surface());
0174               const StripTopology* topology(&(roll->specificTopology()));
0175
0176               // base_bottom, base_top, height, strips, pads (all half length)
0177               auto& parameters(roll->specs()->parameters());
0178               float bottomEdge(parameters[0]);
0179               float topEdge(parameters[1]);
0180               float height(parameters[2]);
0181               float nStrips(parameters[3]);
0182               float nPads(parameters[4]);
0183
0184               LocalPoint lCentre(0., 0., 0.);
0185               GlobalPoint gCentre(bSurface.toGlobal(lCentre));
0186
0187               LocalPoint lTop(0., height, 0.);
0188               GlobalPoint gTop(bSurface.toGlobal(lTop));
0189
0190               LocalPoint lBottom(0., -height, 0.);
0191               GlobalPoint gBottom(bSurface.toGlobal(lBottom));
0192
0193               //   gx, gy, gz, geta, gphi (center)
0194               double cx(gCentre.x());
0195               double cy(gCentre.y());
0196               double cz(gCentre.z());
0197               double ceta(gCentre.eta());
0198               int cphi(static_cast<int>(gCentre.phi().degrees()));
0199               if (cphi < 0)
0200                 cphi += 360;
0201
0202               double tx(gTop.x());
0203               double ty(gTop.y());
0204               double tz(gTop.z());
0205               double teta(gTop.eta());
0206               int tphi(static_cast<int>(gTop.phi().degrees()));
0207               if (tphi < 0)
0208                 tphi += 360;
0209
0210               double bx(gBottom.x());
0211               double by(gBottom.y());
0212               double bz(gBottom.z());
0213               double beta(gBottom.eta());
0214               int bphi(static_cast<int>(gBottom.phi().degrees()));
0215               if (bphi < 0)
0216                 bphi += 360;
0217
0218               // pitch bottom, pitch top, pitch centre
0219               float pitch(roll->pitch());
0220               float topPitch(roll->localPitch(lTop));
0221               float bottomPitch(roll->localPitch(lBottom));
0222
0223               // Type - should be GHA[1-nRolls]
0224               string type(roll->type().name());
0225
0226               // print info about edges
0227               LocalPoint lEdge1(topology->localPosition(0.));
0228               LocalPoint lEdgeN(topology->localPosition((float)nStrips));
0229
0230               double cstrip1(roll->toGlobal(lEdge1).phi().degrees());
0231               double cstripN(roll->toGlobal(lEdgeN).phi().degrees());
0232               double dphi(cstripN - cstrip1);
0233               if (dphi < 0.)
0234                 dphi += 360.;
0235               double deta(abs(beta - teta));
0236               const bool printDetails(true);
0237               if (printDetails) {
0238                 ofos << "    \t\tType: " << type << endl
0239                      << "    \t\tDimensions[cm]: b = " << bottomEdge * 2 << ", B = " << topEdge * 2
0240                      << ", h  = " << height * 2 << endl
0241                      << "    \t\tnStrips = " << nStrips << ", nPads =  " << nPads << endl
0242                      << "    \t\tfirst strip pos = (" << roll->toGlobal(lEdge1).x() << ", "
0243                      << roll->toGlobal(lEdge1).y() << ", " << roll->toGlobal(lEdge1).z() << ")" << endl
0244                      << "    \t\tlast  strip pos = (" << roll->toGlobal(lEdgeN).x() << ", "
0245                      << roll->toGlobal(lEdgeN).y() << ", " << roll->toGlobal(lEdgeN).z() << ")" << endl
0246                      << "    \t\ttop(x,y,z)[cm] = (" << tx << ", " << ty << ", " << tz << "), top(eta,phi) = (" << teta
0247                      << ", " << tphi << ")" << endl
0248                      << "    \t\tcenter(x,y,z)[cm] = (" << cx << ", " << cy << ", " << cz << "), center(eta,phi) = ("
0249                      << ceta << ", " << cphi << ")" << endl
0250                      << "    \t\tbottom(x,y,z)[cm] = (" << bx << ", " << by << ", " << bz << "), bottom(eta,phi) = ("
0251                      << beta << ", " << bphi << ")" << endl
0252                      << "    \t\tpitch (top,center,bottom) = (" << topPitch << ", " << pitch << ", " << bottomPitch
0253                      << "), dEta = " << deta << ", dPhi = " << dphi << endl;
0254               }
0255               ++k;
0256             }
0257             ++j;
0258           }
0259           ++i;
0260         }
0261       }
0262     }
0263   }
0264
0265   ofos << dashedLine_ << " end" << std::endl;
0266 }
0267
0268 //define this as a plug-in
0269 #include "FWCore/Framework/interface/MakerMacros.h"
0270 DEFINE_FWK_MODULE(GEMGeometryAnalyzer);