Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-05-18 03:27:06

0001 #include <string>
0002 #include <vector>
0003 #include <iostream>
0004 #include <map>
0005 
0006 #include "DQM/L1TMonitor/interface/L1TStage2EMTF.h"
0007 
0008 L1TStage2EMTF::L1TStage2EMTF(const edm::ParameterSet& ps)
0009     : daqToken(consumes<l1t::EMTFDaqOutCollection>(ps.getParameter<edm::InputTag>("emtfSource"))),
0010       hitToken(consumes<l1t::EMTFHitCollection>(ps.getParameter<edm::InputTag>("emtfSource"))),
0011       trackToken(consumes<l1t::EMTFTrackCollection>(ps.getParameter<edm::InputTag>("emtfSource"))),
0012       muonToken(consumes<l1t::RegionalMuonCandBxCollection>(ps.getParameter<edm::InputTag>("emtfSource"))),
0013       monitorDir(ps.getUntrackedParameter<std::string>("monitorDir", "")),
0014       verbose(ps.getUntrackedParameter<bool>("verbose", false)) {}
0015 
0016 void L1TStage2EMTF::bookHistograms(DQMStore::IBooker& ibooker, const edm::Run&, const edm::EventSetup&) {
0017   // Monitor Dir
0018   ibooker.setCurrentFolder(monitorDir);
0019 
0020   const std::array<std::string, 6> binNamesErrors{
0021       {"Corruptions", "Synch. Err.", "Synch. Mod.", "BX Mismatch", "Time Misalign", "FMM != Ready"}};
0022 
0023   // DAQ Output Monitor Elements
0024   emtfErrors = ibooker.book1D("emtfErrors", "EMTF Errors", 6, 0, 6);
0025   emtfErrors->setAxisTitle("Error Type (Corruptions Not Implemented)", 1);
0026   emtfErrors->setAxisTitle("Number of Errors", 2);
0027   for (unsigned int bin = 0; bin < binNamesErrors.size(); ++bin) {
0028     emtfErrors->setBinLabel(bin + 1, binNamesErrors[bin], 1);
0029   }
0030 
0031   // CSC LCT Monitor Elements
0032   int nChambs, nWires, nStrips;  // Number of chambers, wiregroups, and halfstrips in each station/ring pair
0033   std::string name, label;
0034   const std::array<std::string, 10> suffix_name{{"42", "41", "32", "31", "22", "21", "13", "12", "11b", "11a"}};
0035   const std::array<std::string, 10> suffix_label{
0036       {"4/2", "4/1", "3/2", "3/1", " 2/2", "2/1", "1/3", "1/2", "1/1b", "1/1a"}};
0037   const std::array<std::string, 12> binNames{
0038       {"ME-N", "ME-4", "ME-3", "ME-2", "ME-1b", "ME-1a", "ME+1a", "ME+1b", "ME+2", "ME+3", "ME+4", "ME+N"}};
0039 
0040   cscLCTBX = ibooker.book2D("cscLCTBX", "CSC LCT BX", 7, -3, 4, 20, 0, 20);
0041   cscLCTBX->setAxisTitle("BX", 1);
0042   for (int xbin = 1, xbin_label = -3; xbin <= 7; ++xbin, ++xbin_label) {
0043     cscLCTBX->setBinLabel(xbin, std::to_string(xbin_label), 1);
0044   }
0045   for (int ybin = 1; ybin <= 10; ++ybin) {
0046     cscLCTBX->setBinLabel(ybin, "ME-" + suffix_label[ybin - 1], 2);
0047     cscLCTBX->setBinLabel(21 - ybin, "ME+" + suffix_label[ybin - 1], 2);
0048   }
0049 
0050   cscLCTOccupancy = ibooker.book2D("cscLCTOccupancy", "CSC Chamber Occupancy", 54, 1, 55, 12, -6, 6);
0051   cscLCTOccupancy->setAxisTitle("Sector (CSCID 1-9 Unlabelled)", 1);
0052   for (int xbin = 1; xbin < 7; ++xbin) {
0053     cscLCTOccupancy->setBinLabel(xbin * 9 - 8, std::to_string(xbin), 1);
0054   }
0055   for (unsigned int ybin = 0; ybin < binNames.size(); ++ybin) {
0056     cscLCTOccupancy->setBinLabel(ybin + 1, binNames[ybin], 2);
0057   }
0058 
0059   //cscOccupancy designed to match the cscDQM plot
0060   cscDQMOccupancy = ibooker.book2D("cscDQMOccupancy", "CSC Chamber Occupancy", 42, 1, 43, 20, 0, 20);
0061   cscDQMOccupancy->setAxisTitle("10#circ Chamber (N=neighbor)", 1);
0062   int count = 0;
0063   for (int xbin = 1; xbin < 43; ++xbin) {
0064     cscDQMOccupancy->setBinLabel(xbin, std::to_string(xbin - count), 1);
0065     if (xbin == 2 || xbin == 9 || xbin == 16 || xbin == 23 || xbin == 30 || xbin == 37) {
0066       ++xbin;
0067       ++count;
0068       cscDQMOccupancy->setBinLabel(xbin, "N", 1);
0069     }
0070   }
0071   for (int ybin = 1; ybin <= 10; ++ybin) {
0072     cscDQMOccupancy->setBinLabel(ybin, "ME-" + suffix_label[ybin - 1], 2);
0073     cscDQMOccupancy->setBinLabel(21 - ybin, "ME+" + suffix_label[ybin - 1], 2);
0074   }
0075   cscDQMOccupancy->getTH2F()->GetXaxis()->SetCanExtend(false);  // Needed to stop multi-thread summing
0076 
0077   mpcLinkErrors = ibooker.book2D("mpcLinkErrors", "MPC Link Errors", 54, 1, 55, 12, -6, 6);
0078   mpcLinkErrors->setAxisTitle("Sector (CSCID 1-9 Unlabelled)", 1);
0079   for (int xbin = 1; xbin < 7; ++xbin) {
0080     mpcLinkErrors->setBinLabel(xbin * 9 - 8, std::to_string(xbin), 1);
0081   }
0082   for (unsigned int ybin = 0; ybin < binNames.size(); ++ybin) {
0083     mpcLinkErrors->setBinLabel(ybin + 1, binNames[ybin], 2);
0084   }
0085 
0086   mpcLinkGood = ibooker.book2D("mpcLinkGood", "MPC Good Links", 54, 1, 55, 12, -6, 6);
0087   mpcLinkGood->setAxisTitle("Sector (CSCID 1-9 Unlabelled)", 1);
0088   for (int xbin = 1; xbin < 7; ++xbin) {
0089     mpcLinkGood->setBinLabel(xbin * 9 - 8, std::to_string(xbin), 1);
0090   }
0091   for (unsigned int ybin = 0; ybin < binNames.size(); ++ybin) {
0092     mpcLinkGood->setBinLabel(ybin + 1, binNames[ybin], 2);
0093   }
0094 
0095   // RPC Monitor Elements
0096   const std::array<std::string, 6> rpc_name{{"43", "42", "33", "32", "22", "12"}};
0097   const std::array<std::string, 6> rpc_label{{"4/3", "4/2", "3/3", "3/2", "2/2", "1/2"}};
0098 
0099   rpcHitBX = ibooker.book2D("rpcHitBX", "RPC Hit BX", 7, -3, 4, 12, 0, 12);
0100   rpcHitBX->setAxisTitle("BX", 1);
0101   for (int xbin = 1, xbin_label = -3; xbin <= 7; ++xbin, ++xbin_label) {
0102     rpcHitBX->setBinLabel(xbin, std::to_string(xbin_label), 1);
0103   }
0104   for (int ybin = 1; ybin <= 6; ++ybin) {
0105     rpcHitBX->setBinLabel(ybin, "RE-" + rpc_label[ybin - 1], 2);
0106     rpcHitBX->setBinLabel(13 - ybin, "RE+" + rpc_label[ybin - 1], 2);
0107   }
0108 
0109   rpcHitOccupancy = ibooker.book2D("rpcHitOccupancy", "RPC Chamber Occupancy", 42, 1, 43, 12, 0, 12);
0110   rpcHitOccupancy->setAxisTitle("Sector (N=neighbor)", 1);
0111   for (int bin = 1; bin <= 6; ++bin) {
0112     rpcHitOccupancy->setBinLabel(bin * 7 - 6, std::to_string(bin), 1);
0113     rpcHitOccupancy->setBinLabel(bin * 7, "N", 1);
0114     rpcHitOccupancy->setBinLabel(bin, "RE-" + rpc_label[bin - 1], 2);
0115     rpcHitOccupancy->setBinLabel(13 - bin, "RE+" + rpc_label[bin - 1], 2);
0116   }
0117   rpcHitOccupancy->getTH2F()->GetXaxis()->SetCanExtend(false);  // Needed to stop multi-thread summing
0118 
0119   // GEM Monitor Elements
0120   // Add GEM Oct 27 2020
0121   hitTypeBX = ibooker.book2D("hitTypeBX", "Hit Type BX", 4, 0.5, 4.5, 7, -3, 4);
0122   hitTypeBX->setBinLabel(1, "CSC", 1);
0123   hitTypeBX->setBinLabel(2, "RPC", 1);
0124   hitTypeBX->setBinLabel(3, "GEM", 1);
0125   hitTypeBX->setBinLabel(4, "Tot", 1);
0126   for (int ybin = 1; ybin < 8; ybin++)
0127     hitTypeBX->setBinLabel(ybin, std::to_string(ybin - 4), 2);
0128 
0129   gemHitBX = ibooker.book2D("gemHitBX", "GEM Hit BX", 7, -3, 4, 2, 0, 2);
0130   gemHitBX->setAxisTitle("BX", 1);
0131   for (int xbin = 1, xbin_label = -3; xbin <= 7; ++xbin, ++xbin_label) {
0132     gemHitBX->setBinLabel(xbin, std::to_string(xbin_label), 1);
0133   }
0134   gemHitBX->setBinLabel(1, "GE-1/1", 2);
0135   gemHitBX->setBinLabel(2, "GE+1/1", 2);
0136 
0137   gemHitOccupancy = ibooker.book2D("gemHitOccupancy", "GEM Chamber Occupancy", 42, 1, 43, 2, 0, 2);
0138   gemHitOccupancy->setAxisTitle("10#circ Chambers (N=neighbor)", 1);
0139   count = 0;
0140   for (int xbin = 1; xbin < 43; ++xbin) {
0141     gemHitOccupancy->setBinLabel(xbin, std::to_string(xbin - count), 1);
0142     if (xbin == 2 || xbin == 9 || xbin == 16 || xbin == 23 || xbin == 30 || xbin == 37) {
0143       ++xbin;
0144       ++count;
0145       gemHitOccupancy->setBinLabel(xbin, "N", 1);
0146     }
0147   }
0148 
0149   gemHitOccupancy->setBinLabel(1, "GE-1/1", 2);
0150   gemHitOccupancy->setBinLabel(2, "GE+1/1", 2);
0151   gemHitOccupancy->getTH2F()->GetXaxis()->SetCanExtend(false);  // Needed to stop multi-thread summing
0152 
0153   // Track Monitor Elements
0154   emtfnTracks = ibooker.book1D("emtfnTracks", "Number of EMTF Tracks per Event", 11, 0, 11);
0155   for (int xbin = 1; xbin <= 10; ++xbin) {
0156     emtfnTracks->setBinLabel(xbin, std::to_string(xbin - 1), 1);
0157   }
0158   emtfnTracks->setBinLabel(11, "Overflow", 1);
0159 
0160   emtfTracknHits = ibooker.book1D("emtfTracknHits", "Number of Hits per EMTF Track", 5, 0, 5);
0161   for (int xbin = 1; xbin <= 5; ++xbin) {
0162     emtfTracknHits->setBinLabel(xbin, std::to_string(xbin - 1), 1);
0163   }
0164 
0165   emtfTrackBX = ibooker.book2D("emtfTrackBX", "EMTF Track Bunch Crossing", 12, -6, 6, 7, -3, 4);
0166   emtfTrackBX->setAxisTitle("Sector (Endcap)", 1);
0167   for (int xbin = 0; xbin < 6; ++xbin) {
0168     emtfTrackBX->setBinLabel(xbin + 1, std::to_string(6 - xbin) + " (-)", 1);
0169     emtfTrackBX->setBinLabel(12 - xbin, std::to_string(6 - xbin) + " (+)", 1);
0170   }
0171   emtfTrackBX->setAxisTitle("Track BX", 2);
0172   for (int ybin = 1, i = -3; ybin <= 7; ++ybin, ++i) {
0173     emtfTrackBX->setBinLabel(ybin, std::to_string(i), 2);
0174   }
0175 
0176   emtfTrackPt = ibooker.book1D("emtfTrackPt", "EMTF Track p_{T}", 256, 1, 257);
0177   emtfTrackPt->setAxisTitle("Track p_{T} [GeV]", 1);
0178 
0179   emtfTrackEta = ibooker.book1D("emtfTrackEta", "EMTF Track #eta", 100, -2.5, 2.5);
0180   emtfTrackEta->setAxisTitle("Track #eta", 1);
0181 
0182   emtfTrackPhi = ibooker.book1D("emtfTrackPhi", "EMTF Track #phi", 126, -3.15, 3.15);
0183   emtfTrackPhi->setAxisTitle("Track #phi", 1);
0184 
0185   emtfTrackOccupancy = ibooker.book2D("emtfTrackOccupancy", "EMTF Track Occupancy", 100, -2.5, 2.5, 126, -3.15, 3.15);
0186   emtfTrackOccupancy->setAxisTitle("#eta", 1);
0187   emtfTrackOccupancy->setAxisTitle("#phi", 2);
0188 
0189   emtfTrackMode = ibooker.book1D("emtfTrackMode", "EMTF Track Mode", 16, 0, 16);
0190   emtfTrackMode->setAxisTitle("Mode", 1);
0191 
0192   emtfTrackQuality = ibooker.book1D("emtfTrackQuality", "EMTF Track Quality", 16, 0, 16);
0193   emtfTrackQuality->setAxisTitle("Quality", 1);
0194 
0195   emtfTrackQualityVsMode = ibooker.book2D("emtfTrackQualityVsMode", "EMTF Track Quality vs Mode", 16, 0, 16, 16, 0, 16);
0196   emtfTrackQualityVsMode->setAxisTitle("Mode", 1);
0197   emtfTrackQualityVsMode->setAxisTitle("Quality", 2);
0198 
0199   RPCvsEMTFTrackMode = ibooker.book2D("RPCvsEMTFTrackMode", "RPC Mode vs EMTF TrackMode", 16, 0, 16, 16, 0, 16);
0200   RPCvsEMTFTrackMode->setAxisTitle("EMTF Mode", 1);
0201   RPCvsEMTFTrackMode->setAxisTitle("RPC Mode", 2);
0202 
0203   for (int bin = 1; bin <= 16; ++bin) {
0204     emtfTrackMode->setBinLabel(bin, std::to_string(bin - 1), 1);
0205     emtfTrackQuality->setBinLabel(bin, std::to_string(bin - 1), 1);
0206     emtfTrackQualityVsMode->setBinLabel(bin, std::to_string(bin - 1), 1);
0207     emtfTrackQualityVsMode->setBinLabel(bin, std::to_string(bin - 1), 2);
0208     RPCvsEMTFTrackMode->setBinLabel(bin, std::to_string(bin - 1), 1);
0209     RPCvsEMTFTrackMode->setBinLabel(bin, std::to_string(bin - 1), 2);
0210   }
0211 
0212   //Chad Freer May 8 2018 (Selected Tracks)
0213   ibooker.setCurrentFolder(monitorDir + "/SelectedTracks");
0214 
0215   //Chad Freer May 8 2018 (High Quality Track Plots)
0216   emtfTrackPtHighQuality = ibooker.book1D("emtfTrackPtHighQuality", "EMTF High Quality Track p_{T}", 256, 1, 257);
0217   emtfTrackPtHighQuality->setAxisTitle("Track p_{T} [GeV] (Quality #geq 12)", 1);
0218 
0219   emtfTrackEtaHighQuality = ibooker.book1D("emtfTrackEtaHighQuality", "EMTF High Quality Track #eta", 100, -2.5, 2.5);
0220   emtfTrackEtaHighQuality->setAxisTitle("Track #eta (Quality #geq 12)", 1);
0221 
0222   emtfTrackPhiHighQuality = ibooker.book1D("emtfTrackPhiHighQuality", "EMTF High Quality #phi", 126, -3.15, 3.15);
0223   emtfTrackPhiHighQuality->setAxisTitle("Track #phi (Quality #geq 12)", 1);
0224 
0225   emtfTrackOccupancyHighQuality = ibooker.book2D(
0226       "emtfTrackOccupancyHighQuality", "EMTF High Quality Track Occupancy", 100, -2.5, 2.5, 126, -3.15, 3.15);
0227   emtfTrackOccupancyHighQuality->setAxisTitle("#eta", 1);
0228   emtfTrackOccupancyHighQuality->setAxisTitle("#phi", 2);
0229 
0230   //Chad Freer may 8 2018 (High Quality and High PT [22 GeV] Track Plots)
0231   emtfTrackPtHighQualityHighPT =
0232       ibooker.book1D("emtfTrackPtHighQualityHighPT", "EMTF High Quality High PT Track p_{T}", 256, 1, 257);
0233   emtfTrackPtHighQualityHighPT->setAxisTitle("Track p_{T} [GeV] (Quality #geq 12 and pT>22)", 1);
0234 
0235   emtfTrackEtaHighQualityHighPT =
0236       ibooker.book1D("emtfTrackEtaHighQualityHighPT", "EMTF High Quality High PT Track #eta", 100, -2.5, 2.5);
0237   emtfTrackEtaHighQualityHighPT->setAxisTitle("Track #eta (Quality #geq 12 and pT>22)", 1);
0238 
0239   emtfTrackPhiHighQualityHighPT =
0240       ibooker.book1D("emtfTrackPhiHighQualityHighPT", "EMTF High Quality High PT #phi", 126, -3.15, 3.15);
0241   emtfTrackPhiHighQualityHighPT->setAxisTitle("Track #phi (Quality #geq 12 and pT>22)", 1);
0242 
0243   emtfTrackOccupancyHighQualityHighPT = ibooker.book2D("emtfTrackOccupancyHighQualityHighPT",
0244                                                        "EMTF High Quality High PT Track Occupancy",
0245                                                        100,
0246                                                        -2.5,
0247                                                        2.5,
0248                                                        126,
0249                                                        -3.15,
0250                                                        3.15);
0251   emtfTrackOccupancyHighQualityHighPT->setAxisTitle("#eta", 1);
0252   emtfTrackOccupancyHighQualityHighPT->setAxisTitle("#phi", 2);
0253   //Chad Freer May 8 2018 (END new plots)
0254 
0255   // CSC Input
0256   ibooker.setCurrentFolder(monitorDir + "/CSCInput");
0257 
0258   for (int hist = 0, i = 0; hist < 20; ++hist, i = hist % 10) {
0259     if (hist < 10) {
0260       name = "MENeg" + suffix_name[i];
0261       label = "ME-" + suffix_label[i];
0262     } else {
0263       name = "MEPos" + suffix_name[9 - i];
0264       label = "ME+" + suffix_label[9 - i];
0265     }
0266 
0267     if (hist < 6) {
0268       nChambs = (i % 2) ? 18 : 36;
0269     } else if (hist > 13) {
0270       nChambs = (i % 2) ? 36 : 18;
0271     } else {
0272       nChambs = 36;
0273     }
0274 
0275     const std::array<int, 10> wiregroups{{64, 96, 64, 96, 64, 112, 32, 64, 48, 48}};
0276     const std::array<int, 10> halfstrips{{160, 160, 160, 160, 160, 160, 128, 160, 128, 96}};
0277 
0278     if (hist < 10) {
0279       nWires = wiregroups[hist];
0280       nStrips = halfstrips[hist];
0281     } else {
0282       nWires = wiregroups[19 - hist];
0283       nStrips = halfstrips[19 - hist];
0284     }
0285 
0286     cscLCTStrip[hist] = ibooker.book1D("cscLCTStrip" + name, "CSC Halfstrip " + label, nStrips, 0, nStrips);
0287     cscLCTStrip[hist]->setAxisTitle("Cathode Halfstrip, " + label, 1);
0288 
0289     cscLCTWire[hist] = ibooker.book1D("cscLCTWire" + name, "CSC Wiregroup " + label, nWires, 0, nWires);
0290     cscLCTWire[hist]->setAxisTitle("Anode Wiregroup, " + label, 1);
0291 
0292     cscChamberStrip[hist] = ibooker.book2D(
0293         "cscChamberStrip" + name, "CSC Halfstrip " + label, nChambs, 1, 1 + nChambs, nStrips, 0, nStrips);
0294     cscChamberStrip[hist]->setAxisTitle("Chamber, " + label, 1);
0295     cscChamberStrip[hist]->setAxisTitle("Cathode Halfstrip", 2);
0296 
0297     cscChamberWire[hist] =
0298         ibooker.book2D("cscChamberWire" + name, "CSC Wiregroup " + label, nChambs, 1, 1 + nChambs, nWires, 0, nWires);
0299     cscChamberWire[hist]->setAxisTitle("Chamber, " + label, 1);
0300     cscChamberWire[hist]->setAxisTitle("Anode Wiregroup", 2);
0301 
0302     for (int bin = 1; bin <= nChambs; ++bin) {
0303       cscChamberStrip[hist]->setBinLabel(bin, std::to_string(bin), 1);
0304       cscChamberWire[hist]->setBinLabel(bin, std::to_string(bin), 1);
0305     }
0306   }
0307 
0308   // RPC Input
0309   ibooker.setCurrentFolder(monitorDir + "/RPCInput");
0310 
0311   for (int hist = 0, i = 0; hist < 12; ++hist, i = hist % 6) {
0312     if (hist < 6) {
0313       name = "RENeg" + rpc_name[i];
0314       label = "RE-" + rpc_label[i];
0315     } else {
0316       name = "REPos" + rpc_name[5 - i];
0317       label = "RE+" + rpc_label[5 - i];
0318     }
0319     rpcHitPhi[hist] = ibooker.book1D("rpcHitPhi" + name, "RPC Hit Phi " + label, 1250, 0, 1250);
0320     rpcHitPhi[hist]->setAxisTitle("#phi", 1);
0321     rpcHitTheta[hist] = ibooker.book1D("rpcHitTheta" + name, "RPC Hit Theta " + label, 32, 0, 32);
0322     rpcHitTheta[hist]->setAxisTitle("#theta", 1);
0323     rpcChamberPhi[hist] = ibooker.book2D("rpcChamberPhi" + name, "RPC Chamber Phi " + label, 36, 1, 37, 1250, 0, 1250);
0324     rpcChamberPhi[hist]->setAxisTitle("Chamber", 1);
0325     rpcChamberPhi[hist]->setAxisTitle("#phi", 2);
0326     rpcChamberTheta[hist] =
0327         ibooker.book2D("rpcChamberTheta" + name, "RPC Chamber Theta " + label, 36, 1, 37, 32, 0, 32);
0328     rpcChamberTheta[hist]->setAxisTitle("Chamber", 1);
0329     rpcChamberTheta[hist]->setAxisTitle("#theta", 2);
0330     for (int xbin = 1; xbin < 37; ++xbin) {
0331       rpcChamberPhi[hist]->setBinLabel(xbin, std::to_string(xbin), 1);
0332       rpcChamberTheta[hist]->setBinLabel(xbin, std::to_string(xbin), 1);
0333     }
0334   }
0335 
0336   // GEM Input Nov 03 2020
0337   ibooker.setCurrentFolder(monitorDir + "/GEMInput");
0338   for (int hist = 0; hist < 2; hist++) {
0339     if (hist == 0) {
0340       name = "GENeg11";
0341       label = "GE-1/1";
0342     }
0343     if (hist == 1) {
0344       name = "GEPos11";
0345       label = "GE+1/1";
0346     }
0347     nChambs = 36;
0348     nStrips = 192;  //use nStrips for number of pads
0349     gemChamberPad[hist] = ibooker.book2D(
0350         "gemChamberPad" + name, "GEM Chamber Pad " + label, nChambs, 1, 1 + nChambs, nStrips, 0, nStrips);  // pads 0-191
0351     gemChamberPad[hist]->setAxisTitle("Chamber, " + label, 1);
0352     gemChamberPad[hist]->setAxisTitle("Pad", 2);
0353     gemChamberPartition[hist] =
0354         ibooker.book2D("gemChamberPartition" + name,
0355                        "GEM Chamber Partition " + label,
0356                        nChambs,
0357                        1,
0358                        1 + nChambs,
0359                        9,
0360                        0,
0361                        9);  // partitions 1-8 or 0-7. There have been changes in different firmware/unpacker verisions.
0362     gemChamberPartition[hist]->setAxisTitle("Chamber, " + label, 1);
0363     gemChamberPartition[hist]->setAxisTitle("Partition", 2);
0364     for (int bin = 1; bin <= nChambs; ++bin) {
0365       gemChamberPad[hist]->setBinLabel(bin, std::to_string(bin), 1);
0366       gemChamberPartition[hist]->setBinLabel(bin, std::to_string(bin), 1);
0367     }
0368   }
0369 
0370   // CSC LCT and RPC Hit Timing
0371   ibooker.setCurrentFolder(monitorDir + "/Timing");
0372 
0373   cscTimingTot = ibooker.book2D("cscTimingTotal", "CSC Total BX ", 42, 1, 43, 20, 0, 20);
0374   cscTimingTot->setAxisTitle("10#circ Chamber (N=neighbor)", 1);
0375 
0376   rpcHitTimingTot = ibooker.book2D("rpcHitTimingTot", "RPC Chamber Occupancy ", 42, 1, 43, 12, 0, 12);
0377   rpcHitTimingTot->setAxisTitle("Sector (N=neighbor)", 1);
0378 
0379   gemHitTimingTot =
0380       ibooker.book2D("gemHitTimingTot", "GEM Chamber Occupancy ", 42, 1, 43, 2, 0, 2);  // Add GEM Timing Oct 27 2020
0381   gemHitTimingTot->setAxisTitle("10#circ Chamber (N=neighbor)", 1);
0382   const std::array<std::string, 5> nameBX{{"BXNeg1", "BXPos1", "BXNeg2", "BXPos2", "BX0"}};
0383   const std::array<std::string, 5> labelBX{{"BX -1", "BX +1", "BX -2", "BX +2", "BX 0"}};
0384 
0385   for (int hist = 0; hist < 5; ++hist) {
0386     count = 0;
0387     cscLCTTiming[hist] =
0388         ibooker.book2D("cscLCTTiming" + nameBX[hist], "CSC Chamber Occupancy " + labelBX[hist], 42, 1, 43, 20, 0, 20);
0389     cscLCTTiming[hist]->setAxisTitle("10#circ Chamber", 1);
0390 
0391     for (int xbin = 1; xbin < 43; ++xbin) {
0392       cscLCTTiming[hist]->setBinLabel(xbin, std::to_string(xbin - count), 1);
0393       if (hist == 0)
0394         cscTimingTot->setBinLabel(xbin, std::to_string(xbin - count), 1);  //only fill once in the loop
0395       if (xbin == 2 || xbin == 9 || xbin == 16 || xbin == 23 || xbin == 30 || xbin == 37) {
0396         ++xbin;
0397         ++count;
0398         cscLCTTiming[hist]->setBinLabel(xbin, "N", 1);
0399         if (hist == 0)
0400           cscTimingTot->setBinLabel(xbin, "N", 1);
0401       }
0402     }
0403 
0404     for (int ybin = 1; ybin <= 10; ++ybin) {
0405       cscLCTTiming[hist]->setBinLabel(ybin, "ME-" + suffix_label[ybin - 1], 2);
0406       cscLCTTiming[hist]->setBinLabel(21 - ybin, "ME+" + suffix_label[ybin - 1], 2);
0407       if (hist == 0)
0408         cscTimingTot->setBinLabel(ybin, "ME-" + suffix_label[ybin - 1], 2);
0409       if (hist == 0)
0410         cscTimingTot->setBinLabel(21 - ybin, "ME+" + suffix_label[ybin - 1], 2);
0411     }
0412     if (hist == 0)
0413       cscTimingTot->getTH2F()->GetXaxis()->SetCanExtend(false);      // Needed to stop multi-thread summing
0414     cscLCTTiming[hist]->getTH2F()->GetXaxis()->SetCanExtend(false);  // Needed to stop multi-thread summing
0415 
0416     rpcHitTiming[hist] =
0417         ibooker.book2D("rpcHitTiming" + nameBX[hist], "RPC Chamber Occupancy " + labelBX[hist], 42, 1, 43, 12, 0, 12);
0418     rpcHitTiming[hist]->setAxisTitle("Sector (N=neighbor)", 1);
0419     for (int bin = 1; bin < 7; ++bin) {
0420       rpcHitTiming[hist]->setBinLabel(bin * 7 - 6, std::to_string(bin), 1);
0421       rpcHitTiming[hist]->setBinLabel(bin * 7, "N", 1);
0422       rpcHitTiming[hist]->setBinLabel(bin, "RE-" + rpc_label[bin - 1], 2);
0423       rpcHitTiming[hist]->setBinLabel(13 - bin, "RE+" + rpc_label[bin - 1], 2);
0424     }
0425     rpcHitTiming[hist]->getTH2F()->GetXaxis()->SetCanExtend(false);  // Needed to stop multi-thread summing
0426     if (hist == 0) {
0427       for (int bin = 1; bin < 7; ++bin) {
0428         rpcHitTimingTot->setBinLabel(bin * 7 - 6, std::to_string(bin), 1);
0429         rpcHitTimingTot->setBinLabel(bin * 7, "N", 1);
0430         rpcHitTimingTot->setBinLabel(bin, "RE-" + rpc_label[bin - 1], 2);
0431         rpcHitTimingTot->setBinLabel(13 - bin, "RE+" + rpc_label[bin - 1], 2);
0432       }
0433       rpcHitTimingTot->getTH2F()->GetXaxis()->SetCanExtend(false);  // Needed to stop multi-thread summing
0434     }
0435     //if (hist == 4) continue; // Don't book for BX = 0
0436 
0437     // Add GEM Timing Oct 27 2020
0438     gemHitTiming[hist] =
0439         ibooker.book2D("gemHitTiming" + nameBX[hist], "GEM Chamber Occupancy " + labelBX[hist], 42, 1, 43, 2, 0, 2);
0440     gemHitTiming[hist]->setAxisTitle("10#circ Chamber", 1);
0441     count = 0;
0442     for (int xbin = 1; xbin < 43; ++xbin) {
0443       gemHitTiming[hist]->setBinLabel(xbin, std::to_string(xbin - count), 1);
0444       if (hist == 0)
0445         gemHitTimingTot->setBinLabel(xbin, std::to_string(xbin - count), 1);  //only fill once in the loop
0446       if (xbin == 2 || xbin == 9 || xbin == 16 || xbin == 23 || xbin == 30 || xbin == 37) {
0447         ++xbin;
0448         ++count;
0449         gemHitTiming[hist]->setBinLabel(xbin, "N", 1);
0450         if (hist == 0)
0451           gemHitTimingTot->setBinLabel(xbin, "N", 1);
0452       }
0453     }
0454 
0455     gemHitTiming[hist]->setBinLabel(1, "GE-1/1", 2);
0456     gemHitTiming[hist]->setBinLabel(2, "GE+1/1", 2);
0457     if (hist == 0) {
0458       gemHitTimingTot->setBinLabel(1, "GE-1/1", 2);
0459       gemHitTimingTot->setBinLabel(2, "GE+1/1", 2);
0460       gemHitTimingTot->getTH2F()->GetXaxis()->SetCanExtend(false);  // Needed to stop multi-thread summing
0461     }
0462     gemHitTiming[hist]->getTH2F()->GetXaxis()->SetCanExtend(false);  // Needed to stop multi-thread summing
0463 
0464     count = 0;
0465     cscLCTTimingFrac[hist] = ibooker.book2D(
0466         "cscLCTTimingFrac" + nameBX[hist], "CSC Chamber Occupancy " + labelBX[hist], 42, 1, 43, 20, 0, 20);
0467     cscLCTTimingFrac[hist]->setAxisTitle("10#circ Chambers", 1);
0468     for (int xbin = 1; xbin < 43; ++xbin) {
0469       cscLCTTimingFrac[hist]->setBinLabel(xbin, std::to_string(xbin - count), 1);
0470       if (xbin == 2 || xbin == 9 || xbin == 16 || xbin == 23 || xbin == 30 || xbin == 37) {
0471         ++xbin;
0472         ++count;
0473         cscLCTTimingFrac[hist]->setBinLabel(xbin, "N", 1);
0474       }
0475     }
0476     for (int ybin = 1; ybin <= 10; ++ybin) {
0477       cscLCTTimingFrac[hist]->setBinLabel(ybin, "ME-" + suffix_label[ybin - 1], 2);
0478       cscLCTTimingFrac[hist]->setBinLabel(21 - ybin, "ME+" + suffix_label[ybin - 1], 2);
0479     }
0480     cscLCTTimingFrac[hist]->getTH2F()->GetXaxis()->SetCanExtend(false);  // Needed to stop multi-thread summing
0481 
0482     rpcHitTimingFrac[hist] = ibooker.book2D(
0483         "rpcHitTimingFrac" + nameBX[hist], "RPC Chamber Fraction in " + labelBX[hist], 42, 1, 43, 12, 0, 12);
0484     rpcHitTimingFrac[hist]->setAxisTitle("Sector (N=neighbor)", 1);
0485     for (int bin = 1; bin < 7; ++bin) {
0486       rpcHitTimingFrac[hist]->setBinLabel(bin * 7 - 6, std::to_string(bin), 1);
0487       rpcHitTimingFrac[hist]->setBinLabel(bin * 7, "N", 1);
0488       rpcHitTimingFrac[hist]->setBinLabel(bin, "RE-" + rpc_label[bin - 1], 2);
0489       rpcHitTimingFrac[hist]->setBinLabel(13 - bin, "RE+" + rpc_label[bin - 1], 2);
0490     }
0491 
0492     // Add GEM Timing Oct 27 2020
0493     gemHitTimingFrac[hist] =
0494         ibooker.book2D("gemHitTimingFrac" + nameBX[hist], "GEM Chamber Occupancy " + labelBX[hist], 42, 1, 43, 2, 0, 2);
0495     gemHitTimingFrac[hist]->setAxisTitle("10#circ Chambers", 1);
0496     count = 0;
0497     for (int xbin = 1; xbin < 43; ++xbin) {
0498       gemHitTimingFrac[hist]->setBinLabel(xbin, std::to_string(xbin - count), 1);
0499       if (xbin == 2 || xbin == 9 || xbin == 16 || xbin == 23 || xbin == 30 || xbin == 37) {
0500         ++xbin;
0501         ++count;
0502         gemHitTimingFrac[hist]->setBinLabel(xbin, "N", 1);
0503       }
0504     }
0505     gemHitTimingFrac[hist]->setBinLabel(1, "GE-1/1", 2);
0506     gemHitTimingFrac[hist]->setBinLabel(2, "GE+1/1", 2);
0507     gemHitTimingFrac[hist]->getTH2F()->GetXaxis()->SetCanExtend(false);  // Needed to stop multi-thread summing
0508   }
0509 
0510   rpcHitTimingInTrack =
0511       ibooker.book2D("rpcHitTimingInTrack", "RPC Hit Timing (matched to track in BX 0)", 7, -3, 4, 12, 0, 12);
0512   rpcHitTimingInTrack->setAxisTitle("BX", 1);
0513   for (int xbin = 1, xbin_label = -3; xbin <= 7; ++xbin, ++xbin_label) {
0514     rpcHitTimingInTrack->setBinLabel(xbin, std::to_string(xbin_label), 1);
0515   }
0516   for (int ybin = 1; ybin <= 6; ++ybin) {
0517     rpcHitTimingInTrack->setBinLabel(ybin, "RE-" + rpc_label[ybin - 1], 2);
0518     rpcHitTimingInTrack->setBinLabel(13 - ybin, "RE+" + rpc_label[ybin - 1], 2);
0519   }
0520 
0521   const std::array<std::string, 3> nameNumStation{{"4Station", "3Station", "2Station"}};
0522   const std::array<std::string, 3> labelNumStation{{"4 Station Track", "3 Station Track", "2 Station Track"}};
0523 
0524   for (int hist = 0; hist < 3; ++hist) {
0525     emtfTrackBXVsCSCLCT[hist] = ibooker.book2D("emtfTrackBXVsCSCLCT" + nameNumStation[hist],
0526                                                "EMTF " + labelNumStation[hist] + " BX vs CSC LCT BX",
0527                                                7,
0528                                                -3,
0529                                                4,
0530                                                7,
0531                                                -3,
0532                                                4);
0533     emtfTrackBXVsCSCLCT[hist]->setAxisTitle("LCT BX", 1);
0534     emtfTrackBXVsCSCLCT[hist]->setAxisTitle("Track BX", 2);
0535     for (int bin = 1, bin_label = -3; bin <= 7; ++bin, ++bin_label) {
0536       emtfTrackBXVsCSCLCT[hist]->setBinLabel(bin, std::to_string(bin_label), 1);
0537       emtfTrackBXVsCSCLCT[hist]->setBinLabel(bin, std::to_string(bin_label), 2);
0538     }
0539     emtfTrackBXVsRPCHit[hist] = ibooker.book2D("emtfTrackBXVsRPCHit" + nameNumStation[hist],
0540                                                "EMTF " + labelNumStation[hist] + " BX vs RPC Hit BX",
0541                                                7,
0542                                                -3,
0543                                                4,
0544                                                7,
0545                                                -3,
0546                                                4);
0547     emtfTrackBXVsRPCHit[hist]->setAxisTitle("Hit BX", 1);
0548     emtfTrackBXVsRPCHit[hist]->setAxisTitle("Track BX", 2);
0549     for (int bin = 1, bin_label = -3; bin <= 7; ++bin, ++bin_label) {
0550       emtfTrackBXVsRPCHit[hist]->setBinLabel(bin, std::to_string(bin_label), 1);
0551       emtfTrackBXVsRPCHit[hist]->setBinLabel(bin, std::to_string(bin_label), 2);
0552     }
0553     // Add GEM vs track BX Dec 05 2020
0554     emtfTrackBXVsGEMHit[hist] = ibooker.book2D("emtfTrackBXVsGEMHit" + nameNumStation[hist],
0555                                                "EMTF " + labelNumStation[hist] + " BX vs GEM Hit BX",
0556                                                7,
0557                                                -3,
0558                                                4,
0559                                                7,
0560                                                -3,
0561                                                4);
0562     emtfTrackBXVsGEMHit[hist]->setAxisTitle("Hit BX", 1);
0563     emtfTrackBXVsGEMHit[hist]->setAxisTitle("Track BX", 2);
0564     for (int bin = 1, bin_label = -3; bin <= 7; ++bin, ++bin_label) {
0565       emtfTrackBXVsGEMHit[hist]->setBinLabel(bin, std::to_string(bin_label), 1);
0566       emtfTrackBXVsGEMHit[hist]->setBinLabel(bin, std::to_string(bin_label), 2);
0567     }
0568   }  // End loop: for (int hist = 0; hist < 3; ++hist)
0569 
0570   // Add mode vs BXdiff comparison Dec 07 2020
0571   const std::array<std::string, 2> nameGEMStation{{"GENeg11", "GEPos11"}};
0572   const std::array<std::string, 2> labelGEMStation{{"GE-1/1", "GE+1/1"}};
0573   const std::array<std::string, 8> nameCSCStation{
0574       {"MENeg4", "MENeg3", "MENeg2", "MENeg1", "MEPos1", "MEPos2", "MEPos3", "MEPos4"}};
0575   const std::array<std::string, 8> labelCSCStation{{"ME-4", "ME-3", "ME-2", "ME-1", "ME+1", "ME+2", "ME+3", "ME+4"}};
0576   const std::array<std::string, 8> nameRPCStation{
0577       {"RENeg4", "RENeg3", "RENeg2", "RENeg1", "REPos1", "REPos2", "REPos3", "REPos4"}};
0578   const std::array<std::string, 8> labelRPCStation{{"RE-4", "RE-3", "RE-2", "RE-1", "RE+1", "RE+2", "RE+3", "RE+4"}};
0579 
0580   for (int iGEM = 0; iGEM < 2; iGEM++) {
0581     emtfTrackModeVsGEMBXDiff[iGEM] = ibooker.book2D("emtfTrackModeVsGEMBXDiff" + nameGEMStation[iGEM],
0582                                                     "EMTF Track Mode vs (Track BX - GEM BX) " + labelGEMStation[iGEM],
0583                                                     9,
0584                                                     -4,
0585                                                     5,
0586                                                     16,
0587                                                     0,
0588                                                     16);
0589     emtfTrackModeVsGEMBXDiff[iGEM]->setAxisTitle("Track BX - GEM BX", 1);
0590     emtfTrackModeVsGEMBXDiff[iGEM]->setAxisTitle("Track Mode", 2);
0591   }
0592   for (int iCSC = 0; iCSC < 8; iCSC++) {
0593     emtfTrackModeVsCSCBXDiff[iCSC] = ibooker.book2D("emtfTrackModeVsCSCBXDiff" + nameCSCStation[iCSC],
0594                                                     "EMTF Track Mode vs (Track BX - LCT BX) " + labelCSCStation[iCSC],
0595                                                     9,
0596                                                     -4,
0597                                                     5,
0598                                                     16,
0599                                                     0,
0600                                                     16);
0601     emtfTrackModeVsCSCBXDiff[iCSC]->setAxisTitle("Track BX - LCT BX", 1);
0602     emtfTrackModeVsCSCBXDiff[iCSC]->setAxisTitle("Track Mode", 2);
0603   }
0604   for (int iRPC = 0; iRPC < 8; iRPC++) {
0605     emtfTrackModeVsRPCBXDiff[iRPC] = ibooker.book2D("emtfTrackModeVsRPCBXDiff" + nameRPCStation[iRPC],
0606                                                     "EMTF Track Mode vs (Track BX - RPC BX) " + labelRPCStation[iRPC],
0607                                                     9,
0608                                                     -4,
0609                                                     5,
0610                                                     16,
0611                                                     0,
0612                                                     16);
0613     emtfTrackModeVsRPCBXDiff[iRPC]->setAxisTitle("Track BX - RPC BX", 1);
0614     emtfTrackModeVsRPCBXDiff[iRPC]->setAxisTitle("Track Mode", 2);
0615   }
0616 
0617   // GEM vs CSC Dec 06 2020
0618   ibooker.setCurrentFolder(monitorDir + "/GEMVsCSC");
0619   for (int hist = 0; hist < 2; hist++) {
0620     gemHitPhi[hist] = ibooker.book2D(
0621         "gemHitPhi" + nameGEMStation[hist], "GEM Hit Phi " + labelGEMStation[hist], 4921, 0, 4921, 6, 1, 7);
0622     gemHitTheta[hist] = ibooker.book2D(
0623         "gemHitTheta" + nameGEMStation[hist], "GEM Hit Theta " + labelGEMStation[hist], 128, 0, 128, 6, 1, 7);
0624     gemHitVScscLCTPhi[hist] = ibooker.book2D("gemHitVScscLCTPhi" + nameGEMStation[hist],
0625                                              "GEM Hit Phi - CSC LCT Phi " + labelGEMStation[hist],
0626                                              1200,
0627                                              -600,
0628                                              600,
0629                                              36,
0630                                              1,
0631                                              37);  // one chamber is 10 degrees, 60 integer phi per degree
0632     gemHitVScscLCTTheta[hist] = ibooker.book2D("gemHitVScscLCTTheta" + nameGEMStation[hist],
0633                                                "GEM Hit Theta - CSC LCT Theta " + labelGEMStation[hist],
0634                                                20,
0635                                                -10,
0636                                                10,
0637                                                36,
0638                                                1,
0639                                                37);  // 0.1 eta is at most 9.5 integer theta (between eta 1.5 and 1.6)
0640     gemHitVScscLCTBX[hist] = ibooker.book2D("gemHitVScscLCTBX" + nameGEMStation[hist],
0641                                             "GEM Hit BX - CSC LCT BX " + labelGEMStation[hist],
0642                                             9,
0643                                             -4,
0644                                             5,
0645                                             36,
0646                                             1,
0647                                             37);
0648 
0649     gemHitPhi[hist]->setAxisTitle("Integer #phi", 1);
0650     gemHitTheta[hist]->setAxisTitle("Integer #theta", 1);
0651     gemHitVScscLCTPhi[hist]->setAxisTitle("Integer #phi", 1);
0652     gemHitVScscLCTTheta[hist]->setAxisTitle("Integer #theta", 1);
0653     gemHitVScscLCTBX[hist]->setAxisTitle("GEM BX - CSC BX", 1);
0654 
0655     gemHitPhi[hist]->setAxisTitle("Sector", 2);
0656     gemHitTheta[hist]->setAxisTitle("Sector", 2);
0657     gemHitVScscLCTPhi[hist]->setAxisTitle("Chamber", 2);
0658     gemHitVScscLCTTheta[hist]->setAxisTitle("Chamber", 2);
0659     gemHitVScscLCTBX[hist]->setAxisTitle("Chamber", 2);
0660   }
0661 
0662   // Muon Cand
0663   ibooker.setCurrentFolder(monitorDir + "/MuonCand");
0664 
0665   // Regional Muon Candidate Monitor Elements
0666   emtfMuonBX = ibooker.book1D("emtfMuonBX", "EMTF Muon Cand BX", 7, -3, 4);
0667   emtfMuonBX->setAxisTitle("BX", 1);
0668   for (int xbin = 1, bin_label = -3; xbin <= 7; ++xbin, ++bin_label) {
0669     emtfMuonBX->setBinLabel(xbin, std::to_string(bin_label), 1);
0670   }
0671 
0672   emtfMuonhwPt = ibooker.book1D("emtfMuonhwPt", "EMTF Muon Cand p_{T}", 512, 0, 512);
0673   emtfMuonhwPt->setAxisTitle("Hardware p_{T}", 1);
0674 
0675   emtfMuonhwEta = ibooker.book1D("emtfMuonhwEta", "EMTF Muon Cand #eta", 460, -230, 230);
0676   emtfMuonhwEta->setAxisTitle("Hardware #eta", 1);
0677 
0678   emtfMuonhwPhi = ibooker.book1D("emtfMuonhwPhi", "EMTF Muon Cand #phi", 145, -40, 105);
0679   emtfMuonhwPhi->setAxisTitle("Hardware #phi", 1);
0680 
0681   emtfMuonhwQual = ibooker.book1D("emtfMuonhwQual", "EMTF Muon Cand Quality", 16, 0, 16);
0682   emtfMuonhwQual->setAxisTitle("Quality", 1);
0683   for (int xbin = 1; xbin <= 16; ++xbin) {
0684     emtfMuonhwQual->setBinLabel(xbin, std::to_string(xbin - 1), 1);
0685   }
0686 }
0687 
0688 // CSCOccupancy chamber mapping for neighbor inclusive plots
0689 int chamber_bin(int station, int ring, int chamber) {
0690   int chamber_bin_index = 0;
0691   if (station > 1 && (ring % 2) == 1) {
0692     chamber_bin_index = (chamber * 2) + ((chamber + 1) / 3);
0693   } else {
0694     chamber_bin_index = chamber + ((chamber + 3) / 6);
0695   }
0696   return chamber_bin_index;
0697 };
0698 
0699 void L1TStage2EMTF::analyze(const edm::Event& e, const edm::EventSetup& c) {
0700   if (verbose)
0701     edm::LogInfo("L1TStage2EMTF") << "L1TStage2EMTF: analyze..." << std::endl;
0702 
0703   // DAQ Output
0704   edm::Handle<l1t::EMTFDaqOutCollection> DaqOutCollection;
0705   e.getByToken(daqToken, DaqOutCollection);
0706 
0707   for (auto DaqOut = DaqOutCollection->begin(); DaqOut != DaqOutCollection->end(); ++DaqOut) {
0708     const l1t::emtf::MECollection* MECollection = DaqOut->PtrMECollection();
0709     for (auto ME = MECollection->begin(); ME != MECollection->end(); ++ME) {
0710       if (ME->SE())
0711         emtfErrors->Fill(1);
0712       if (ME->SM())
0713         emtfErrors->Fill(2);
0714       if (ME->BXE())
0715         emtfErrors->Fill(3);
0716       if (ME->AF())
0717         emtfErrors->Fill(4);
0718     }
0719 
0720     const l1t::emtf::EventHeader* EventHeader = DaqOut->PtrEventHeader();
0721     if (!EventHeader->Rdy())
0722       emtfErrors->Fill(5);
0723 
0724     // Fill MPC input link errors
0725     int offset = (EventHeader->Sector() - 1) * 9;
0726     int endcap = EventHeader->Endcap();
0727     l1t::emtf::Counters CO = DaqOut->GetCounters();
0728     const std::array<std::array<int, 9>, 5> counters{
0729         {{{CO.ME1a_1(),
0730            CO.ME1a_2(),
0731            CO.ME1a_3(),
0732            CO.ME1a_4(),
0733            CO.ME1a_5(),
0734            CO.ME1a_6(),
0735            CO.ME1a_7(),
0736            CO.ME1a_8(),
0737            CO.ME1a_9()}},
0738          {{CO.ME1b_1(),
0739            CO.ME1b_2(),
0740            CO.ME1b_3(),
0741            CO.ME1b_4(),
0742            CO.ME1b_5(),
0743            CO.ME1b_6(),
0744            CO.ME1b_7(),
0745            CO.ME1b_8(),
0746            CO.ME1b_9()}},
0747          {{CO.ME2_1(), CO.ME2_2(), CO.ME2_3(), CO.ME2_4(), CO.ME2_5(), CO.ME2_6(), CO.ME2_7(), CO.ME2_8(), CO.ME2_9()}},
0748          {{CO.ME3_1(), CO.ME3_2(), CO.ME3_3(), CO.ME3_4(), CO.ME3_5(), CO.ME3_6(), CO.ME3_7(), CO.ME3_8(), CO.ME3_9()}},
0749          {{CO.ME4_1(), CO.ME4_2(), CO.ME4_3(), CO.ME4_4(), CO.ME4_5(), CO.ME4_6(), CO.ME4_7(), CO.ME4_8(), CO.ME4_9()}}}};
0750     for (int i = 0; i < 5; i++) {
0751       for (int j = 0; j < 9; j++) {
0752         if (counters.at(i).at(j) != 0)
0753           mpcLinkErrors->Fill(j + 1 + offset, endcap * (i + 0.5), counters.at(i).at(j));
0754         else
0755           mpcLinkGood->Fill(j + 1 + offset, endcap * (i + 0.5));
0756       }
0757     }
0758     if (CO.ME1n_3() == 1)
0759       mpcLinkErrors->Fill(1 + offset, endcap * 5.5);
0760     if (CO.ME1n_6() == 1)
0761       mpcLinkErrors->Fill(2 + offset, endcap * 5.5);
0762     if (CO.ME1n_9() == 1)
0763       mpcLinkErrors->Fill(3 + offset, endcap * 5.5);
0764     if (CO.ME2n_3() == 1)
0765       mpcLinkErrors->Fill(4 + offset, endcap * 5.5);
0766     if (CO.ME2n_9() == 1)
0767       mpcLinkErrors->Fill(5 + offset, endcap * 5.5);
0768     if (CO.ME3n_3() == 1)
0769       mpcLinkErrors->Fill(6 + offset, endcap * 5.5);
0770     if (CO.ME3n_9() == 1)
0771       mpcLinkErrors->Fill(7 + offset, endcap * 5.5);
0772     if (CO.ME4n_3() == 1)
0773       mpcLinkErrors->Fill(8 + offset, endcap * 5.5);
0774     if (CO.ME4n_9() == 1)
0775       mpcLinkErrors->Fill(9 + offset, endcap * 5.5);
0776     if (CO.ME1n_3() == 0)
0777       mpcLinkGood->Fill(1 + offset, endcap * 5.5);
0778     if (CO.ME1n_6() == 0)
0779       mpcLinkGood->Fill(2 + offset, endcap * 5.5);
0780     if (CO.ME1n_9() == 0)
0781       mpcLinkGood->Fill(3 + offset, endcap * 5.5);
0782     if (CO.ME2n_3() == 0)
0783       mpcLinkGood->Fill(4 + offset, endcap * 5.5);
0784     if (CO.ME2n_9() == 0)
0785       mpcLinkGood->Fill(5 + offset, endcap * 5.5);
0786     if (CO.ME3n_3() == 0)
0787       mpcLinkGood->Fill(6 + offset, endcap * 5.5);
0788     if (CO.ME3n_9() == 0)
0789       mpcLinkGood->Fill(7 + offset, endcap * 5.5);
0790     if (CO.ME4n_3() == 0)
0791       mpcLinkGood->Fill(8 + offset, endcap * 5.5);
0792     if (CO.ME4n_9() == 0)
0793       mpcLinkGood->Fill(9 + offset, endcap * 5.5);
0794   }
0795 
0796   // Hits (CSC LCTs and RPC hits)
0797   edm::Handle<l1t::EMTFHitCollection> HitCollection;
0798   e.getByToken(hitToken, HitCollection);
0799 
0800   // Maps CSC station and ring to the monitor element index and uses symmetry of the endcaps
0801   const std::map<std::pair<int, int>, int> histIndexCSC = {{{1, 4}, 9},
0802                                                            {{1, 1}, 8},
0803                                                            {{1, 2}, 7},
0804                                                            {{1, 3}, 6},
0805                                                            {{2, 1}, 5},
0806                                                            {{2, 2}, 4},
0807                                                            {{3, 1}, 3},
0808                                                            {{3, 2}, 2},
0809                                                            {{4, 1}, 1},
0810                                                            {{4, 2}, 0}};
0811 
0812   // Maps CSC BX from -2 to 2 to monitor element cscLCTTIming
0813   const std::map<int, int> histIndexBX = {{0, 4}, {-1, 0}, {1, 1}, {-2, 2}, {2, 3}};
0814 
0815   // Maps RPC station and ring to the monitor element index and uses symmetry of the endcaps
0816   const std::map<std::pair<int, int>, int> histIndexRPC = {
0817       {{4, 3}, 0}, {{4, 2}, 1}, {{3, 3}, 2}, {{3, 2}, 3}, {{2, 2}, 4}, {{1, 2}, 5}};
0818 
0819   // Reverse 'rotate by 2' for RPC subsector
0820   auto get_subsector_rpc_cppf = [](int subsector_rpc) { return ((subsector_rpc + 3) % 6) + 1; };
0821 
0822   for (auto Hit = HitCollection->begin(); Hit != HitCollection->end(); ++Hit) {
0823     int endcap = Hit->Endcap();
0824     int sector = Hit->Sector();
0825     int station = Hit->Station();
0826     int ring = Hit->Ring();
0827     int cscid = Hit->CSC_ID();
0828     int chamber = Hit->Chamber();
0829     int strip = Hit->Strip();
0830     int wire = Hit->Wire();
0831     int cscid_offset = (sector - 1) * 9;
0832 
0833     int hist_index = 0;
0834     if (ring == 4 && strip >= 128)
0835       strip -= 128;
0836 
0837     if (Hit->Is_CSC() == true) {
0838       hist_index = histIndexCSC.at({station, ring});
0839       if (endcap > 0)
0840         hist_index = 19 - hist_index;
0841       cscLCTBX->Fill(Hit->BX(), hist_index);
0842       float evt_wgt = (Hit->Station() > 1 && Hit->Ring() == 1) ? 0.5 : 1.0;
0843       if (Hit->Neighbor() == false) {
0844         //Map for cscDQMOccupancy plot
0845         cscDQMOccupancy->Fill(chamber_bin(station, ring, chamber), hist_index, evt_wgt);
0846         if (station > 1 && (ring % 2) == 1) {
0847           cscDQMOccupancy->Fill(chamber_bin(station, ring, chamber) - 1, hist_index, evt_wgt);
0848         }
0849         cscLCTStrip[hist_index]->Fill(strip);
0850         cscLCTWire[hist_index]->Fill(wire);
0851         cscChamberStrip[hist_index]->Fill(chamber, strip);
0852         cscChamberWire[hist_index]->Fill(chamber, wire);
0853         if (Hit->Subsector() == 1) {
0854           cscLCTOccupancy->Fill(cscid + cscid_offset, endcap * (station - 0.5));
0855         } else {
0856           cscLCTOccupancy->Fill(cscid + cscid_offset, endcap * (station + 0.5));
0857         }
0858       } else {
0859         // Map neighbor chambers to "fake" CSC IDs: 1/3 --> 1, 1/6 --> 2, 1/9 --> 3, 2/3 --> 4, 2/9 --> 5, etc.
0860         int cscid_n = (station == 1 ? (cscid / 3) : (station * 2) + ((cscid - 3) / 6));
0861         cscLCTOccupancy->Fill(cscid_n + cscid_offset, endcap * 5.5);
0862       }
0863       if (Hit->Neighbor() == true) {
0864         cscDQMOccupancy->Fill((sector % 6 + 1) * 7 - 4, hist_index, evt_wgt);
0865       }
0866     }
0867 
0868     if (Hit->Is_RPC() == true) {
0869       hist_index = histIndexRPC.at({station, ring});
0870       if (endcap > 0)
0871         hist_index = 11 - hist_index;
0872 
0873       rpcHitBX->Fill(Hit->BX(), hist_index);
0874 
0875       if (Hit->Neighbor() == false) {
0876         rpcHitPhi[hist_index]->Fill(Hit->Phi_fp() / 4);
0877         rpcHitTheta[hist_index]->Fill(Hit->Theta_fp() / 4);
0878         rpcChamberPhi[hist_index]->Fill(chamber, Hit->Phi_fp() / 4);
0879         rpcChamberTheta[hist_index]->Fill(chamber, Hit->Theta_fp() / 4);
0880         rpcHitOccupancy->Fill((Hit->Sector_RPC() - 1) * 7 + get_subsector_rpc_cppf(Hit->Subsector_RPC()),
0881                               hist_index + 0.5);
0882       } else if (Hit->Neighbor() == true) {
0883         rpcHitOccupancy->Fill((Hit->Sector_RPC() - 1) * 7 + 7, hist_index + 0.5);
0884       }
0885     }
0886 
0887     // Add GEM Oct 27 2020
0888     hitTypeBX->Fill(4, Hit->BX());
0889     if (Hit->Is_CSC() == true)
0890       hitTypeBX->Fill(1, Hit->BX());
0891     else if (Hit->Is_RPC() == true)
0892       hitTypeBX->Fill(2, Hit->BX());
0893     else if (Hit->Is_GEM() == true)
0894       hitTypeBX->Fill(3, Hit->BX());
0895 
0896     if (Hit->Is_GEM() == true) {
0897       gemHitBX->Fill(Hit->BX(), (endcap > 0) ? 1.5 : 0.5);
0898       hist_index = (endcap > 0) ? 1 : 0;
0899       if (Hit->Neighbor() == false) {
0900         gemChamberPad[hist_index]->Fill(chamber, Hit->Pad());
0901         gemChamberPartition[hist_index]->Fill(chamber, Hit->Partition());
0902         gemHitOccupancy->Fill(chamber_bin(1, 1, chamber), (endcap > 0) ? 1.5 : 0.5);  // follow CSC convention
0903       } else {
0904         gemChamberPad[hist_index]->Fill((Hit->Sector() % 6) * 6 + 2, Hit->Pad());
0905         gemChamberPartition[hist_index]->Fill((Hit->Sector() % 6) * 6 + 2, Hit->Partition());
0906         gemHitOccupancy->Fill((Hit->Sector() % 6 + 1) * 7 - 4, (endcap > 0) ? 1.5 : 0.5);  // follow CSC convention
0907       }
0908     }  // End of if (Hit->Is_GEM() == true)
0909   }    // End of for (auto Hit = HitCollection->begin(); Hit != HitCollection->end(); ++Hit)
0910 
0911   // Tracks
0912   edm::Handle<l1t::EMTFTrackCollection> TrackCollection;
0913   e.getByToken(trackToken, TrackCollection);
0914 
0915   int nTracks = TrackCollection->size();
0916 
0917   emtfnTracks->Fill(std::min(nTracks, emtfnTracks->getTH1F()->GetNbinsX() - 1));
0918 
0919   for (auto Track = TrackCollection->begin(); Track != TrackCollection->end(); ++Track) {
0920     int endcap = Track->Endcap();
0921     int sector = Track->Sector();
0922     float eta = Track->Eta();
0923     float phi_glob_rad = Track->Phi_glob() * M_PI / 180.;
0924     int mode = Track->Mode();
0925     int quality = Track->GMT_quality();
0926     int numHits = Track->NumHits();
0927     int modeNeighbor = Track->Mode_neighbor();
0928     int modeRPC = Track->Mode_RPC();
0929     int singleMuQuality = 12;
0930     int singleMuPT = 22;
0931 
0932     // Only plot if there are <= 1 neighbor hits in the track to avoid spikes at sector boundaries
0933     if (modeNeighbor >= 2 && modeNeighbor != 4 && modeNeighbor != 8)
0934       continue;
0935 
0936     emtfTracknHits->Fill(numHits);
0937     emtfTrackBX->Fill(endcap * (sector - 0.5), Track->BX());
0938     emtfTrackPt->Fill(Track->Pt());
0939     emtfTrackEta->Fill(eta);
0940 
0941     emtfTrackOccupancy->Fill(eta, phi_glob_rad);
0942     emtfTrackMode->Fill(mode);
0943     emtfTrackQuality->Fill(quality);
0944     emtfTrackQualityVsMode->Fill(mode, quality);
0945     RPCvsEMTFTrackMode->Fill(mode, modeRPC);
0946     emtfTrackPhi->Fill(phi_glob_rad);
0947 
0948     if (quality >= singleMuQuality) {
0949       emtfTrackPtHighQuality->Fill(Track->Pt());
0950       emtfTrackEtaHighQuality->Fill(eta);
0951       emtfTrackPhiHighQuality->Fill(phi_glob_rad);
0952       emtfTrackOccupancyHighQuality->Fill(eta, phi_glob_rad);
0953       if (Track->Pt() >= singleMuPT) {
0954         emtfTrackPtHighQualityHighPT->Fill(Track->Pt());
0955         emtfTrackEtaHighQualityHighPT->Fill(eta);
0956         emtfTrackPhiHighQualityHighPT->Fill(phi_glob_rad);
0957         emtfTrackOccupancyHighQualityHighPT->Fill(eta, phi_glob_rad);
0958       }
0959     }
0960 
0961     ////////////////////////////////////////////////////
0962     ///  Begin block for CSC LCT and RPC hit timing  ///
0963     ////////////////////////////////////////////////////
0964     {
0965       // LCT and RPC Timing
0966       if (numHits < 2 || numHits > 4)
0967         continue;
0968       int numHitsInTrack_BX0 = 0;
0969       unsigned int hist_index2 = 4 - numHits;
0970 
0971       for (const auto& iTrkHit : Track->Hits()) {
0972         if (iTrkHit.Is_CSC()) {
0973           emtfTrackBXVsCSCLCT[hist_index2]->Fill(iTrkHit.BX(), Track->BX());
0974           int iCSC = (endcap > 0) ? (iTrkHit.Station() + 3) : (4 - iTrkHit.Station());
0975           emtfTrackModeVsCSCBXDiff[iCSC]->Fill(Track->BX() - iTrkHit.BX(),
0976                                                mode);  // Add mode vs BXdiff comparison Dec 07 2020
0977         } else if (iTrkHit.Is_RPC()) {
0978           emtfTrackBXVsRPCHit[hist_index2]->Fill(iTrkHit.BX(), Track->BX());
0979           int iRPC = (endcap > 0) ? (iTrkHit.Station() + 3) : (4 - iTrkHit.Station());
0980           emtfTrackModeVsRPCBXDiff[iRPC]->Fill(Track->BX() - iTrkHit.BX(),
0981                                                mode);  // Add mode vs BXdiff comparison Dec 07 2020
0982         } else if (iTrkHit.Is_GEM()) {
0983           emtfTrackBXVsGEMHit[hist_index2]->Fill(iTrkHit.BX(), Track->BX());
0984           int iGEM = (endcap > 0) ? 1 : 0;
0985           emtfTrackModeVsGEMBXDiff[iGEM]->Fill(Track->BX() - iTrkHit.BX(),
0986                                                mode);  // Add mode vs BXdiff comparison Dec 07 2020
0987         }
0988       }
0989 
0990       // Select well-timed tracks: >= 3 hits, with <= 1 in BX != 0
0991       if (numHits < 3)
0992         continue;
0993       for (const auto& jTrkHit : Track->Hits()) {
0994         if (jTrkHit.BX() == 0)
0995           numHitsInTrack_BX0++;
0996       }
0997       if (numHitsInTrack_BX0 < numHits - 1)
0998         continue;
0999 
1000       for (const auto& TrkHit : Track->Hits()) {
1001         int trackHitBX = TrkHit.BX();
1002         if (std::abs(trackHitBX) > 2)
1003           continue;  // Should never happen, but just to be safe ...
1004         //int cscid        = TrkHit.CSC_ID();
1005         int ring = TrkHit.Ring();
1006         int station = TrkHit.Station();
1007         int sector = TrkHit.Sector();
1008         //int cscid_offset = (sector - 1) * 9;//no longer needed after new time plots (maybe useful for future plots)
1009         int neighbor = TrkHit.Neighbor();
1010         int endcap = TrkHit.Endcap();
1011         int chamber = TrkHit.Chamber();
1012 
1013         int hist_index = 0;
1014         float evt_wgt = (TrkHit.Station() > 1 && TrkHit.Ring() == 1) ? 0.5 : 1.0;
1015 
1016         if (TrkHit.Is_CSC() == true) {
1017           hist_index = histIndexCSC.at({station, ring});
1018           if (endcap > 0)
1019             hist_index = 19 - hist_index;
1020           if (neighbor == false) {
1021             cscLCTTiming[histIndexBX.at(trackHitBX)]->Fill(chamber_bin(station, ring, chamber), hist_index, evt_wgt);
1022             cscTimingTot->Fill(chamber_bin(station, ring, chamber), hist_index, evt_wgt);
1023             if (station > 1 && (ring % 2) == 1) {
1024               cscLCTTiming[histIndexBX.at(trackHitBX)]->Fill(
1025                   chamber_bin(station, ring, chamber) - 1, hist_index, evt_wgt);
1026               cscTimingTot->Fill(chamber_bin(station, ring, chamber) - 1, hist_index, evt_wgt);
1027             }
1028           } else {
1029             // Map neighbor chambers to "fake" CSC IDs: 1/3 --> 1, 1/6 --> 2, 1/9 --> 3, 2/3 --> 4, 2/9 --> 5, etc.
1030             //int cscid_n = (station == 1 ? (cscid / 3) : (station * 2) + ((cscid - 3) / 6) );
1031             cscLCTTiming[histIndexBX.at(trackHitBX)]->Fill((sector % 6 + 1) * 7 - 4, hist_index, evt_wgt);
1032             cscTimingTot->Fill((sector % 6 + 1) * 7 - 4, hist_index, evt_wgt);
1033           }
1034 
1035           // Fill RPC timing with matched CSC LCTs
1036           if (trackHitBX == 0 && ring == 2) {
1037             for (auto Hit = HitCollection->begin(); Hit != HitCollection->end(); ++Hit) {
1038               if (Hit->Is_RPC() == false || neighbor == true)
1039                 continue;
1040               if (std::abs(Track->Eta() - Hit->Eta()) > 0.1)
1041                 continue;
1042               if (Hit->Endcap() != endcap || Hit->Station() != station || Hit->Chamber() != chamber)
1043                 continue;
1044               if (std::abs(Hit->BX()) > 2)
1045                 continue;
1046 
1047               hist_index = histIndexRPC.at({Hit->Station(), Hit->Ring()});
1048               if (Hit->Endcap() > 0)
1049                 hist_index = 11 - hist_index;
1050               rpcHitTimingInTrack->Fill(Hit->BX(), hist_index + 0.5);
1051               rpcHitTiming[histIndexBX.at(Hit->BX())]->Fill(
1052                   (Hit->Sector_RPC() - 1) * 7 + get_subsector_rpc_cppf(Hit->Subsector_RPC()), hist_index + 0.5);
1053               rpcHitTimingTot->Fill((Hit->Sector_RPC() - 1) * 7 + get_subsector_rpc_cppf(Hit->Subsector_RPC()),
1054                                     hist_index + 0.5);
1055             }  // End loop: for (auto Hit = HitCollection->begin(); Hit != HitCollection->end(); ++Hit)
1056           }    // End conditional: if (trackHitBX == 0 && ring == 2)
1057 
1058           // Fill GEM timing with matched CSC LCTs
1059           if (trackHitBX == 0 && station == 1 && ring == 1) {  // GEM only in station 1
1060             for (auto Hit = HitCollection->begin(); Hit != HitCollection->end(); ++Hit) {
1061               if (Hit->Is_GEM() == false)
1062                 continue;
1063               if (std::abs(Track->Eta() - Hit->Eta()) > 0.1)
1064                 continue;
1065               if (Hit->Endcap() != endcap || Hit->Station() != station || Hit->Chamber() != chamber ||
1066                   Hit->Neighbor() != neighbor)  //different neighbor requirement from RPC
1067                 continue;
1068               if (std::abs(Hit->BX()) > 2)
1069                 continue;
1070 
1071               if (neighbor == false) {
1072                 gemHitTiming[histIndexBX.at(Hit->BX())]->Fill(chamber_bin(1, 1, chamber), (endcap > 0) ? 1.5 : 0.5);
1073                 gemHitTimingTot->Fill(chamber_bin(1, 1, chamber), (endcap > 0) ? 1.5 : 0.5);
1074                 int ihist = (endcap > 0) ? 1 : 0;
1075                 gemHitPhi[ihist]->Fill(Hit->Phi_fp(), sector);
1076                 gemHitTheta[ihist]->Fill(Hit->Theta_fp(), sector);
1077                 gemHitVScscLCTPhi[ihist]->Fill(Hit->Phi_fp() - TrkHit.Phi_fp(), chamber);  // GEM vs CSC Dec 06 2020
1078                 gemHitVScscLCTTheta[ihist]->Fill(Hit->Theta_fp() - TrkHit.Theta_fp(), chamber);
1079                 gemHitVScscLCTBX[ihist]->Fill(Hit->BX() - TrkHit.BX(), chamber);
1080               } else {
1081                 gemHitTiming[histIndexBX.at(Hit->BX())]->Fill((sector % 6 + 1) * 7 - 4, (endcap > 0) ? 1.5 : 0.5);
1082                 gemHitTimingTot->Fill((sector % 6 + 1) * 7 - 4, (endcap > 0) ? 1.5 : 0.5);
1083               }
1084 
1085             }  // End loop: for (auto Hit = HitCollection->begin(); Hit != HitCollection->end(); ++Hit)
1086           }    // End conditional: if (trackHitBX == 0 && station == 1 && ring == 1)
1087         }      // End conditional: if (TrkHit.Is_CSC() == true)
1088 
1089         if (TrkHit.Is_RPC() == true && neighbor == false) {
1090           hist_index = histIndexRPC.at({station, ring});
1091           if (endcap > 0)
1092             hist_index = 11 - hist_index;
1093 
1094           rpcHitTimingInTrack->Fill(trackHitBX, hist_index + 0.5);
1095           rpcHitTiming[histIndexBX.at(trackHitBX)]->Fill(
1096               (TrkHit.Sector_RPC() - 1) * 7 + get_subsector_rpc_cppf(TrkHit.Subsector_RPC()), hist_index + 0.5);
1097           rpcHitTimingTot->Fill((TrkHit.Sector_RPC() - 1) * 7 + get_subsector_rpc_cppf(TrkHit.Subsector_RPC()),
1098                                 hist_index + 0.5);
1099         }  // End conditional: if (TrkHit.Is_RPC() == true && neighbor == false)
1100         if (TrkHit.Is_RPC() == true && neighbor == true) {
1101           hist_index = histIndexRPC.at({station, ring});
1102           if (endcap > 0)
1103             hist_index = 11 - hist_index;
1104           rpcHitTiming[histIndexBX.at(trackHitBX)]->Fill((TrkHit.Sector_RPC() - 1) * 7, hist_index + 0.5);
1105         }
1106         // Add GEM Timing Oct 27 2020
1107         if (TrkHit.Is_GEM() == true) {
1108           if (neighbor == false) {
1109             gemHitTiming[histIndexBX.at(trackHitBX)]->Fill(chamber_bin(1, 1, chamber), (endcap > 0) ? 1.5 : 0.5);
1110             gemHitTimingTot->Fill(chamber_bin(1, 1, chamber), (endcap > 0) ? 1.5 : 0.5);
1111           } else {
1112             gemHitTiming[histIndexBX.at(trackHitBX)]->Fill((sector % 6 + 1) * 7 - 4, (endcap > 0) ? 1.5 : 0.5);
1113             gemHitTimingTot->Fill((sector % 6 + 1) * 7 - 4, (endcap > 0) ? 1.5 : 0.5);
1114           }
1115         }  // End condition: if (TrkHit.Is_GEM() == true)
1116       }    // End loop: for (int iHit = 0; iHit < numHits; ++iHit)
1117     }
1118     //////////////////////////////////////////////////
1119     ///  End block for CSC LCT and RPC hit timing  ///
1120     //////////////////////////////////////////////////
1121 
1122   }  // End loop: for (auto Track = TrackCollection->begin(); Track != TrackCollection->end(); ++Track)
1123 
1124   // CSC LCT and RPC Hit Timing Efficieny
1125   for (int hist_index = 0; hist_index < 5; ++hist_index) {
1126     cscLCTTimingFrac[hist_index]->getTH2F()->Divide(cscLCTTiming[hist_index]->getTH2F(), cscTimingTot->getTH2F());
1127     rpcHitTimingFrac[hist_index]->getTH2F()->Divide(rpcHitTiming[hist_index]->getTH2F(), rpcHitTimingTot->getTH2F());
1128     gemHitTimingFrac[hist_index]->getTH2F()->Divide(gemHitTiming[hist_index]->getTH2F(), gemHitTimingTot->getTH2F());
1129   }
1130 
1131   // Regional Muon Candidates
1132   edm::Handle<l1t::RegionalMuonCandBxCollection> MuonBxCollection;
1133   e.getByToken(muonToken, MuonBxCollection);
1134 
1135   for (int itBX = MuonBxCollection->getFirstBX(); itBX <= MuonBxCollection->getLastBX(); ++itBX) {
1136     for (l1t::RegionalMuonCandBxCollection::const_iterator Muon = MuonBxCollection->begin(itBX);
1137          Muon != MuonBxCollection->end(itBX);
1138          ++Muon) {
1139       emtfMuonBX->Fill(itBX);
1140       emtfMuonhwPt->Fill(Muon->hwPt());
1141       emtfMuonhwEta->Fill(Muon->hwEta());
1142       emtfMuonhwPhi->Fill(Muon->hwPhi());
1143       emtfMuonhwQual->Fill(Muon->hwQual());
1144     }
1145   }
1146 }