L1TMuonEndCap/src/EMTFSubsystemCollector.cc

0001 #include "L1Trigger/L1TMuonEndCap/interface/EMTFSubsystemCollector.h"
0002 #include "DataFormats/L1TMuon/interface/L1TMuonSubsystems.h"
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "DataFormats/Common/interface/Handle.h"
0005
0006 #include "Geometry/RPCGeometry/interface/RPCGeometry.h"  // needed to handle RPCRecHit
0007
0008 #include "helper.h"  // adjacent_cluster
0009
0010 // _____________________________________________________________________________
0011 // Specialized for DT
0012 template <>
0013 void EMTFSubsystemCollector::extractPrimitives(emtf::DTTag tag,
0014                                                const GeometryTranslator* tp_geom,
0015                                                const edm::Event& iEvent,
0016                                                const edm::EDGetToken& token1,
0017                                                const edm::EDGetToken& token2,
0018                                                TriggerPrimitiveCollection& out) const {
0019   edm::Handle<emtf::DTTag::digi_collection> phiContainer;
0020   iEvent.getByToken(token1, phiContainer);
0021
0022   edm::Handle<emtf::DTTag::theta_digi_collection> thetaContainer;
0023   iEvent.getByToken(token2, thetaContainer);
0024
0025   TriggerPrimitiveCollection muon_primitives;
0026
0027   // Adapted from L1Trigger/L1TMuonBarrel/src/L1TMuonBarrelKalmanStubProcessor.cc
0028   constexpr int minPhiQuality = 0;
0029   constexpr int minBX = -3;
0030   constexpr int maxBX = 3;
0031
0032   for (int bx = minBX; bx <= maxBX; bx++) {
0033     for (int wheel = -2; wheel <= 2; wheel++) {
0034       for (int sector = 0; sector < 12; sector++) {
0035         for (int station = 1; station < 5; station++) {
0036           if (wheel == -1 || wheel == 0 || wheel == 1)
0037             continue;  // do not include wheels -1, 0, +1
0038           if (station == 4)
0039             continue;  // do not include MB4
0040
0041           // According to Michalis, in legacy BMTF, the second stub was coming as BXNUM=-1.
0042           // This is a code convention now, but you need bx-1 to get the proper second stub.
0043           emtf::DTTag::theta_digi_type const* theta_segm = thetaContainer->chThetaSegm(wheel, station, sector, bx);
0044           emtf::DTTag::digi_type const* phi_segm_high = phiContainer->chPhiSegm1(wheel, station, sector, bx);
0045           emtf::DTTag::digi_type const* phi_segm_low = phiContainer->chPhiSegm2(wheel, station, sector, bx - 1);
0046
0047           // Find theta BTI group(s)
0048           bool has_theta_segm = false;
0049           int bti_group1 = -1;
0050           int bti_group2 = -1;
0051
0052           // Case with theta segment
0053           if (theta_segm != nullptr) {
0054             has_theta_segm = true;
0055
0056             for (unsigned int i = 0; i < 7; ++i) {
0057               if (theta_segm->position(i) != 0) {
0058                 if (bti_group1 < 0) {
0059                   bti_group1 = i;
0060                   bti_group2 = i;
0061                 } else {
0062                   bti_group2 = i;
0063                 }
0064               }
0065             }
0066             emtf_assert(bti_group1 != -1 && bti_group2 != -1);
0067           }
0068
0069           // 1st phi segment
0070           if (phi_segm_high != nullptr) {
0071             if (phi_segm_high->code() >= minPhiQuality) {
0072               DTChamberId detid(phi_segm_high->whNum(), phi_segm_high->stNum(), phi_segm_high->scNum() + 1);
0073               if (has_theta_segm) {
0074                 muon_primitives.emplace_back(detid, *phi_segm_high, *theta_segm, bti_group1);
0075               } else {
0076                 muon_primitives.emplace_back(detid, *phi_segm_high, 1);
0077               }
0078             }
0079           }
0080
0081           // 2nd phi segment
0082           if (phi_segm_low != nullptr) {
0083             if (phi_segm_low->code() >= minPhiQuality) {
0084               DTChamberId detid(phi_segm_low->whNum(), phi_segm_low->stNum(), phi_segm_low->scNum() + 1);
0085               if (has_theta_segm) {
0086                 muon_primitives.emplace_back(detid, *phi_segm_low, *theta_segm, bti_group2);
0087               } else {
0088                 muon_primitives.emplace_back(detid, *phi_segm_low, 2);
0089               }
0090             }
0091           }
0092
0093           // Duplicate DT muon primitives, if more than one theta segment, but only one phi segment
0094           if (phi_segm_high != nullptr && phi_segm_low == nullptr && bti_group1 != bti_group2) {
0095             DTChamberId detid(phi_segm_high->whNum(), phi_segm_high->stNum(), phi_segm_high->scNum() + 1);
0096             muon_primitives.emplace_back(detid, *phi_segm_high, *theta_segm, bti_group2);
0097           }
0098
0099         }  // end loop over station
0100       }    // end loop over sector
0101     }      // end loop over wheel
0102   }        // end loop over bx
0103
0104   // Remove duplicates using erase-remove idiom,
0105   // assuming the vector is already sorted
0106   muon_primitives.erase(std::unique(muon_primitives.begin(), muon_primitives.end()), muon_primitives.end());
0107   std::copy(muon_primitives.begin(), muon_primitives.end(), std::back_inserter(out));
0108   return;
0109 }
0110
0111 // _____________________________________________________________________________
0112 // Specialized for CSC
0113 template <>
0114 void EMTFSubsystemCollector::extractPrimitives(emtf::CSCTag tag,
0115                                                const GeometryTranslator* tp_geom,
0116                                                const edm::Event& iEvent,
0117                                                const edm::EDGetToken& token,
0118                                                TriggerPrimitiveCollection& out) const {
0119   edm::Handle<emtf::CSCTag::digi_collection> cscDigis;
0120   iEvent.getByToken(token, cscDigis);
0121
0122   auto chamber = cscDigis->begin();
0123   auto chend = cscDigis->end();
0124   for (; chamber != chend; ++chamber) {
0125     auto digi = (*chamber).second.first;
0126     auto dend = (*chamber).second.second;
0127     for (; digi != dend; ++digi) {
0128       out.emplace_back((*chamber).first, *digi);
0129     }
0130   }
0131   return;
0132 }
0133
0134 // _____________________________________________________________________________
0135 // Specialized for RPC
0136 template <>
0137 void EMTFSubsystemCollector::extractPrimitives(emtf::RPCTag tag,
0138                                                const GeometryTranslator* tp_geom,
0139                                                const edm::Event& iEvent,
0140                                                const edm::EDGetToken& token,
0141                                                TriggerPrimitiveCollection& out) const {
0142   edm::Handle<emtf::RPCTag::digi_collection> rpcDigis;
0143   iEvent.getByToken(token, rpcDigis);
0144
0145   TriggerPrimitiveCollection muon_primitives;
0146
0147   auto chamber = rpcDigis->begin();
0148   auto chend = rpcDigis->end();
0149   for (; chamber != chend; ++chamber) {
0150     auto digi = (*chamber).second.first;
0151     auto dend = (*chamber).second.second;
0152     for (; digi != dend; ++digi) {
0153       if ((*chamber).first.region() != 0) {  // 0 is barrel
0154         if ((*chamber).first.station() <= 2 && (*chamber).first.ring() == 3)
0155           continue;  // do not include RE1/3, RE2/3
0156         if ((*chamber).first.station() >= 3 && (*chamber).first.ring() == 1)
0157           continue;  // do not include RE3/1, RE4/1 (iRPC)
0158
0159         muon_primitives.emplace_back((*chamber).first, *digi);
0160       }
0161     }
0162   }
0163
0164   // Cluster the RPC digis
0165   TriggerPrimitiveCollection clus_muon_primitives;
0166   cluster_rpc(muon_primitives, clus_muon_primitives);
0167
0168   // Output
0169   std::copy(clus_muon_primitives.begin(), clus_muon_primitives.end(), std::back_inserter(out));
0170   return;
0171 }
0172
0173 // Specialized for RPC (using RecHits)
0174 template <>
0175 void EMTFSubsystemCollector::extractPrimitives(emtf::RPCTag tag,
0176                                                const GeometryTranslator* tp_geom,
0177                                                const edm::Event& iEvent,
0178                                                const edm::EDGetToken& token1,
0179                                                const edm::EDGetToken& token2,
0180                                                TriggerPrimitiveCollection& out) const {
0181   constexpr int maxClusterSize = 3;
0182
0183   //edm::Handle<RPCTag::digi_collection> rpcDigis;
0184   //iEvent.getByToken(token1, rpcDigis);
0185
0186   edm::Handle<emtf::RPCTag::rechit_collection> rpcRecHits;
0187   iEvent.getByToken(token2, rpcRecHits);
0188
0189   auto rechit = rpcRecHits->begin();
0190   auto rhend = rpcRecHits->end();
0191   for (; rechit != rhend; ++rechit) {
0192     const RPCDetId& detid = rechit->rpcId();
0193     const RPCRoll* roll = dynamic_cast<const RPCRoll*>(tp_geom->getRPCGeometry().roll(detid));
0194     if (roll == nullptr)
0195       continue;
0196
0197     if (detid.region() != 0) {  // 0 is barrel
0198       if (detid.station() <= 2 && detid.ring() == 3)
0199         continue;  // do not include RE1/3, RE2/3
0200       if (detid.station() >= 3 && detid.ring() == 1)
0201         continue;  // do not include RE3/1, RE4/1 (iRPC)
0202
0203       if (rechit->clusterSize() <= maxClusterSize) {
0204         out.emplace_back(detid, *rechit);
0205       }
0206     }
0207   }
0208   return;
0209 }
0210
0211 // _____________________________________________________________________________
0212 // Specialized for iRPC
0213 template <>
0214 void EMTFSubsystemCollector::extractPrimitives(emtf::IRPCTag tag,
0215                                                const GeometryTranslator* tp_geom,
0216                                                const edm::Event& iEvent,
0217                                                const edm::EDGetToken& token,
0218                                                TriggerPrimitiveCollection& out) const {
0219   edm::Handle<emtf::IRPCTag::digi_collection> irpcDigis;
0220   iEvent.getByToken(token, irpcDigis);
0221
0222   TriggerPrimitiveCollection muon_primitives;
0223
0224   auto chamber = irpcDigis->begin();
0225   auto chend = irpcDigis->end();
0226   for (; chamber != chend; ++chamber) {
0227     auto digi = (*chamber).second.first;
0228     auto dend = (*chamber).second.second;
0229     for (; digi != dend; ++digi) {
0230       if ((*chamber).first.region() != 0) {  // 0 is barrel
0231         if (!((*chamber).first.station() >= 3 && (*chamber).first.ring() == 1))
0232           continue;  // only RE3/1, RE4/1 (iRPC)
0233
0234         muon_primitives.emplace_back((*chamber).first, *digi);
0235       }
0236     }
0237   }
0238
0239   // Cluster the iRPC digis
0240   TriggerPrimitiveCollection clus_muon_primitives;
0241   cluster_rpc(muon_primitives, clus_muon_primitives);
0242
0243   // Output
0244   std::copy(clus_muon_primitives.begin(), clus_muon_primitives.end(), std::back_inserter(out));
0245   return;
0246 }
0247
0248 // Specialized for iRPC (using RecHits)
0249 template <>
0250 void EMTFSubsystemCollector::extractPrimitives(emtf::IRPCTag tag,
0251                                                const GeometryTranslator* tp_geom,
0252                                                const edm::Event& iEvent,
0253                                                const edm::EDGetToken& token1,
0254                                                const edm::EDGetToken& token2,
0255                                                TriggerPrimitiveCollection& out) const {
0256   constexpr int maxClusterSize = 6;
0257
0258   //edm::Handle<emtf::IRPCTag::digi_collection> irpcDigis;
0259   //iEvent.getByToken(token1, irpcDigis);
0260
0261   edm::Handle<emtf::IRPCTag::rechit_collection> irpcRecHits;
0262   iEvent.getByToken(token2, irpcRecHits);
0263
0264   auto rechit = irpcRecHits->begin();
0265   auto rhend = irpcRecHits->end();
0266   for (; rechit != rhend; ++rechit) {
0267     const RPCDetId& detid = rechit->rpcId();
0268     const RPCRoll* roll = dynamic_cast<const RPCRoll*>(tp_geom->getRPCGeometry().roll(detid));
0269     if (roll == nullptr)
0270       continue;
0271
0272     if (detid.region() != 0) {  // 0 is barrel
0273       if (!(detid.station() >= 3 && detid.ring() == 1))
0274         continue;  // only RE3/1, RE4/1 (iRPC)
0275
0276       if (rechit->clusterSize() <= maxClusterSize) {
0277         out.emplace_back(detid, *rechit);
0278       }
0279     }
0280   }
0281   return;
0282 }
0283
0284 // _____________________________________________________________________________
0285 // Specialized for CPPF
0286 template <>
0287 void EMTFSubsystemCollector::extractPrimitives(emtf::CPPFTag tag,
0288                                                const GeometryTranslator* tp_geom,
0289                                                const edm::Event& iEvent,
0290                                                const edm::EDGetToken& token,
0291                                                TriggerPrimitiveCollection& out) const {
0292   edm::Handle<emtf::CPPFTag::digi_collection> cppfDigis;
0293   iEvent.getByToken(token, cppfDigis);
0294
0295   for (const auto& digi : *cppfDigis) {
0296     out.emplace_back(digi.rpcId(), digi);
0297   }
0298   return;
0299 }
0300
0301 // _____________________________________________________________________________
0302 // Specialized for GEM
0303 template <>
0304 void EMTFSubsystemCollector::extractPrimitives(emtf::GEMTag tag,
0305                                                const GeometryTranslator* tp_geom,
0306                                                const edm::Event& iEvent,
0307                                                const edm::EDGetToken& token,
0308                                                TriggerPrimitiveCollection& out) const {
0309   edm::Handle<emtf::GEMTag::digi_collection> gemDigis;
0310   iEvent.getByToken(token, gemDigis);
0311
0312   TriggerPrimitiveCollection muon_primitives;
0313
0314   auto chamber = gemDigis->begin();
0315   auto chend = gemDigis->end();
0316   for (; chamber != chend; ++chamber) {
0317     auto digi = (*chamber).second.first;
0318     auto dend = (*chamber).second.second;
0319     for (; digi != dend; ++digi) {
0320       auto detid = (*chamber).first;
0321       // temporarily ignore 16-partition GE2/1 clusters, because the EMTF
0322       // is not yet adapted to handle these objects
0323       if (detid.isGE21() and digi->nPartitions() == GEMPadDigi::GE21SplitStrip)
0324         continue;
0325       muon_primitives.emplace_back((*chamber).first, *digi);
0326     }
0327   }
0328
0329   // Make GEM coincidence pads
0330   TriggerPrimitiveCollection copad_muon_primitives;
0331   make_copad_gem(muon_primitives, copad_muon_primitives);
0332
0333   // Output
0334   std::copy(copad_muon_primitives.begin(), copad_muon_primitives.end(), std::back_inserter(out));
0335   return;
0336 }
0337
0338 // _____________________________________________________________________________
0339 // Specialized for ME0
0340 template <>
0341 void EMTFSubsystemCollector::extractPrimitives(emtf::ME0Tag tag,
0342                                                const GeometryTranslator* tp_geom,
0343                                                const edm::Event& iEvent,
0344                                                const edm::EDGetToken& token,
0345                                                TriggerPrimitiveCollection& out) const {
0346   edm::Handle<emtf::ME0Tag::digi_collection> me0Digis;
0347   iEvent.getByToken(token, me0Digis);
0348
0349   auto chamber = me0Digis->begin();
0350   auto chend = me0Digis->end();
0351   for (; chamber != chend; ++chamber) {
0352     auto digi = (*chamber).second.first;
0353     auto dend = (*chamber).second.second;
0354     for (; digi != dend; ++digi) {
0355       out.emplace_back((*chamber).first, *digi);
0356     }
0357   }
0358   return;
0359 }
0360
0361 // _____________________________________________________________________________
0362 // RPC functions
0363 void EMTFSubsystemCollector::cluster_rpc(const TriggerPrimitiveCollection& muon_primitives,
0364                                          TriggerPrimitiveCollection& clus_muon_primitives) const {
0365   // Define operator to select RPC digis
0366   struct {
0367     typedef TriggerPrimitive value_type;
0368     bool operator()(const value_type& x) const { return (x.subsystem() == L1TMuon::kRPC); }
0369   } rpc_digi_select;
0370
0371   // Define operator to sort the RPC digis prior to clustering.
0372   // Use rawId, bx and strip as the sorting id. RPC rawId fully specifies
0373   // sector, subsector, endcap, station, ring, layer, roll. Strip is used as
0374   // the least significant sorting id.
0375   struct {
0376     typedef TriggerPrimitive value_type;
0377     bool operator()(const value_type& lhs, const value_type& rhs) const {
0378       bool cmp = (std::make_pair(std::make_pair(lhs.rawId(), lhs.getRPCData().bx), lhs.getRPCData().strip) <
0379                   std::make_pair(std::make_pair(rhs.rawId(), rhs.getRPCData().bx), rhs.getRPCData().strip));
0380       return cmp;
0381     }
0382   } rpc_digi_less;
0383
0384   struct {
0385     typedef TriggerPrimitive value_type;
0386     bool operator()(const value_type& lhs, const value_type& rhs) const {
0387       bool cmp = (std::make_pair(std::make_pair(lhs.rawId(), lhs.getRPCData().bx), lhs.getRPCData().strip) ==
0388                   std::make_pair(std::make_pair(rhs.rawId(), rhs.getRPCData().bx), rhs.getRPCData().strip));
0389       return cmp;
0390     }
0391   } rpc_digi_equal;
0392
0393   // Define operators for the nearest-neighbor clustering algorithm.
0394   // If two digis are next to each other (check strip_hi on the 'left', and
0395   // strip_low on the 'right'), cluster them (increment strip_hi on the 'left')
0396   struct {
0397     typedef TriggerPrimitive value_type;
0398     bool operator()(const value_type& lhs, const value_type& rhs) const {
0399       bool cmp = ((lhs.rawId() == rhs.rawId()) && (lhs.getRPCData().bx == rhs.getRPCData().bx) &&
0400                   (lhs.getRPCData().strip_hi + 1 == rhs.getRPCData().strip_low));
0401       return cmp;
0402     }
0403   } rpc_digi_adjacent;
0404
0405   struct {
0406     typedef TriggerPrimitive value_type;
0407     void operator()(value_type& lhs, value_type& rhs) {  // pass by reference
0408       lhs.accessRPCData().strip_hi += 1;
0409     }
0410   } rpc_digi_cluster;
0411
0412   // ___________________________________________________________________________
0413   // Do clustering using C++ <algorithm> functions
0414
0415   // 1. Select RPC digis
0416   clus_muon_primitives.clear();
0417   std::copy_if(
0418       muon_primitives.begin(), muon_primitives.end(), std::back_inserter(clus_muon_primitives), rpc_digi_select);
0419
0420   // 2. Sort
0421   std::stable_sort(clus_muon_primitives.begin(), clus_muon_primitives.end(), rpc_digi_less);
0422
0423   // 3. Remove duplicates
0424   clus_muon_primitives.erase(std::unique(clus_muon_primitives.begin(), clus_muon_primitives.end(), rpc_digi_equal),
0425                              clus_muon_primitives.end());
0426
0427   // 4. Cluster adjacent digis
0428   clus_muon_primitives.erase(
0429       adjacent_cluster(clus_muon_primitives.begin(), clus_muon_primitives.end(), rpc_digi_adjacent, rpc_digi_cluster),
0430       clus_muon_primitives.end());
0431 }
0432
0433 // _____________________________________________________________________________
0434 // GEM functions
0435 void EMTFSubsystemCollector::make_copad_gem(const TriggerPrimitiveCollection& muon_primitives,
0436                                             TriggerPrimitiveCollection& copad_muon_primitives) const {
0437   // Use the inner layer (layer 1) hit coordinates as output, and the outer
0438   // layer (layer 2) as coincidence
0439
0440   // Adapted from L1Trigger/CSCTriggerPrimitives/src/GEMCoPadProcessor.cc
0441   constexpr unsigned int maxDeltaBX = 1;
0442   constexpr unsigned int maxDeltaRoll = 1;
0443   constexpr unsigned int maxDeltaPadGE11 = 3;  // it was 2
0444   constexpr unsigned int maxDeltaPadGE21 = 2;
0445
0446   // Make sure that the difference is calculated using signed integer, and
0447   // output the absolute difference (as unsigned integer)
0448   auto calculate_delta = [](int a, int b) -> unsigned int { return std::abs(a - b); };
0449
0450   // Create maps of GEM pads (key = detid), split by layer
0451   std::map<uint32_t, TriggerPrimitiveCollection> in_pads_layer1, in_pads_layer2;
0452
0453   auto tp_it = muon_primitives.begin();
0454   auto tp_end = muon_primitives.end();
0455
0456   for (; tp_it != tp_end; ++tp_it) {
0457     GEMDetId detid = tp_it->detId<GEMDetId>();
0458     emtf_assert(detid.layer() == 1 || detid.layer() == 2);
0459     emtf_assert(1 <= detid.roll() && detid.roll() <= 8);
0460     uint32_t layer = detid.layer();
0461
0462     // Remove layer number and roll number from detid
0463     detid = GEMDetId(detid.region(), detid.ring(), detid.station(), 0, detid.chamber(), 0);
0464
0465     if (layer == 1) {
0466       in_pads_layer1[detid.rawId()].push_back(*tp_it);
0467     } else {
0468       in_pads_layer2[detid.rawId()].push_back(*tp_it);
0469     }
0470   }
0471
0472   // Build coincidences
0473   copad_muon_primitives.clear();
0474
0475   auto map_tp_it = in_pads_layer1.begin();
0476   auto map_tp_end = in_pads_layer1.end();
0477
0478   for (; map_tp_it != map_tp_end; ++map_tp_it) {
0479     const GEMDetId& detid = map_tp_it->first;
0480     const TriggerPrimitiveCollection& pads = map_tp_it->second;
0481
0482     // find all corresponding ids with layer 2
0483     auto found = in_pads_layer2.find(detid);
0484
0485     // empty range = no possible coincidence pads
0486     if (found == in_pads_layer2.end())
0487       continue;
0488
0489     // now let's correlate the pads in two layers of this partition
0490     const TriggerPrimitiveCollection& co_pads = found->second;
0491     for (auto p = pads.begin(); p != pads.end(); ++p) {
0492       bool has_copad = false;
0493
0494       for (auto co_p = co_pads.begin(); co_p != co_pads.end() && !has_copad; ++co_p) {
0495         // check the match in pad
0496         unsigned int deltaPad = calculate_delta(p->getGEMData().pad, co_p->getGEMData().pad);
0497         if ((detid.station() == 1 && deltaPad > maxDeltaPadGE11) ||
0498             (detid.station() == 2 && deltaPad > maxDeltaPadGE21))
0499           continue;
0500         // check the match in BX
0501         unsigned int deltaBX = calculate_delta(p->getGEMData().bx, co_p->getGEMData().bx);
0502         if (deltaBX > maxDeltaBX)
0503           continue;
0504         // check the match in roll
0505         unsigned int deltaRoll = calculate_delta(p->detId<GEMDetId>().roll(), co_p->detId<GEMDetId>().roll());
0506         if (deltaRoll > maxDeltaRoll)
0507           continue;
0508         has_copad = true;
0509       }  // end loop over co_pads
0510
0511       // make a new coincidence pad digi
0512       if (has_copad) {
0513         copad_muon_primitives.push_back(*p);
0514       }
0515     }  // end loop over pads
0516   }    // end loop over in_pads_layer1
0517 }