All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Classes | Public Types | Public Member Functions | Private Member Functions | Private Attributes | List of all members
trkf::KalmanFilterTrajectoryFitter Class Reference
Inheritance diagram for trkf::KalmanFilterTrajectoryFitter:

Classes

struct  Config
 
struct  Inputs
 
struct  Options
 

Public Types

using Parameters = art::EDProducer::Table< Config >
 

Public Member Functions

 KalmanFilterTrajectoryFitter (Parameters const &p)
 
 KalmanFilterTrajectoryFitter (KalmanFilterTrajectoryFitter const &)=delete
 
 KalmanFilterTrajectoryFitter (KalmanFilterTrajectoryFitter &&)=delete
 
KalmanFilterTrajectoryFitteroperator= (KalmanFilterTrajectoryFitter const &)=delete
 
KalmanFilterTrajectoryFitteroperator= (KalmanFilterTrajectoryFitter &&)=delete
 

Private Member Functions

void produce (art::Event &e) override
 
double setMomValue (const recob::TrackTrajectory *ptraj, const double pMC, const int pId) const
 
int setPId () const
 
bool setDirFlip (const recob::TrackTrajectory *ptraj, TVector3 &mcdir) const
 
void restoreInputPoints (const recob::TrackTrajectory &track, const std::vector< art::Ptr< recob::Hit >> &inHits, recob::Track &outTrack, std::vector< art::Ptr< recob::Hit >> &outHits) const
 

Private Attributes

Parameters p_
 
TrackStatePropagator prop
 
trkf::TrackKalmanFitter kalmanFitter
 
trkf::TrackMomentumCalculator tmc {}
 
art::InputTag trajectoryInputTag
 
art::InputTag simTrackInputTag
 
bool isTT
 

Detailed Description

Definition at line 34 of file KalmanFilterTrajectoryFitter_module.cc.

Member Typedef Documentation

using trkf::KalmanFilterTrajectoryFitter::Parameters = art::EDProducer::Table<Config>

Definition at line 104 of file KalmanFilterTrajectoryFitter_module.cc.

Constructor & Destructor Documentation

KalmanFilterTrajectoryFitter::KalmanFilterTrajectoryFitter ( Parameters const &  p)
explicit

Definition at line 138 of file KalmanFilterTrajectoryFitter_module.cc.

