File indexing completed on 2021-10-20 10:41:39
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021 #include <memory>
0022
0023
0024 #include "FWCore/Framework/interface/Frameworkfwd.h"
0025 #include "FWCore/Framework/interface/stream/EDProducer.h"
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/EventSetup.h"
0028 #include "FWCore/Framework/interface/MakerMacros.h"
0029 #include "FWCore/Framework/interface/ConsumesCollector.h"
0030 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0031
0032
0033 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0034 #include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h"
0035
0036
0037 #include "DataFormats/MuonReco/interface/Muon.h"
0038 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0039 #include "DataFormats/TrackReco/interface/Track.h"
0040 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0041 #include "TrackingTools/GeomPropagators/interface/Propagator.h"
0042 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0043 #include "Geometry/DTGeometry/interface/DTGeometry.h"
0044 #include "Geometry/CSCGeometry/interface/CSCGeometry.h"
0045 #include "DataFormats/MuonDetId/interface/MuonSubdetId.h"
0046 #include "TrackingTools/Records/interface/TrackingComponentsRecord.h"
0047 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0048 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0049 #include "MagneticField/Engine/interface/MagneticField.h"
0050 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0051 #include "DataFormats/TrackingRecHit/interface/TrackingRecHitFwd.h"
0052 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
0053 #include "RecoMuon/TransientTrackingRecHit/interface/MuonTransientTrackingRecHitBuilder.h"
0054 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
0055 #include "Geometry/CommonDetUnit/interface/GlobalTrackingGeometry.h"
0056 #include "Geometry/Records/interface/GlobalTrackingGeometryRecord.h"
0057 #include "TrackingTools/TrackRefitter/interface/TrackTransformer.h"
0058
0059
0060
0061
0062
0063 class TrackerToMuonPropagator : public edm::stream::EDProducer<> {
0064 public:
0065 explicit TrackerToMuonPropagator(const edm::ParameterSet&);
0066 ~TrackerToMuonPropagator() override;
0067
0068 private:
0069 void produce(edm::Event&, const edm::EventSetup&) override;
0070
0071
0072
0073
0074 const edm::ESGetToken<Propagator, TrackingComponentsRecord> m_esTokenProp;
0075 const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> m_esTokenTk;
0076 const edm::ESGetToken<DTGeometry, MuonGeometryRecord> m_esTokenDT;
0077 const edm::ESGetToken<CSCGeometry, MuonGeometryRecord> m_esTokenCSC;
0078 const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> m_esTokenMF;
0079 const edm::ESGetToken<GlobalTrackingGeometry, GlobalTrackingGeometryRecord> m_esTokenGTGeo;
0080
0081 edm::InputTag m_globalMuons, m_globalMuonTracks;
0082
0083 bool m_refitTracker;
0084 TrackTransformer* m_trackTransformer;
0085 };
0086
0087
0088
0089
0090
0091
0092
0093
0094
0095
0096
0097
0098 TrackerToMuonPropagator::TrackerToMuonPropagator(const edm::ParameterSet& iConfig)
0099 : m_esTokenProp(esConsumes(edm::ESInputTag("", iConfig.getParameter<std::string>("propagator")))),
0100 m_esTokenTk(esConsumes()),
0101 m_esTokenDT(esConsumes()),
0102 m_esTokenCSC(esConsumes()),
0103 m_esTokenMF(esConsumes()),
0104 m_esTokenGTGeo(esConsumes()) {
0105 m_globalMuons = iConfig.getParameter<edm::InputTag>("globalMuons");
0106 m_globalMuonTracks = iConfig.getParameter<edm::InputTag>("globalMuonTracks");
0107 m_refitTracker = iConfig.getParameter<bool>("refitTrackerTrack");
0108 if (m_refitTracker) {
0109 m_trackTransformer =
0110 new TrackTransformer(iConfig.getParameter<edm::ParameterSet>("trackerTrackTransformer"), consumesCollector());
0111 } else
0112 m_trackTransformer = nullptr;
0113
0114 produces<std::vector<Trajectory>>();
0115 produces<TrajTrackAssociationCollection>();
0116 }
0117
0118 TrackerToMuonPropagator::~TrackerToMuonPropagator() {
0119
0120
0121 }
0122
0123
0124
0125
0126
0127
0128 void TrackerToMuonPropagator::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0129 if (m_trackTransformer)
0130 m_trackTransformer->setServices(iSetup);
0131
0132 edm::Handle<reco::MuonCollection> globalMuons;
0133 iEvent.getByLabel(m_globalMuons, globalMuons);
0134
0135 edm::Handle<reco::TrackCollection> globalMuonTracks;
0136 iEvent.getByLabel(m_globalMuonTracks, globalMuonTracks);
0137
0138 const Propagator* propagator = &iSetup.getData(m_esTokenProp);
0139 const TrackerGeometry* trackerGeometry = &iSetup.getData(m_esTokenTk);
0140 const DTGeometry* dtGeometry = &iSetup.getData(m_esTokenDT);
0141 const CSCGeometry* cscGeometry = &iSetup.getData(m_esTokenCSC);
0142 const MagneticField* magneticField = &iSetup.getData(m_esTokenMF);
0143 const GlobalTrackingGeometry* globalGeometry = &iSetup.getData(m_esTokenGTGeo);
0144
0145
0146
0147 MuonTransientTrackingRecHitBuilder muonTransBuilder;
0148
0149
0150 auto trajectoryCollection = std::make_unique<std::vector<Trajectory>>();
0151
0152
0153 std::map<edm::Ref<std::vector<Trajectory>>::key_type, edm::Ref<reco::TrackCollection>::key_type> reference_map;
0154 edm::Ref<std::vector<Trajectory>>::key_type trajCounter = 0;
0155
0156 for (reco::MuonCollection::const_iterator globalMuon = globalMuons->begin(); globalMuon != globalMuons->end();
0157 ++globalMuon) {
0158
0159 edm::Ref<reco::TrackCollection>::key_type trackCounter = 0;
0160 reco::TrackCollection::const_iterator globalMuonTrack = globalMuonTracks->begin();
0161 for (; globalMuonTrack != globalMuonTracks->end(); ++globalMuonTrack) {
0162 trackCounter++;
0163 if (fabs(globalMuon->combinedMuon()->phi() - globalMuonTrack->phi()) < 1e-10 &&
0164 fabs(globalMuon->combinedMuon()->eta() - globalMuonTrack->eta()) < 1e-10)
0165 break;
0166 }
0167 if (globalMuonTrack == globalMuonTracks->end()) {
0168 throw cms::Exception("BadConfig") << "The tracks label doesn't correspond to the same objects as the muons label"
0169 << std::endl;
0170 }
0171
0172 TrajectoryStateOnSurface tracker_tsos;
0173 DetId outerDetId;
0174 if (m_refitTracker) {
0175 std::vector<Trajectory> trackerTrajectories = m_trackTransformer->transform(*globalMuon->track());
0176 if (trackerTrajectories.size() == 1) {
0177 const Trajectory trackerTrajectory = *(trackerTrajectories.begin());
0178
0179
0180 tracker_tsos = trackerTrajectory.firstMeasurement().forwardPredictedState();
0181 outerDetId = trackerTrajectory.firstMeasurement().recHit()->geographicalId();
0182 } else
0183 continue;
0184 } else {
0185
0186 GlobalPoint outerPosition(globalMuon->track()->outerPosition().x(),
0187 globalMuon->track()->outerPosition().y(),
0188 globalMuon->track()->outerPosition().z());
0189 GlobalVector outerMomentum(globalMuon->track()->outerMomentum().x(),
0190 globalMuon->track()->outerMomentum().y(),
0191 globalMuon->track()->outerMomentum().z());
0192 int charge = globalMuon->track()->charge();
0193 const reco::Track::CovarianceMatrix outerStateCovariance = globalMuon->track()->outerStateCovariance();
0194 outerDetId = DetId(globalMuon->track()->outerDetId());
0195
0196
0197 GlobalTrajectoryParameters globalTrajParams(outerPosition, outerMomentum, charge, magneticField);
0198 CurvilinearTrajectoryError curviError(outerStateCovariance);
0199 FreeTrajectoryState tracker_state(globalTrajParams, curviError);
0200
0201
0202 tracker_tsos =
0203 TrajectoryStateOnSurface(globalTrajParams, curviError, trackerGeometry->idToDet(outerDetId)->surface());
0204 }
0205
0206 TrajectoryStateOnSurface last_tsos = tracker_tsos;
0207
0208
0209 edm::OwnVector<TrackingRecHit> muonHits;
0210 std::vector<TrajectoryStateOnSurface> TSOSes;
0211 for (auto const& hit : globalMuon->combinedMuon()->recHits()) {
0212 DetId id = hit->geographicalId();
0213
0214 TrajectoryStateOnSurface extrapolation;
0215 bool extrapolated = false;
0216 if (id.det() == DetId::Muon && id.subdetId() == MuonSubdetId::DT) {
0217 extrapolation = propagator->propagate(last_tsos, dtGeometry->idToDet(id)->surface());
0218 extrapolated = true;
0219 } else if (id.det() == DetId::Muon && id.subdetId() == MuonSubdetId::CSC) {
0220 extrapolation = propagator->propagate(last_tsos, cscGeometry->idToDet(id)->surface());
0221 extrapolated = true;
0222 }
0223
0224 if (extrapolated && extrapolation.isValid()) {
0225 muonHits.push_back(hit->clone());
0226 TSOSes.push_back(extrapolation);
0227 }
0228 }
0229
0230
0231 if (!muonHits.empty()) {
0232 PTrajectoryStateOnDet const& PTraj = trajectoryStateTransform::persistentState(tracker_tsos, outerDetId.rawId());
0233 TrajectorySeed trajectorySeed(PTraj, muonHits, alongMomentum);
0234 Trajectory trajectory(trajectorySeed, alongMomentum);
0235
0236 for (unsigned int i = 0; i < muonHits.size(); i++) {
0237 TrajectoryMeasurement::ConstRecHitPointer hitPtr(muonTransBuilder.build(&(muonHits[i]), globalGeometry));
0238 TrajectoryStateOnSurface TSOS = TSOSes[i];
0239 trajectory.push(TrajectoryMeasurement(TSOS, TSOS, TSOS, hitPtr));
0240 }
0241
0242 trajectoryCollection->push_back(trajectory);
0243
0244
0245 trajCounter++;
0246 reference_map[trajCounter] = trackCounter;
0247
0248 }
0249
0250 }
0251
0252 unsigned int numTrajectories = trajectoryCollection->size();
0253
0254
0255 edm::OrphanHandle<std::vector<Trajectory>> ohTrajs = iEvent.put(std::move(trajectoryCollection));
0256
0257
0258 auto trajTrackMap = std::make_unique<TrajTrackAssociationCollection>();
0259
0260 for (trajCounter = 0; trajCounter < numTrajectories; trajCounter++) {
0261 edm::Ref<reco::TrackCollection>::key_type trackCounter = reference_map[trajCounter];
0262
0263 trajTrackMap->insert(edm::Ref<std::vector<Trajectory>>(ohTrajs, trajCounter),
0264 edm::Ref<reco::TrackCollection>(globalMuonTracks, trackCounter));
0265 }
0266
0267 iEvent.put(std::move(trajTrackMap));
0268 }
0269
0270
0271 DEFINE_FWK_MODULE(TrackerToMuonPropagator);