Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-06-27 22:37:19

0001 /** Derived from DTRadiiAnalyzer by Nicola Amapane
0002  *
0003  *  \author M. Maggi - INFN Bari
0004  */
0005 
0006 #include <memory>
0007 #include <fstream>
0008 #include <FWCore/Framework/interface/Frameworkfwd.h>
0009 
0010 #include <FWCore/Framework/interface/one/EDAnalyzer.h>
0011 #include <FWCore/Framework/interface/Event.h>
0012 #include <FWCore/Framework/interface/EventSetup.h>
0013 #include <FWCore/ParameterSet/interface/ParameterSet.h>
0014 
0015 #include "Geometry/RPCGeometry/interface/RPCGeometry.h"
0016 #include "Geometry/RPCGeometry/interface/RPCRollSpecs.h"
0017 #include "Geometry/RPCGeometry/interface/RPCGeomServ.h"
0018 #include <Geometry/Records/interface/MuonGeometryRecord.h>
0019 #include <Geometry/CommonTopologies/interface/RectangularStripTopology.h>
0020 #include <Geometry/CommonTopologies/interface/TrapezoidalStripTopology.h>
0021 
0022 #include <string>
0023 #include <cmath>
0024 #include <vector>
0025 #include <iomanip>
0026 #include <set>
0027 
0028 using namespace std;
0029 
0030 class RPCRadiiAnalyzer : public edm::one::EDAnalyzer<> {
0031 public:
0032   RPCRadiiAnalyzer(const edm::ParameterSet& pset);
0033 
0034   ~RPCRadiiAnalyzer() override;
0035 
0036   void beginJob() override {}
0037   void analyze(edm::Event const& iEvent, edm::EventSetup const&) override;
0038   void endJob() override {}
0039 
0040   const std::string& myName() { return myName_; }
0041 
0042 private:
0043   const edm::ESGetToken<RPCGeometry, MuonGeometryRecord> tokRPC_;
0044   const int dashedLineWidth_;
0045   const std::string dashedLine_;
0046   const std::string myName_;
0047   std::ofstream ofos;
0048 };
0049 
0050 RPCRadiiAnalyzer::RPCRadiiAnalyzer(const edm::ParameterSet& /*iConfig*/)
0051     : tokRPC_{esConsumes<RPCGeometry, MuonGeometryRecord>(edm::ESInputTag{})},
0052       dashedLineWidth_(104),
0053       dashedLine_(std::string(dashedLineWidth_, '-')),
0054       myName_("RPCRadiiAnalyzer") {
0055   ofos.open("MytestOutput.out");
0056   std::cout << "======================== Opening output file" << std::endl;
0057 }
0058 
0059 RPCRadiiAnalyzer::~RPCRadiiAnalyzer() {
0060   ofos.close();
0061   std::cout << "======================== Closing output file" << std::endl;
0062 }
0063 
0064 void RPCRadiiAnalyzer::analyze(const edm::Event& /*iEvent*/, const edm::EventSetup& iSetup) {
0065   const RPCGeometry* pDD = &iSetup.getData(tokRPC_);
0066 
0067   std::cout << myName() << ": Analyzer..." << std::endl;
0068   std::cout << "start " << dashedLine_ << std::endl;
0069 
0070   std::cout << " Geometry node for RPCGeom is  " << &(*pDD) << std::endl;
0071   cout << " I have " << pDD->detTypes().size() << " detTypes" << endl;
0072   cout << " I have " << pDD->detUnits().size() << " detUnits" << endl;
0073   cout << " I have " << pDD->dets().size() << " dets" << endl;
0074   cout << " I have " << pDD->rolls().size() << " rolls" << endl;
0075   cout << " I have " << pDD->chambers().size() << " chambers" << endl;
0076 
0077   std::cout << myName() << ": Begin iteration over geometry..." << std::endl;
0078   std::cout << "iter " << dashedLine_ << std::endl;
0079 
0080   const double dPi = 3.14159265358;
0081   const double radToDeg = 180. / dPi;  //@@ Where to get pi from?
0082 
0083   for (auto it : pDD->dets()) {
0084     //      //----------------------- RPCCHAMBER TEST -------------------------------------------------------
0085 
0086     if (dynamic_cast<const RPCChamber*>(it) != nullptr) {
0087       const RPCChamber* ch = dynamic_cast<const RPCChamber*>(it);
0088 
0089       //RPCDetId detId=ch->id();
0090 
0091       std::vector<const RPCRoll*> rolls = (ch->rolls());
0092       for (auto& roll : rolls) {
0093         if (roll->id().region() == -1 && roll->id().station() > 0)  // &&
0094                                                                     //     (*r)->id().ring() == 2)
0095         {
0096           //        std::cout<<"RPCDetId = "<<(*r)->id()<<std::endl;
0097           RPCGeomServ geosvc(roll->id());
0098           LocalPoint centre(0., 0., 0.);
0099           GlobalPoint gc = roll->toGlobal(centre);
0100           double phic = double(gc.phi()) * radToDeg;
0101           double radii = double(gc.perp());
0102           std::cout << geosvc.name() << " phi=" << phic << " r=" << radii << " detName " << roll->specs()->detName()
0103                     << " s=" << roll->id().sector() << " subs=" << roll->id().subsector() << std::endl;
0104         }
0105       }
0106     }
0107   }
0108   std::cout << dashedLine_ << " end" << std::endl;
0109 }
0110 
0111 //define this as a plug-in
0112 #include <FWCore/Framework/interface/MakerMacros.h>
0113 DEFINE_FWK_MODULE(RPCRadiiAnalyzer);