140  : EDProducer{p}
141  , p_(p)
142  , prop{p_().propagator}
143  , kalmanFitter{&prop, p_().fitter}
144  , trajectoryInputTag{p_().inputs().inputTrajectoryLabel()}
145 {
146  if (p_().options().pFromMC() || p_().options().dirFromMC())
147  simTrackInputTag = art::InputTag{p_().inputs().inputMCLabel()};
148 
149  isTT = p_().inputs().isTrackTrajectory();
150 
151  produces<std::vector<recob::Track>>();
152  produces<art::Assns<recob::Track, recob::Hit>>();
153  produces<art::Assns<recob::Track, recob::Hit, recob::TrackHitMeta>>();
154  if (isTT) { produces<art::Assns<recob::TrackTrajectory, recob::Track>>(); }
155  else {
156  produces<art::Assns<recob::Trajectory, recob::Track>>();
157  }
158  if (p_().options().produceTrackFitHitInfo()) {
159  produces<std::vector<std::vector<recob::TrackFitHitInfo>>>();
160  }
161  if (p_().options().produceSpacePoints()) {
162  produces<std::vector<recob::SpacePoint>>();
163  produces<art::Assns<recob::Hit, recob::SpacePoint>>();
164  }
165 
166  //throw expections to avoid possible silent failures due to incompatible configuration options
167 
168  unsigned int nDirs = 0;
169  if (p_().options().dirFromMC()) nDirs++;
170  if (p_().options().dirFromVec()) nDirs++;
171  if (p_().options().alwaysInvertDir()) nDirs++;
172  if (nDirs > 1) {
173  throw cet::exception("KalmanFilterTrajectoryFitter")
174  << "Incompatible configuration parameters: only at most one can be set to true among "
175  "dirFromMC, dirFromVec, and alwaysInvertDir."
176  << "\n";
177  }
178 
179  unsigned int nPFroms = 0;
180  if (p_().options().pFromLength()) nPFroms++;
181  if (p_().options().pFromMC()) nPFroms++;
182  if (nPFroms > 1) {
183  throw cet::exception("KalmanFilterTrajectoryFitter")
184  << "Incompatible configuration parameters: only at most one can be set to true among "
185  "pFromLength, and pFromMC."
186  << "\n"; //pFromMSChi2,
187  }
188 
189  if (p_().options().keepInputTrajectoryPoints()) {
190  if (p_().fitter().sortHitsByPlane() || p_().fitter().sortOutputHitsMinLength() ||
191  p_().fitter().skipNegProp()) {
192  throw cet::exception("KalmanFilterTrajectoryFitter")
193  << "Incompatible configuration parameters: keepInputTrajectoryPoints needs the following "
194  "fitter options all set to false: sortHitsByPlane, sortOutputHitsMinLength, skipNegProp."
195  << "\n";
196  }
197  }
198 }
pdgs p
Definition: selectors.fcl:22
then echo echo For and will not be changed by echo further linking echo echo B echo The symbol is in the uninitialized data multiple common symbols may appear with the echo same name If the symbol is defined the common echo symbols are treated as undefined references For more echo details on common see the discussion of warn common echo in *Note Linker options
trkf::KalmanFilterTrajectoryFitter::KalmanFilterTrajectoryFitter ( KalmanFilterTrajectoryFitter const &  )
delete
trkf::KalmanFilterTrajectoryFitter::KalmanFilterTrajectoryFitter ( KalmanFilterTrajectoryFitter &&  )
delete

Member Function Documentation

KalmanFilterTrajectoryFitter& trkf::KalmanFilterTrajectoryFitter::operator= ( KalmanFilterTrajectoryFitter const &  )
delete
KalmanFilterTrajectoryFitter& trkf::KalmanFilterTrajectoryFitter::operator= ( KalmanFilterTrajectoryFitter &&  )
delete
void KalmanFilterTrajectoryFitter::produce ( art::Event &  e)
overrideprivate

Definition at line 201 of file KalmanFilterTrajectoryFitter_module.cc.

