Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 09:55:50

0001 #include "L1Trigger/Phase2L1GMT/interface/L1TPhase2GMTEndcapStubProcessor.h"
0002 #include "L1Trigger/L1TMuon/interface/MuonTriggerPrimitive.h"
0003 #include <cmath>
0004 #include <iostream>
0005 #include <string>
0006 #include <sstream>
0007 
0008 L1TPhase2GMTEndcapStubProcessor::L1TPhase2GMTEndcapStubProcessor() : minBX_(-3), maxBX_(3) {}
0009 
0010 L1TPhase2GMTEndcapStubProcessor::L1TPhase2GMTEndcapStubProcessor(const edm::ParameterSet& iConfig)
0011     : minBX_(iConfig.getParameter<int>("minBX")),
0012       maxBX_(iConfig.getParameter<int>("maxBX")),
0013       coord1LSB_(iConfig.getParameter<double>("coord1LSB")),
0014       coord2LSB_(iConfig.getParameter<double>("coord2LSB")),
0015       eta1LSB_(iConfig.getParameter<double>("eta1LSB")),
0016       eta2LSB_(iConfig.getParameter<double>("eta2LSB")),
0017       etaMatch_(iConfig.getParameter<double>("etaMatch")),
0018       phiMatch_(iConfig.getParameter<double>("phiMatch")),
0019       verbose_(iConfig.getParameter<unsigned int>("verbose")) {}
0020 
0021 L1TPhase2GMTEndcapStubProcessor::~L1TPhase2GMTEndcapStubProcessor() {}
0022 
0023 l1t::MuonStub L1TPhase2GMTEndcapStubProcessor::buildCSCOnlyStub(const CSCDetId& detid,
0024                                                                 const CSCCorrelatedLCTDigi& digi,
0025                                                                 const L1TMuon::GeometryTranslator* translator) {
0026   int endcap = detid.endcap();
0027   int station = detid.station();
0028   int chamber = detid.chamber();
0029   int ring = detid.ring();
0030 
0031   L1TMuon::TriggerPrimitive primitive(detid, digi);
0032 
0033   const GlobalPoint& gp = translator->getGlobalPoint(primitive);
0034 
0035   int phi = int(gp.phi().value() / coord1LSB_);
0036   int eta1 = int(gp.eta() / eta1LSB_);
0037 
0038   int wheel = 0;
0039   int sign = endcap == 1 ? -1 : 1;
0040 
0041   if (ring == 3)
0042     wheel = sign * 3;
0043   else if (ring == 2)
0044     wheel = sign * 4;
0045   else if (ring == 1)
0046     wheel = sign * 5;
0047 
0048   int sector = fabs(chamber);
0049 
0050   int bx = digi.getBX() - 8;
0051   int quality = 1;
0052 
0053   uint tfLayer = 0;
0054   if ((ring == 3 || ring == 2) && station == 1)  //ME1/3
0055     tfLayer = 4;
0056   else if (ring == 1 && station == 1)  //ME1/3
0057     tfLayer = 0;
0058   else if (station == 2)  //ME2/2
0059     tfLayer = 2;
0060   else if (station == 3)  //ME3/2
0061     tfLayer = 1;
0062   else if (station == 4)  //ME4/2
0063     tfLayer = 3;
0064 
0065   l1t::MuonStub stub(wheel, sector, station, tfLayer, phi, 0, 0, bx, quality, eta1, 0, 1, 0);
0066 
0067   stub.setOfflineQuantities(gp.phi().value(), 0.0, gp.eta(), 0.0);
0068   return stub;
0069 }
0070 
0071 l1t::MuonStub L1TPhase2GMTEndcapStubProcessor::buildRPCOnlyStub(const RPCDetId& detid,
0072                                                                 const RPCDigi& digi,
0073                                                                 const L1TMuon::GeometryTranslator* translator) {
0074   L1TMuon::TriggerPrimitive primitive(detid, digi);
0075   const GlobalPoint& gp = translator->getGlobalPoint(primitive);
0076 
0077   int phi2 = int(gp.phi().value() / coord2LSB_);
0078   int eta2 = int(gp.eta() / eta2LSB_);
0079 
0080   int wheel = -(6 - detid.ring()) * detid.region();
0081   int sector = (detid.sector() - 1) * 6 + detid.subsector();
0082   int station = detid.station();
0083   bool tag = detid.trIndex();
0084   int bx = digi.bx();
0085   int quality = 2;
0086 
0087   int ring = detid.ring();
0088 
0089   uint tfLayer = 0;
0090 
0091   if ((ring == 3 || ring == 2) && station == 1)  //ME1/3
0092     tfLayer = 4;
0093   //  else if (ring==1 && station==1) //ME1/3
0094   //    tfLayer=0;
0095   else if (station == 2)  //ME2/2
0096     tfLayer = 2;
0097   else if (station == 3)  //ME3/2
0098     tfLayer = 1;
0099   else if (station == 4)  //ME4/2
0100     tfLayer = 3;
0101 
0102   l1t::MuonStub stub(wheel, sector, station, tfLayer, 0, phi2, tag, bx, quality, 0, eta2, 2, 0);
0103   stub.setOfflineQuantities(gp.phi().value(), gp.phi().value(), gp.eta(), gp.eta());
0104   return stub;
0105 }
0106 
0107 l1t::MuonStubCollection L1TPhase2GMTEndcapStubProcessor::combineStubs(const l1t::MuonStubCollection& cscStubs,
0108                                                                       const l1t::MuonStubCollection& rpcStubs) {
0109   l1t::MuonStubCollection out;
0110   l1t::MuonStubCollection usedRPC;
0111   l1t::MuonStubCollection cleanedRPC;
0112 
0113   l1t::MuonStubCollection cleanedCSC;
0114 
0115   //clean ME11 ambiguities
0116   l1t::MuonStubCollection allCSC = cscStubs;
0117 
0118   while (!allCSC.empty()) {
0119     l1t::MuonStub stub = allCSC[0];
0120     l1t::MuonStubCollection freeCSC;
0121     for (uint i = 1; i < allCSC.size(); ++i) {
0122       if ((stub.etaRegion() == allCSC[i].etaRegion()) && (stub.depthRegion() == allCSC[i].depthRegion()) &&
0123           (fabs(deltaPhi(stub.offline_coord1(), allCSC[i].offline_coord1())) < 0.001)) {
0124         if (fabs(stub.offline_eta1() - allCSC[i].offline_eta1()) > 0.001) {
0125           stub.setEta(stub.eta1(), allCSC[i].eta1(), 3);
0126           stub.setOfflineQuantities(stub.offline_coord1(), 0.0, stub.offline_eta1(), allCSC[i].offline_eta1());
0127         }
0128       } else {
0129         freeCSC.push_back(allCSC[i]);
0130       }
0131     }
0132     cleanedCSC.push_back(stub);
0133     allCSC = freeCSC;
0134   }
0135 
0136   for (const auto& csc : cleanedCSC) {
0137     int nRPC = 0;
0138     float phiF = 0.0;
0139     float etaF = 0.0;
0140     int phi = 0;
0141     int eta = 0;
0142     for (const auto& rpc : rpcStubs) {
0143       if (csc.depthRegion() != rpc.depthRegion())
0144         continue;
0145       if (fabs(deltaPhi(csc.offline_coord1(), rpc.offline_coord2())) < phiMatch_ &&
0146           fabs(csc.offline_eta1() - rpc.offline_eta2()) < etaMatch_ && csc.bxNum() == rpc.bxNum()) {
0147         phiF += rpc.offline_coord2();
0148         etaF += rpc.offline_eta2();
0149         phi += rpc.coord2();
0150         eta += rpc.eta2();
0151         nRPC++;
0152         usedRPC.push_back(rpc);
0153       }
0154     }
0155 
0156     int finalRPCPhi = 0;
0157     int finalRPCEta = 0;
0158     double offline_finalRPCPhi = 0;
0159     double offline_finalRPCEta = 0;
0160     if (nRPC != 0) {
0161       finalRPCPhi = phi / nRPC;
0162       finalRPCEta = eta / nRPC;
0163       offline_finalRPCPhi = phiF / nRPC;
0164       offline_finalRPCEta = etaF / nRPC;
0165       l1t::MuonStub stub(csc.etaRegion(),
0166                          csc.phiRegion(),
0167                          csc.depthRegion(),
0168                          csc.tfLayer(),
0169                          csc.coord1(),
0170                          finalRPCPhi,
0171                          0,
0172                          csc.bxNum(),
0173                          3,
0174                          csc.eta1(),
0175                          finalRPCEta,
0176                          3,
0177                          0);
0178       stub.setOfflineQuantities(csc.offline_coord1(), offline_finalRPCPhi, csc.offline_eta1(), offline_finalRPCEta);
0179       out.push_back(stub);
0180     } else {
0181       out.push_back(csc);
0182     }
0183   }
0184 
0185   //clean the RPC from the used ones
0186 
0187   for (const auto& rpc : rpcStubs) {
0188     bool keep = true;
0189     for (const auto& rpc2 : usedRPC) {
0190       if (rpc == rpc2) {
0191         keep = false;
0192         break;
0193       }
0194     }
0195     if (keep)
0196       cleanedRPC.push_back(rpc);
0197   }
0198 
0199   while (!cleanedRPC.empty()) {
0200     l1t::MuonStubCollection freeRPC;
0201 
0202     int nRPC = 1;
0203     float phiF = cleanedRPC[0].offline_coord2();
0204     float etaF = cleanedRPC[0].offline_eta2();
0205     int phi = cleanedRPC[0].coord2();
0206     int eta = cleanedRPC[0].eta2();
0207 
0208     for (unsigned i = 1; i < cleanedRPC.size(); ++i) {
0209       if (fabs(deltaPhi(cleanedRPC[0].offline_coord2(), cleanedRPC[i].offline_coord2())) < phiMatch_ &&
0210           cleanedRPC[0].depthRegion() == cleanedRPC[i].depthRegion() &&
0211           fabs(cleanedRPC[0].offline_eta2() - cleanedRPC[i].offline_eta2()) < etaMatch_ &&
0212           cleanedRPC[0].bxNum() == cleanedRPC[i].bxNum()) {
0213         phiF += cleanedRPC[i].offline_coord2();
0214         etaF += cleanedRPC[i].offline_eta2();
0215         phi += cleanedRPC[i].coord2();
0216         eta += cleanedRPC[i].eta2();
0217         nRPC++;
0218       } else {
0219         freeRPC.push_back(cleanedRPC[i]);
0220       }
0221     }
0222     l1t::MuonStub stub(cleanedRPC[0].etaRegion(),
0223                        cleanedRPC[0].phiRegion(),
0224                        cleanedRPC[0].depthRegion(),
0225                        cleanedRPC[0].tfLayer(),
0226                        0,
0227                        phi / nRPC,
0228                        0,
0229                        cleanedRPC[0].bxNum(),
0230                        2,
0231                        0,
0232                        eta / nRPC,
0233                        2,
0234                        0);
0235     stub.setOfflineQuantities(phiF / nRPC, phiF / nRPC, etaF / nRPC, etaF / nRPC);
0236     out.push_back(stub);
0237     cleanedRPC = freeRPC;
0238   };
0239   return out;
0240 }
0241 
0242 l1t::MuonStubCollection L1TPhase2GMTEndcapStubProcessor::makeStubs(
0243     const MuonDigiCollection<CSCDetId, CSCCorrelatedLCTDigi>& csc,
0244     const MuonDigiCollection<RPCDetId, RPCDigi>& cleaned,
0245     const L1TMuon::GeometryTranslator* t,
0246     const edm::EventSetup& iSetup) {
0247   l1t::MuonStubCollection cscStubs;
0248   auto chamber = csc.begin();
0249   auto chend = csc.end();
0250   for (; chamber != chend; ++chamber) {
0251     auto digi = (*chamber).second.first;
0252     auto dend = (*chamber).second.second;
0253     for (; digi != dend; ++digi) {
0254       l1t::MuonStub stub = buildCSCOnlyStub((*chamber).first, *digi, t);
0255       if (stub.bxNum() >= minBX_ && stub.bxNum() <= maxBX_)
0256         cscStubs.push_back(stub);
0257     }
0258   }
0259 
0260   l1t::MuonStubCollection rpcStubs;
0261 
0262   auto rpcchamber = cleaned.begin();
0263   auto rpcchend = cleaned.end();
0264   for (; rpcchamber != rpcchend; ++rpcchamber) {
0265     if ((*rpcchamber).first.region() == 0)
0266       continue;
0267     auto digi = (*rpcchamber).second.first;
0268     auto dend = (*rpcchamber).second.second;
0269     for (; digi != dend; ++digi) {
0270       l1t::MuonStub stub = buildRPCOnlyStub((*rpcchamber).first, *digi, t);
0271       if (stub.bxNum() >= minBX_ && stub.bxNum() <= maxBX_)
0272         rpcStubs.push_back(stub);
0273     }
0274   }
0275 
0276   l1t::MuonStubCollection combinedStubs = combineStubs(cscStubs, rpcStubs);
0277 
0278   if (verbose_) {
0279     printf("CSC Stubs\n");
0280     for (const auto& stub : cscStubs)
0281       printf(
0282           "CSC Stub bx=%d TF=%d etaRegion=%d phiRegion=%d depthRegion=%d  coord1=%f,%d coord2=%f,%d eta1=%f,%d "
0283           "eta2=%f,%d quality=%d etaQuality=%d\n",
0284           stub.bxNum(),
0285           stub.tfLayer(),
0286           stub.etaRegion(),
0287           stub.phiRegion(),
0288           stub.depthRegion(),
0289           stub.offline_coord1(),
0290           stub.coord1(),
0291           stub.offline_coord2(),
0292           stub.coord2(),
0293           stub.offline_eta1(),
0294           stub.eta1(),
0295           stub.offline_eta2(),
0296           stub.eta2(),
0297           stub.quality(),
0298           stub.etaQuality());
0299     printf("RPC Stubs\n");
0300     for (const auto& stub : rpcStubs)
0301       printf(
0302           "RPC Stub bx=%d TF=%d etaRegion=%d phiRegion=%d depthRegion=%d  coord1=%f,%d coord2=%f,%d eta1=%f,%d "
0303           "eta2=%f,%d quality=%d etaQuality=%d\n",
0304           stub.bxNum(),
0305           stub.tfLayer(),
0306           stub.etaRegion(),
0307           stub.phiRegion(),
0308           stub.depthRegion(),
0309           stub.offline_coord1(),
0310           stub.coord1(),
0311           stub.offline_coord2(),
0312           stub.coord2(),
0313           stub.offline_eta1(),
0314           stub.eta1(),
0315           stub.offline_eta2(),
0316           stub.eta2(),
0317           stub.quality(),
0318           stub.etaQuality());
0319     for (const auto& stub : combinedStubs)
0320       printf(
0321           "Combined Stub bx=%d TF=%d etaRegion=%d phiRegion=%d depthRegion=%d  coord1=%f,%d coord2=%f,%d eta1=%f,%d "
0322           "eta2=%f,%d quality=%d etaQuality=%d\n",
0323           stub.bxNum(),
0324           stub.tfLayer(),
0325           stub.etaRegion(),
0326           stub.phiRegion(),
0327           stub.depthRegion(),
0328           stub.offline_coord1(),
0329           stub.coord1(),
0330           stub.offline_coord2(),
0331           stub.coord2(),
0332           stub.offline_eta1(),
0333           stub.eta1(),
0334           stub.offline_eta2(),
0335           stub.eta2(),
0336           stub.quality(),
0337           stub.etaQuality());
0338   }
0339 
0340   return combinedStubs;
0341 }