Back to home page

Project CMSSW displayed by LXR

 
 

    


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

0001 // -*- C++ -*-
0002 //
0003 // Package:    RPCCSC
0004 // Class:      RPCCSC
0005 //
0006 /**\class RPCCSC RPCCSC.cc TESTCSCRPC/RPCCSC/src/RPCCSC.cc
0007 
0008 Description: <one line class summary>
0009 
0010 Implementation:
0011 <Notes on implementation>
0012  */
0013 //
0014 // Original Author:  Haiyun Teng
0015 //         Created:  Wed Feb 25 18:09:15 CET 2009
0016 //
0017 //
0018 
0019 // system include files
0020 #include <memory>
0021 
0022 // user include files
0023 #include "FWCore/Framework/interface/Frameworkfwd.h"
0024 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0025 
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/MakerMacros.h"
0028 
0029 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0030 
0031 #include <memory>
0032 #include "FWCore/Framework/interface/MakerMacros.h"
0033 #include <DataFormats/MuonDetId/interface/RPCDetId.h>
0034 #include <Geometry/RPCGeometry/interface/RPCGeomServ.h>
0035 #include <Geometry/CommonDetUnit/interface/GeomDet.h>
0036 #include <Geometry/Records/interface/MuonGeometryRecord.h>
0037 #include <Geometry/CommonTopologies/interface/RectangularStripTopology.h>
0038 #include <Geometry/CommonTopologies/interface/TrapezoidalStripTopology.h>
0039 #include <Geometry/CSCGeometry/interface/CSCGeometry.h>
0040 #include <Geometry/RPCGeometry/interface/RPCGeometry.h>
0041 #include <Geometry/CSCGeometry/interface/CSCChamber.h>
0042 #include <DataFormats/MuonDetId/interface/CSCDetId.h>
0043 
0044 //
0045 // class decleration
0046 //
0047 class CSCStationIndex {
0048 public:
0049   CSCStationIndex() : _region(0), _station(0), _ring(0), _chamber(0) {}
0050   CSCStationIndex(int region, int station, int ring, int chamber)
0051       : _region(region), _station(station), _ring(ring), _chamber(chamber) {}
0052   ~CSCStationIndex() {}
0053   int region() const { return _region; }
0054   int station() const { return _station; }
0055   int ring() const { return _ring; }
0056   int chamber() const { return _chamber; }
0057   bool operator<(const CSCStationIndex& cscind) const {
0058     if (cscind.region() != this->region())
0059       return cscind.region() < this->region();
0060     else if (cscind.station() != this->station())
0061       return cscind.station() < this->station();
0062     else if (cscind.ring() != this->ring())
0063       return cscind.ring() < this->ring();
0064     else if (cscind.chamber() != this->chamber())
0065       return cscind.chamber() < this->chamber();
0066     return false;
0067   }
0068 
0069 private:
0070   int _region;
0071   int _station;
0072   int _ring;
0073   int _chamber;
0074 };
0075 
0076 class RPCCSC : public edm::one::EDAnalyzer<> {
0077 public:
0078   explicit RPCCSC(const edm::ParameterSet&);
0079   ~RPCCSC() override;
0080 
0081   void beginJob() override {}
0082   void analyze(edm::Event const& iEvent, edm::EventSetup const&) override;
0083   void endJob() override {}
0084 
0085 private:
0086   const edm::ESGetToken<RPCGeometry, MuonGeometryRecord> tokRPC_;
0087   const edm::ESGetToken<CSCGeometry, MuonGeometryRecord> tokCSC_;
0088   std::map<CSCStationIndex, std::set<RPCDetId> > rollstoreCSC;
0089 };
0090 
0091 RPCCSC::RPCCSC(const edm::ParameterSet& /*iConfig*/)
0092     : tokRPC_{esConsumes<RPCGeometry, MuonGeometryRecord>(edm::ESInputTag{})},
0093       tokCSC_{esConsumes<CSCGeometry, MuonGeometryRecord>(edm::ESInputTag{})} {}
0094 
0095 RPCCSC::~RPCCSC() {}
0096 
0097 // ------------ method called once each job just before starting event loop  ------------
0098 void RPCCSC::analyze(const edm::Event& /*iEvent*/, const edm::EventSetup& iSetup) {
0099   using namespace std;
0100   const RPCGeometry* rpcGeometry = &iSetup.getData(tokRPC_);
0101   const CSCGeometry* cscGeometry = &iSetup.getData(tokCSC_);
0102 
0103   for (TrackingGeometry::DetContainer::const_iterator it = rpcGeometry->dets().begin(); it < rpcGeometry->dets().end();
0104        it++) {
0105     if (dynamic_cast<const RPCChamber*>(*it) != nullptr) {
0106       const RPCChamber* ch = dynamic_cast<const RPCChamber*>(*it);
0107       std::vector<const RPCRoll*> roles = (ch->rolls());
0108       for (std::vector<const RPCRoll*>::const_iterator r = roles.begin(); r != roles.end(); ++r) {
0109         RPCDetId rpcId = (*r)->id();
0110         int region = rpcId.region();
0111         //booking all histograms
0112         RPCGeomServ rpcsrv(rpcId);
0113         std::string nameRoll = rpcsrv.name();
0114         //std::cout<<"Booking for "<<nameRoll<<std::endl;
0115 
0116         if (region != 0) {
0117           //      const TrapezoidalStripTopology* topE_=dynamic_cast<const TrapezoidalStripTopology*>(&((*r)->topology()));
0118           //      float stripl = topE_->stripLength();
0119           //      float stripw = topE_->pitch();
0120           int region = rpcId.region();
0121           int station = rpcId.station();
0122           int ring = rpcId.ring();
0123           int cscring = ring;
0124           int cscstation = station;
0125           RPCGeomServ rpcsrv(rpcId);
0126           int rpcsegment = rpcsrv.segment();  //This replace rpcsrv.segment();
0127           //std::cout<<"My segment="<<mySegment(rpcId)<<" GeomServ="<<rpcsrv.segment()<<std::endl;
0128           int cscchamber = rpcsegment;                        //FIX THIS ACCORDING TO RPCGeomServ::segment()Definition
0129           if ((station == 2 || station == 3) && ring == 3) {  //Adding Ring 3 of RPC to the CSC Ring 2
0130             cscring = 2;
0131           }
0132 
0133           CSCStationIndex ind(region, cscstation, cscring, cscchamber);
0134           std::set<RPCDetId> myrolls;
0135           if (rollstoreCSC.find(ind) != rollstoreCSC.end()) {
0136             myrolls = rollstoreCSC[ind];
0137           }
0138           myrolls.insert(rpcId);
0139           rollstoreCSC[ind] = myrolls;
0140         }
0141       }
0142     }
0143   }
0144   for (TrackingGeometry::DetContainer::const_iterator it = rpcGeometry->dets().begin(); it < rpcGeometry->dets().end();
0145        it++) {
0146     if (dynamic_cast<const RPCChamber*>(*it) != nullptr) {
0147       const RPCChamber* ch = dynamic_cast<const RPCChamber*>(*it);
0148       std::vector<const RPCRoll*> roles = (ch->rolls());
0149       for (std::vector<const RPCRoll*>::const_iterator r = roles.begin(); r != roles.end(); ++r) {
0150         RPCDetId rpcId = (*r)->id();
0151         int region = rpcId.region();
0152         if (region != 0 && (rpcId.ring() == 2 || rpcId.ring() == 3)) {
0153           int region = rpcId.region();
0154           int station = rpcId.station();
0155           int ring = rpcId.ring();
0156           int cscring = ring;
0157 
0158           if ((station == 2 || station == 3) && ring == 3)
0159             cscring = 2;  //CSC Ring 2 covers rpc ring 2 & 3
0160 
0161           int cscstation = station;
0162           RPCGeomServ rpcsrv(rpcId);
0163           int rpcsegment = rpcsrv.segment();
0164 
0165           int cscchamber = rpcsegment + 1;
0166           if (cscchamber == 37)
0167             cscchamber = 1;
0168           CSCStationIndex ind(region, cscstation, cscring, cscchamber);
0169           std::set<RPCDetId> myrolls;
0170           if (rollstoreCSC.find(ind) != rollstoreCSC.end())
0171             myrolls = rollstoreCSC[ind];
0172           myrolls.insert(rpcId);
0173           rollstoreCSC[ind] = myrolls;
0174 
0175           cscchamber = rpcsegment - 1;
0176           if (cscchamber == 0)
0177             cscchamber = 36;
0178           CSCStationIndex indDos(region, cscstation, cscring, cscchamber);
0179           std::set<RPCDetId> myrollsDos;
0180           if (rollstoreCSC.find(indDos) != rollstoreCSC.end())
0181             myrollsDos = rollstoreCSC[indDos];
0182           myrollsDos.insert(rpcId);
0183           rollstoreCSC[indDos] = myrollsDos;
0184         }
0185       }
0186     }
0187   }
0188 
0189   //adding more rpcs
0190 
0191   // Now check binding
0192   const CSCGeometry::ChamberContainer& cscChambers = cscGeometry->chambers();
0193   for (auto cscChamber : cscChambers) {
0194     CSCDetId CSCId = cscChamber->id();
0195 
0196     int cscEndCap = CSCId.endcap();
0197     int cscStation = CSCId.station();
0198     int cscRing = CSCId.ring();
0199     //     int cscChamber = CSCId.chamber();
0200     int rpcRegion = 1;
0201     if (cscEndCap == 2)
0202       rpcRegion = -1;  //Relacion entre las endcaps
0203     int rpcRing = cscRing;
0204     if (cscRing == 4)
0205       rpcRing = 1;
0206     int rpcStation = cscStation;
0207     int rpcSegment = CSCId.chamber();
0208 
0209     //std::cout<<"CSC \t \t Getting chamber from Geometry"<<std::endl;
0210     const CSCChamber* TheChamber = cscGeometry->chamber(CSCId);
0211     //std::cout<<"CSC \t \t Getting ID from Chamber"<<std::endl;
0212 
0213     std::set<RPCDetId> rollsForThisCSC = rollstoreCSC[CSCStationIndex(rpcRegion, rpcStation, rpcRing, rpcSegment)];
0214     if (CSCId.ring() != 1)
0215       std::cout << "CSC for" << CSCId << " " << rollsForThisCSC.size() << " rolls." << std::endl;
0216 
0217     for (auto iteraRoll : rollsForThisCSC) {
0218       const RPCRoll* rollasociated = rpcGeometry->roll(iteraRoll);
0219       RPCDetId rpcId = rollasociated->id();
0220       RPCGeomServ rpcsrv(rpcId);
0221 
0222       const BoundPlane& RPCSurface = rollasociated->surface();
0223 
0224       GlobalPoint CenterPointRollGlobal = RPCSurface.toGlobal(LocalPoint(0, 0, 0));
0225       GlobalPoint CenterPointCSCGlobal = TheChamber->toGlobal(LocalPoint(0, 0, 0));
0226 
0227       //LocalPoint CenterRollinCSCFrame = TheChamber->toLocal(CenterPointRollGlobal);
0228 
0229       float rpcphi = 0;
0230       float cscphi = 0;
0231 
0232       (CenterPointRollGlobal.barePhi() < 0) ? rpcphi = 2 * 3.141592 + CenterPointRollGlobal.barePhi()
0233                                             : rpcphi = CenterPointRollGlobal.barePhi();
0234 
0235       (CenterPointCSCGlobal.barePhi() < 0) ? cscphi = 2 * 3.1415926536 + CenterPointCSCGlobal.barePhi()
0236                                            : cscphi = CenterPointCSCGlobal.barePhi();
0237 
0238       float df = fabs(cscphi - rpcphi);
0239       float dr = fabs(CenterPointRollGlobal.perp() - CenterPointCSCGlobal.perp());
0240       float diffz = CenterPointRollGlobal.z() - CenterPointCSCGlobal.z();
0241       float dfg = df * 180. / 3.14159265;
0242 
0243       std::cout << "CSC \t " << rpcsrv.segment() << rpcsrv.name() << " dr=" << dr << " dz=" << diffz << " dfg=" << dfg
0244                 << std::endl;
0245 
0246       bool print = false;
0247 
0248       if ((dr > 200. || fabs(diffz) > 55. || dfg > 1.) && print) {
0249         std::cout << "\t \t problem CSC Station= " << CSCId.station() << " Ring= " << CSCId.ring()
0250                   << " Chamber= " << CSCId.chamber() << " cscphi=" << cscphi * 180 / 3.14159265
0251                   << "\t RPC Station= " << rpcId.station() << " ring= " << rpcId.ring() << " segment =-> "
0252                   << rpcsrv.segment() << " rollphi=" << rpcphi * 180 / 3.14159265 << "\t dfg=" << dfg << " dz=" << diffz
0253                   << " dr=" << dr << std::endl;
0254       }
0255     }
0256   }
0257 }
0258 
0259 //define this as a plug-in
0260 DEFINE_FWK_MODULE(RPCCSC);