202 {
203 
204  auto outputTracks = std::make_unique<std::vector<recob::Track>>();
205  auto outputHitsMeta =
206  std::make_unique<art::Assns<recob::Track, recob::Hit, recob::TrackHitMeta>>();
207  auto outputHits = std::make_unique<art::Assns<recob::Track, recob::Hit>>();
208  auto outputHitInfo = std::make_unique<std::vector<std::vector<recob::TrackFitHitInfo>>>();
209 
210  //only one will be filled and pushed into the event:
211  auto outputTTjTAssn = std::make_unique<art::Assns<recob::TrackTrajectory, recob::Track>>();
212  auto outputTjTAssn = std::make_unique<art::Assns<recob::Trajectory, recob::Track>>();
213 
214  auto const tid = e.getProductID<std::vector<recob::Track>>();
215  auto const tidgetter = e.productGetter(tid);
216 
217  auto outputSpacePoints = std::make_unique<std::vector<recob::SpacePoint>>();
218  auto outputHitSpacePointAssn = std::make_unique<art::Assns<recob::Hit, recob::SpacePoint>>();
219  auto const spid = e.getProductID<std::vector<recob::SpacePoint>>();
220  auto const spidgetter = e.productGetter(spid);
221 
222  //FIXME, eventually remove this (ok only for single particle MC)
223  double pMC = -1.;
224  TVector3 mcdir;
225  if (p_().options().pFromMC() || p_().options().dirFromMC()) {
226  art::ValidHandle<std::vector<sim::MCTrack>> simTracks =
227  e.getValidHandle<std::vector<sim::MCTrack>>(simTrackInputTag);
228  for (unsigned int iMC = 0; iMC < simTracks->size(); ++iMC) {
229  const sim::MCTrack& mctrack = simTracks->at(iMC);
230  //fiducial cuts on MC tracks
231  if (mctrack.PdgCode() != 13) continue;
232  if (mctrack.Process() != "primary") continue;
233  pMC = mctrack.Start().Momentum().P() * 0.001;
234  mcdir = TVector3(mctrack.Start().Momentum().X() * 0.001 / pMC,
235  mctrack.Start().Momentum().Y() * 0.001 / pMC,
236  mctrack.Start().Momentum().Z() * 0.001 / pMC);
237  break;
238  }
239  //std::cout << "mc momentum value = " << pval << " GeV" << std::endl;
240  }
241 
242  unsigned int nTrajs = 0;
243 
244  art::Handle<std::vector<recob::TrackTrajectory>> inputTrackTrajectoryH;
245  art::Handle<std::vector<recob::Trajectory>> inputTrajectoryH;
246  const std::vector<recob::TrackTrajectory>* trackTrajectoryVec = nullptr;
247  const std::vector<recob::Trajectory>* trajectoryVec = nullptr;
248  const art::Assns<recob::TrackTrajectory, recob::Hit>* trackTrajectoryHitsAssn = nullptr;
249  const art::Assns<recob::Trajectory, recob::Hit>* trajectoryHitsAssn = nullptr;
250  if (isTT) {
251  bool ok = e.getByLabel(trajectoryInputTag, inputTrackTrajectoryH);
252  if (!ok)
253  throw cet::exception("KalmanFilterTrajectoryFitter")
254  << "Cannot find recob::TrackTrajectory art::Handle with inputTag " << trajectoryInputTag;
255  trackTrajectoryVec = inputTrackTrajectoryH.product();
256  trackTrajectoryHitsAssn =
257  e.getValidHandle<art::Assns<recob::TrackTrajectory, recob::Hit>>(trajectoryInputTag)
258  .product();
259  nTrajs = trackTrajectoryVec->size();
260  }
261  else {
262  bool ok = e.getByLabel(trajectoryInputTag, inputTrajectoryH);
263  if (!ok)
264  throw cet::exception("KalmanFilterTrajectoryFitter")
265  << "Cannot find recob::Trajectory art::Handle with inputTag " << trajectoryInputTag;
266  trajectoryVec = inputTrajectoryH.product();
267  trajectoryHitsAssn =
268  e.getValidHandle<art::Assns<recob::Trajectory, recob::Hit>>(trajectoryInputTag).product();
269  nTrajs = trajectoryVec->size();
270  }
271 
272  auto const detProp = art::ServiceHandle<detinfo::DetectorPropertiesService const>()->DataFor(e);
273 
274  for (unsigned int iTraj = 0; iTraj < nTrajs; ++iTraj) {
275 
276  const recob::TrackTrajectory& inTraj =
277  (isTT ? trackTrajectoryVec->at(iTraj) :
278  recob::TrackTrajectory(trajectoryVec->at(iTraj),
279  std::vector<recob::TrajectoryPointFlags>()));
280  // this is not computationally optimal, but at least preserves the order
281  // unlike FindManyP
282  std::vector<art::Ptr<recob::Hit>> inHits;
283  if (isTT) {
284  for (auto it = trackTrajectoryHitsAssn->begin(); it != trackTrajectoryHitsAssn->end(); ++it) {
285  if (it->first.key() == iTraj)
286  inHits.push_back(it->second);
287  else if (inHits.size() > 0)
288  break;
289  }
290  }
291  else {
292  for (auto it = trajectoryHitsAssn->begin(); it != trajectoryHitsAssn->end(); ++it) {
293  if (it->first.key() == iTraj)
294  inHits.push_back(it->second);
295  else if (inHits.size() > 0)
296  break;
297  }
298  }
299  const int pId = setPId();
300  const double mom = setMomValue(&inTraj, pMC, pId);
301  const bool flipDir = setDirFlip(&inTraj, mcdir);
302 
303  recob::Track outTrack;
304  std::vector<art::Ptr<recob::Hit>> outHits;
305  trkmkr::OptionalOutputs optionals;
306  if (p_().options().produceTrackFitHitInfo()) optionals.initTrackFitInfos();
307  bool fitok = kalmanFitter.fitTrack(detProp,
308  inTraj,
309  iTraj,
310  SMatrixSym55(),
311  SMatrixSym55(),
312  inHits, // inFlags,
313  mom,
314  pId,
315  flipDir,
316  outTrack,
317  outHits,
318  optionals);
319  if (!fitok) continue;
320 
321  if (p_().options().keepInputTrajectoryPoints()) {
322  restoreInputPoints(inTraj, inHits, outTrack, outHits);
323  }
324 
325  outputTracks->emplace_back(std::move(outTrack));
326  art::Ptr<recob::Track> aptr(tid, outputTracks->size() - 1, tidgetter);
327  unsigned int ip = 0;
328  for (auto const& trhit : outHits) {
329  //the fitter produces collections with 1-1 match between hits and point
330  recob::TrackHitMeta metadata(ip, -1);
331  outputHitsMeta->addSingle(aptr, trhit, metadata);
332  outputHits->addSingle(aptr, trhit);
333  if (p_().options().produceSpacePoints() && outputTracks->back().HasValidPoint(ip)) {
334  auto& tp = outputTracks->back().Trajectory().LocationAtPoint(ip);
335  double fXYZ[3] = {tp.X(), tp.Y(), tp.Z()};
336  double fErrXYZ[6] = {0};
337  recob::SpacePoint sp(fXYZ, fErrXYZ, -1.);
338  outputSpacePoints->emplace_back(std::move(sp));
339  art::Ptr<recob::SpacePoint> apsp(spid, outputSpacePoints->size() - 1, spidgetter);
340  outputHitSpacePointAssn->addSingle(trhit, apsp);
341  }
342  ip++;
343  }
344  outputHitInfo->emplace_back(optionals.trackFitHitInfos());
345  if (isTT) {
346  outputTTjTAssn->addSingle(art::Ptr<recob::TrackTrajectory>(inputTrackTrajectoryH, iTraj),
347  aptr);
348  }
349  else {
350  outputTjTAssn->addSingle(art::Ptr<recob::Trajectory>(inputTrajectoryH, iTraj), aptr);
351  }
352  }
353  e.put(std::move(outputTracks));
354  e.put(std::move(outputHitsMeta));
355  e.put(std::move(outputHits));
356  if (p_().options().produceTrackFitHitInfo()) { e.put(std::move(outputHitInfo)); }
357  if (p_().options().produceSpacePoints()) {
358  e.put(std::move(outputSpacePoints));
359  e.put(std::move(outputHitSpacePointAssn));
360  }
361  if (isTT)
362  e.put(std::move(outputTTjTAssn));
363  else
364  e.put(std::move(outputTjTAssn));
365 }
void restoreInputPoints(const recob::TrackTrajectory &track, const std::vector< art::Ptr< recob::Hit >> &inHits, recob::Track &outTrack, std::vector< art::Ptr< recob::Hit >> &outHits) const
void initTrackFitInfos()
initialize the output vector of TrackFitHitInfos
Definition: TrackMaker.h:161
recob::tracking::SMatrixSym55 SMatrixSym55
Definition: TrackState.h:15
bool fitTrack(detinfo::DetectorPropertiesData const &detProp, const recob::TrackTrajectory &traj, int tkID, const SMatrixSym55 &covVtx, const SMatrixSym55 &covEnd, const std::vector< art::Ptr< recob::Hit >> &hits, const double pval, const int pdgid, const bool flipDirection, recob::Track &outTrack, std::vector< art::Ptr< recob::Hit >> &outHits, trkmkr::OptionalOutputs &optionals) const
Fit track starting from TrackTrajectory.
Data related to recob::Hit associated with recob::Track.The purpose is to collect several variables t...
Definition: TrackHitMeta.h:43
A trajectory in space reconstructed from hits.
double setMomValue(const recob::TrackTrajectory *ptraj, const double pMC, const int pId) const
int PdgCode() const
Definition: MCTrack.h:41
const TLorentzVector & Momentum() const
Definition: MCStep.h:38
std::vector< recob::TrackFitHitInfo > trackFitHitInfos()
get the output vector of TrackFitHitInfos by releasing and moving
Definition: TrackMaker.h:185
const std::string & Process() const
Definition: MCTrack.h:43
then echo echo For and will not be changed by echo further linking echo echo B echo The symbol is in the uninitialized data multiple common symbols may appear with the echo same name If the symbol is defined the common echo symbols are treated as undefined references For more echo details on common see the discussion of warn common echo in *Note Linker options
do i e
const MCStep & Start() const
Definition: MCTrack.h:44
bool setDirFlip(const recob::TrackTrajectory *ptraj, TVector3 &mcdir) const
Struct holding optional TrackMaker outputs.
Definition: TrackMaker.h:125
auto const detProp
Track from a non-cascading particle.A recob::Track consists of a recob::TrackTrajectory, plus additional members relevant for a &quot;fitted&quot; track:
void KalmanFilterTrajectoryFitter::restoreInputPoints ( const recob::TrackTrajectory track,
const std::vector< art::Ptr< recob::Hit >> &  inHits,
recob::Track outTrack,
std::vector< art::Ptr< recob::Hit >> &  outHits 
) const
private

