File indexing completed on 2021-10-01 22:40:33
0001 #include "FWCore/Framework/interface/Event.h"
0002
0003 #include "DataFormats/MuonSeed/interface/L3MuonTrajectorySeedCollection.h"
0004 #include "DataFormats/TrackerRecHit2D/interface/FastTrackerRecHit.h"
0005 #include "DataFormats/TrackReco/interface/Track.h"
0006 #include "DataFormats/Math/interface/deltaPhi.h"
0007
0008 #include "SimDataFormats/Track/interface/SimTrackContainer.h"
0009
0010 #include "RecoTracker/TkTrackingRegions/interface/RectangularEtaPhiTrackingRegion.h"
0011
0012 #include "RecoMuon/GlobalTrackingTools/interface/MuonTrackingRegionBuilder.h"
0013
0014 #include "FastSimulation/Muons/plugins/FastTSGFromL2Muon.h"
0015
0016
0017
0018 #include <set>
0019
0020 FastTSGFromL2Muon::FastTSGFromL2Muon(const edm::ParameterSet& cfg) {
0021 produces<L3MuonTrajectorySeedCollection>();
0022
0023 edm::ParameterSet serviceParameters = cfg.getParameter<edm::ParameterSet>("ServiceParameters");
0024
0025 thePtCut = cfg.getParameter<double>("PtCut");
0026
0027 theL2CollectionLabel = cfg.getParameter<edm::InputTag>("MuonCollectionLabel");
0028 theSeedCollectionLabels = cfg.getParameter<std::vector<edm::InputTag> >("SeedCollectionLabels");
0029 theSimTrackCollectionLabel = cfg.getParameter<edm::InputTag>("SimTrackCollectionLabel");
0030
0031 edm::ParameterSet regionBuilderPSet = cfg.getParameter<edm::ParameterSet>("MuonTrackingRegionBuilder");
0032 theRegionBuilder = std::make_unique<MuonTrackingRegionBuilder>(regionBuilderPSet, consumesCollector());
0033 }
0034
0035 FastTSGFromL2Muon::~FastTSGFromL2Muon() {}
0036
0037 void FastTSGFromL2Muon::beginRun(edm::Run const& run, edm::EventSetup const& es) {
0038
0039
0040
0041
0042
0043
0044
0045
0046
0047
0048
0049
0050
0051
0052 }
0053
0054 void FastTSGFromL2Muon::produce(edm::Event& ev, const edm::EventSetup& es) {
0055
0056 std::unique_ptr<L3MuonTrajectorySeedCollection> result(new L3MuonTrajectorySeedCollection());
0057
0058
0059 theRegionBuilder->setEvent(ev, es);
0060
0061
0062 edm::Handle<edm::SimTrackContainer> theSimTracks;
0063 ev.getByLabel(theSimTrackCollectionLabel, theSimTracks);
0064
0065
0066 edm::Handle<reco::TrackCollection> l2muonH;
0067 ev.getByLabel(theL2CollectionLabel, l2muonH);
0068
0069
0070 unsigned seedCollections = theSeedCollectionLabels.size();
0071 std::vector<edm::Handle<edm::View<TrajectorySeed> > > theSeeds;
0072 theSeeds.resize(seedCollections);
0073 unsigned seed_size = 0;
0074 for (unsigned iseed = 0; iseed < seedCollections; ++iseed) {
0075 ev.getByLabel(theSeedCollectionLabels[iseed], theSeeds[iseed]);
0076 seed_size += theSeeds[iseed]->size();
0077 }
0078
0079
0080 unsigned int imu = 0;
0081 unsigned int imuMax = l2muonH->size();
0082
0083 for (; imu != imuMax; ++imu) {
0084
0085 reco::TrackRef muRef(l2muonH, imu);
0086
0087
0088 if (muRef->pt() < thePtCut || muRef->innerMomentum().Rho() < thePtCut || muRef->innerMomentum().R() < 2.5)
0089 continue;
0090
0091
0092 std::unique_ptr<RectangularEtaPhiTrackingRegion> region = theRegionBuilder->region(muRef);
0093
0094
0095 std::vector<TrajectorySeed> tkSeeds;
0096 std::set<unsigned> tkIds;
0097 tkSeeds.reserve(seed_size);
0098 for (unsigned iseed = 0; iseed < seedCollections; ++iseed) {
0099 edm::Handle<edm::View<TrajectorySeed> > aSeedCollection = theSeeds[iseed];
0100 unsigned nSeeds = aSeedCollection->size();
0101 for (unsigned seednr = 0; seednr < nSeeds; ++seednr) {
0102
0103 const BasicTrajectorySeed* aSeed = &((*aSeedCollection)[seednr]);
0104
0105
0106 int simTrackId = static_cast<FastTrackerRecHit const&>(*aSeed->recHits().begin()).simTrackId(0);
0107
0108
0109 std::set<unsigned>::iterator tkId = tkIds.find(simTrackId);
0110 if (tkId != tkIds.end())
0111 continue;
0112
0113 const SimTrack& theSimTrack = (*theSimTracks)[simTrackId];
0114
0115 if (clean(muRef, region.get(), aSeed, theSimTrack))
0116 tkSeeds.push_back(*aSeed);
0117 tkIds.insert(simTrackId);
0118
0119 }
0120
0121 }
0122
0123
0124
0125
0126
0127
0128
0129
0130 unsigned int is = 0;
0131 unsigned int isMax = tkSeeds.size();
0132 for (; is != isMax; ++is) {
0133 result->push_back(L3MuonTrajectorySeed(tkSeeds[is], muRef));
0134 }
0135
0136 }
0137
0138
0139
0140
0141
0142
0143
0144 ev.put(std::move(result));
0145 }
0146
0147 bool FastTSGFromL2Muon::clean(reco::TrackRef muRef,
0148 RectangularEtaPhiTrackingRegion* region,
0149 const BasicTrajectorySeed* aSeed,
0150 const SimTrack& theSimTrack) {
0151
0152 const PixelRecoRange<float>& etaRange = region->etaRange();
0153 double etaSeed = theSimTrack.momentum().Eta();
0154 double etaLimit = (fabs(fabs(etaRange.max()) - fabs(etaRange.mean())) < 0.05)
0155 ? 0.05
0156 : fabs(fabs(etaRange.max()) - fabs(etaRange.mean()));
0157 bool inEtaRange = etaSeed >= (etaRange.mean() - etaLimit) && etaSeed <= (etaRange.mean() + etaLimit);
0158 if (!inEtaRange)
0159 return false;
0160
0161
0162 const TkTrackingRegionsMargin<float>& phiMargin = region->phiMargin();
0163 double phiSeed = theSimTrack.momentum().Phi();
0164 double phiLimit = (phiMargin.right() < 0.05) ? 0.05 : phiMargin.right();
0165 bool inPhiRange = (fabs(deltaPhi(phiSeed, double(region->direction().phi()))) < phiLimit);
0166 if (!inPhiRange)
0167 return false;
0168
0169
0170 double ptSeed = std::sqrt(theSimTrack.momentum().Perp2());
0171 double ptMin = (region->ptMin() > 3.5) ? 3.5 : region->ptMin();
0172 bool inPtRange = ptSeed >= ptMin && ptSeed <= 2 * (muRef->pt());
0173 return inPtRange;
0174 }