Muon/plugins/HLTMuonL3PreFilter.cc

0001 /** \class HLTMuonL3PreFilter
0002  *
0003  * See header file for documentation
0004  *
0005  *  \author J. Alcaraz, J-R Vlimant
0006  *
0007  */
0008
0009 #include "HLTMuonL3PreFilter.h"
0010
0011 #include "DataFormats/Common/interface/Handle.h"
0012 #include "DataFormats/Common/interface/RefToBase.h"
0013
0014 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0015 #include "DataFormats/HLTReco/interface/TriggerRefsCollections.h"
0016
0017 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0018 #include "DataFormats/TrackReco/interface/Track.h"
0019 #include "DataFormats/MuonReco/interface/MuonTrackLinks.h"
0020 #include "DataFormats/MuonReco/interface/Muon.h"
0021
0022 #include "DataFormats/MuonSeed/interface/L3MuonTrajectorySeed.h"
0023 #include "DataFormats/MuonSeed/interface/L3MuonTrajectorySeedCollection.h"
0024 #include "DataFormats/TrajectorySeed/interface/TrajectorySeedCollection.h"
0025
0026 #include "FWCore/Utilities/interface/InputTag.h"
0027
0028 //
0029 // constructors and destructor
0030 //
0031 using namespace std;
0032 using namespace edm;
0033 using namespace reco;
0034 using namespace trigger;
0035
0036 HLTMuonL3PreFilter::HLTMuonL3PreFilter(const ParameterSet& iConfig)
0037     : HLTFilter(iConfig),
0038       propSetup_(iConfig, consumesCollector()),
0039       beamspotTag_(iConfig.getParameter<edm::InputTag>("BeamSpotTag")),
0040       beamspotToken_(consumes<reco::BeamSpot>(beamspotTag_)),
0041       candTag_(iConfig.getParameter<InputTag>("CandTag")),
0042       candToken_(consumes<reco::RecoChargedCandidateCollection>(candTag_)),
0043       previousCandTag_(iConfig.getParameter<InputTag>("PreviousCandTag")),
0044       previousCandToken_(consumes<trigger::TriggerFilterObjectWithRefs>(previousCandTag_)),
0045       requireL3MuonTrajectorySeed_(iConfig.getParameter<bool>("requireL3MuonTrajectorySeed")),
0046       l1CandTag_(iConfig.getParameter<InputTag>("L1CandTag")),
0047       l1CandToken_(requireL3MuonTrajectorySeed_ ? edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs>()
0048                                                 : consumes<trigger::TriggerFilterObjectWithRefs>(l1CandTag_)),
0049       recoMuTag_(iConfig.getParameter<InputTag>("inputMuonCollection")),
0050       recoMuToken_(requireL3MuonTrajectorySeed_ ? edm::EDGetTokenT<reco::MuonCollection>()
0051                                                 : consumes<reco::MuonCollection>(recoMuTag_)),
0052       min_N_(iConfig.getParameter<int>("MinN")),
0053       max_Eta_(iConfig.getParameter<double>("MaxEta")),
0054       min_Nhits_(iConfig.getParameter<int>("MinNhits")),
0055       max_Dr_(iConfig.getParameter<double>("MaxDr")),
0056       min_Dr_(iConfig.getParameter<double>("MinDr")),
0057       max_Dz_(iConfig.getParameter<double>("MaxDz")),
0058       min_DxySig_(iConfig.getParameter<double>("MinDxySig")),
0059       min_Pt_(iConfig.getParameter<double>("MinPt")),
0060       nsigma_Pt_(iConfig.getParameter<double>("NSigmaPt")),
0061       max_NormalizedChi2_(iConfig.getParameter<double>("MaxNormalizedChi2")),
0062       max_DXYBeamSpot_(iConfig.getParameter<double>("MaxDXYBeamSpot")),
0063       min_DXYBeamSpot_(iConfig.getParameter<double>("MinDXYBeamSpot")),
0064       min_NmuonHits_(iConfig.getParameter<int>("MinNmuonHits")),
0065       max_PtDifference_(iConfig.getParameter<double>("MaxPtDifference")),
0066       min_TrackPt_(iConfig.getParameter<double>("MinTrackPt")),
0067       min_MuonStations_L3fromL1_(iConfig.getParameter<int>("minMuonStations")),
0068       allowedTypeMask_L3fromL1_(iConfig.getParameter<unsigned int>("allowedTypeMask")),
0069       requiredTypeMask_L3fromL1_(iConfig.getParameter<unsigned int>("requiredTypeMask")),
0070       maxNormalizedChi2_L3fromL1_(iConfig.getParameter<double>("MaxNormalizedChi2_L3FromL1")),
0071       trkMuonId_(muon::SelectionType(iConfig.getParameter<unsigned int>("trkMuonId"))),
0072       L1MatchingdR_(iConfig.getParameter<double>("L1MatchingdR")),
0073       L1MatchingdR2_(L1MatchingdR_ * L1MatchingdR_),
0074       matchPreviousCand_(iConfig.getParameter<bool>("MatchToPreviousCand")),
0075
0076       devDebug_(false),
0077       theL3LinksLabel(iConfig.getParameter<InputTag>("InputLinks")),
0078       linkToken_(requireL3MuonTrajectorySeed_ ? edm::EDGetTokenT<reco::MuonTrackLinksCollection>()
0079                                               : consumes<reco::MuonTrackLinksCollection>(theL3LinksLabel)) {
0080   if (L1MatchingdR_ <= 0.) {
0081     throw cms::Exception("HLTMuonL3PreFilterConfiguration")
0082         << "invalid value for parameter \"L1MatchingdR\" (must be > 0): " << L1MatchingdR_;
0083   }
0084   LogDebug("HLTMuonL3PreFilter") << " CandTag/MinN/MaxEta/MinNhits/MaxDr/MinDr/MaxDz/MinDxySig/MinPt/NSigmaPt : "
0085                                  << candTag_.encode() << " " << min_N_ << " " << max_Eta_ << " " << min_Nhits_ << " "
0086                                  << max_Dr_ << " " << min_Dr_ << " " << max_Dz_ << " " << min_DxySig_ << " " << min_Pt_
0087                                  << " " << nsigma_Pt_;
0088 }
0089
0090 HLTMuonL3PreFilter::~HLTMuonL3PreFilter() = default;
0091
0092 void HLTMuonL3PreFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0093   edm::ParameterSetDescription desc;
0094   makeHLTFilterDescription(desc);
0095   desc.add<edm::InputTag>("BeamSpotTag", edm::InputTag("hltOfflineBeamSpot"));
0096   desc.add<edm::InputTag>("CandTag", edm::InputTag("hltL3MuonCandidates"));
0097   //  desc.add<edm::InputTag>("PreviousCandTag",edm::InputTag("hltDiMuonL2PreFiltered0"));
0098   desc.add<edm::InputTag>("PreviousCandTag", edm::InputTag(""));
0099   desc.add<edm::InputTag>("L1CandTag", edm::InputTag(""));
0100   desc.add<edm::InputTag>("inputMuonCollection", edm::InputTag(""));
0101   desc.add<int>("MinN", 1);
0102   desc.add<double>("MaxEta", 2.5);
0103   desc.add<int>("MinNhits", 0);
0104   desc.add<double>("MaxDr", 2.0);
0105   desc.add<double>("MinDr", -1.0);
0106   desc.add<double>("MaxDz", 9999.0);
0107   desc.add<double>("MinDxySig", -1.0);
0108   desc.add<double>("MinPt", 3.0);
0109   desc.add<double>("NSigmaPt", 0.0);
0110   desc.add<double>("MaxNormalizedChi2", 9999.0);
0111   desc.add<double>("MaxDXYBeamSpot", 9999.0);
0112   desc.add<double>("MinDXYBeamSpot", -1.0);
0113   desc.add<int>("MinNmuonHits", 0);
0114   desc.add<double>("MaxPtDifference", 9999.0);
0115   desc.add<double>("MinTrackPt", 0.0);
0116   desc.add<int>("minMuonStations", -1);
0117   desc.add<int>("minTrkHits", -1);
0118   desc.add<int>("minMuonHits", -1);
0119   desc.add<unsigned int>("allowedTypeMask", 255);
0120   desc.add<unsigned int>("requiredTypeMask", 0);
0121   desc.add<double>("MaxNormalizedChi2_L3FromL1", 9999.);
0122   desc.add<unsigned int>("trkMuonId", 0);
0123   desc.add<double>("L1MatchingdR", 0.3);
0124   desc.add<bool>("MatchToPreviousCand", true);
0125   desc.add<edm::InputTag>("InputLinks", edm::InputTag(""));
0126   desc.add<bool>("requireL3MuonTrajectorySeed", false);
0127   PropagateToMuonSetup::fillPSetDescription(desc);
0128   descriptions.add("hltMuonL3PreFilter", desc);
0129 }
0130
0131 //
0132 // member functions
0133 //
0134
0135 // ------------ method called to produce the data  ------------
0136 bool HLTMuonL3PreFilter::hltFilter(Event& iEvent,
0137                                    const EventSetup& iSetup,
0138                                    trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0139   // All HLT filters must create and fill an HLT filter object,
0140   // recording any reconstructed physics objects satisfying (or not)
0141   // this HLT filter, and place it in the Event.
0142
0143   auto const prop = propSetup_.init(iSetup);
0144
0145   if (saveTags())
0146     filterproduct.addCollectionTag(candTag_);
0147
0148   // Read RecoChargedCandidates from L3MuonCandidateProducer:
0149   Handle<RecoChargedCandidateCollection> mucands;
0150   iEvent.getByToken(candToken_, mucands);
0151
0152   // Read L2 triggered objects:
0153   Handle<TriggerFilterObjectWithRefs> previousLevelCands;
0154   iEvent.getByToken(previousCandToken_, previousLevelCands);
0155   vector<RecoChargedCandidateRef> vl2cands;
0156   previousLevelCands->getObjects(TriggerMuon, vl2cands);
0157
0158   // Read BeamSpot information:
0159   Handle<BeamSpot> recoBeamSpotHandle;
0160   iEvent.getByToken(beamspotToken_, recoBeamSpotHandle);
0161   const BeamSpot& beamSpot = *recoBeamSpotHandle;
0162
0163   // Number of objects passing the L3 Trigger:
0164   int n = 0;
0165
0166   // sort them by L2Track
0167   std::map<reco::TrackRef, std::vector<RecoChargedCandidateRef> > L2toL3s;
0168   // map the L3 cands matched to a L1 to their position in the recoMuon collection
0169   std::map<unsigned int, RecoChargedCandidateRef> MuonToL3s;
0170
0171   // Test to see if we can use L3MuonTrajectorySeeds:
0172   if (mucands->empty())
0173     return false;
0174   auto const& tk = (*mucands)[0].track();
0175   bool useL3MTS = false;
0176
0177   if (tk->seedRef().isNonnull()) {
0178     auto a = dynamic_cast<const L3MuonTrajectorySeed*>(tk->seedRef().get());
0179
0180     if (requireL3MuonTrajectorySeed_ && not a) {
0181       throw cms::Exception("HLTMuonL3PreFilter")
0182           << "requireL3MuonTrajectorySeed is required to be true, but cannot cast the input "
0183              "RecoChargedCandidateCollection to a L3MuonTrajectorySeed./n Please cross-check the configuration of this "
0184              "module";
0185     }
0186
0187     useL3MTS = a != nullptr;
0188   }
0189
0190   // If we can use L3MuonTrajectory seeds run the older code:
0191   if (useL3MTS) {
0192     LogDebug("HLTMuonL3PreFilter") << "HLTMuonL3PreFilter::hltFilter is in mode: useL3MTS";
0193
0194     unsigned int maxI = mucands->size();
0195     for (unsigned int i = 0; i != maxI; ++i) {
0196       const TrackRef& tk = (*mucands)[i].track();
0197       edm::Ref<L3MuonTrajectorySeedCollection> l3seedRef =
0198           tk->seedRef().castTo<edm::Ref<L3MuonTrajectorySeedCollection> >();
0199       TrackRef staTrack = l3seedRef->l2Track();
0200       LogDebug("HLTMuonL3PreFilter") << "L2 from: " << iEvent.getStableProvenance(staTrack.id()).moduleLabel()
0201                                      << " index: " << staTrack.key();
0202       L2toL3s[staTrack].push_back(RecoChargedCandidateRef(mucands, i));
0203     }
0204   }  //end of useL3MTS
0205
0206   // Using normal TrajectorySeeds:
0207   else {
0208     LogDebug("HLTMuonL3PreFilter") << "HLTMuonL3PreFilter::hltFilter is in mode: not useL3MTS";
0209
0210     // Read Links collection:
0211     edm::Handle<reco::MuonTrackLinksCollection> links;
0212     iEvent.getByToken(linkToken_, links);
0213
0214     edm::Handle<trigger::TriggerFilterObjectWithRefs> level1Cands;
0215     std::vector<l1t::MuonRef> vl1cands;
0216
0217     bool check_l1match = true;
0218
0219     // Loop over RecoChargedCandidates:
0220     for (unsigned int i(0); i < mucands->size(); ++i) {
0221       RecoChargedCandidateRef cand(mucands, i);
0222       TrackRef tk = cand->track();
0223
0224       if (!matchPreviousCand_) {
0225         MuonToL3s[i] = RecoChargedCandidateRef(cand);
0226       } else {
0227         check_l1match = true;
0228         for (auto const& link : *links) {
0229           // Using the same method that was used to create the links between L3 and L2
0230           // ToDo: there should be a better way than dR,dPt matching
0231           const reco::Track& trackerTrack = *link.trackerTrack();
0232
0233           float dR2 = deltaR2(tk->eta(), tk->phi(), trackerTrack.eta(), trackerTrack.phi());
0234           float dPt = std::abs(tk->pt() - trackerTrack.pt());
0235           if (tk->pt() != 0)
0236             dPt = dPt / tk->pt();
0237
0238           if (dR2 < 0.02 * 0.02 and dPt < 0.001) {
0239             const TrackRef staTrack = link.standAloneTrack();
0240             L2toL3s[staTrack].push_back(RecoChargedCandidateRef(cand));
0241             check_l1match = false;
0242           }
0243         }  //MTL loop
0244
0245         if (!l1CandTag_.label().empty() && check_l1match) {
0246           auto const propagated = prop.extrapolate(*tk);
0247           auto const etaForMatch = propagated.isValid() ? propagated.globalPosition().eta() : cand->eta();
0248           auto const phiForMatch = propagated.isValid() ? (double)propagated.globalPosition().phi() : cand->phi();
0249           iEvent.getByToken(l1CandToken_, level1Cands);
0250           level1Cands->getObjects(trigger::TriggerL1Mu, vl1cands);
0251           const unsigned int nL1Muons(vl1cands.size());
0252           for (unsigned int il1 = 0; il1 != nL1Muons; ++il1) {
0253             if (deltaR2(etaForMatch, phiForMatch, vl1cands[il1]->eta(), vl1cands[il1]->phi()) < L1MatchingdR2_) {
0254               MuonToL3s[i] = RecoChargedCandidateRef(cand);
0255             }
0256           }
0257         }
0258       }
0259     }  //RCC loop
0260   }    //end of using normal TrajectorySeeds
0261
0262   // look at all mucands,  check cuts and add to filter object
0263   LogDebug("HLTMuonL3PreFilter") << "looking at: " << L2toL3s.size() << " L2->L3s from: " << mucands->size();
0264   for (const auto& L2toL3s_it : L2toL3s) {
0265     if (!triggeredByLevel2(L2toL3s_it.first, vl2cands))
0266       continue;
0267
0268     //loop over the L3Tk reconstructed for this L2.
0269     unsigned int iTk = 0;
0270     unsigned int maxItk = L2toL3s_it.second.size();
0271     for (; iTk != maxItk; iTk++) {
0272       const RecoChargedCandidateRef& cand = L2toL3s_it.second[iTk];
0273       if (!applySelection(cand, beamSpot))
0274         continue;
0275
0276       filterproduct.addObject(TriggerMuon, cand);
0277       n++;
0278       break;  // and go on with the next L2 association
0279     }
0280   }  ////loop over L2s from L3 grouping
0281
0282   // now loop on L3 from L1
0283   edm::Handle<reco::MuonCollection> recomuons;
0284
0285   // do the consume only if the MuonToL3s is not empty
0286   // this happens certainly if requireL3MuonTrajectorySeed_ is true
0287   if (!MuonToL3s.empty()) {
0288     iEvent.getByToken(recoMuToken_, recomuons);
0289   }
0290
0291   for (const auto& MuonToL3s_it : MuonToL3s) {
0292     const reco::Muon& muon(recomuons->at(MuonToL3s_it.first));
0293
0294     // applys specific cuts for TkMu
0295     if ((muon.type() & allowedTypeMask_L3fromL1_) == 0)
0296       continue;
0297     if ((muon.type() & requiredTypeMask_L3fromL1_) != requiredTypeMask_L3fromL1_)
0298       continue;
0299     if (muon.numberOfMatchedStations() < min_MuonStations_L3fromL1_)
0300       continue;
0301     if (!muon.globalTrack().isNull()) {
0302       if (muon.globalTrack()->normalizedChi2() > maxNormalizedChi2_L3fromL1_)
0303         continue;
0304     }
0305     if (muon.isTrackerMuon() && !muon::isGoodMuon(muon, trkMuonId_))
0306       continue;
0307
0308     const RecoChargedCandidateRef& cand = MuonToL3s_it.second;
0309     // apply common selection
0310     if (!applySelection(cand, beamSpot))
0311       continue;
0312     filterproduct.addObject(TriggerMuon, cand);
0313     n++;
0314
0315     break;  // and go on with the next L3 from L1
0316   }
0317
0318   vector<RecoChargedCandidateRef> vref;
0319   filterproduct.getObjects(TriggerMuon, vref);
0320   for (auto& i : vref) {
0321     RecoChargedCandidateRef candref = RecoChargedCandidateRef(i);
0322     TrackRef tk = candref->get<TrackRef>();
0323     LogDebug("HLTMuonL3PreFilter") << " Track passing filter: trackRef pt= " << tk->pt() << " (" << candref->pt()
0324                                    << ") "
0325                                    << ", eta: " << tk->eta() << " (" << candref->eta() << ") ";
0326   }
0327
0328   // filter decision
0329   const bool accept(n >= min_N_);
0330
0331   LogDebug("HLTMuonL3PreFilter") << " >>>>> Result of HLTMuonL3PreFilter is " << accept
0332                                  << ", number of muons passing thresholds= " << n;
0333
0334   return accept;
0335 }
0336
0337 bool HLTMuonL3PreFilter::triggeredByLevel2(const TrackRef& staTrack, vector<RecoChargedCandidateRef>& vcands) const {
0338   bool ok = false;
0339   for (auto& vcand : vcands) {
0340     if (vcand->get<TrackRef>() == staTrack) {
0341       ok = true;
0342       LogDebug("HLTMuonL3PreFilter") << "The L2 track triggered";
0343       break;
0344     }
0345   }
0346   return ok;
0347 }
0348
0349 bool HLTMuonL3PreFilter::applySelection(const RecoChargedCandidateRef& cand, const BeamSpot& beamSpot) const {
0350   // eta cut
0351   if (std::abs(cand->eta()) > max_Eta_)
0352     return false;
0353
0354   TrackRef tk = cand->track();
0355   LogDebug("HLTMuonL3PreFilter") << " Muon in loop, q*pt= " << tk->charge() * tk->pt() << " ("
0356                                  << cand->charge() * cand->pt() << ") "
0357                                  << ", eta= " << tk->eta() << " (" << cand->eta() << ") "
0358                                  << ", hits= " << tk->numberOfValidHits() << ", d0= " << tk->d0()
0359                                  << ", dz= " << tk->dz();
0360
0361   // cut on number of hits
0362   if (tk->numberOfValidHits() < min_Nhits_)
0363     return false;
0364
0365   //max dr cut
0366   auto dr =
0367       std::abs((-(cand->vx() - beamSpot.x0()) * cand->py() + (cand->vy() - beamSpot.y0()) * cand->px()) / cand->pt());
0368   if (dr > max_Dr_)
0369     return false;
0370
0371   //min dr cut
0372   if (dr < min_Dr_)
0373     return false;
0374
0375   //dz cut
0376   if (std::abs((cand->vz() - beamSpot.z0()) -
0377                ((cand->vx() - beamSpot.x0()) * cand->px() + (cand->vy() - beamSpot.y0()) * cand->py()) / cand->pt() *
0378                    cand->pz() / cand->pt()) > max_Dz_)
0379     return false;
0380
0381   // dxy significance cut (safeguard against bizarre values)
0382   if (min_DxySig_ > 0 && (tk->dxyError() <= 0 || std::abs(tk->dxy(beamSpot.position()) / tk->dxyError()) < min_DxySig_))
0383     return false;
0384
0385   //normalizedChi2 cut
0386   if (tk->normalizedChi2() > max_NormalizedChi2_)
0387     return false;
0388
0389   //dxy beamspot cut
0390   float absDxy = std::abs(tk->dxy(beamSpot.position()));
0391   if (absDxy > max_DXYBeamSpot_ || absDxy < min_DXYBeamSpot_)
0392     return false;
0393
0394   //min muon hits cut
0395   const reco::HitPattern& trackHits = tk->hitPattern();
0396   if (trackHits.numberOfValidMuonHits() < min_NmuonHits_)
0397     return false;
0398
0399   //pt difference cut
0400   double candPt = cand->pt();
0401   double trackPt = tk->pt();
0402
0403   if (std::abs(candPt - trackPt) > max_PtDifference_)
0404     return false;
0405
0406   //track pt cut
0407   if (trackPt < min_TrackPt_)
0408     return false;
0409
0410   // Pt threshold cut
0411   double pt = cand->pt();
0412   double err0 = tk->error(0);
0413   double abspar0 = std::abs(tk->parameter(0));
0414   double ptLx = pt;
0415   // convert 50% efficiency threshold to 90% efficiency threshold
0416   if (abspar0 > 0)
0417     ptLx += nsigma_Pt_ * err0 / abspar0 * pt;
0418   LogTrace("HLTMuonL3PreFilter") << " ...Muon in loop, trackkRef pt= " << tk->pt() << ", ptLx= " << ptLx << " cand pT "
0419                                  << cand->pt();
0420   if (ptLx < min_Pt_)
0421     return false;
0422
0423   return true;
0424 }
0425
0426 // declare this class as a framework plugin
0427 #include "FWCore/Framework/interface/MakerMacros.h"
0428 DEFINE_FWK_MODULE(HLTMuonL3PreFilter);