Definition at line 368 of file KalmanFilterTrajectoryFitter_module.cc.

373 {
374  const auto np = outTrack.NumberTrajectoryPoints();
375  std::vector<Point_t> positions(np);
376  std::vector<Vector_t> momenta(np);
377  std::vector<recob::TrajectoryPointFlags> outFlags(np);
378  //
379  for (unsigned int p = 0; p < np; ++p) {
380  auto flag = outTrack.FlagsAtPoint(p);
381  auto mom = outTrack.VertexMomentum();
382  auto op = flag.fromHit();
383  positions[op] = track.LocationAtPoint(op);
384  momenta[op] = mom * track.DirectionAtPoint(op);
385  auto mask = flag.mask();
388  outFlags[op] = recob::TrajectoryPointFlags(op, mask);
389  }
390  auto covs = outTrack.Covariances();
391  outTrack = recob::Track(
392  recob::TrackTrajectory(std::move(positions), std::move(momenta), std::move(outFlags), true),
393  outTrack.ParticleId(),
394  outTrack.Chi2(),
395  outTrack.Ndof(),
396  std::move(covs.first),
397  std::move(covs.second),
398  outTrack.ID());
399  //
400  outHits.clear();
401  for (auto h : inHits)
402  outHits.push_back(h);
403 }
double VertexMomentum() const
static constexpr Flag_t NoPoint
The trajectory point is not defined.
T DirectionAtPoint(unsigned int p) const
Direction at point p. Use e.g. as:
then echo unknown compiler flag
pdgs p
Definition: selectors.fcl:22
size_t NumberTrajectoryPoints() const
Various functions related to the presence and the number of (valid) points.
std::pair< SMatrixSym55, SMatrixSym55 > Covariances() const
while getopts h
constexpr mask_t< EnumType > mask(EnumType bit, OtherBits...otherBits)
Returns a mask with all specified bits set.
A trajectory in space reconstructed from hits.
T LocationAtPoint(unsigned int p) const
Position at point p. Use e.g. as:
PointFlags_t const & FlagsAtPoint(size_t i) const
TrackCollectionProxyElement< TrackCollProxy > Track
Proxy to an element of a proxy collection of recob::Track objects.
Set of flags pertaining a point of the track.
bool KalmanFilterTrajectoryFitter::setDirFlip ( const recob::TrackTrajectory ptraj,
TVector3 &  mcdir 
) const
private

