test/stubs/CSCDetIdAnalyzer.cc

0001 // Test CSCDetId & CSCIndexer 13.11.2007 ptc
0002
0003 #include <FWCore/Framework/interface/one/EDAnalyzer.h>
0004 #include <FWCore/Framework/interface/ESHandle.h>
0005 #include <FWCore/Framework/interface/EventSetup.h>
0006 #include <FWCore/Framework/interface/MakerMacros.h>
0007
0008 #include <Geometry/Records/interface/MuonGeometryRecord.h>
0009 #include <Geometry/CSCGeometry/interface/CSCGeometry.h>
0010 #include <Geometry/CSCGeometry/interface/CSCLayer.h>
0011 #include <DataFormats/MuonDetId/interface/CSCIndexer.h>
0012 #include <DataFormats/GeometryVector/interface/Pi.h>
0013
0014 #include <cmath>
0015 #include <iomanip>  // for setw() etc.
0016 #include <string>
0017 #include <vector>
0018
0019 class CSCDetIdAnalyzer : public edm::one::EDAnalyzer<> {
0020 public:
0021   explicit CSCDetIdAnalyzer(const edm::ParameterSet&);
0022   ~CSCDetIdAnalyzer() override = default;
0023
0024   virtual void analyze(const edm::Event&, const edm::EventSetup&);
0025
0026   const std::string& myName() { return myName_; }
0027
0028 private:
0029   const int dashedLineWidth_;
0030   const std::string dashedLine_;
0031   const std::string myName_;
0032   const edm::ESGetToken<CSCGeometry, MuonGeometryRecord> tokGeom_;
0033 };
0034
0035 CSCDetIdAnalyzer::CSCDetIdAnalyzer(const edm::ParameterSet& iConfig)
0036     : dashedLineWidth_(140),
0037       dashedLine_(std::string(dashedLineWidth_, '-')),
0038       myName_("CSCDetIdAnalyzer"),
0039       tokGeom_(esConsumes()) {
0040   std::cout << dashedLine_ << std::endl;
0041   std::cout << "Welcome to " << myName_ << std::endl;
0042   std::cout << dashedLine_ << std::endl;
0043   std::cout << "I will build the CSC geometry, then iterate over all layers." << std::endl;
0044   std::cout << "From each CSCDetId I will build the associated linear index, skipping ME1a layers." << std::endl;
0045   std::cout << "I will build this index once from the layer labels and once from the CSCDetId, and check they agree."
0046             << std::endl;
0047   std::cout << "I will output one line per layer, listing the CSCDetId, the labels, and these two indexes."
0048             << std::endl;
0049   std::cout << "I will append the strip-channel indexes for the two edge strips and the central strip." << std::endl;
0050   std::cout << "Finally, I will rebuild a CSCDetId from the layer index, and check it is the same as the original,"
0051             << std::endl;
0052   std::cout << "and rebuild a CSCDetId from the final strip-channel index, and check it is the same as the original."
0053             << std::endl;
0054   std::cout << "If any of these tests fail, you will see assert failure messages in the output." << std::endl;
0055   std::cout << "If there are no such failures then the tests passed." << std::endl;
0056   std::cout << dashedLine_ << std::endl;
0057 }
0058
0059 void CSCDetIdAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0060   const double dPi = Geom::pi();
0061   const double radToDeg = 180. / dPi;
0062
0063   std::cout << myName() << ": Analyzer..." << std::endl;
0064   std::cout << "start " << dashedLine_ << std::endl;
0065   std::cout << "pi = " << dPi << ", radToDeg = " << radToDeg << std::endl;
0066
0067   const edm::ESHandle<CSCGeometry>& pDD = iSetup.getHandle(tokGeom_);
0068   std::cout << " Geometry node for CSCGeom is  " << &(*pDD) << std::endl;
0069   std::cout << " I have " << pDD->detTypes().size() << " detTypes" << std::endl;
0070   std::cout << " I have " << pDD->detUnits().size() << " detUnits" << std::endl;
0071   std::cout << " I have " << pDD->dets().size() << " dets" << std::endl;
0072   std::cout << " I have " << pDD->layers().size() << " layers" << std::endl;
0073   std::cout << " I have " << pDD->chambers().size() << " chambers" << std::endl;
0074
0075   std::cout << myName() << ": Begin iteration over geometry..." << std::endl;
0076   std::cout << "iter " << dashedLine_ << std::endl;
0077
0078   std::cout << "\n  #     id(dec)      id(oct)                        "
0079                "lindex     lindex2      cindex       label       strip  sindex   strip  sindex   strip  sindex"
0080             << std::endl;
0081
0082   // Construct an indexer object
0083   CSCIndexer* theIndexer = new CSCIndexer;
0084
0085   int icount = 0;
0086   int icountAll = 0;
0087
0088   // Iterate over the DetUnits in the CSCGeometry
0089   for (const auto& it : pDD->detUnits()) {
0090     // Check each DetUnit really is a CSC layer
0091     auto layer = dynamic_cast<CSCLayer const*>(it);
0092
0093     if (layer) {
0094       ++icountAll;  // how many layers we see
0095
0096       // What's its DetId?
0097
0098       DetId detId = layer->geographicalId();
0099       int id = detId();  // or detId.rawId()
0100       CSCDetId cscDetId = layer->id();
0101
0102       // There's going to be a lot of messing with field width (and precision) so
0103       // save input values...
0104       int iw = std::cout.width();      // save current width
0105       int ip = std::cout.precision();  // save current precision
0106
0107       short ie = CSCDetId::endcap(id);
0108       short is = CSCDetId::station(id);
0109       short ir = CSCDetId::ring(id);
0110       short ic = CSCDetId::chamber(id);
0111       short il = CSCDetId::layer(id);
0112
0113       if (ir == 4)
0114         continue;  // DOES NOT HANDLE ME1a SEPARATELY FROM ME11
0115       ++icount;    // count ONLY non-ME1a layers!
0116
0117       std::cout << std::setw(4) << icount << std::setw(12) << id << std::oct << std::setw(12) << id << std::dec
0118                 << std::setw(iw) << "   E" << ie << " S" << is << " R" << ir << " C" << std::setw(2) << ic
0119                 << std::setw(iw) << " L" << il;
0120
0121       unsigned lind = theIndexer->layerIndex(ie, is, ir, ic, il);
0122       unsigned cind = theIndexer->startChamberIndexInEndcap(ie, is, ir) + ic - 1;
0123       unsigned lind2 = theIndexer->layerIndex(cscDetId);
0124
0125       //       std::cout << std::setw(12) << std::setw(12) << lind << std::setw(12) << lind2 << "     " << std::endl;
0126       std::cout << std::setw(12) << lind << std::setw(12) << lind2 << std::setw(12) << cind << std::setw(12)
0127                 << theIndexer->chamberLabelFromChamberIndex(cind) << "     ";
0128
0129       // Index a few strips
0130       unsigned short nstrips = theIndexer->stripChannelsPerLayer(is, ir);
0131       unsigned int sc1 = theIndexer->stripChannelIndex(ie, is, ir, ic, il, 1);
0132       unsigned int scm = theIndexer->stripChannelIndex(ie, is, ir, ic, il, nstrips / 2);
0133       unsigned int scn = theIndexer->stripChannelIndex(ie, is, ir, ic, il, nstrips);
0134
0135       std::cout << "      1  " << std::setw(6) << sc1 << "      " << nstrips / 2 << "  " << std::setw(6) << scm
0136                 << "      " << nstrips << "  " << std::setw(6) << scn << std::endl;
0137
0138       // Reset the values we changed
0139       std::cout << std::setprecision(ip) << std::setw(iw);
0140
0141       // ASSERTS
0142       // =======
0143
0144       // Check layer indices are consistent
0145       //       std::cout << "lind2 = " << lind2 << ", lind=" << lind << std::endl;
0146       assert(lind2 == lind);
0147
0148       // Build CSCDetId from this index and check it's same as original
0149       CSCDetId cscDetId2 = theIndexer->detIdFromLayerIndex(lind);
0150       // std::cout << "cscDetId2 = E" << cscDetId2.endcap() << " S" << cscDetId2.station() << " R" << cscDetId2.ring() << " C" << cscDetId2.chamber() << " L" << cscDetId2.layer() << std::endl;
0151       assert(cscDetId2 == cscDetId);
0152
0153       // Build CSCDetId from the strip-channel index for strip 'nstrips' and check it matches
0154       std::pair<CSCDetId, unsigned short int> p = theIndexer->detIdFromStripChannelIndex(scn);
0155       CSCDetId cscDetId3 = p.first;
0156       unsigned short iscn = p.second;
0157       // std::cout << "scn=" << scn << "  iscn=" << iscn << std::endl;
0158       // std::cout << "cscDetId3 = E" << cscDetId3.endcap() << " S" << cscDetId3.station() << " R" << cscDetId3.ring() << " C" << cscDetId3.chamber() << " L" << cscDetId3.layer() << std::endl;
0159       assert(iscn == nstrips);
0160       assert(cscDetId3 == cscDetId);
0161
0162       // Check idToDetUnit
0163       const GeomDetUnit* gdu = pDD->idToDetUnit(detId);
0164       assert(gdu == layer);
0165       // Check idToDet
0166       const GeomDet* gd = pDD->idToDet(detId);
0167       assert(gd == layer);
0168     } else {
0169       std::cout << "Something wrong ... could not dynamic_cast Det* to CSCLayer* " << std::endl;
0170     }
0171   }
0172
0173   delete theIndexer;
0174
0175   std::cout << dashedLine_ << " end" << std::endl;
0176 }
0177
0178 //define this as a plug-in
0179 DEFINE_FWK_MODULE(CSCDetIdAnalyzer);