MuonTools/plugins/MuDTMuonExtTableProducer.cc

0001 /** \class MuDTMuonExtTableProducer MuDTMuonExtTableProducer.cc DPGAnalysis/MuonTools/plugins/MuDTMuonExtTableProducer.cc
0002  *
0003  * Helper class : the muon filler
0004  *
0005  * \author L. Lunerti (INFN BO)
0006  *
0007  *
0008  */
0009
0010 #include "DataFormats/DetId/interface/DetId.h"
0011 #include "DataFormats/MuonDetId/interface/DTChamberId.h"
0012 #include "DataFormats/MuonReco/interface/MuonChamberMatch.h"
0013 #include "DataFormats/MuonReco/interface/MuonSegmentMatch.h"
0014 #include "DataFormats/MuonReco/interface/MuonSelectors.h"
0015
0016 #include "DataFormats/DTRecHit/interface/DTRecSegment4D.h"
0017 #include "DataFormats/DTRecHit/interface/DTRecSegment4DCollection.h"
0018 #include "DataFormats/TrackingRecHit/interface/TrackingRecHitFwd.h"
0019 #include "DataFormats/MuonDetId/interface/DTChamberId.h"
0020
0021 #include "DataFormats/Math/interface/deltaR.h"
0022
0023 #include "TString.h"
0024 #include "TRegexp.h"
0025
0026 #include <numeric>
0027 #include <vector>
0028
0029 #include "DPGAnalysis/MuonTools/interface/MuBaseFlatTableProducer.h"
0030
0031 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0032 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0033
0034 #include "DataFormats/Common/interface/View.h"
0035 #include "DataFormats/MuonReco/interface/Muon.h"
0036 #include "DataFormats/VertexReco/interface/Vertex.h"
0037 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0038 #include "DataFormats/MuonReco/interface/MuonIsolation.h"
0039 #include "DataFormats/MuonReco/interface/MuonPFIsolation.h"
0040 #include "DataFormats/DTRecHit/interface/DTRecSegment4D.h"
0041
0042 #include "DataFormats/Common/interface/TriggerResults.h"
0043 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0044 #include "DataFormats/HLTReco/interface/TriggerObject.h"
0045 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0046
0047 class MuDTMuonExtTableProducer : public MuBaseFlatTableProducer {
0048 public:
0049   /// Constructor
0050   MuDTMuonExtTableProducer(const edm::ParameterSet&);
0051
0052   /// Fill descriptors
0053   static void fillDescriptions(edm::ConfigurationDescriptions&);
0054
0055 protected:
0056   /// Fill tree branches for a given events
0057   void fillTable(edm::Event&) final;
0058
0059   /// Get info from the ES by run
0060   void getFromES(const edm::Run&, const edm::EventSetup&) final;
0061
0062 private:
0063   /// Tokens
0064   nano_mu::EDTokenHandle<edm::View<reco::Muon>> m_muToken;
0065   nano_mu::EDTokenHandle<DTRecSegment4DCollection> m_dtSegmentToken;
0066
0067   nano_mu::EDTokenHandle<edm::TriggerResults> m_trigResultsToken;
0068   nano_mu::EDTokenHandle<trigger::TriggerEvent> m_trigEventToken;
0069
0070   /// Fill matches table
0071   bool m_fillMatches;
0072
0073   /// Name of the triggers used by muon filler for trigger matching
0074   std::string m_trigName;
0075   std::string m_isoTrigName;
0076
0077   /// DT Geometry
0078   nano_mu::ESTokenHandle<DTGeometry, MuonGeometryRecord, edm::Transition::BeginRun> m_dtGeometry;
0079
0080   /// HLT config provider
0081   HLTConfigProvider m_hltConfig;
0082
0083   /// Indices of the triggers used by muon filler for trigger matching
0084   std::vector<int> m_trigIndices;
0085   std::vector<int> m_isoTrigIndices;
0086
0087   bool hasTrigger(std::vector<int>&,
0088                   const trigger::TriggerObjectCollection&,
0089                   edm::Handle<trigger::TriggerEvent>&,
0090                   const reco::Muon&);
0091 };
0092
0093 MuDTMuonExtTableProducer::MuDTMuonExtTableProducer(const edm::ParameterSet& config)
0094     : MuBaseFlatTableProducer(config),
0095       m_muToken{config, consumesCollector(), "src"},
0096       m_dtSegmentToken{config, consumesCollector(), "dtSegmentSrc"},
0097       m_trigResultsToken{config, consumesCollector(), "trigResultsSrc"},
0098       m_trigEventToken{config, consumesCollector(), "trigEventSrc"},
0099       m_fillMatches{config.getParameter<bool>("fillMatches")},
0100       m_trigName{config.getParameter<std::string>("trigName")},
0101       m_isoTrigName{config.getParameter<std::string>("isoTrigName")},
0102       m_dtGeometry{consumesCollector()} {
0103   produces<nanoaod::FlatTable>();
0104   if (m_fillMatches) {
0105     produces<nanoaod::FlatTable>("matches");
0106     produces<nanoaod::FlatTable>("staMatches");
0107   }
0108 }
0109
0110 void MuDTMuonExtTableProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0111   edm::ParameterSetDescription desc;
0112
0113   desc.add<std::string>("name", "muon");
0114   desc.add<edm::InputTag>("src", edm::InputTag{"patMuons"});
0115   desc.add<edm::InputTag>("dtSegmentSrc", edm::InputTag{"dt4DSegments"});
0116   desc.add<edm::InputTag>("trigEventSrc", edm::InputTag{"hltTriggerSummaryAOD::HLT"});
0117   desc.add<edm::InputTag>("trigResultsSrc", edm::InputTag{"TriggerResults::HLT"});
0118
0119   desc.add<bool>("fillMatches", true);
0120
0121   desc.add<std::string>("trigName", "HLT_Mu55*");
0122   desc.add<std::string>("isoTrigName", "HLT_IsoMu2*");
0123
0124   descriptions.addWithDefaultLabel(desc);
0125 }
0126
0127 void MuDTMuonExtTableProducer::getFromES(const edm::Run& run, const edm::EventSetup& environment) {
0128   m_dtGeometry.getFromES(environment);
0129
0130   bool changed{true};
0131   m_hltConfig.init(run, environment, "HLT", changed);
0132
0133   const bool enableWildcard{true};
0134
0135   TString tName = TString(m_trigName);
0136   TRegexp tNamePattern = TRegexp(tName, enableWildcard);
0137
0138   for (unsigned iPath = 0; iPath < m_hltConfig.size(); ++iPath) {
0139     TString pathName = TString(m_hltConfig.triggerName(iPath));
0140     if (pathName.Contains(tNamePattern))
0141       m_trigIndices.push_back(static_cast<int>(iPath));
0142   }
0143
0144   tName = TString(m_isoTrigName);
0145   tNamePattern = TRegexp(tName, enableWildcard);
0146
0147   for (unsigned iPath = 0; iPath < m_hltConfig.size(); ++iPath) {
0148     TString pathName = TString(m_hltConfig.triggerName(iPath));
0149     if (pathName.Contains(tNamePattern))
0150       m_isoTrigIndices.push_back(static_cast<int>(iPath));
0151   }
0152 }
0153
0154 void MuDTMuonExtTableProducer::fillTable(edm::Event& ev) {
0155   unsigned int nMuons{0};
0156
0157   std::vector<bool> firesIsoTrig;
0158   std::vector<bool> firesTrig;
0159
0160   std::vector<int> nMatches;
0161   std::vector<int> staMu_nMatchSeg;
0162
0163   std::vector<uint32_t> matches_begin;
0164   std::vector<uint32_t> matches_end;
0165
0166   std::vector<uint32_t> staMatches_begin;
0167   std::vector<uint32_t> staMatches_end;
0168
0169   std::vector<int16_t> matches_wheel;
0170   std::vector<int16_t> matches_sector;
0171   std::vector<int16_t> matches_station;
0172
0173   std::vector<float> matches_x;
0174   std::vector<float> matches_y;
0175
0176   std::vector<float> matches_phi;
0177   std::vector<float> matches_eta;
0178   std::vector<float> matches_edgeX;
0179   std::vector<float> matches_edgeY;
0180
0181   std::vector<float> matches_dXdZ;
0182   std::vector<float> matches_dYdZ;
0183
0184   std::vector<uint32_t> staMatches_MuSegIdx;
0185
0186   auto muons = m_muToken.conditionalGet(ev);
0187   auto segments = m_dtSegmentToken.conditionalGet(ev);
0188
0189   auto triggerResults = m_trigResultsToken.conditionalGet(ev);
0190   auto triggerEvent = m_trigEventToken.conditionalGet(ev);
0191
0192   if (muons.isValid() && segments.isValid()) {
0193     for (const auto& muon : (*muons)) {
0194       if (triggerResults.isValid() && triggerEvent.isValid()) {
0195         const auto& triggerObjects = triggerEvent->getObjects();
0196
0197         bool hasIsoTrig = hasTrigger(m_isoTrigIndices, triggerObjects, triggerEvent, muon);
0198         bool hasTrig = hasTrigger(m_trigIndices, triggerObjects, triggerEvent, muon);
0199
0200         firesIsoTrig.push_back(hasIsoTrig);
0201         firesTrig.push_back(hasTrig);
0202
0203       } else {
0204         firesIsoTrig.push_back(false);
0205         firesTrig.push_back(false);
0206       }
0207
0208       size_t iMatches = 0;
0209       size_t iSegMatches = 0;
0210
0211       if (m_fillMatches) {
0212         matches_begin.push_back(matches_wheel.size());
0213
0214         if (muon.isMatchesValid()) {
0215           for (const auto& match : muon.matches()) {
0216             if (match.id.det() == DetId::Muon && match.id.subdetId() == MuonSubdetId::DT) {
0217               DTChamberId dtId(match.id.rawId());
0218               const auto chamb = m_dtGeometry->chamber(static_cast<DTChamberId>(match.id));
0219
0220               matches_wheel.push_back(dtId.wheel());
0221               matches_sector.push_back(dtId.sector());
0222               matches_station.push_back(dtId.station());
0223
0224               matches_x.push_back(match.x);
0225               matches_y.push_back(match.y);
0226
0227               matches_phi.push_back(chamb->toGlobal(LocalPoint(match.x, match.y, 0.)).phi());
0228               matches_eta.push_back(chamb->toGlobal(LocalPoint(match.x, match.y, 0.)).eta());
0229
0230               matches_edgeX.push_back(match.edgeX);
0231               matches_edgeY.push_back(match.edgeY);
0232
0233               matches_dXdZ.push_back(match.dXdZ);
0234               matches_dYdZ.push_back(match.dYdZ);
0235
0236               ++iMatches;
0237             }
0238           }
0239         }
0240
0241         matches_end.push_back(matches_wheel.size());
0242
0243         //SEGMENT MATCHING VARIABLES
0244
0245         staMatches_begin.push_back(staMatches_MuSegIdx.size());
0246
0247         if (!muon.outerTrack().isNull()) {
0248           reco::TrackRef outerTrackRef = muon.outerTrack();
0249
0250           auto recHitIt = outerTrackRef->recHitsBegin();
0251           auto recHitEnd = outerTrackRef->recHitsEnd();
0252
0253           for (; recHitIt != recHitEnd; ++recHitIt) {
0254             DetId detId = (*recHitIt)->geographicalId();
0255
0256             if (detId.det() == DetId::Muon && detId.subdetId() == MuonSubdetId::DT) {
0257               const auto dtSegmentSta = dynamic_cast<const DTRecSegment4D*>((*recHitIt));
0258               int iSeg = 0;
0259
0260               for (const auto& segment : (*segments)) {
0261                 if (dtSegmentSta && dtSegmentSta->chamberId().station() == segment.chamberId().station() &&
0262                     dtSegmentSta->chamberId().wheel() == segment.chamberId().wheel() &&
0263                     dtSegmentSta->chamberId().sector() == segment.chamberId().sector() &&
0264                     std::abs(dtSegmentSta->localPosition().x() - segment.localPosition().x()) < 0.001 &&
0265                     std::abs(dtSegmentSta->localPosition().y() - segment.localPosition().y()) < 0.001 &&
0266                     std::abs(dtSegmentSta->localDirection().x() - segment.localDirection().x()) < 0.001 &&
0267                     std::abs(dtSegmentSta->localDirection().y() - segment.localDirection().y()) < 0.001) {
0268                   staMatches_MuSegIdx.push_back(iSeg);
0269                   ++iSegMatches;
0270                 }
0271
0272                 ++iSeg;
0273               }  //loop over segments
0274             }
0275
0276           }  //loop over recHits
0277         }
0278
0279         staMatches_end.push_back(staMatches_MuSegIdx.size());
0280       }
0281
0282       nMatches.push_back(iMatches);
0283       staMu_nMatchSeg.push_back(iSegMatches);
0284
0285       ++nMuons;
0286     }
0287   }
0288
0289   auto table = std::make_unique<nanoaod::FlatTable>(nMuons, m_name, false, true);
0290
0291   addColumn(table,
0292             "firesIsoTrig",
0293             firesIsoTrig,
0294             "True if the muon is matched to an isolated trigger"
0295             "<br />specified in the ntuple config file");
0296
0297   addColumn(table,
0298             "firesTrig",
0299             firesTrig,
0300             "True if the muon is matched to a  (non isolated)trigger"
0301             "<br />specified in the ntuple config file");
0302
0303   addColumn(table, "nMatches", nMatches, "Number of muon chamber matches (DT only)");
0304   addColumn(table, "staMu_nMatchSeg", staMu_nMatchSeg, "Number of segments used in the standalone track (DT only)");
0305
0306   addColumn(table,
0307             "matches_begin",
0308             matches_begin,
0309             "begin() of range of quantities for a given muon in the *_matches_* vectors");
0310   addColumn(
0311       table, "matches_end", matches_end, "end() of range of quantities for a given muon in the *_matches_* vectors");
0312
0313   addColumn(table,
0314             "staMatches_begin",
0315             staMatches_begin,
0316             "begin() of range of quantities for a given muon in the matches_staMuSegIdx vector");
0317   addColumn(table,
0318             "staMatches_end",
0319             staMatches_end,
0320             "end() of range of quantities for a given muon in the matches_staMuSegIdx vector");
0321
0322   ev.put(std::move(table));
0323
0324   if (m_fillMatches) {
0325     auto sum = [](std::vector<int> v) { return std::accumulate(v.begin(), v.end(), 0); };
0326
0327     auto tabMatches = std::make_unique<nanoaod::FlatTable>(sum(nMatches), m_name + "_matches", false, false);
0328
0329     tabMatches->setDoc("RECO muon matches_* vectors");
0330
0331     addColumn(tabMatches, "x", matches_x, "x position of the extrapolated track on the matched DT chamber");
0332     addColumn(tabMatches, "y", matches_y, "y position of the extrapolated track on the matched DT chamber");
0333
0334     addColumn(tabMatches,
0335               "edgeX",
0336               matches_edgeX,
0337               "distance in x of the extrapolated track to the matched DT chamber border (<0 == inside the chamber)");
0338     addColumn(tabMatches,
0339               "edgeY",
0340               matches_edgeY,
0341               "distance in y of the extrapolated track to the matched DT chamber border (<0 == inside the chamber)");
0342
0343     addColumn(tabMatches, "wheel", matches_wheel, "matched DT chamber wheel");
0344     addColumn(tabMatches, "sector", matches_sector, "matched DT chamber sector");
0345     addColumn(tabMatches, "station", matches_station, "matched DT chamber station");
0346
0347     addColumn(
0348         tabMatches, "phi", matches_phi, "phi of the (x,y) position on the matched DT chamber (global reference frame)");
0349     addColumn(tabMatches,
0350               "eta",
0351               matches_eta,
0352               " eta of the (x,y) position on the matched cDT hamber (global reference frame)");
0353
0354     addColumn(tabMatches, "dXdZ", matches_dXdZ, "dXdZ of the extrapolated track on the matched DT chamber");
0355     addColumn(tabMatches, "dYdZ", matches_dYdZ, "dYdZ of the extrapolated track on the matched DT chamber");
0356
0357     ev.put(std::move(tabMatches), "matches");
0358
0359     auto tabStaMatches =
0360         std::make_unique<nanoaod::FlatTable>(sum(staMu_nMatchSeg), m_name + "_staMatches", false, false);
0361
0362     tabStaMatches->setDoc("RECO muon staMatches_* vector");
0363
0364     addColumn(tabStaMatches,
0365               "MuSegIdx",
0366               staMatches_MuSegIdx,
0367               "Index of DT segments used in the standalone track it corresponds"
0368               "<br />to the index of a given segment in the ntuple seg_* collection");
0369
0370     ev.put(std::move(tabStaMatches), "staMatches");
0371   }
0372 }
0373
0374 bool MuDTMuonExtTableProducer::hasTrigger(std::vector<int>& trigIndices,
0375                                           const trigger::TriggerObjectCollection& trigObjs,
0376                                           edm::Handle<trigger::TriggerEvent>& trigEvent,
0377                                           const reco::Muon& muon) {
0378   float dRMatch = 999.;
0379   for (int trigIdx : trigIndices) {
0380     const std::vector<std::string> trigModuleLabels = m_hltConfig.moduleLabels(trigIdx);
0381
0382     const unsigned trigModuleIndex =
0383         std::find(trigModuleLabels.begin(), trigModuleLabels.end(), "hltBoolEnd") - trigModuleLabels.begin() - 1;
0384     const unsigned hltFilterIndex = trigEvent->filterIndex(edm::InputTag(trigModuleLabels[trigModuleIndex], "", "HLT"));
0385     if (hltFilterIndex < trigEvent->sizeFilters()) {
0386       const trigger::Keys keys = trigEvent->filterKeys(hltFilterIndex);
0387       const trigger::Vids vids = trigEvent->filterIds(hltFilterIndex);
0388       const unsigned nTriggers = vids.size();
0389
0390       for (unsigned iTrig = 0; iTrig < nTriggers; ++iTrig) {
0391         trigger::TriggerObject trigObj = trigObjs[keys[iTrig]];
0392         float dR = deltaR(muon, trigObj);
0393         if (dR < dRMatch)
0394           dRMatch = dR;
0395       }
0396     }
0397   }
0398
0399   return dRMatch < 0.1;  //CB should become programmable
0400 }
0401
0402 #include "FWCore/PluginManager/interface/ModuleDef.h"
0403 #include "FWCore/Framework/interface/MakerMacros.h"
0404
0405 DEFINE_FWK_MODULE(MuDTMuonExtTableProducer);