Definition at line 426 of file KalmanFilterTrajectoryFitter_module.cc.

428 {
429  bool result = false;
430  if (p_().options().alwaysInvertDir()) { return true; }
431  else if (p_().options().dirFromMC()) {
432  auto tdir = ptraj->VertexDirection();
433  if ((mcdir.X() * tdir.X() + mcdir.Y() * tdir.Y() + mcdir.Z() * tdir.Z()) < 0.) result = true;
434  }
435  else if (p_().options().dirFromVec()) {
436  std::array<float, 3> dir = p_().options().dirVec();
437  auto tdir = ptraj->VertexDirection();
438  if ((dir[0] * tdir.X() + dir[1] * tdir.Y() + dir[2] * tdir.Z()) < 0.) result = true;
439  }
440  return result;
441 }
Vector_t VertexDirection() const
Returns the direction of the trajectory at the first point.
tuple dir
Definition: dropbox.py:28
then echo echo For and will not be changed by echo further linking echo echo B echo The symbol is in the uninitialized data multiple common symbols may appear with the echo same name If the symbol is defined the common echo symbols are treated as undefined references For more echo details on common see the discussion of warn common echo in *Note Linker options
double KalmanFilterTrajectoryFitter::setMomValue ( const recob::TrackTrajectory ptraj,
const double  pMC,
const int  pId 
) const
private

