File indexing completed on 2025-01-09 23:33:54
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024 #include "DataFormats/L1TMuonPhase2/interface/TrackerMuon.h"
0025 #include "DataFormats/Math/interface/deltaR.h"
0026 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0027 #include "DataFormats/MuonSeed/interface/L2MuonTrajectorySeedCollection.h"
0028 #include "DataFormats/TrackReco/interface/Track.h"
0029 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0030 #include "FWCore/Framework/interface/ConsumesCollector.h"
0031 #include "FWCore/Framework/interface/ESHandle.h"
0032 #include "FWCore/Framework/interface/Event.h"
0033 #include "FWCore/Framework/interface/EventSetup.h"
0034 #include "FWCore/Framework/interface/Frameworkfwd.h"
0035 #include "FWCore/Framework/interface/MakerMacros.h"
0036 #include "FWCore/Framework/interface/stream/EDProducer.h"
0037 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0038 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0039 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0040 #include "FWCore/Utilities/interface/ESGetToken.h"
0041 #include "FWCore/Utilities/interface/InputTag.h"
0042
0043 #include <unordered_set>
0044
0045 class Phase2HLTMuonSelectorForL3 : public edm::stream::EDProducer<> {
0046 public:
0047
0048 Phase2HLTMuonSelectorForL3(const edm::ParameterSet&);
0049
0050
0051 ~Phase2HLTMuonSelectorForL3() override = default;
0052
0053
0054 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0055
0056
0057 void produce(edm::Event&, const edm::EventSetup&) override;
0058
0059 private:
0060 const edm::EDGetTokenT<l1t::TrackerMuonCollection> l1TkMuCollToken_;
0061 const edm::EDGetTokenT<reco::TrackCollection> l2MuCollectionToken_;
0062 const edm::EDGetTokenT<reco::TrackCollection> l3TrackCollectionToken_;
0063
0064 const bool IOFirst_;
0065 const double matchingDr_;
0066 const bool applyL3Filters_;
0067 const double maxNormalizedChi2_, maxPtDifference_;
0068 const int minNhits_, minNhitsMuons_, minNhitsPixel_, minNhitsTracker_;
0069
0070
0071 const bool rejectL3Track(l1t::TrackerMuonRef l1TkMuRef, reco::TrackRef l3TrackRef) const;
0072 };
0073
0074
0075 Phase2HLTMuonSelectorForL3::Phase2HLTMuonSelectorForL3(const edm::ParameterSet& iConfig)
0076 : l1TkMuCollToken_(consumes<l1t::TrackerMuonCollection>(iConfig.getParameter<edm::InputTag>("l1TkMuons"))),
0077 l2MuCollectionToken_(consumes<reco::TrackCollection>(iConfig.getParameter<edm::InputTag>("l2MuonsUpdVtx"))),
0078 l3TrackCollectionToken_(consumes<reco::TrackCollection>(iConfig.getParameter<edm::InputTag>("l3Tracks"))),
0079 IOFirst_(iConfig.getParameter<bool>("IOFirst")),
0080 matchingDr_(iConfig.getParameter<double>("matchingDr")),
0081 applyL3Filters_(iConfig.getParameter<bool>("applyL3Filters")),
0082 maxNormalizedChi2_(iConfig.getParameter<double>("MaxNormalizedChi2")),
0083 maxPtDifference_(iConfig.getParameter<double>("MaxPtDifference")),
0084 minNhits_(iConfig.getParameter<int>("MinNhits")),
0085 minNhitsMuons_(iConfig.getParameter<int>("MinNhitsMuons")),
0086 minNhitsPixel_(iConfig.getParameter<int>("MinNhitsPixel")),
0087 minNhitsTracker_(iConfig.getParameter<int>("MinNhitsTracker")) {
0088 if (IOFirst_) {
0089 produces<reco::TrackCollection>("L2MuToReuse");
0090 produces<reco::TrackCollection>("L3IOTracksFiltered");
0091 } else {
0092 produces<l1t::TrackerMuonCollection>("L1TkMuToReuse");
0093 produces<reco::TrackCollection>("L3OITracksFiltered");
0094 }
0095 }
0096
0097 void Phase2HLTMuonSelectorForL3::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0098 edm::ParameterSetDescription desc;
0099 desc.add<edm::InputTag>("l1TkMuons", edm::InputTag("l1tTkMuonsGmt"));
0100 desc.add<edm::InputTag>("l2MuonsUpdVtx", edm::InputTag("hltL2MuonsFromL1TkMuon", "UpdatedAtVtx"));
0101 desc.add<edm::InputTag>("l3Tracks", edm::InputTag("hltIter2Phase2L3FromL1TkMuonMerged"));
0102 desc.add<bool>("IOFirst", true);
0103 desc.add<double>("matchingDr", 0.02);
0104 desc.add<bool>("applyL3Filters", true);
0105 desc.add<int>("MinNhits", 1);
0106 desc.add<double>("MaxNormalizedChi2", 5.0);
0107 desc.add<int>("MinNhitsMuons", 0);
0108 desc.add<int>("MinNhitsPixel", 1);
0109 desc.add<int>("MinNhitsTracker", 6);
0110 desc.add<double>("MaxPtDifference", 999.0);
0111 descriptions.add("Phase2HLTMuonSelectorForL3", desc);
0112 }
0113
0114
0115
0116 void Phase2HLTMuonSelectorForL3::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0117 const std::string metname = "RecoMuon|Phase2HLTMuonSelectorForL3";
0118
0119
0120 auto l3TracksCollectionH = iEvent.getHandle(l3TrackCollectionToken_);
0121
0122 if (IOFirst_) {
0123 LogDebug(metname) << "Inside-Out reconstruction done first, looping over L2 muons";
0124
0125
0126 auto const l2MuonsCollectionH = iEvent.getHandle(l2MuCollectionToken_);
0127
0128
0129 std::unique_ptr<reco::TrackCollection> L2MuToReuse = std::make_unique<reco::TrackCollection>();
0130 std::unique_ptr<reco::TrackCollection> L3IOTracksFiltered = std::make_unique<reco::TrackCollection>();
0131
0132
0133 std::unordered_set<size_t> goodL3Indexes;
0134
0135
0136 for (size_t l2MuIndex = 0; l2MuIndex != l2MuonsCollectionH->size(); ++l2MuIndex) {
0137 reco::TrackRef l2MuRef(l2MuonsCollectionH, l2MuIndex);
0138 bool reuseL2 = true;
0139
0140
0141 edm::RefToBase<TrajectorySeed> seedRef = l2MuRef->seedRef();
0142 edm::Ref<L2MuonTrajectorySeedCollection> l2Seed = seedRef.castTo<edm::Ref<L2MuonTrajectorySeedCollection>>();
0143 l1t::TrackerMuonRef l1TkMuRef = l2Seed->l1TkMu();
0144
0145
0146 if (l1TkMuRef.isNonnull()) {
0147
0148 LogDebug(metname) << "Looping over L3 tracks";
0149 for (size_t l3MuIndex = 0; l3MuIndex != l3TracksCollectionH->size(); ++l3MuIndex) {
0150 reco::TrackRef l3TrackRef(l3TracksCollectionH, l3MuIndex);
0151 bool rejectL3 = true;
0152
0153 if (applyL3Filters_) {
0154 LogDebug(metname) << "Checking L3 Track quality";
0155 rejectL3 = rejectL3Track(l1TkMuRef, l3TrackRef);
0156 if (!rejectL3) {
0157 LogDebug(metname) << "Adding good quality L3 IO track to filtered collection";
0158 goodL3Indexes.insert(l3MuIndex);
0159 }
0160 }
0161
0162 float dR2 = deltaR2(l1TkMuRef->phEta(), l1TkMuRef->phPhi(), l3TrackRef->eta(), l3TrackRef->phi());
0163 LogDebug(metname) << "deltaR2: " << dR2;
0164 if (dR2 < matchingDr_ * matchingDr_) {
0165 LogDebug(metname) << "Found L2 muon that matches the L3 track";
0166 reuseL2 = applyL3Filters_ ? rejectL3 : false;
0167 LogDebug(metname) << "Reuse L2: " << reuseL2;
0168 }
0169 }
0170 } else {
0171 LogDebug(metname) << "Found L2 muon without an associated L1TkMu";
0172 }
0173 if (reuseL2) {
0174 LogDebug(metname) << "Found a L2 muon to be reused";
0175 L2MuToReuse->push_back(*l2MuRef);
0176 }
0177 }
0178
0179
0180 for (const size_t index : goodL3Indexes) {
0181 L3IOTracksFiltered->push_back(*(reco::TrackRef(l3TracksCollectionH, index)));
0182 }
0183
0184 LogDebug(metname) << "Placing L2 Muons to be reused in the event";
0185 iEvent.put(std::move(L2MuToReuse), "L2MuToReuse");
0186 LogDebug(metname) << "Placing good quality L3 IO Tracks in the event";
0187 iEvent.put(std::move(L3IOTracksFiltered), "L3IOTracksFiltered");
0188 } else {
0189 LogDebug(metname) << "Outside-In reconstruction done first, looping over L1Tk muons";
0190
0191
0192 auto const l1TkMuonsCollectionH = iEvent.getHandle(l1TkMuCollToken_);
0193
0194
0195 std::unique_ptr<l1t::TrackerMuonCollection> L1TkMuToReuse = std::make_unique<l1t::TrackerMuonCollection>();
0196 std::unique_ptr<reco::TrackCollection> L3OITracksFiltered = std::make_unique<reco::TrackCollection>();
0197
0198
0199 std::unordered_set<size_t> goodL3Indexes;
0200
0201
0202 for (size_t l1TkMuIndex = 0; l1TkMuIndex != l1TkMuonsCollectionH->size(); ++l1TkMuIndex) {
0203 l1t::TrackerMuonRef l1TkMuRef(l1TkMuonsCollectionH, l1TkMuIndex);
0204 bool reuseL1TkMu = true;
0205
0206
0207 LogDebug(metname) << "Looping over L3 tracks";
0208 for (size_t l3MuIndex = 0; l3MuIndex != l3TracksCollectionH->size(); ++l3MuIndex) {
0209 reco::TrackRef l3TrackRef(l3TracksCollectionH, l3MuIndex);
0210 bool rejectL3 = true;
0211
0212 if (applyL3Filters_) {
0213 LogDebug(metname) << "Checking L3 Track quality";
0214 rejectL3 = rejectL3Track(l1TkMuRef, l3TrackRef);
0215 if (!rejectL3) {
0216 LogDebug(metname) << "Adding good quality L3 OI track to filtered collection";
0217 goodL3Indexes.insert(l3MuIndex);
0218 }
0219 }
0220
0221 float dR2 = deltaR2(l1TkMuRef->phEta(), l1TkMuRef->phPhi(), l3TrackRef->eta(), l3TrackRef->phi());
0222 LogDebug(metname) << "deltaR2: " << dR2;
0223 if (dR2 < matchingDr_ * matchingDr_) {
0224 LogDebug(metname) << "Found L1TkMu that matches the L3 track";
0225 reuseL1TkMu = applyL3Filters_ ? rejectL3 : false;
0226 LogDebug(metname) << "Reuse L1TkMu: " << reuseL1TkMu;
0227 }
0228 }
0229 if (reuseL1TkMu) {
0230 LogDebug(metname) << "Found a L1TkMu to be reused";
0231 L1TkMuToReuse->push_back(*l1TkMuRef);
0232 }
0233 }
0234
0235
0236 for (const size_t index : goodL3Indexes) {
0237 L3OITracksFiltered->push_back(*(reco::TrackRef(l3TracksCollectionH, index)));
0238 }
0239
0240 LogDebug(metname) << "Placing L1Tk Muons to be reused in the event";
0241 iEvent.put(std::move(L1TkMuToReuse), "L1TkMuToReuse");
0242 LogDebug(metname) << "Placing good quality L3 OI Tracks in the event";
0243 iEvent.put(std::move(L3OITracksFiltered), "L3OITracksFiltered");
0244 }
0245 }
0246
0247 const bool Phase2HLTMuonSelectorForL3::rejectL3Track(l1t::TrackerMuonRef l1TkMuRef, reco::TrackRef l3TrackRef) const {
0248 const std::string metname = "RecoMuon|Phase2HLTMuonSelectorForL3";
0249
0250 bool nHitsCut = l3TrackRef->numberOfValidHits() < minNhits_;
0251 bool chi2Cut = l3TrackRef->normalizedChi2() > maxNormalizedChi2_;
0252 bool nHitsMuonsCut = l3TrackRef->hitPattern().numberOfValidMuonHits() < minNhitsMuons_;
0253 bool nHitsPixelCut = l3TrackRef->hitPattern().numberOfValidPixelHits() < minNhitsPixel_;
0254 bool nHitsTrackerCut = l3TrackRef->hitPattern().trackerLayersWithMeasurement() < minNhitsTracker_;
0255 bool ptCut = std::abs(l3TrackRef->pt() - l1TkMuRef->phPt()) > maxPtDifference_ * l3TrackRef->pt();
0256
0257 bool reject = nHitsCut or chi2Cut or nHitsMuonsCut or nHitsPixelCut or nHitsTrackerCut or ptCut;
0258
0259 LogDebug(metname) << "nHits: " << l3TrackRef->numberOfValidHits() << " | chi2: " << l3TrackRef->normalizedChi2()
0260 << " | nHitsMuon: " << l3TrackRef->hitPattern().numberOfValidMuonHits()
0261 << " | nHitsPixel: " << l3TrackRef->hitPattern().numberOfValidPixelHits()
0262 << " | nHitsTracker: " << l3TrackRef->hitPattern().trackerLayersWithMeasurement();
0263 LogDebug(metname) << "Reject L3 Track: " << reject;
0264 return reject;
0265 }
0266
0267 #include "FWCore/Framework/interface/MakerMacros.h"
0268 DEFINE_FWK_MODULE(Phase2HLTMuonSelectorForL3);