File indexing completed on 2024-04-06 12:11:23
0001
0002 #include <memory>
0003 #include <vector>
0004 #include <map>
0005
0006
0007 #include "FWCore/Framework/interface/stream/EDProducer.h"
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/Framework/interface/EventSetup.h"
0010 #include "FWCore/Framework/interface/ESHandle.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012 #include "FWCore/Utilities/interface/InputTag.h"
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014
0015
0016 #include "DataFormats/Common/interface/Handle.h"
0017 #include "DataFormats/Common/interface/OwnVector.h"
0018 #include "DataFormats/TrackCandidate/interface/TrackCandidateCollection.h"
0019 #include "DataFormats/TrajectorySeed/interface/TrajectorySeedCollection.h"
0020 #include "DataFormats/TrackerRecHit2D/interface/FastTrackerRecHit.h"
0021 #include "DataFormats/TrackerRecHit2D/interface/FastTrackerRecHitCollection.h"
0022 #include "SimDataFormats/Track/interface/SimTrack.h"
0023 #include "SimDataFormats/Track/interface/SimTrackContainer.h"
0024
0025
0026 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0027 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0028 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0029 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0030 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
0031 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0032 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h"
0033 #include "TrackingTools/GeomPropagators/interface/Propagator.h"
0034 #include "TrackingTools/Records/interface/TrackingComponentsRecord.h"
0035
0036
0037 #include "FastSimulation/Tracking/interface/TrackingLayer.h"
0038 #include "FastSimulation/Tracking/interface/FastTrackingUtilities.h"
0039 #include "FastSimulation/Tracking/interface/FastTrackerRecHitSplitter.h"
0040 #include "FastSimulation/Tracking/interface/SeedMatcher.h"
0041
0042 class TrackCandidateProducer : public edm::stream::EDProducer<> {
0043 public:
0044 explicit TrackCandidateProducer(const edm::ParameterSet& conf);
0045
0046 void produce(edm::Event&, const edm::EventSetup&) override;
0047
0048 private:
0049
0050 edm::EDGetTokenT<edm::View<TrajectorySeed> > seedToken;
0051 edm::EDGetTokenT<FastTrackerRecHitCombinationCollection> recHitCombinationsToken;
0052 edm::EDGetTokenT<std::vector<bool> > hitMasksToken;
0053 edm::EDGetTokenT<edm::SimTrackContainer> simTrackToken;
0054 std::string propagatorLabel;
0055 const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> magneticFieldESToken_;
0056 const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> trackerGeometryESToken_;
0057 const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> trackerTopologyESToken_;
0058 const edm::ESGetToken<Propagator, TrackingComponentsRecord> propagatorESToken_;
0059
0060
0061 bool rejectOverlaps;
0062 bool splitHits;
0063 FastTrackerRecHitSplitter hitSplitter;
0064 double maxSeedMatchEstimator;
0065 };
0066
0067 TrackCandidateProducer::TrackCandidateProducer(const edm::ParameterSet& conf)
0068 : propagatorLabel(conf.getParameter<std::string>("propagator")),
0069 magneticFieldESToken_(esConsumes()),
0070 trackerGeometryESToken_(esConsumes()),
0071 trackerTopologyESToken_(esConsumes()),
0072 propagatorESToken_(esConsumes(edm::ESInputTag("", propagatorLabel))),
0073 hitSplitter() {
0074
0075 produces<TrackCandidateCollection>();
0076
0077
0078 if (conf.exists("hitMasks")) {
0079 hitMasksToken = consumes<std::vector<bool> >(conf.getParameter<edm::InputTag>("hitMasks"));
0080 }
0081 seedToken = consumes<edm::View<TrajectorySeed> >(conf.getParameter<edm::InputTag>("src"));
0082 recHitCombinationsToken =
0083 consumes<FastTrackerRecHitCombinationCollection>(conf.getParameter<edm::InputTag>("recHitCombinations"));
0084 simTrackToken = consumes<edm::SimTrackContainer>(conf.getParameter<edm::InputTag>("simTracks"));
0085
0086
0087 maxSeedMatchEstimator = conf.getUntrackedParameter<double>("maxSeedMatchEstimator", 0);
0088 rejectOverlaps = conf.getParameter<bool>("OverlapCleaning");
0089 splitHits = conf.getParameter<bool>("SplitHits");
0090 }
0091
0092 void TrackCandidateProducer::produce(edm::Event& e, const edm::EventSetup& es) {
0093
0094 auto const& magneticField = es.getData(magneticFieldESToken_);
0095 auto const& trackerGeometry = es.getData(trackerGeometryESToken_);
0096 auto const& trackerTopology = es.getData(trackerTopologyESToken_);
0097 auto const& propagator = es.getData(propagatorESToken_);
0098
0099
0100 edm::Handle<edm::View<TrajectorySeed> > seeds;
0101 e.getByToken(seedToken, seeds);
0102
0103 edm::Handle<FastTrackerRecHitCombinationCollection> recHitCombinations;
0104 e.getByToken(recHitCombinationsToken, recHitCombinations);
0105
0106 edm::Handle<edm::SimTrackContainer> simTracks;
0107 e.getByToken(simTrackToken, simTracks);
0108
0109 edm::Handle<std::vector<bool> > hitMasks;
0110 if (!hitMasksToken.isUninitialized()) {
0111 e.getByToken(hitMasksToken, hitMasks);
0112 }
0113
0114
0115 std::unique_ptr<TrackCandidateCollection> output(new TrackCandidateCollection);
0116
0117
0118 for (unsigned seedIndex = 0; seedIndex < seeds->size(); ++seedIndex) {
0119 const TrajectorySeed seed = (*seeds)[seedIndex];
0120 std::vector<int32_t> recHitCombinationIndices;
0121
0122
0123 if (seed.nHits() == 0) {
0124 recHitCombinationIndices = SeedMatcher::matchRecHitCombinations(
0125 seed, *recHitCombinations, *simTracks, maxSeedMatchEstimator, propagator, magneticField, trackerGeometry);
0126 }
0127
0128 else {
0129 int32_t icomb = fastTrackingUtilities::getRecHitCombinationIndex(seed);
0130 recHitCombinationIndices.push_back(icomb);
0131 }
0132
0133
0134 for (auto icomb : recHitCombinationIndices) {
0135 if (icomb < 0 || unsigned(icomb) >= recHitCombinations->size()) {
0136 throw cms::Exception("TrackCandidateProducer")
0137 << " found seed with recHitCombination out or range: " << icomb << std::endl;
0138 }
0139 const FastTrackerRecHitCombination& recHitCombination = (*recHitCombinations)[icomb];
0140
0141
0142 std::vector<const FastTrackerRecHit*> selectedRecHits;
0143
0144
0145 for (auto const& recHit : seed.recHits()) {
0146 selectedRecHits.push_back(static_cast<const FastTrackerRecHit*>(&recHit));
0147 }
0148
0149
0150 const FastTrackerRecHit* lastHitToSkip = nullptr;
0151 if (!selectedRecHits.empty()) {
0152 lastHitToSkip = selectedRecHits.back();
0153 }
0154
0155
0156 bool hitsAlongMomentum = (seed.direction() == alongMomentum);
0157
0158
0159 for (unsigned hitIndex = hitsAlongMomentum ? 0 : recHitCombination.size() - 1;
0160 hitIndex < recHitCombination.size();
0161 hitsAlongMomentum ? ++hitIndex : --hitIndex) {
0162 const FastTrackerRecHit* selectedRecHit = recHitCombination[hitIndex].get();
0163
0164
0165 if (lastHitToSkip) {
0166 if (lastHitToSkip->sameId(selectedRecHit)) {
0167 lastHitToSkip = nullptr;
0168 }
0169 continue;
0170 }
0171
0172
0173 if (hitMasks.isValid() && fastTrackingUtilities::hitIsMasked(selectedRecHit, *hitMasks)) {
0174 continue;
0175 }
0176
0177
0178
0179
0180 if (!rejectOverlaps || selectedRecHits.empty() ||
0181 (TrackingLayer::createFromDetId(selectedRecHits.back()->geographicalId(), trackerTopology) !=
0182 TrackingLayer::createFromDetId(selectedRecHit->geographicalId(), trackerTopology))) {
0183 selectedRecHits.push_back(selectedRecHit);
0184 }
0185
0186
0187
0188
0189 else if (fastTrackingUtilities::hitLocalError(selectedRecHit) <
0190 fastTrackingUtilities::hitLocalError(selectedRecHits.back())) {
0191 selectedRecHits.back() = selectedRecHit;
0192 }
0193 }
0194
0195
0196 edm::OwnVector<TrackingRecHit> hitsForTrackCandidate;
0197 for (unsigned index = 0; index < selectedRecHits.size(); ++index) {
0198 if (splitHits) {
0199
0200 hitSplitter.split(*selectedRecHits[index], hitsForTrackCandidate, hitsAlongMomentum);
0201 } else {
0202 hitsForTrackCandidate.push_back(selectedRecHits[index]->clone());
0203 }
0204 }
0205
0206
0207 fastTrackingUtilities::setRecHitCombinationIndex(hitsForTrackCandidate, icomb);
0208
0209
0210
0211 DetId seedDetId(seed.startingState().detId());
0212 const GeomDet* gdet = trackerGeometry.idToDet(seedDetId);
0213 TrajectoryStateOnSurface seedTSOS =
0214 trajectoryStateTransform::transientState(seed.startingState(), &(gdet->surface()), &magneticField);
0215
0216 const GeomDet* initialLayer = trackerGeometry.idToDet(hitsForTrackCandidate.front().geographicalId());
0217 const TrajectoryStateOnSurface initialTSOS = propagator.propagate(seedTSOS, initialLayer->surface());
0218
0219 if (!initialTSOS.isValid())
0220 continue;
0221 PTrajectoryStateOnDet PTSOD = trajectoryStateTransform::persistentState(
0222 initialTSOS, hitsForTrackCandidate.front().geographicalId().rawId());
0223
0224
0225 output->push_back(
0226 TrackCandidate(hitsForTrackCandidate, seed, PTSOD, edm::RefToBase<TrajectorySeed>(seeds, seedIndex)));
0227 }
0228 }
0229
0230
0231 e.put(std::move(output));
0232 }
0233
0234 #include "FWCore/Framework/interface/MakerMacros.h"
0235 DEFINE_FWK_MODULE(TrackCandidateProducer);