File indexing completed on 2023-03-17 11:22:19
0001 #include "RecoTracker/FinalTrackSelectors/interface/TrackCollectionCloner.h"
0002
0003 void TrackCollectionCloner::fill(edm::ParameterSetDescription& desc) {
0004 desc.addUntracked<bool>("copyExtras", true);
0005 desc.addUntracked<bool>("copyTrajectories", false);
0006 }
0007
0008 TrackCollectionCloner::Producer::Producer(edm::Event& ievt, TrackCollectionCloner const& cloner)
0009 : copyExtras_(cloner.copyExtras_), copyTrajectories_(cloner.copyTrajectories_), evt(ievt) {
0010 selTracks_ = std::make_unique<reco::TrackCollection>();
0011 if (copyExtras_) {
0012 selTrackExtras_ = std::make_unique<reco::TrackExtraCollection>();
0013 selHits_ = std::make_unique<TrackingRecHitCollection>();
0014 }
0015 if (copyTrajectories_) {
0016 selTrajs_ = std::make_unique<std::vector<Trajectory> >();
0017 selTTAss_ = std::make_unique<TrajTrackAssociationCollection>(evt.getRefBeforePut<std::vector<Trajectory> >(),
0018 evt.getRefBeforePut<reco::TrackCollection>());
0019 }
0020 }
0021
0022
0023 void TrackCollectionCloner::Producer::operator()(Tokens const& tokens, std::vector<unsigned int> const& selected) {
0024 auto rTracks = evt.getRefBeforePut<reco::TrackCollection>();
0025
0026 TrackingRecHitRefProd rHits;
0027 reco::TrackExtraRefProd rTrackExtras;
0028 if (copyExtras_) {
0029 rHits = evt.getRefBeforePut<TrackingRecHitCollection>();
0030 rTrackExtras = evt.getRefBeforePut<reco::TrackExtraCollection>();
0031 }
0032
0033 edm::RefProd<std::vector<Trajectory> > trajRefProd;
0034 if (copyTrajectories_) {
0035 trajRefProd = evt.getRefBeforePut<std::vector<Trajectory> >();
0036 }
0037
0038 std::vector<Trajectory> dummy;
0039
0040 auto const& trajIn = copyTrajectories_ ? tokens.trajectories(evt) : dummy;
0041
0042 auto const& tracksIn = tokens.tracks(evt);
0043 for (auto k : selected) {
0044 auto const& trk = tracksIn[k];
0045 selTracks_->emplace_back(trk);
0046 if (copyTrajectories_) {
0047
0048 selTrajs_->emplace_back(trajIn[k]);
0049 assert(selTrajs_->back().measurements().size() == trk.recHitsSize());
0050 selTTAss_->insert(edm::Ref<std::vector<Trajectory> >(trajRefProd, selTrajs_->size() - 1),
0051 reco::TrackRef(rTracks, selTracks_->size() - 1));
0052 }
0053
0054 if (!copyExtras_)
0055 continue;
0056
0057
0058 selTrackExtras_->emplace_back(trk.outerPosition(),
0059 trk.outerMomentum(),
0060 trk.outerOk(),
0061 trk.innerPosition(),
0062 trk.innerMomentum(),
0063 trk.innerOk(),
0064 trk.outerStateCovariance(),
0065 trk.outerDetId(),
0066 trk.innerStateCovariance(),
0067 trk.innerDetId(),
0068 trk.seedDirection(),
0069 trk.seedRef());
0070 selTracks_->back().setExtra(reco::TrackExtraRef(rTrackExtras, selTrackExtras_->size() - 1));
0071 auto& tx = selTrackExtras_->back();
0072 tx.setResiduals(trk.residuals());
0073 auto nh1 = trk.recHitsSize();
0074 tx.setHits(rHits, selHits_->size(), nh1);
0075 tx.setTrajParams(trk.extra()->trajParams(), trk.extra()->chi2sX5());
0076 assert(tx.trajParams().size() == tx.recHitsSize());
0077
0078 for (auto hit = trk.recHitsBegin(); hit != trk.recHitsEnd(); ++hit) {
0079 selHits_->push_back((*hit)->clone());
0080 }
0081 }
0082 }
0083
0084 TrackCollectionCloner::Producer::~Producer() {
0085 selTracks_->shrink_to_fit();
0086 auto tsize = selTracks_->size();
0087 evt.put(std::move(selTracks_));
0088 if (copyExtras_) {
0089 selTrackExtras_->shrink_to_fit();
0090 assert(selTrackExtras_->size() == tsize);
0091 selHits_->shrink_to_fit();
0092 evt.put(std::move(selTrackExtras_));
0093 evt.put(std::move(selHits_));
0094 }
0095 if (copyTrajectories_) {
0096 selTrajs_->shrink_to_fit();
0097 assert(selTrajs_->size() == tsize);
0098 assert(selTTAss_->size() == tsize);
0099 evt.put(std::move(selTrajs_));
0100 evt.put(std::move(selTTAss_));
0101 }
0102 }