Definition at line 406 of file KalmanFilterTrajectoryFitter_module.cc.

409 {
410  double result = p_().options().pval();
411  if (p_().options().pFromLength()) { result = tmc.GetTrackMomentum(ptraj->Length(), pId); }
412  else if (p_().options().pFromMC() && pMC > 0.) {
413  result = pMC;
414  }
415  return result;
416 }
double Length(size_t startAt=0) const
Returns the approximate length of the trajectory.
then echo echo For and will not be changed by echo further linking echo echo B echo The symbol is in the uninitialized data multiple common symbols may appear with the echo same name If the symbol is defined the common echo symbols are treated as undefined references For more echo details on common see the discussion of warn common echo in *Note Linker options
double GetTrackMomentum(double trkrange, int pdg) const
int KalmanFilterTrajectoryFitter::setPId ( ) const
private

Definition at line 419 of file KalmanFilterTrajectoryFitter_module.cc.

420 {
421  int result = p_().options().pdgId();
422  return result;
423 }

Member Data Documentation

bool trkf::KalmanFilterTrajectoryFitter::isTT
private

Definition at line 125 of file KalmanFilterTrajectoryFitter_module.cc.

trkf::TrackKalmanFitter trkf::KalmanFilterTrajectoryFitter::kalmanFitter
private

Definition at line 119 of file KalmanFilterTrajectoryFitter_module.cc.

Parameters trkf::KalmanFilterTrajectoryFitter::p_
private

Definition at line 117 of file KalmanFilterTrajectoryFitter_module.cc.

TrackStatePropagator trkf::KalmanFilterTrajectoryFitter::prop
private

Definition at line 118 of file KalmanFilterTrajectoryFitter_module.cc.

art::InputTag trkf::KalmanFilterTrajectoryFitter::simTrackInputTag
private

Definition at line 123 of file KalmanFilterTrajectoryFitter_module.cc.

trkf::TrackMomentumCalculator trkf::KalmanFilterTrajectoryFitter::tmc {}
private

Definition at line 120 of file KalmanFilterTrajectoryFitter_module.cc.

art::InputTag trkf::KalmanFilterTrajectoryFitter::trajectoryInputTag
private

Definition at line 122 of file KalmanFilterTrajectoryFitter_module.cc.


The documentation for this class was generated from the following file: