Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 10:55:14

0001 #include <vector>
0002 #include <string>
0003 
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 
0008 #include "DQMServices/Core/interface/MonitorElement.h"
0009 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0010 #include "FWCore/Framework/interface/MakerMacros.h"
0011 
0012 #include "Geometry/RPCGeometry/interface/RPCRoll.h"
0013 #include "Geometry/RPCGeometry/interface/RPCGeometry.h"
0014 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0015 #include "DataFormats/L1TMuon/interface/CPPFDigi.h"
0016 
0017 // class decleration
0018 
0019 class L1TStage2CPPF : public DQMEDAnalyzer {
0020 public:
0021   // class constructor
0022   L1TStage2CPPF(const edm::ParameterSet& ps);
0023   // class destructor
0024   ~L1TStage2CPPF() override;
0025 
0026   // member functions
0027   edm::ESHandle<RPCGeometry> rpcGeom;
0028 
0029 protected:
0030   void analyze(const edm::Event&, const edm::EventSetup&) override;
0031   void bookHistograms(DQMStore::IBooker&, const edm::Run&, const edm::EventSetup&) override;
0032 
0033   // data members
0034 private:
0035   std::string monitorDir;
0036   bool verbose;
0037   float global_phi;
0038   const edm::EDGetTokenT<l1t::CPPFDigiCollection> cppfDigiToken_;
0039   int EMTF_sector;
0040   int EMTF_subsector;
0041   int EMTF_bx;
0042 
0043   std::vector<int> EMTFsector1bins;
0044   std::vector<int> EMTFsector2bins;
0045   std::vector<int> EMTFsector3bins;
0046   std::vector<int> EMTFsector4bins;
0047   std::vector<int> EMTFsector5bins;
0048   std::vector<int> EMTFsector6bins;
0049 
0050   std::map<int, std::vector<int>> fill_info;
0051 
0052   MonitorElement* Occupancy_EMTFSector;
0053   MonitorElement* Track_Bx;
0054 };
0055 
0056 L1TStage2CPPF::L1TStage2CPPF(const edm::ParameterSet& ps)
0057     : monitorDir(ps.getUntrackedParameter<std::string>("monitorDir", "")),
0058       verbose(ps.getUntrackedParameter<bool>("verbose", false)),
0059       global_phi(-1000),
0060       cppfDigiToken_(consumes<l1t::CPPFDigiCollection>(ps.getParameter<edm::InputTag>("cppfSource"))) {}
0061 
0062 L1TStage2CPPF::~L1TStage2CPPF() {}
0063 
0064 void L1TStage2CPPF::bookHistograms(DQMStore::IBooker& ibooker, const edm::Run& iRun, const edm::EventSetup& eveSetup) {
0065   ibooker.setCurrentFolder(monitorDir);
0066 
0067   Occupancy_EMTFSector = ibooker.book2D("Occupancy_EMTFSector", "Occupancy_EMTFSector", 36, 1., 37., 12, 1., 13.);
0068   Track_Bx = ibooker.book2D("Track_Bx", "Track_Bx", 12, 1., 13., 7, -3., 4.);
0069 }
0070 
0071 void L1TStage2CPPF::analyze(const edm::Event& eve, const edm::EventSetup& eveSetup) {
0072   if (verbose) {
0073     edm::LogInfo("L1TStage2CPPF") << "L1TStage2CPPF: analyze....";
0074   }
0075 
0076   edm::Handle<l1t::CPPFDigiCollection> CppfDigis;
0077   eve.getByToken(cppfDigiToken_, CppfDigis);
0078 
0079   //Fill the specific bin for each EMTF sector
0080   EMTFsector1bins.clear();
0081   EMTFsector2bins.clear();
0082   EMTFsector3bins.clear();
0083   EMTFsector4bins.clear();
0084   EMTFsector5bins.clear();
0085   EMTFsector6bins.clear();
0086   for (int i = 1; i < 7; i++) {
0087     EMTFsector1bins.push_back(i);
0088     EMTFsector2bins.push_back(i + 6);
0089     EMTFsector3bins.push_back(i + 12);
0090     EMTFsector4bins.push_back(i + 18);
0091     EMTFsector5bins.push_back(i + 24);
0092     EMTFsector6bins.push_back(i + 30);
0093   }
0094   //FIll the map for each EMTF sector
0095   fill_info[1] = EMTFsector1bins;
0096   fill_info[2] = EMTFsector2bins;
0097   fill_info[3] = EMTFsector3bins;
0098   fill_info[4] = EMTFsector4bins;
0099   fill_info[5] = EMTFsector5bins;
0100   fill_info[6] = EMTFsector6bins;
0101 
0102   for (auto& cppf_digis : *CppfDigis) {
0103     RPCDetId rpcId = cppf_digis.rpcId();
0104     int ring = rpcId.ring();
0105     int station = rpcId.station();
0106     int region = rpcId.region();
0107     int subsector = rpcId.subsector();
0108 
0109     //Region -
0110     if (region == -1) {
0111       //for Occupancy
0112       EMTF_sector = rpcId.sector();
0113       EMTF_subsector = fill_info[EMTF_sector][subsector - 1];
0114 
0115       if ((station == 4) && (ring == 3)) {
0116         Occupancy_EMTFSector->Fill(EMTF_subsector, 1);
0117       } else if ((station == 4) && (ring == 2)) {
0118         Occupancy_EMTFSector->Fill(EMTF_subsector, 2);
0119       } else if ((station == 3) && (ring == 3)) {
0120         Occupancy_EMTFSector->Fill(EMTF_subsector, 3);
0121       } else if ((station == 3) && (ring == 2)) {
0122         Occupancy_EMTFSector->Fill(EMTF_subsector, 4);
0123       } else if ((station == 2) && (ring == 2)) {
0124         Occupancy_EMTFSector->Fill(EMTF_subsector, 5);
0125       } else if ((station == 1) && (ring == 2)) {
0126         Occupancy_EMTFSector->Fill(EMTF_subsector, 6);
0127       }
0128 
0129       //for Track_Bx
0130       if (EMTF_sector >= 1 && EMTF_sector <= 6) {
0131         EMTF_bx = cppf_digis.bx();
0132         Track_Bx->Fill(7 - EMTF_sector, EMTF_bx);
0133       }
0134     }
0135     //Region +
0136     if (region == 1) {
0137       //for Occupancy
0138       EMTF_sector = rpcId.sector();
0139       EMTF_subsector = fill_info[EMTF_sector][subsector - 1];
0140 
0141       if ((station == 1) && (ring == 2)) {
0142         Occupancy_EMTFSector->Fill(EMTF_subsector, 7);
0143       } else if ((station == 2) && (ring == 2)) {
0144         Occupancy_EMTFSector->Fill(EMTF_subsector, 8);
0145       } else if ((station == 3) && (ring == 2)) {
0146         Occupancy_EMTFSector->Fill(EMTF_subsector, 9);
0147       } else if ((station == 3) && (ring == 3)) {
0148         Occupancy_EMTFSector->Fill(EMTF_subsector, 10);
0149       } else if ((station == 4) && (ring == 2)) {
0150         Occupancy_EMTFSector->Fill(EMTF_subsector, 11);
0151       } else if ((station == 4) && (ring == 3)) {
0152         Occupancy_EMTFSector->Fill(EMTF_subsector, 12);
0153       }
0154 
0155       //for Track_Bx
0156       if (EMTF_sector >= 1 && EMTF_sector <= 6) {
0157         EMTF_bx = cppf_digis.bx();
0158         Track_Bx->Fill(6 + EMTF_sector, EMTF_bx);
0159       }
0160     }
0161   }  //loop over CPPFDigis
0162 }
0163 DEFINE_FWK_MODULE(L1TStage2CPPF);