Bugfixes in simtrack to GEM/CSC digi and stub matching

DataFormats/GEMDigi/interface/GEMDigi.h

@@ -14,7 +14,7 @@
 
 class GEMDigi {
 public:
-  explicit GEMDigi(int strip, int bx);
+  explicit GEMDigi(uint16_t strip, int16_t bx);
   GEMDigi();
 
   bool operator==(const GEMDigi& digi) const;
@@ -23,8 +23,8 @@ class GEMDigi {
   bool isValid() const;
 
   // return the strip number. counts from 0.
-  int strip() const { return strip_; }
-  int bx() const { return bx_; }
+  uint16_t strip() const { return strip_; }
+  int16_t bx() const { return bx_; }
 
   void print() const;
 

DataFormats/GEMDigi/src/GEMDigi.cc

@@ -1,7 +1,7 @@
 #include "DataFormats/GEMDigi/interface/GEMDigi.h"
 #include <iostream>
 
-GEMDigi::GEMDigi(int strip, int bx) : strip_(strip), bx_(bx) {}
+GEMDigi::GEMDigi(uint16_t strip, int16_t bx) : strip_(strip), bx_(bx) {}
 
 GEMDigi::GEMDigi() : strip_(65535), bx_(-99) {}
 

SimDataFormats/GEMDigiSimLink/interface/GEMDigiSimLink.h

@@ -35,7 +35,7 @@ class GEMDigiSimLink {
   ~GEMDigiSimLink() { ; }
 
   unsigned int getStrip() const { return _digi.first; }
-  unsigned int getBx() const { return _digi.second; }
+  int getBx() const { return _digi.second; }
   Local3DPoint getEntryPoint() const { return _entryPoint; }
   LocalVector getMomentumAtEntry() const { return _momentumAtEntry; }
   float getTimeOfFlight() const { return _timeOfFlight; }

Validation/MuonCSCDigis/src/CSCDigiMatcher.cc

@@ -57,16 +57,6 @@ void CSCDigiMatcher::match(const SimTrack& t, const SimVertex& v) {
 }
 
 void CSCDigiMatcher::matchComparatorsToSimTrack(const CSCComparatorDigiCollection& comparators) {
-  for (auto detUnitIt = comparators.begin(); detUnitIt != comparators.end(); ++detUnitIt) {
-    const CSCDetId& id = (*detUnitIt).first;
-    const auto& range = (*detUnitIt).second;
-    for (auto digiIt = range.first; digiIt != range.second; ++digiIt) {
-      if (id.station() == 1 and (id.ring() == 1 or id.ring() == 4))
-        if (verboseComparator_)
-          cout << "CSCid " << id << " Comparator digi (comparator, comparator, Tbin ) " << (*digiIt) << endl;
-    }
-  }
-
   const auto& det_ids = muonSimHitMatcher_->detIds(0);
   for (const auto& id : det_ids) {
     CSCDetId layer_id(id);

Validation/MuonCSCDigis/src/CSCStubMatcher.cc

@@ -1,6 +1,5 @@
 #include "Validation/MuonCSCDigis/interface/CSCStubMatcher.h"
 #include "L1Trigger/CSCCommonTrigger/interface/CSCConstants.h"
-
 #include <algorithm>
 
 using namespace std;
@@ -87,12 +86,9 @@ void CSCStubMatcher::matchCLCTsToSimTrack(const CSCCLCTDigiCollection& clcts) {
     }
 
     int ring = ch_id.ring();
-    if (ring == 4)
-      ring = 1;  //use ME1b id to get CLCTs
-    CSCDetId ch_id2(ch_id.endcap(), ch_id.station(), ring, ch_id.chamber(), 0);
 
     // do not consider CSCs with too few hits
-    if (cscDigiMatcher_->nLayersWithStripInChamber(ch_id2) < minNHitsChamberCLCT_)
+    if (cscDigiMatcher_->nLayersWithStripInChamber(ch_id) < minNHitsChamberCLCT_)
       continue;
 
     // get the comparator digis in this chamber
@@ -104,7 +100,7 @@ void CSCStubMatcher::matchCLCTsToSimTrack(const CSCCLCTDigiCollection& clcts) {
 
     // print out the digis
     if (verboseCLCT_) {
-      cout << "clct: digi_strips " << ch_id2 << endl;
+      cout << "clct: comparators " << ch_id << endl;
       int layer = 0;
       for (const auto& p : comps) {
         layer++;
@@ -115,6 +111,12 @@ void CSCStubMatcher::matchCLCTsToSimTrack(const CSCCLCTDigiCollection& clcts) {
       }
     }
 
+    //use ME1b id to get CLCTs
+    const bool isME1a(ch_id.station() == 1 and ch_id.ring() == 4);
+    if (isME1a)
+      ring = 1;
+    CSCDetId ch_id2(ch_id.endcap(), ch_id.station(), ring, ch_id.chamber(), 0);
+
     const auto& clcts_in_det = clcts.get(ch_id2);
 
     for (auto c = clcts_in_det.first; c != clcts_in_det.second; ++c) {
@@ -137,15 +139,28 @@ void CSCStubMatcher::matchCLCTsToSimTrack(const CSCCLCTDigiCollection& clcts) {
       for (const auto& p : comps) {
         layer++;
         for (const auto& q : p) {
+          if (verboseCLCT_)
+            cout << "L" << layer << " " << q << " " << q.getHalfStrip() << " " << std::endl;
           for (const auto& clctComp : (*c).getHits()[layer - 1]) {
-            if (q.getHalfStrip() == clctComp) {
+            if (clctComp == 65535)
+              continue;
+            if (verboseCLCT_) {
+              std::cout << "\t" << clctComp << " " << endl;
+            }
+            if (q.getHalfStrip() == clctComp or (isME1a and q.getHalfStrip() + 128 == clctComp)) {
               nMatches++;
+              if (verboseCLCT_) {
+                cout << "\t\tnMatches " << nMatches << std::endl;
+              }
             }
           }
         }
+        if (verboseCLCT_) {
+          cout << endl;
+        }
       }
 
-      // require at least 3 good half-strip matches
+      // require at least 3 good matches
       if (nMatches < 3)
         continue;
 
@@ -257,19 +272,23 @@ void CSCStubMatcher::matchLCTsToSimTrack(const CSCCorrelatedLCTDigiCollection& l
     for (const auto& lct : lcts_tmp) {
       iLct++;
 
-      bool lct_matched(false);
       bool lct_clct_match(false);
       bool lct_alct_match(false);
       bool lct_gem1_match(false);
       bool lct_gem2_match(false);
 
-      if (verboseLCT_)
-        cout << "in LCT, getCLCT " << lct.getCLCT() << " getALCT " << lct.getALCT() << endl;
-
+      if (verboseLCT_) {
+        cout << ch_id << " " << ch_id2 << endl;
+        cout << lct << endl;
+        cout << "getCLCT " << lct.getCLCT() << "\ngetALCT " << lct.getALCT() << "\ngetGEM1 " << lct.getGEM1()
+             << "\ngetGEM2 " << lct.getGEM2() << endl;
+      }
       // Check if matched to an CLCT
       for (const auto& p : clctsInChamber(id)) {
         if (p == lct.getCLCT()) {
           lct_clct_match = true;
+          if (verboseLCT_)
+            cout << "\t...lct_clct_match" << endl;
           break;
         }
       }
@@ -278,6 +297,8 @@ void CSCStubMatcher::matchLCTsToSimTrack(const CSCCorrelatedLCTDigiCollection& l
       for (const auto& p : alctsInChamber(id)) {
         if (p == lct.getALCT()) {
           lct_alct_match = true;
+          if (verboseLCT_)
+            cout << "\t...lct_alct_match" << endl;
           break;
         }
       }
@@ -289,6 +310,8 @@ void CSCStubMatcher::matchLCTsToSimTrack(const CSCCorrelatedLCTDigiCollection& l
         for (const auto& p : gemDigiMatcher_->padsInChamber(gemDetIdL1.rawId())) {
           if (p == lct.getGEM1()) {
             lct_gem1_match = true;
+            if (verboseLCT_)
+              cout << "\t...lct_gem1_match" << endl;
             break;
           }
         }
@@ -298,17 +321,24 @@ void CSCStubMatcher::matchLCTsToSimTrack(const CSCCorrelatedLCTDigiCollection& l
         for (const auto& p : gemDigiMatcher_->padsInChamber(gemDetIdL2.rawId())) {
           if (p == lct.getGEM2()) {
             lct_gem2_match = true;
+            if (verboseLCT_)
+              cout << "\t...lct_gem2_match" << endl;
             break;
           }
         }
       }
 
-      lct_matched = ((lct_clct_match and lct_alct_match) or (lct_alct_match and lct_gem1_match and lct_gem2_match) or
-                     (lct_clct_match and lct_gem1_match and lct_gem2_match));
+      const bool alct_clct = lct_clct_match and lct_alct_match;
+      const bool alct_gem = lct_alct_match and lct_gem1_match and lct_gem2_match;
+      const bool clct_gem = lct_clct_match and lct_gem1_match and lct_gem2_match;
+
+      bool lct_tight_matched = alct_clct or alct_gem or clct_gem;
+      bool lct_loose_matched = lct_clct_match or lct_alct_match;
+      bool lct_matched = lct_loose_matched or lct_tight_matched;
 
       if (lct_matched) {
         if (verboseLCT_)
-          cout << "this LCT matched to simtrack in chamber " << ch_id << endl;
+          cout << "...was matched" << endl;
         if (std::find(chamber_to_lcts_[id].begin(), chamber_to_lcts_[id].end(), lct) == chamber_to_lcts_[id].end()) {
           chamber_to_lcts_[id].emplace_back(lct);
         }
@@ -331,7 +361,6 @@ void CSCStubMatcher::matchMPLCTsToSimTrack(const CSCCorrelatedLCTDigiCollection&
       if (!lct->isValid())
         continue;
 
-      // std::cout << "MPC Stub ALL " << *lct << std::endl;
       chamber_to_mplcts_all_[id].emplace_back(*lct);
 
       // check if this stub corresponds with a previously matched stub

Validation/MuonGEMDigis/BuildFile.xml

@@ -1,6 +1,7 @@
 <use name="rootcore"/>
 <use name="Validation/MuonHits"/>
 <use name="DataFormats/GEMDigi"/>
+<use name="SimDataFormats/GEMDigiSimLink"/>
 <export>
   <lib name="1"/>
 </export>

Validation/MuonGEMDigis/interface/GEMDigiMatcher.h

@@ -10,16 +10,19 @@
 
 #include "Validation/MuonHits/interface/GEMSimHitMatcher.h"
 #include "Geometry/GEMGeometry/interface/GEMGeometry.h"
+#include "SimDataFormats/GEMDigiSimLink/interface/GEMDigiSimLink.h"
 #include "DataFormats/GEMDigi/interface/GEMDigiCollection.h"
 #include "DataFormats/GEMDigi/interface/GEMPadDigiCollection.h"
 #include "DataFormats/GEMDigi/interface/GEMPadDigiClusterCollection.h"
 #include "DataFormats/GEMDigi/interface/GEMCoPadDigiCollection.h"
+#include "DataFormats/Common/interface/DetSetVector.h"
 
 #include <vector>
 #include <map>
 #include <set>
 #include <unordered_set>
 
+typedef std::vector<GEMDigiSimLink> GEMDigiSimLinkContainer;
 typedef std::vector<GEMDigi> GEMDigiContainer;
 typedef std::vector<GEMPadDigi> GEMPadDigiContainer;
 typedef std::vector<GEMPadDigiCluster> GEMPadDigiClusterContainer;
@@ -39,6 +42,8 @@ class GEMDigiMatcher {
   // do the matching
   void match(const SimTrack& t, const SimVertex& v);
 
+  std::set<unsigned int> detIdsSimLink(int gem_type = MuonHitHelper::GEM_ALL) const;
+
   // partition GEM detIds with digis
   std::set<unsigned int> detIdsDigi(int gem_type = MuonHitHelper::GEM_ALL) const;
   std::set<unsigned int> detIdsPad(int gem_type = MuonHitHelper::GEM_ALL) const;
@@ -55,6 +60,8 @@ class GEMDigiMatcher {
   std::set<unsigned int> superChamberIdsCluster(int gem_type = MuonHitHelper::GEM_ALL) const;
   std::set<unsigned int> superChamberIdsCoPad(int gem_type = MuonHitHelper::GEM_ALL) const;
 
+  const GEMDigiSimLinkContainer& simLinksInDetId(unsigned int) const;
+
   // GEM digis from a particular partition, chamber or superchamber
   const GEMDigiContainer& digisInDetId(unsigned int) const;
   const GEMDigiContainer& digisInChamber(unsigned int) const;
@@ -78,6 +85,7 @@ class GEMDigiMatcher {
   // #layers with digis from this simtrack
   int nLayersWithDigisInSuperChamber(unsigned int) const;
   int nLayersWithPadsInSuperChamber(unsigned int) const;
+  int nLayersWithClustersInSuperChamber(unsigned int) const;
 
   /// How many pads in GEM did this simtrack get in total?
   int nPads() const;
@@ -102,18 +110,21 @@ class GEMDigiMatcher {
   std::shared_ptr<GEMSimHitMatcher> muonSimHitMatcher() { return muonSimHitMatcher_; }
 
 private:
+  void matchDigisSLToSimTrack(const edm::DetSetVector<GEMDigiSimLink>&);
   void matchDigisToSimTrack(const GEMDigiCollection&);
   void matchPadsToSimTrack(const GEMPadDigiCollection&);
   void matchClustersToSimTrack(const GEMPadDigiClusterCollection&);
   void matchCoPadsToSimTrack(const GEMCoPadDigiCollection&);
 
   void clear();
 
+  edm::EDGetTokenT<edm::DetSetVector<GEMDigiSimLink>> gemSimLinkToken_;
   edm::EDGetTokenT<GEMDigiCollection> gemDigiToken_;
   edm::EDGetTokenT<GEMPadDigiCollection> gemPadToken_;
   edm::EDGetTokenT<GEMPadDigiClusterCollection> gemClusterToken_;
   edm::EDGetTokenT<GEMCoPadDigiCollection> gemCoPadToken_;
 
+  edm::Handle<edm::DetSetVector<GEMDigiSimLink>> gemDigisSLH_;
   edm::Handle<GEMDigiCollection> gemDigisH_;
   edm::Handle<GEMPadDigiCollection> gemPadsH_;
   edm::Handle<GEMPadDigiClusterCollection> gemClustersH_;
@@ -127,18 +138,24 @@ class GEMDigiMatcher {
   template <class T>
   std::set<unsigned int> selectDetIds(const T&, int) const;
 
+  bool simMuOnly_;
+  bool discardEleHits_;
+
   int minBXDigi_, maxBXDigi_;
   int minBXPad_, maxBXPad_;
   int minBXCluster_, maxBXCluster_;
   int minBXCoPad_, maxBXCoPad_;
 
   int matchDeltaStrip_;
 
+  bool verboseSimLink_;
   bool verboseDigi_;
   bool verbosePad_;
   bool verboseCluster_;
   bool verboseCoPad_;
 
+  std::map<unsigned int, GEMDigiSimLinkContainer> detid_to_simLinks_;
+
   std::map<unsigned int, GEMDigiContainer> detid_to_digis_;
   std::map<unsigned int, GEMDigiContainer> chamber_to_digis_;
   std::map<unsigned int, GEMDigiContainer> superchamber_to_digis_;
@@ -153,6 +170,7 @@ class GEMDigiMatcher {
 
   std::map<unsigned int, GEMCoPadDigiContainer> superchamber_to_copads_;
 
+  GEMDigiSimLinkContainer no_gem_simLinks_;
   GEMDigiContainer no_gem_digis_;
   GEMPadDigiContainer no_gem_pads_;
   GEMPadDigiClusterContainer no_gem_clusters_;

Validation/MuonGEMDigis/python/muonGEMDigiPSet.py

@@ -1,6 +1,12 @@
 import FWCore.ParameterSet.Config as cms
 
 muonGEMDigiPSet = cms.PSet(
+    gemSimLink = cms.PSet(
+        verbose = cms.int32(0),
+        inputTag = cms.InputTag("simMuonGEMDigis","GEM"),
+        simMuOnly = cms.bool(True),
+        discardEleHits = cms.bool(True),
+    ),
     gemStripDigi = cms.PSet(
         verbose = cms.int32(0),
         inputTag = cms.InputTag("simMuonGEMDigis"),
@@ -30,7 +36,7 @@
     gemCoPadDigi = cms.PSet(
         verbose = cms.int32(0),
         inputTag = cms.InputTag("simCscTriggerPrimitiveDigis"),
-        minBX = cms.int32(-1),
-        maxBX = cms.int32(1),
+        minBX = cms.int32(0),
+        maxBX = cms.int32(0),
     ),
 )