CMS 3D CMS Logo

TrackCollectionCloner.cc
Go to the documentation of this file.
2 
4  desc.addUntracked<bool>("copyExtras", true);
5  desc.addUntracked<bool>("copyTrajectories", false);
6 }
7 
9  : copyExtras_(cloner.copyExtras_), copyTrajectories_(cloner.copyTrajectories_), evt(ievt) {
10  selTracks_ = std::make_unique<reco::TrackCollection>();
11  if (copyExtras_) {
12  selTrackExtras_ = std::make_unique<reco::TrackExtraCollection>();
13  selHits_ = std::make_unique<TrackingRecHitCollection>();
14  }
15  if (copyTrajectories_) {
16  selTrajs_ = std::make_unique<std::vector<Trajectory> >();
17  selTTAss_ = std::make_unique<TrajTrackAssociationCollection>(evt.getRefBeforePut<std::vector<Trajectory> >(),
19  }
20 }
21 
23 void TrackCollectionCloner::Producer::operator()(Tokens const& tokens, std::vector<unsigned int> const& selected) {
24  auto rTracks = evt.getRefBeforePut<reco::TrackCollection>();
25 
27  reco::TrackExtraRefProd rTrackExtras;
28  if (copyExtras_) {
29  rHits = evt.getRefBeforePut<TrackingRecHitCollection>();
30  rTrackExtras = evt.getRefBeforePut<reco::TrackExtraCollection>();
31  }
32 
34  if (copyTrajectories_) {
35  trajRefProd = evt.getRefBeforePut<std::vector<Trajectory> >();
36  }
37 
38  std::vector<Trajectory> dummy;
39 
40  auto const& trajIn = copyTrajectories_ ? tokens.trajectories(evt) : dummy;
41 
42  auto const& tracksIn = tokens.tracks(evt);
43  for (auto k : selected) {
44  auto const& trk = tracksIn[k];
45  selTracks_->emplace_back(trk); // clone and store
46  if (copyTrajectories_) {
47  // we assume tracks and trajectories are one-to-one and the assocMap is useless
48  selTrajs_->emplace_back(trajIn[k]);
49  assert(selTrajs_->back().measurements().size() == trk.recHitsSize());
50  selTTAss_->insert(edm::Ref<std::vector<Trajectory> >(trajRefProd, selTrajs_->size() - 1),
51  reco::TrackRef(rTracks, selTracks_->size() - 1));
52  }
53 
54  if (!copyExtras_)
55  continue;
56 
57  // TrackExtras
58  selTrackExtras_->emplace_back(trk.outerPosition(),
59  trk.outerMomentum(),
60  trk.outerOk(),
61  trk.innerPosition(),
62  trk.innerMomentum(),
63  trk.innerOk(),
64  trk.outerStateCovariance(),
65  trk.outerDetId(),
66  trk.innerStateCovariance(),
67  trk.innerDetId(),
68  trk.seedDirection(),
69  trk.seedRef());
70  selTracks_->back().setExtra(reco::TrackExtraRef(rTrackExtras, selTrackExtras_->size() - 1));
71  auto& tx = selTrackExtras_->back();
72  tx.setResiduals(trk.residuals());
73  auto nh1 = trk.recHitsSize();
74  tx.setHits(rHits, selHits_->size(), nh1);
75  tx.setTrajParams(trk.extra()->trajParams(), trk.extra()->chi2sX5());
76  assert(tx.trajParams().size() == tx.recHitsSize());
77  // TrackingRecHits
78  for (auto hit = trk.recHitsBegin(); hit != trk.recHitsEnd(); ++hit) {
79  selHits_->push_back((*hit)->clone());
80  }
81  }
82 }
83 
85  selTracks_->shrink_to_fit();
86  auto tsize = selTracks_->size();
87  evt.put(std::move(selTracks_));
88  if (copyExtras_) {
89  selTrackExtras_->shrink_to_fit();
90  assert(selTrackExtras_->size() == tsize);
91  selHits_->shrink_to_fit();
92  evt.put(std::move(selTrackExtras_));
93  evt.put(std::move(selHits_));
94  }
95  if (copyTrajectories_) {
96  selTrajs_->shrink_to_fit();
97  assert(selTrajs_->size() == tsize);
98  assert(selTTAss_->size() == tsize);
99  evt.put(std::move(selTrajs_));
100  evt.put(std::move(selTTAss_));
101  }
102 }
std::unique_ptr< std::vector< Trajectory > > selTrajs_
bool copyExtras_
copy only the tracks, not extras and rechits (for AOD)
std::vector< Track > TrackCollection
collection of Tracks
Definition: TrackFwd.h:14
Trktree trk
Definition: Trktree.cc:2
static void fill(edm::ParameterSetDescription &desc)
assert(be >=bs)
std::unique_ptr< reco::TrackExtraCollection > selTrackExtras_
std::unique_ptr< reco::TrackCollection > selTracks_
std::unique_ptr< TrackingRecHitCollection > selHits_
bool copyTrajectories_
copy also trajectories and trajectory->track associations
RefProd< PROD > getRefBeforePut()
Definition: Event.h:158
std::vector< TrackExtra > TrackExtraCollection
collection of TrackExtra objects
Definition: TrackExtraFwd.h:10
std::vector< Trajectory > const & trajectories(edm::Event &evt) const
bool copyExtras_
copy only the tracks, not extras and rechits (for AOD)
std::unique_ptr< TrajTrackAssociationCollection > selTTAss_
Producer(edm::Event &ievt, TrackCollectionCloner const &cloner)
reco::TrackCollection const & tracks(edm::Event &evt) const
def move(src, dest)
Definition: eostools.py:511
bool copyTrajectories_
copy also trajectories and trajectory->track associations
void operator()(Tokens const &tokens, std::vector< unsigned int > const &selected)
process one event