LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
trkf::KalmanFilterAlg Class Reference

#include "KalmanFilterAlg.h"

Public Member Functions

 KalmanFilterAlg (const fhicl::ParameterSet &pset)
 Constructor. More...
 
bool getTrace () const
 Trace config parameters. More...
 
int getPlane () const
 Preferred view plane. More...
 
void setTrace (bool trace)
 Set trace config parameter. More...
 
void setPlane (int plane)
 Set preferred view plane. More...
 
bool buildTrack (const KTrack &trk, KGTrack &trg, const Propagator &prop, const Propagator::PropDirection dir, KHitContainer &hits, bool linear) const
 Make a new track. More...
 
bool smoothTrack (KGTrack &trg, KGTrack *trg1, const Propagator &prop) const
 Smooth track. More...
 
bool extendTrack (KGTrack &trg, const Propagator &prop, KHitContainer &hits) const
 Add hits to existing track. More...
 
bool fitMomentumRange (const KGTrack &trg, KETrack &tremom) const
 Estimate track momentum using range. More...
 
bool fitMomentumMS (const KGTrack &trg, const Propagator &prop, KETrack &tremom) const
 Estimate track momentum using multiple scattering. More...
 
bool fitMomentum (const KGTrack &trg, const Propagator &prop, KETrack &tremom) const
 Estimate track momentum using either range or multiple scattering. More...
 
bool updateMomentum (const KETrack &tremom, const Propagator &prop, KGTrack &trg) const
 Set track momentum at each track surface. More...
 
bool smoothTrackIter (int niter, KGTrack &trg, const Propagator &prop) const
 Iteratively smooth a track. More...
 
void cleanTrack (KGTrack &trg) const
 Clean track by removing noise hits near endpoints. More...
 

Private Attributes

bool fTrace
 Trace flag. More...
 
double fMaxPErr
 Maximum pointing error for free propagation. More...
 
double fGoodPErr
 Pointing error threshold for switching to free propagation. More...
 
double fMaxIncChisq
 Maximum incremental chisquare to accept a hit. More...
 
double fMaxSeedIncChisq
 Maximum incremental chisquare to accept a hit in seed phase. More...
 
double fMaxSmoothIncChisq
 Maximum incremental chisquare to accept a hit in smooth phase. More...
 
double fMaxEndChisq
 Maximum incremental chisquare for endpoint hit. More...
 
int fMinLHits
 Minimum number of hits to turn off linearized propagation. More...
 
double fMaxLDist
 Maximum distance for linearized propagation. More...
 
double fMaxPredDist
 Maximum prediciton distance to accept a hit. More...
 
double fMaxSeedPredDist
 Maximum prediciton distance to accept a hit in seed phase. More...
 
double fMaxPropDist
 Maximum propagation distance to candidate surface. More...
 
double fMinSortDist
 Sort low distance threshold. More...
 
double fMaxSortDist
 Sort high distance threshold. More...
 
int fMaxSamePlane
 Maximum consecutive hits in same plane. More...
 
double fGapDist
 Minimum gap distance. More...
 
int fMaxNoiseHits
 Maximum number of hits in noise cluster. More...
 
double fMinSampleDist
 Minimum sample distance (for momentum measurement). More...
 
bool fFitMomRange
 Fit momentum using range. More...
 
bool fFitMomMS
 Fit momentum using multiple scattering. More...
 
bool fGTrace
 Graphical trace flag. More...
 
double fGTraceWW
 Window width. More...
 
double fGTraceWH
 Window height. More...
 
double fGTraceXMin
 Graphical trace minimum x. More...
 
double fGTraceXMax
 Graphical trace maximum x. More...
 
std::vector< double > fGTraceZMin
 Graphical trace minimum z for each view. More...
 
std::vector< double > fGTraceZMax
 Graphical trace maximum z for each view. More...
 
int fPlane
 Preferred view plane. More...
 
std::vector< std::unique_ptr< TCanvas > > fCanvases
 Graphical trace canvases. More...
 
std::vector< TVirtualPad * > fPads
 View pads in current canvas. More...
 
TVirtualPad * fInfoPad
 Information pad. More...
 
TPaveText * fMessages
 Message box. More...
 
std::map< int, TMarker * > fMarkerMap
 Markers in current canvas. More...
 

Detailed Description

Definition at line 68 of file KalmanFilterAlg.h.

Constructor & Destructor Documentation

trkf::KalmanFilterAlg::KalmanFilterAlg ( const fhicl::ParameterSet pset)

Constructor.

Definition at line 223 of file KalmanFilterAlg.cxx.

References fFitMomMS, fFitMomRange, fGapDist, fGoodPErr, fGTrace, fGTraceWH, fGTraceWW, fGTraceXMax, fGTraceXMin, fGTraceZMax, fGTraceZMin, fMaxEndChisq, fMaxIncChisq, fMaxLDist, fMaxNoiseHits, fMaxPErr, fMaxPredDist, fMaxPropDist, fMaxSamePlane, fMaxSeedIncChisq, fMaxSeedPredDist, fMaxSmoothIncChisq, fMaxSortDist, fMinLHits, fMinSampleDist, fMinSortDist, fTrace, and fhicl::ParameterSet::get().

224  : fTrace(false)
225  , fMaxPErr(0.)
226  , fGoodPErr(0.)
227  , fMaxIncChisq(0.)
228  , fMaxSeedIncChisq(0.)
229  , fMaxSmoothIncChisq(0.)
230  , fMaxEndChisq(0.)
231  , fMinLHits(0)
232  , fMaxLDist(0.)
233  , fMaxPredDist(0.)
234  , fMaxSeedPredDist(0.)
235  , fMaxPropDist(0.)
236  , fMinSortDist(0.)
237  , fMaxSortDist(0.)
238  , fMaxSamePlane(0)
239  , fGapDist(0.)
240  , fMaxNoiseHits(0)
241  , fMinSampleDist(0.)
242  , fFitMomRange(false)
243  , fFitMomMS(false)
244  , fGTrace(false)
245  , fGTraceWW(0)
246  , fGTraceWH(0)
247  , fPlane(-1)
248  , fInfoPad(0)
249  , fMessages(0)
250 {
251  mf::LogInfo("KalmanFilterAlg") << "KalmanFilterAlg instantiated.";
252 
253  fTrace = pset.get<bool>("Trace");
254  fMaxPErr = pset.get<double>("MaxPErr");
255  fGoodPErr = pset.get<double>("GoodPErr");
256  fMaxIncChisq = pset.get<double>("MaxIncChisq");
257  fMaxSeedIncChisq = pset.get<double>("MaxSeedIncChisq");
258  fMaxSmoothIncChisq = pset.get<double>("MaxSmoothIncChisq");
259  fMaxEndChisq = pset.get<double>("MaxEndChisq");
260  fMinLHits = pset.get<int>("MinLHits");
261  fMaxLDist = pset.get<double>("MaxLDist");
262  fMaxPredDist = pset.get<double>("MaxPredDist");
263  fMaxSeedPredDist = pset.get<double>("MaxSeedPredDist");
264  fMaxPropDist = pset.get<double>("MaxPropDist");
265  fMinSortDist = pset.get<double>("MinSortDist");
266  fMaxSortDist = pset.get<double>("MaxSortDist");
267  fMaxSamePlane = pset.get<int>("MaxSamePlane");
268  fGapDist = pset.get<double>("GapDist");
269  fMaxNoiseHits = pset.get<double>("MaxNoiseHits");
270  fMinSampleDist = pset.get<double>("MinSampleDist");
271  fFitMomRange = pset.get<bool>("FitMomRange");
272  fFitMomMS = pset.get<bool>("FitMomMS");
273  fGTrace = pset.get<bool>("GTrace");
274  if (fGTrace) {
275  fGTraceWW = pset.get<double>("GTraceWW");
276  fGTraceWH = pset.get<double>("GTraceWH");
277  fGTraceXMin = pset.get<double>("GTraceXMin");
278  fGTraceXMax = pset.get<double>("GTraceXMax");
279  fGTraceZMin = pset.get<std::vector<double>>("GTraceZMin");
280  fGTraceZMax = pset.get<std::vector<double>>("GTraceZMax");
281  }
282 }
bool fGTrace
Graphical trace flag.
bool fFitMomRange
Fit momentum using range.
double fMaxPropDist
Maximum propagation distance to candidate surface.
double fGapDist
Minimum gap distance.
bool fFitMomMS
Fit momentum using multiple scattering.
double fMaxEndChisq
Maximum incremental chisquare for endpoint hit.
bool fTrace
Trace flag.
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
std::vector< double > fGTraceZMin
Graphical trace minimum z for each view.
double fMaxLDist
Maximum distance for linearized propagation.
double fMaxSortDist
Sort high distance threshold.
double fMaxSeedPredDist
Maximum prediciton distance to accept a hit in seed phase.
double fMaxPredDist
Maximum prediciton distance to accept a hit.
double fGTraceWH
Window height.
TVirtualPad * fInfoPad
Information pad.
T get(std::string const &key) const
Definition: ParameterSet.h:314
double fGTraceXMin
Graphical trace minimum x.
double fGoodPErr
Pointing error threshold for switching to free propagation.
TPaveText * fMessages
Message box.
double fMaxSeedIncChisq
Maximum incremental chisquare to accept a hit in seed phase.
double fMaxSmoothIncChisq
Maximum incremental chisquare to accept a hit in smooth phase.
double fMaxIncChisq
Maximum incremental chisquare to accept a hit.
double fMinSampleDist
Minimum sample distance (for momentum measurement).
double fMinSortDist
Sort low distance threshold.
double fGTraceXMax
Graphical trace maximum x.
double fGTraceWW
Window width.
std::vector< double > fGTraceZMax
Graphical trace maximum z for each view.
int fMaxSamePlane
Maximum consecutive hits in same plane.
int fMaxNoiseHits
Maximum number of hits in noise cluster.
double fMaxPErr
Maximum pointing error for free propagation.
int fMinLHits
Minimum number of hits to turn off linearized propagation.
int fPlane
Preferred view plane.

Member Function Documentation

bool trkf::KalmanFilterAlg::buildTrack ( const KTrack trk,
KGTrack trg,
const Propagator prop,
const Propagator::PropDirection  dir,
KHitContainer hits,
bool  linear 
) const

Make a new track.

Add hits to track.

Arguments:

trk - Starting track. trg - Result global track. prop - Propagator. dir - Direction. hits - Candidate hits.

Returns: True if success.

This method makes a unidirectional Kalman fit in the specified direction, visiting each surface of the passed KHitContainer and updating the track. In case of multiple measurements on the same surface, keep (at most) the one with the smallest incremental chisquare. Any measurements that fail the incremental chisquare cut are rejected. Resolved hits (accepted or rejected) are moved to the unused list in KHitContainer. The container is sorted at the start of the method, and may be resorted during the progress of the fit.

Definition at line 306 of file KalmanFilterAlg.cxx.

References util::abs(), trkf::KGTrack::addTrack(), trkf::KFitTrack::BACKWARD, trkf::Propagator::BACKWARD, trkf::KFitTrack::BACKWARD_PREDICTED, cleanTrack(), larg4::dist(), trkf::KGTrack::endTrack(), fCanvases, fGoodPErr, fGTrace, fGTraceWH, fGTraceWW, fGTraceXMax, fGTraceXMin, fGTraceZMax, fGTraceZMin, fInfoPad, fMarkerMap, fMaxLDist, fMaxPErr, fMaxPropDist, fMaxSamePlane, fMaxSeedIncChisq, fMaxSeedPredDist, fMaxSortDist, fMessages, fMinLHits, trkf::KFitTrack::FORWARD, trkf::Propagator::FORWARD, trkf::KFitTrack::FORWARD_PREDICTED, fPads, fPlane, fTrace, trkf::KFitTrack::getChisq(), trkf::KHitBase::getChisq(), trkf::KHitGroup::getHits(), trkf::KHitBase::getID(), trkf::KHitBase::getMeasPlane(), trkf::KHitGroup::getPath(), trkf::KFitTrack::getPath(), trkf::KHitGroup::getPlane(), trkf::KHitBase::getPredDistance(), trkf::KHitContainer::getSorted(), trkf::KHitGroup::getSurface(), trkf::KTrack::getSurface(), trkf::KHitContainer::getUnsorted(), trkf::KHitContainer::getUnused(), trkf::KGTrack::isValid(), trkf::KTrack::isValid(), leg, trkf::Propagator::noise_prop(), trkf::KGTrack::numHits(), trkf::KFitTrack::OPTIMAL, trkf::KETrack::PointingError(), trkf::KHitBase::predict(), trkf::KFitTrack::setChisq(), trkf::KFitTrack::setPath(), trkf::KFitTrack::setStat(), trkf::KHitContainer::sort(), trkf::KGTrack::startTrack(), trkf::Propagator::UNKNOWN, and lar::dump::vector().

Referenced by trkf::Track3DKalmanHitAlg::makeKalmanTracks(), and trkf::TrackKalmanCheater::produce().

312 {
313  // Direction must be forward or backward (unknown is not allowed).
314 
316  throw cet::exception("KalmanFilterAlg") << "No direction for Kalman fit.\n";
317 
318  // Set up canvas for graphical trace.
319 
320  if (fGTrace) {
321 
322  // Make a new canvas with a unique name.
323 
324  static int cnum = 0;
325  ++cnum;
326  std::ostringstream ostr;
327  ostr << "khit" << cnum;
328  fCanvases.emplace_back(
329  new TCanvas(ostr.str().c_str(), ostr.str().c_str(), fGTraceWW, fGTraceWH));
330  fPads.clear();
331  fMarkerMap.clear();
332  int nview = 3;
333  fCanvases.back()->Divide(2, nview / 2 + 1);
334 
335  // Make subpad for each view.
336 
337  for (int iview = 0; iview < nview; ++iview) {
338 
339  std::ostringstream ostr;
340  ostr << "Plane " << iview;
341 
342  fCanvases.back()->cd(iview + 1);
343  double zmin = 0.06;
344  double zmax = 0.94;
345  double xmin = 0.04;
346  double xmax = 0.95;
347  TPad* p = new TPad(ostr.str().c_str(), ostr.str().c_str(), zmin, xmin, zmax, xmax);
348  p->SetBit(kCanDelete); // Give away ownership.
349  p->Range(fGTraceZMin[iview], fGTraceXMin, fGTraceZMax[iview], fGTraceXMax);
350  p->SetFillStyle(4000); // Transparent.
351  p->Draw();
352  fPads.push_back(p);
353 
354  // Draw label.
355 
356  TText* t = new TText(zmax - 0.02, xmax - 0.03, ostr.str().c_str());
357  t->SetBit(kCanDelete); // Give away ownership.
358  t->SetTextAlign(33);
359  t->SetTextFont(42);
360  t->SetTextSize(0.04);
361  t->Draw();
362 
363  // Draw axes.
364 
365  TGaxis* pz1 =
366  new TGaxis(zmin, xmin, zmax, xmin, fGTraceZMin[iview], fGTraceZMax[iview], 510, "");
367  pz1->SetBit(kCanDelete); // Give away ownership.
368  pz1->Draw();
369 
370  TGaxis* px1 = new TGaxis(zmin, xmin, zmin, xmax, fGTraceXMin, fGTraceXMax, 510, "");
371  px1->SetBit(kCanDelete); // Give away ownership.
372  px1->Draw();
373 
374  TGaxis* pz2 =
375  new TGaxis(zmin, xmax, zmax, xmax, fGTraceZMin[iview], fGTraceZMax[iview], 510, "-");
376  pz2->SetBit(kCanDelete); // Give away ownership.
377  pz2->Draw();
378 
379  TGaxis* px2 = new TGaxis(zmax, xmin, zmax, xmax, fGTraceXMin, fGTraceXMax, 510, "L+");
380  px2->SetBit(kCanDelete); // Give away ownership.
381  px2->Draw();
382  }
383 
384  // Draw legend.
385 
386  fCanvases.back()->cd(nview + 1);
387  fInfoPad = gPad;
388  TLegend* leg = new TLegend(0.6, 0.5, 0.99, 0.99);
389  leg->SetBit(kCanDelete); // Give away ownership.
390 
391  TLegendEntry* entry = 0;
392  TMarker* m = 0;
393 
394  m = new TMarker(0., 0., 20);
395  m->SetBit(kCanDelete);
396  m->SetMarkerSize(1.2);
397  m->SetMarkerColor(kRed);
398  entry = leg->AddEntry(m, "Hits on Track", "P");
399  entry->SetBit(kCanDelete);
400 
401  m = new TMarker(0., 0., 20);
402  m->SetBit(kCanDelete);
403  m->SetMarkerSize(1.2);
404  m->SetMarkerColor(kOrange);
405  entry = leg->AddEntry(m, "Smoothed Hits on Track", "P");
406  entry->SetBit(kCanDelete);
407 
408  m = new TMarker(0., 0., 20);
409  m->SetBit(kCanDelete);
410  m->SetMarkerSize(1.2);
411  m->SetMarkerColor(kBlack);
412  entry = leg->AddEntry(m, "Available Hits", "P");
413  entry->SetBit(kCanDelete);
414 
415  m = new TMarker(0., 0., 20);
416  m->SetBit(kCanDelete);
417  m->SetMarkerSize(1.2);
418  m->SetMarkerColor(kBlue);
419  entry = leg->AddEntry(m, "Rejected Hits", "P");
420  entry->SetBit(kCanDelete);
421 
422  m = new TMarker(0., 0., 20);
423  m->SetBit(kCanDelete);
424  m->SetMarkerSize(1.2);
425  m->SetMarkerColor(kGreen);
426  entry = leg->AddEntry(m, "Removed Hits", "P");
427  entry->SetBit(kCanDelete);
428 
429  m = new TMarker(0., 0., 20);
430  m->SetBit(kCanDelete);
431  m->SetMarkerSize(1.2);
432  m->SetMarkerColor(kMagenta);
433  entry = leg->AddEntry(m, "Seed Hits", "P");
434  entry->SetBit(kCanDelete);
435 
436  m = new TMarker(0., 0., 20);
437  m->SetBit(kCanDelete);
438  m->SetMarkerSize(1.2);
439  m->SetMarkerColor(kCyan);
440  entry = leg->AddEntry(m, "Unreachable Seed Hits", "P");
441  entry->SetBit(kCanDelete);
442 
443  leg->Draw();
444 
445  // Draw message box.
446 
447  fMessages = new TPaveText(0.01, 0.01, 0.55, 0.99, "");
448  fMessages->SetBit(kCanDelete);
449  fMessages->SetTextFont(42);
450  fMessages->SetTextSize(0.04);
451  fMessages->SetTextAlign(12);
452  fMessages->Draw();
453  TText* t = fMessages->AddText("Enter buildTrack");
454  t->SetBit(kCanDelete);
455 
456  fCanvases.back()->Update();
457  }
458 
459  // Sort container using this seed track.
460 
461  hits.sort(trk, true, prop, Propagator::UNKNOWN);
462 
463  // Draw hits and populate hit->marker map.
464 
465  if (fGTrace && fCanvases.size() > 0) {
466 
467  // Loop over sorted KHitGroups.
468  // Paint sorted seed hits magenta.
469 
470  const std::list<KHitGroup>& groups = hits.getSorted();
471  for (auto const& gr : groups) {
472 
473  // Loop over hits in this group.
474 
475  const std::vector<std::shared_ptr<const KHitBase>>& phits = gr.getHits();
476  for (auto const& phit : phits) {
477  const KHitBase& hit = *phit;
478  int pl = hit.getMeasPlane();
479  if (pl >= 0 && pl < int(fPads.size())) {
480  double z = 0.;
481  double x = 0.;
482  hit_position(hit, z, x);
483  TMarker* marker = new TMarker(z, x, 20);
484  fMarkerMap[hit.getID()] = marker;
485  fPads[pl]->cd();
486  marker->SetBit(kCanDelete); // Give away ownership.
487  marker->SetMarkerSize(0.5);
488  marker->SetMarkerColor(kMagenta);
489  marker->Draw();
490  }
491  }
492  }
493 
494  // Loop over unsorted KHitGroups.
495  // Paint unsorted seed hits cyan.
496  // There should be few, if any, unsorted seed hits.
497 
498  const std::list<KHitGroup>& ugroups = hits.getUnsorted();
499  for (auto const& gr : ugroups) {
500 
501  // Loop over hits in this group.
502 
503  const std::vector<std::shared_ptr<const KHitBase>>& phits = gr.getHits();
504  for (auto const& phit : phits) {
505  const KHitBase& hit = *phit;
506  int pl = hit.getMeasPlane();
507  if (pl >= 0 && pl < int(fPads.size())) {
508  double z = 0.;
509  double x = 0.;
510  hit_position(hit, z, x);
511  TMarker* marker = new TMarker(z, x, 20);
512  fMarkerMap[hit.getID()] = marker;
513  fPads[pl]->cd();
514  marker->SetBit(kCanDelete); // Give away ownership.
515  marker->SetMarkerSize(0.5);
516  marker->SetMarkerColor(kCyan);
517  marker->Draw();
518  }
519  }
520  }
521  fCanvases.back()->Update();
522  }
523 
524  // Loop over measurements (KHitGroup) from sorted list.
525 
526  double tchisq = 0.; // Cumulative chisquare.
527  double path = 0.; // Cumulative path distance.
528  int step = 0; // Step count.
529  int nsame = 0; // Number of consecutive measurements in same plane.
530  int last_plane = -1; // Last plane.
531  bool has_pref_plane = false; // Set to true when first preferred plane hit is added to track.
532 
533  // Make a copy of the starting track, in the form of a KFitTrack,
534  // which we will update as we go.
535 
536  TrackError err;
537  trk.getSurface()->getStartingError(err);
538  KETrack tre(trk, err);
539  KFitTrack trf(tre, path, tchisq);
540 
541  // Also use the starting track as the reference track for linearized
542  // propagation, until the track is established with reasonably small
543  // errors.
544 
545  KTrack ref(trk);
546  KTrack* pref = &ref;
547 
548  mf::LogInfo log("KalmanFilterAlg");
549 
550  // Loop over measurement groups (KHitGroups).
551 
552  while (hits.getSorted().size() > 0) {
553  ++step;
554  if (fTrace) {
555  log << "Build Step " << step << "\n";
556  log << "KGTrack has " << trg.numHits() << " hits.\n";
557  log << trf;
558  }
559 
560  // Get an iterator for the next KHitGroup.
561 
563  if (dir == Propagator::FORWARD)
564  it = hits.getSorted().begin();
565  else {
566  it = hits.getSorted().end();
567  --it;
568  }
569  const KHitGroup& gr = *it;
570 
571  if (fTrace) {
572  double path_est = gr.getPath();
573  log << "Next surface: " << *(gr.getSurface()) << "\n";
574  log << " Estimated path length of hit group = " << path_est << "\n";
575  }
576 
577  // Get the next prediction surface. If this KHitGroup is on the
578  // preferred plane, use that as the prediction surface.
579  // Otherwise, use the current track surface as the prediction
580  // surface.
581 
582  std::shared_ptr<const Surface> psurf = trf.getSurface();
583  if (gr.getPlane() < 0) throw cet::exception("KalmanFilterAlg") << "negative plane?\n";
584  if (fPlane < 0 || gr.getPlane() < 0 || fPlane == gr.getPlane()) psurf = gr.getSurface();
585 
586  // Propagate track to the prediction surface.
587 
588  std::optional<double> dist = prop.noise_prop(trf, psurf, Propagator::UNKNOWN, true, pref);
589  if (dist && std::abs(*dist) > fMaxPropDist) dist = std::nullopt;
590  double ds = 0.;
591 
592  if (dist) {
593 
594  // Propagation succeeded.
595  // Update cumulative path distance and track status.
596 
597  ds = *dist;
598  path += ds;
599  trf.setPath(path);
600  if (dir == Propagator::FORWARD)
601  trf.setStat(KFitTrack::FORWARD_PREDICTED);
602  else {
603  trf.setStat(KFitTrack::BACKWARD_PREDICTED);
604  }
605  if (fTrace) {
606  log << "After propagation\n";
607  log << " Incremental propagation distance = " << ds << "\n";
608  log << " Path length of prediction surface = " << path << "\n";
609  log << "KGTrack has " << trg.numHits() << " hits.\n";
610  log << trf;
611  }
612 
613  // Loop over measurements in this group.
614 
615  const std::vector<std::shared_ptr<const KHitBase>>& hits = gr.getHits();
616  double best_chisq = 0.;
617  std::shared_ptr<const KHitBase> best_hit;
618  for (std::vector<std::shared_ptr<const KHitBase>>::const_iterator ihit = hits.begin();
619  ihit != hits.end();
620  ++ihit) {
621  const KHitBase& hit = **ihit;
622 
623  // Turn this hit blue.
624 
625  if (fGTrace && fCanvases.size() > 0) {
626  auto marker_it = fMarkerMap.find(hit.getID());
627  if (marker_it != fMarkerMap.end()) {
628  TMarker* marker = marker_it->second;
629  marker->SetMarkerColor(kBlue);
630  }
631  //fCanvases.back()->Update();
632  }
633 
634  // Update predction using current track hypothesis and get
635  // incremental chisquare.
636 
637  if (hit.predict(trf, prop)) {
638  double chisq = hit.getChisq();
639  double preddist = hit.getPredDistance();
640  if (fTrace) {
641  log << "Trying Hit.\n"
642  << hit << "\nchisq = " << chisq << "\n"
643  << "prediction distance = " << preddist << "\n";
644  }
645  if ((!has_pref_plane || abs(preddist) < fMaxSeedPredDist) &&
646  (best_hit.get() == 0 || chisq < best_chisq)) {
647  best_chisq = chisq;
648  if (chisq < fMaxSeedIncChisq) best_hit = *ihit;
649  }
650  }
651  }
652  if (fTrace) {
653  log << "Best hit incremental chisquare = " << best_chisq << "\n";
654  if (best_hit.get() != 0) {
655  log << "Hit after prediction\n";
656  log << *best_hit;
657  }
658  else
659  log << "No hit passed chisquare cut.\n";
660  }
661  if (fGTrace && fCanvases.size() > 0) fCanvases.back()->Update();
662 
663  // If we found a best measurement, and if the incremental
664  // chisquare passes the cut, add it to the track and update
665  // fit information.
666 
667  bool update_ok = false;
668  if (best_hit.get() != 0) {
669  KFitTrack trf0(trf);
670  best_hit->update(trf);
671  update_ok = trf.isValid();
672  if (!update_ok) trf = trf0;
673  }
674  if (update_ok) {
675  ds += best_hit->getPredDistance();
676  tchisq += best_chisq;
677  trf.setChisq(tchisq);
678  if (dir == Propagator::FORWARD)
679  trf.setStat(KFitTrack::FORWARD);
680  else {
681  trf.setStat(KFitTrack::BACKWARD);
682  }
683 
684  // If the pointing error got too large, and there is no
685  // reference track, quit.
686 
687  if (pref == 0 && trf.PointingError() > fMaxPErr) {
688  if (fTrace) log << "Quitting because pointing error got too large.\n";
689  break;
690  }
691 
692  // Test number of consecutive measurements in the same plane.
693 
694  if (gr.getPlane() >= 0) {
695  if (gr.getPlane() == last_plane)
696  ++nsame;
697  else {
698  nsame = 1;
699  last_plane = gr.getPlane();
700  }
701  }
702  else {
703  nsame = 0;
704  last_plane = -1;
705  }
706  if (nsame <= fMaxSamePlane) {
707 
708  // Turn best hit red.
709 
710  if (fGTrace && fCanvases.size() > 0) {
711  int pl = best_hit->getMeasPlane();
712  if (pl >= 0 && pl < int(fPads.size())) {
713  auto marker_it = fMarkerMap.find(best_hit->getID());
714  if (marker_it != fMarkerMap.end()) {
715  TMarker* marker = marker_it->second;
716  marker->SetMarkerColor(kRed);
717 
718  // Redraw marker so that it will be on top.
719 
720  fPads[pl]->cd();
721  marker->Draw();
722  }
723  }
724  fCanvases.back()->Update();
725  }
726 
727  // Make a KHitTrack and add it to the KGTrack.
728 
729  KHitTrack trh(trf, best_hit);
730  trg.addTrack(trh);
731  if (fPlane == gr.getPlane()) has_pref_plane = true;
732 
733  // Decide if we want to kill the reference track.
734 
735  if (!linear && pref != 0 && int(trg.numHits()) >= fMinLHits &&
736  (trf.PointingError() < fGoodPErr || path > fMaxLDist)) {
737  pref = 0;
738  if (fTrace) log << "Killing reference track.\n";
739  }
740 
741  if (fTrace) {
742  log << "After update\n";
743  log << "KGTrack has " << trg.numHits() << " hits.\n";
744  log << trf;
745  }
746  }
747  }
748  }
749 
750  // The current KHitGroup is now resolved.
751  // Move it to unused list.
752 
753  hits.getUnused().splice(hits.getUnused().end(), hits.getSorted(), it);
754 
755  // If the propagation distance was the wrong direction, resort the measurements.
756 
757  if (pref == 0 && dist &&
758  ((dir == Propagator::FORWARD && (ds < fMinSortDist || ds > fMaxSortDist)) ||
759  (dir == Propagator::BACKWARD && (-ds < fMinSortDist || -ds > fMaxSortDist)))) {
760  if (fTrace) log << "Resorting measurements.\n";
761  hits.sort(trf, true, prop, dir);
762  }
763  }
764 
765  // Clean track.
766 
767  cleanTrack(trg);
768 
769  // Set the fit status of the last added KHitTrack to optimal and get
770  // the final chisquare and path length.
771 
772  double fchisq = 0.;
773  path = 0.;
774  if (trg.isValid()) {
775  path = trg.endTrack().getPath() - trg.startTrack().getPath();
776  if (dir == Propagator::FORWARD) {
777  trg.endTrack().setStat(KFitTrack::OPTIMAL);
778  fchisq = trg.endTrack().getChisq();
779  }
780  else {
781  trg.startTrack().setStat(KFitTrack::OPTIMAL);
782  fchisq = trg.startTrack().getChisq();
783  }
784  }
785 
786  // Summary.
787 
788  log << "KalmanFilterAlg build track summary.\n"
789  << "Build direction = " << (dir == Propagator::FORWARD ? "FORWARD" : "BACKWARD") << "\n"
790  << "Track has " << trg.numHits() << " hits.\n"
791  << "Track length = " << path << "\n"
792  << "Track chisquare = " << fchisq << "\n";
793 
794  // Done. Return success if we added at least one measurement.
795 
796  bool ok = trg.numHits() > 0;
797 
798  if (fGTrace && fCanvases.size() > 0) {
799  TText* t = 0;
800  if (ok)
801  t = fMessages->AddText("Exit buildTrack, status success");
802  else
803  t = fMessages->AddText("Exit buildTrack, status fail");
804  t->SetBit(kCanDelete);
805  fInfoPad->Modified();
806  fCanvases.back()->Update();
807  }
808 
809  return ok;
810 }
Float_t x
Definition: compare.C:6
bool fGTrace
Graphical trace flag.
intermediate_table::iterator iterator
double fMaxPropDist
Maximum propagation distance to candidate surface.
std::map< int, TMarker * > fMarkerMap
Markers in current canvas.
KSymMatrix< 5 >::type TrackError
Track error matrix, dimension 5x5.
bool fTrace
Trace flag.
Double_t z
Definition: plot.C:276
constexpr auto abs(T v)
Returns the absolute value of the argument.
void cleanTrack(KGTrack &trg) const
Clean track by removing noise hits near endpoints.
std::vector< double > fGTraceZMin
Graphical trace minimum z for each view.
double fMaxLDist
Maximum distance for linearized propagation.
double fMaxSortDist
Sort high distance threshold.
double fMaxSeedPredDist
Maximum prediciton distance to accept a hit in seed phase.
auto vector(Vector const &v)
Returns a manipulator which will print the specified array.
Definition: DumpUtils.h:289
void hits()
Definition: readHits.C:15
double fGTraceWH
Window height.
TVirtualPad * fInfoPad
Information pad.
std::vector< std::unique_ptr< TCanvas > > fCanvases
Graphical trace canvases.
double fGTraceXMin
Graphical trace minimum x.
double fGoodPErr
Pointing error threshold for switching to free propagation.
TPaveText * fMessages
Message box.
double fMaxSeedIncChisq
Maximum incremental chisquare to accept a hit in seed phase.
Detector simulation of raw signals on wires.
TDirectory * dir
Definition: macro.C:5
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
double fGTraceXMax
Graphical trace maximum x.
std::vector< TVirtualPad * > fPads
View pads in current canvas.
TLegend * leg
Definition: compare.C:67
double fGTraceWW
Window width.
std::vector< double > fGTraceZMax
Graphical trace maximum z for each view.
int fMaxSamePlane
Maximum consecutive hits in same plane.
double fMaxPErr
Maximum pointing error for free propagation.
int fMinLHits
Minimum number of hits to turn off linearized propagation.
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
int fPlane
Preferred view plane.
void trkf::KalmanFilterAlg::cleanTrack ( KGTrack trg) const

Clean track by removing noise hits near endpoints.

Clean track by removing noise hits near endpoints.

Arguments:

trg - Track.

Definition at line 2097 of file KalmanFilterAlg.cxx.

References fCanvases, fGapDist, fGTrace, fMarkerMap, fMaxEndChisq, fMaxNoiseHits, trkf::KHitBase::getChisq(), trkf::KHitTrack::getHit(), trkf::KHitBase::getID(), trkf::KHitBase::getMeasPlane(), and trkf::KGTrack::getTrackMap().

Referenced by buildTrack(), and extendTrack().

2098 {
2099  // Get hold of a modifiable track map from trg.
2100 
2101  std::multimap<double, KHitTrack>& trackmap = trg.getTrackMap();
2102 
2103  // Do an indefinite loop until we no longer find any dirt.
2104 
2105  bool done = false;
2106  while (!done) {
2107 
2108  // If track map has fewer than fMaxNoiseHits tracks, then this is a
2109  // noise track. Clear the map, making the whole track invalid.
2110 
2111  int ntrack = trackmap.size();
2112  if (ntrack <= fMaxNoiseHits) {
2113  for (auto const& trackmap_entry : trackmap) {
2114  const KHitTrack& trh = trackmap_entry.second;
2115  const KHitBase& hit = *(trh.getHit());
2116  auto marker_it = fMarkerMap.find(hit.getID());
2117  if (marker_it != fMarkerMap.end()) {
2118  TMarker* marker = marker_it->second;
2119  marker->SetMarkerColor(kGreen);
2120  }
2121  }
2122  trackmap.clear();
2123  done = true;
2124  break;
2125  }
2126 
2127  // Make sure the first two and last two tracks belong to different
2128  // views. If not, remove the first or last track.
2129 
2130  if (ntrack >= 2) {
2131 
2132  // Check start.
2133 
2134  std::multimap<double, KHitTrack>::iterator it = trackmap.begin();
2135  const KHitTrack& trh1a = (*it).second;
2136  const KHitBase& hit1a = *(trh1a.getHit());
2137  int plane1 = hit1a.getMeasPlane();
2138  double chisq1 = hit1a.getChisq();
2139  ++it;
2140  const KHitTrack& trh2a = (*it).second;
2141  const KHitBase& hit2a = *(trh2a.getHit());
2142  int plane2 = hit2a.getMeasPlane();
2143  double chisq2 = hit2a.getChisq();
2144  if ((plane1 >= 0 && plane2 >= 0 && plane1 == plane2) || chisq1 > fMaxEndChisq ||
2145  chisq2 > fMaxEndChisq) {
2146  auto marker_it = fMarkerMap.find(hit1a.getID());
2147  if (marker_it != fMarkerMap.end()) {
2148  TMarker* marker = marker_it->second;
2149  marker->SetMarkerColor(kGreen);
2150  }
2151  trackmap.erase(trackmap.begin(), it);
2152  done = false;
2153  continue;
2154  }
2155 
2156  // Check end.
2157 
2158  it = trackmap.end();
2159  --it;
2160  const KHitTrack& trh1b = (*it).second;
2161  const KHitBase& hit1b = *(trh1b.getHit());
2162  plane1 = hit1b.getMeasPlane();
2163  chisq1 = hit1b.getChisq();
2164  --it;
2165  const KHitTrack& trh2b = (*it).second;
2166  const KHitBase& hit2b = *(trh2b.getHit());
2167  plane2 = hit2b.getMeasPlane();
2168  chisq2 = hit2b.getChisq();
2169  if ((plane1 >= 0 && plane2 >= 0 && plane1 == plane2) || chisq1 > fMaxEndChisq ||
2170  chisq2 > fMaxEndChisq) {
2171  ++it;
2172  auto marker_it = fMarkerMap.find(hit1b.getID());
2173  if (marker_it != fMarkerMap.end()) {
2174  TMarker* marker = marker_it->second;
2175  marker->SetMarkerColor(kGreen);
2176  }
2177  trackmap.erase(it, trackmap.end());
2178  done = false;
2179  continue;
2180  }
2181  }
2182 
2183  // Loop over successive pairs of elements of track map. Look for
2184  // adjacent elements with distance separation greater than
2185  // fGapDist.
2186 
2187  std::multimap<double, KHitTrack>::iterator it = trackmap.begin();
2188  std::multimap<double, KHitTrack>::iterator jt = trackmap.end();
2189  int nb = 0; // Number of elements from begin to jt.
2190  int ne = ntrack; // Number of elements it to end.
2191  bool found_noise = false;
2192  for (; it != trackmap.end(); ++it, ++nb, --ne) {
2193  if (jt == trackmap.end())
2194  jt = trackmap.begin();
2195  else {
2196  if (nb < 1)
2197  throw cet::exception("KalmanFilterAlg")
2198  << "KalmanFilterAlg::cleanTrack(): nb not positive\n";
2199  if (ne < 1)
2200  throw cet::exception("KalmanFilterAlg")
2201  << "KalmanFilterAlg::cleanTrack(): ne not positive\n";
2202 
2203  double disti = (*it).first;
2204  double distj = (*jt).first;
2205  double sep = disti - distj;
2206  if (sep < 0.)
2207  throw cet::exception("KalmanFilterAlg")
2208  << "KalmanFilterAlg::fitMomentumRange(): negative separation\n";
2209 
2210  if (sep > fGapDist) {
2211 
2212  // Found a gap. See if we want to trim track.
2213 
2214  if (nb <= fMaxNoiseHits) {
2215 
2216  // Trim front.
2217 
2218  found_noise = true;
2219  for (auto jt = trackmap.begin(); jt != it; ++jt) {
2220  const KHitTrack& trh = jt->second;
2221  const KHitBase& hit = *(trh.getHit());
2222  auto marker_it = fMarkerMap.find(hit.getID());
2223  if (marker_it != fMarkerMap.end()) {
2224  TMarker* marker = marker_it->second;
2225  marker->SetMarkerColor(kGreen);
2226  }
2227  }
2228  trackmap.erase(trackmap.begin(), it);
2229  break;
2230  }
2231  if (ne <= fMaxNoiseHits) {
2232 
2233  // Trim back.
2234 
2235  found_noise = true;
2236  for (auto jt = it; jt != trackmap.end(); ++jt) {
2237  const KHitTrack& trh = jt->second;
2238  const KHitBase& hit = *(trh.getHit());
2239  auto marker_it = fMarkerMap.find(hit.getID());
2240  if (marker_it != fMarkerMap.end()) {
2241  TMarker* marker = marker_it->second;
2242  marker->SetMarkerColor(kGreen);
2243  }
2244  }
2245  trackmap.erase(it, trackmap.end());
2246  break;
2247  }
2248  }
2249  ++jt;
2250  }
2251  }
2252 
2253  // Decide if we are done.
2254 
2255  done = !found_noise;
2256  }
2257  if (fGTrace && fCanvases.size() > 0) fCanvases.back()->Update();
2258 }
bool fGTrace
Graphical trace flag.
intermediate_table::iterator iterator
double fGapDist
Minimum gap distance.
std::map< int, TMarker * > fMarkerMap
Markers in current canvas.
double fMaxEndChisq
Maximum incremental chisquare for endpoint hit.
std::vector< std::unique_ptr< TCanvas > > fCanvases
Graphical trace canvases.
Detector simulation of raw signals on wires.
int fMaxNoiseHits
Maximum number of hits in noise cluster.
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
bool trkf::KalmanFilterAlg::extendTrack ( KGTrack trg,
const Propagator prop,
KHitContainer hits 
) const

Add hits to existing track.

Add hits to existing track.

Arguments:

trg - Track to extend. prop - Propagator. hits - Hit collection to choose hits from.

This method extends a KGTrack by adding hits. The KGTrack must have previously been produced by a unidirectional Kalman fit (it should be optimal at one end). This method finds the optimal end and extends the track in that direction. If any hits are added, the originally optimal end has its status reset to forward or backward, and the new endpoint is optimal. In any case, the final result is unidirectionally fit KGTrack.

Definition at line 1152 of file KalmanFilterAlg.cxx.

References util::abs(), trkf::KGTrack::addTrack(), trkf::KFitTrack::BACKWARD, trkf::Propagator::BACKWARD, trkf::KFitTrack::BACKWARD_PREDICTED, cleanTrack(), dir, larg4::dist(), trkf::KGTrack::endTrack(), fCanvases, fGTrace, fInfoPad, fMarkerMap, fMaxIncChisq, fMaxPErr, fMaxPredDist, fMaxPropDist, fMaxSamePlane, fMaxSortDist, fMessages, trkf::KFitTrack::FORWARD, trkf::Propagator::FORWARD, trkf::KFitTrack::FORWARD_PREDICTED, fPads, fPlane, fTrace, trkf::KFitTrack::getChisq(), trkf::KHitBase::getChisq(), trkf::KETrack::getError(), trkf::KHitGroup::getHits(), trkf::KHitBase::getID(), trkf::KHitBase::getMeasPlane(), trkf::KHitGroup::getPath(), trkf::KFitTrack::getPath(), trkf::KHitGroup::getPlane(), trkf::KHitBase::getPredDistance(), trkf::KHitContainer::getSorted(), trkf::KFitTrack::getStat(), trkf::KHitGroup::getSurface(), trkf::KTrack::getSurface(), trkf::KHitContainer::getUnsorted(), trkf::KHitContainer::getUnused(), trkf::KTrack::getVector(), trkf::KGTrack::isValid(), trkf::KTrack::isValid(), trkf::Propagator::noise_prop(), trkf::KGTrack::numHits(), trkf::KFitTrack::OPTIMAL, trkf::KETrack::PointingError(), trkf::KHitBase::predict(), trkf::KFitTrack::setChisq(), trkf::KFitTrack::setPath(), trkf::KFitTrack::setStat(), trkf::KHitContainer::sort(), trkf::KGTrack::startTrack(), trkf::Propagator::UNKNOWN, and lar::dump::vector().

Referenced by trkf::Track3DKalmanHitAlg::extendandsmoothLoop().

1155 {
1156  // Default result failure.
1157 
1158  bool result = false;
1159 
1160  if (fGTrace && fCanvases.size() > 0) {
1161  TText* t = fMessages->AddText("Enter extendTrack");
1162  t->SetBit(kCanDelete);
1163  fInfoPad->Modified();
1164  fCanvases.back()->Update();
1165  }
1166 
1167  // Remember the original number of measurement.
1168 
1169  unsigned int nhits0 = trg.numHits();
1170 
1171  // It is an error if the KGTrack is not valid.
1172 
1173  if (trg.isValid()) {
1174  mf::LogInfo log("KalmanFilterAlg");
1175 
1176  // Examine the track endpoints and figure out which end of the
1177  // track to extend. The track is always extended from the optimal
1178  // end. It is an error if neither end point is optimal, or both
1179  // endoints are optimal. Reset the status of the optimal, and
1180  // make a copy of the starting track fit. Also get starting path
1181  // and chisquare.
1182 
1183  KHitTrack& trh0 = trg.startTrack();
1184  KHitTrack& trh1 = trg.endTrack();
1185  KFitTrack::FitStatus stat0 = trh0.getStat();
1186  KFitTrack::FitStatus stat1 = trh1.getStat();
1187  bool dofit = false;
1189  KFitTrack trf;
1190  double path = 0.;
1191  double tchisq = 0.;
1192 
1193  if (stat0 == KFitTrack::OPTIMAL) {
1194  if (stat1 == KFitTrack::OPTIMAL) {
1195 
1196  // Both ends optimal (do nothing, return failure).
1197 
1198  dofit = false;
1199  result = false;
1200  return result;
1201  }
1202  else {
1203 
1204  // Start optimal. Extend backward.
1205 
1206  dofit = true;
1207  dir = Propagator::BACKWARD;
1208  trh0.setStat(KFitTrack::BACKWARD);
1209  trf = trh0;
1210  path = trh0.getPath();
1211  tchisq = trh0.getChisq();
1212  }
1213  }
1214  else {
1215  if (stat1 == KFitTrack::OPTIMAL) {
1216 
1217  // End optimal. Extend forward.
1218 
1219  dofit = true;
1220  dir = Propagator::FORWARD;
1221  trh1.setStat(KFitTrack::FORWARD);
1222  trf = trh1;
1223  path = trh1.getPath();
1224  tchisq = trh1.getChisq();
1225 
1226  // Make sure forward extend track momentum is over some
1227  // minimum value.
1228 
1229  if (trf.getVector()(4) > 5.) {
1230  trf.getVector()(4) = 5.;
1231  trf.getError()(4, 4) = 5.;
1232  }
1233  }
1234  else {
1235 
1236  // Neither end optimal (do nothing, return failure).
1237 
1238  dofit = false;
1239  result = false;
1240  return result;
1241  }
1242  }
1243  if (dofit) {
1244 
1245  // Sort hit container using starting track.
1246 
1247  hits.sort(trf, true, prop, dir);
1248 
1249  // Draw and add hits in hit->marker map that are not already there.
1250 
1251  if (fGTrace && fCanvases.size() > 0) {
1252 
1253  // Loop over sorted KHitGroups.
1254  // Paint sorted hits black.
1255 
1256  const std::list<KHitGroup>& groups = hits.getSorted();
1257  for (auto const& gr : groups) {
1258 
1259  // Loop over hits in this group.
1260 
1261  const std::vector<std::shared_ptr<const KHitBase>>& phits = gr.getHits();
1262  for (auto const& phit : phits) {
1263  const KHitBase& hit = *phit;
1264  int pl = hit.getMeasPlane();
1265  if (pl >= 0 && pl < int(fPads.size())) {
1266 
1267  // Is this hit already in map?
1268 
1269  TMarker* marker = 0;
1270  auto marker_it = fMarkerMap.find(hit.getID());
1271  if (marker_it == fMarkerMap.end()) {
1272  double z = 0.;
1273  double x = 0.;
1274  hit_position(hit, z, x);
1275  marker = new TMarker(z, x, 20);
1276  fMarkerMap[hit.getID()] = marker;
1277  fPads[pl]->cd();
1278  marker->SetBit(kCanDelete); // Give away ownership.
1279  marker->SetMarkerSize(0.5);
1280  marker->Draw();
1281  }
1282  else
1283  marker = marker_it->second;
1284  marker->SetMarkerColor(kBlack);
1285  }
1286  }
1287  }
1288 
1289  // Loop over unsorted KHitGroups.
1290  // Paint unsorted hits blue.
1291 
1292  const std::list<KHitGroup>& ugroups = hits.getUnsorted();
1293  for (auto const& gr : ugroups) {
1294 
1295  // Loop over hits in this group.
1296 
1297  const std::vector<std::shared_ptr<const KHitBase>>& phits = gr.getHits();
1298  for (auto const& phit : phits) {
1299  const KHitBase& hit = *phit;
1300  int pl = hit.getMeasPlane();
1301  if (pl >= 0 && pl < int(fPads.size())) {
1302 
1303  // Is this hit already in map?
1304 
1305  TMarker* marker = 0;
1306  auto marker_it = fMarkerMap.find(hit.getID());
1307  if (marker_it == fMarkerMap.end()) {
1308  double z = 0.;
1309  double x = 0.;
1310  hit_position(hit, z, x);
1311  marker = new TMarker(z, x, 20);
1312  fMarkerMap[hit.getID()] = marker;
1313  fPads[pl]->cd();
1314  marker->SetBit(kCanDelete); // Give away ownership.
1315  marker->SetMarkerSize(0.5);
1316  marker->Draw();
1317  }
1318  else
1319  marker = marker_it->second;
1320  marker->SetMarkerColor(kBlue);
1321  }
1322  }
1323  }
1324  fCanvases.back()->Update();
1325  }
1326 
1327  //std::cout << "extendTrack: marker map has " << fMarkerMap.size() << " entries." << std::endl;
1328 
1329  // Extend loop starts here.
1330 
1331  int step = 0;
1332  int nsame = 0;
1333  int last_plane = -1;
1334  while (hits.getSorted().size() > 0) {
1335  ++step;
1336  if (fTrace) {
1337  log << "Extend Step " << step << "\n";
1338  log << "KGTrack has " << trg.numHits() << " hits.\n";
1339  log << trf;
1340  }
1341 
1342  // Get an iterator for the next KHitGroup.
1343 
1345  switch (dir) {
1346  case Propagator::FORWARD: it = hits.getSorted().begin(); break;
1347  case Propagator::BACKWARD:
1348  it = hits.getSorted().end();
1349  --it;
1350  break;
1351  default:
1352  throw cet::exception("KalmanFilterAlg")
1353  << "KalmanFilterAlg::extendTrack(): invalid direction\n";
1354  } // switch
1355  const KHitGroup& gr = *it;
1356 
1357  if (fTrace) {
1358  double path_est = gr.getPath();
1359  log << "Next surface: " << *(gr.getSurface()) << "\n";
1360  log << " Estimated path length of hit group = " << path_est << "\n";
1361  }
1362 
1363  // Get the next prediction surface. If this KHitGroup is on the
1364  // preferred plane, use that as the prediction surface.
1365  // Otherwise, use the current track surface as the prediction
1366  // surface.
1367 
1368  std::shared_ptr<const Surface> psurf = trf.getSurface();
1369  if (gr.getPlane() < 0)
1370  throw cet::exception("KalmanFilterAlg")
1371  << "KalmanFilterAlg::extendTrack(): negative plane?\n";
1372  if (fPlane < 0 || gr.getPlane() < 0 || fPlane == gr.getPlane()) psurf = gr.getSurface();
1373 
1374  // Propagate track to the prediction surface.
1375 
1376  std::optional<double> dist = prop.noise_prop(trf, psurf, Propagator::UNKNOWN, true);
1377  if (dist && std::abs(*dist) > fMaxPropDist) dist = std::nullopt;
1378  double ds = 0.;
1379 
1380  if (dist) {
1381 
1382  // Propagation succeeded.
1383  // Update cumulative path distance and track status.
1384 
1385  ds = *dist;
1386  path += ds;
1387  trf.setPath(path);
1388  if (dir == Propagator::FORWARD)
1389  trf.setStat(KFitTrack::FORWARD_PREDICTED);
1390  else {
1391  trf.setStat(KFitTrack::BACKWARD_PREDICTED);
1392  }
1393  if (fTrace) {
1394  log << "After propagation\n";
1395  log << " Incremental distance = " << ds << "\n";
1396  log << " Track path length = " << path << "\n";
1397  log << "KGTrack has " << trg.numHits() << " hits.\n";
1398  log << trf;
1399  }
1400 
1401  // Loop over measurements in this group.
1402 
1403  const std::vector<std::shared_ptr<const KHitBase>>& hits = gr.getHits();
1404  double best_chisq = 0.;
1405  std::shared_ptr<const KHitBase> best_hit;
1406  for (std::vector<std::shared_ptr<const KHitBase>>::const_iterator ihit = hits.begin();
1407  ihit != hits.end();
1408  ++ihit) {
1409  const KHitBase& hit = **ihit;
1410 
1411  // Turn this hit blue.
1412 
1413  if (fGTrace && fCanvases.size() > 0) {
1414  auto marker_it = fMarkerMap.find(hit.getID());
1415  if (marker_it != fMarkerMap.end()) {
1416  TMarker* marker = marker_it->second;
1417  marker->SetMarkerColor(kBlue);
1418  }
1419  //fCanvases.back()->Update();
1420  }
1421 
1422  // Update predction using current track hypothesis and get
1423  // incremental chisquare.
1424 
1425  bool ok = hit.predict(trf, prop);
1426  if (ok) {
1427  double chisq = hit.getChisq();
1428  double preddist = hit.getPredDistance();
1429  if (abs(preddist) < fMaxPredDist && (best_hit.get() == 0 || chisq < best_chisq)) {
1430  best_chisq = chisq;
1431  if (chisq < fMaxIncChisq) best_hit = *ihit;
1432  }
1433  }
1434  }
1435  if (fTrace) {
1436  log << "Best hit incremental chisquare = " << best_chisq << "\n";
1437  if (best_hit.get() != 0) {
1438  log << "Hit after prediction\n";
1439  log << *best_hit;
1440  }
1441  else
1442  log << "No hit passed chisquare cut.\n";
1443  }
1444  if (fGTrace && fCanvases.size() > 0) fCanvases.back()->Update();
1445 
1446  // If we found a best measurement, and if the incremental
1447  // chisquare passes the cut, add it to the track and update
1448  // fit information.
1449 
1450  bool update_ok = false;
1451  if (best_hit.get() != 0) {
1452  KFitTrack trf0(trf);
1453  best_hit->update(trf);
1454  update_ok = trf.isValid();
1455  if (!update_ok) trf = trf0;
1456  }
1457  if (update_ok) {
1458  ds += best_hit->getPredDistance();
1459  tchisq += best_chisq;
1460  trf.setChisq(tchisq);
1461  if (dir == Propagator::FORWARD)
1462  trf.setStat(KFitTrack::FORWARD);
1463  else {
1464  trf.setStat(KFitTrack::BACKWARD);
1465  }
1466 
1467  // If the pointing error got too large, quit.
1468 
1469  if (trf.PointingError() > fMaxPErr) {
1470  if (fTrace) log << "Quitting because pointing error got too large.\n";
1471  break;
1472  }
1473 
1474  // Test number of consecutive measurements in the same plane.
1475 
1476  if (gr.getPlane() >= 0) {
1477  if (gr.getPlane() == last_plane)
1478  ++nsame;
1479  else {
1480  nsame = 1;
1481  last_plane = gr.getPlane();
1482  }
1483  }
1484  else {
1485  nsame = 0;
1486  last_plane = -1;
1487  }
1488  if (nsame <= fMaxSamePlane) {
1489 
1490  // Turn best hit red.
1491 
1492  if (fGTrace && fCanvases.size() > 0) {
1493  int pl = best_hit->getMeasPlane();
1494  if (pl >= 0 && pl < int(fPads.size())) {
1495  auto marker_it = fMarkerMap.find(best_hit->getID());
1496  if (marker_it != fMarkerMap.end()) {
1497  TMarker* marker = marker_it->second;
1498  marker->SetMarkerColor(kRed);
1499 
1500  // Redraw marker so that it will be on top.
1501 
1502  fPads[pl]->cd();
1503  marker->Draw();
1504  }
1505  }
1506  fCanvases.back()->Update();
1507  }
1508 
1509  // Make a KHitTrack and add it to the KGTrack.
1510 
1511  KHitTrack trh(trf, best_hit);
1512  trg.addTrack(trh);
1513 
1514  if (fTrace) {
1515  log << "After update\n";
1516  log << "KGTrack has " << trg.numHits() << " hits.\n";
1517  log << trf;
1518  }
1519  }
1520  }
1521  }
1522 
1523  // The current KHitGroup is now resolved.
1524  // Move it to unused list.
1525 
1526  hits.getUnused().splice(hits.getUnused().end(), hits.getSorted(), it);
1527 
1528  // If the propagation distance was the wrong direction, resort the measurements.
1529 
1530  if (dist && ((dir == Propagator::FORWARD && (ds < fMinSortDist || ds > fMaxSortDist)) ||
1531  (dir == Propagator::BACKWARD && (-ds < fMinSortDist || -ds > fMaxSortDist)))) {
1532  if (fTrace) log << "Resorting measurements.\n";
1533  hits.sort(trf, true, prop, dir);
1534  }
1535  }
1536  }
1537 
1538  // Clean track.
1539 
1540  cleanTrack(trg);
1541 
1542  // Set the fit status of the last added KHitTrack to optimal and
1543  // get the final chisquare and path length.
1544 
1545  double fchisq = 0.;
1546  path = 0.;
1547  if (trg.isValid()) {
1548  path = trg.endTrack().getPath() - trg.startTrack().getPath();
1549  switch (dir) {
1550  case Propagator::FORWARD:
1551  trg.endTrack().setStat(KFitTrack::OPTIMAL);
1552  fchisq = trg.endTrack().getChisq();
1553  break;
1554  case Propagator::BACKWARD:
1555  trg.startTrack().setStat(KFitTrack::OPTIMAL);
1556  fchisq = trg.startTrack().getChisq();
1557  break;
1558  default:
1559  throw cet::exception("KalmanFilterAlg")
1560  << "KalmanFilterAlg::extendTrack(): invalid direction [II]\n";
1561  } // switch
1562  }
1563 
1564  // Summary.
1565 
1566  log << "KalmanFilterAlg extend track summary.\n"
1567  << "Extend direction = " << (dir == Propagator::FORWARD ? "FORWARD" : "BACKWARD") << "\n"
1568  << "Track has " << trg.numHits() << " hits.\n"
1569  << "Track length = " << path << "\n"
1570  << "Track chisquare = " << fchisq << "\n";
1571  // log << trg << "\n";
1572  }
1573 
1574  // Done.
1575 
1576  result = (trg.numHits() > nhits0);
1577 
1578  if (fGTrace && fCanvases.size() > 0) {
1579  TText* t = 0;
1580  if (result)
1581  t = fMessages->AddText("Exit extendTrack, status success");
1582  else
1583  t = fMessages->AddText("Exit extendTrack, status fail");
1584  t->SetBit(kCanDelete);
1585  fInfoPad->Modified();
1586  fCanvases.back()->Update();
1587  }
1588 
1589  return result;
1590 }
Float_t x
Definition: compare.C:6
bool fGTrace
Graphical trace flag.
intermediate_table::iterator iterator
double fMaxPropDist
Maximum propagation distance to candidate surface.
std::map< int, TMarker * > fMarkerMap
Markers in current canvas.
bool fTrace
Trace flag.
Double_t z
Definition: plot.C:276
constexpr auto abs(T v)
Returns the absolute value of the argument.
void cleanTrack(KGTrack &trg) const
Clean track by removing noise hits near endpoints.
FitStatus
Fit status enum.
Definition: KFitTrack.h:39
double fMaxSortDist
Sort high distance threshold.
double fMaxPredDist
Maximum prediciton distance to accept a hit.
auto vector(Vector const &v)
Returns a manipulator which will print the specified array.
Definition: DumpUtils.h:289
void hits()
Definition: readHits.C:15
TVirtualPad * fInfoPad
Information pad.
std::vector< std::unique_ptr< TCanvas > > fCanvases
Graphical trace canvases.
TPaveText * fMessages
Message box.
double fMaxIncChisq
Maximum incremental chisquare to accept a hit.
Detector simulation of raw signals on wires.
TDirectory * dir
Definition: macro.C:5
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
std::vector< TVirtualPad * > fPads
View pads in current canvas.
int fMaxSamePlane
Maximum consecutive hits in same plane.
double fMaxPErr
Maximum pointing error for free propagation.
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
PropDirection
Propagation direction enum.
Definition: Propagator.h:94
int fPlane
Preferred view plane.
bool trkf::KalmanFilterAlg::fitMomentum ( const KGTrack trg,
const Propagator prop,
KETrack tremom 
) const

Estimate track momentum using either range or multiple scattering.

Estimate track momentum using either range or multiple scattering.

Arguments:

trg - Global track whose momentum is to be updated. prop - Propagator. tremom - Track with momentum estimate.

Returns: True if success.

Definition at line 1853 of file KalmanFilterAlg.cxx.

References fFitMomMS, fFitMomRange, fitMomentumMS(), fitMomentumRange(), trkf::KETrack::getError(), trkf::KTrack::getVector(), trkf::KGTrack::startTrack(), and updateMomentum().

Referenced by trkf::Track3DKalmanHitAlg::fitnupdateMomentum(), and trkf::TrackKalmanCheater::produce().

1856 {
1857  mf::LogInfo log("KalmanFilterAlg");
1858  double invp_range = 0.;
1859  double invp_ms = 0.;
1860 
1861  // Get multiple scattering momentum estimate.
1862 
1863  KETrack tremom_ms;
1864  bool ok_ms = false;
1865  if (fFitMomMS) {
1866  ok_ms = fitMomentumMS(trg, prop, tremom_ms);
1867  if (ok_ms) {
1868  KGTrack trg_ms(trg);
1869  ok_ms = updateMomentum(tremom_ms, prop, trg_ms);
1870  if (ok_ms) {
1871  invp_ms = trg_ms.startTrack().getVector()(4);
1872  double var_invp = trg_ms.startTrack().getError()(4, 4);
1873  double p = 0.;
1874  if (invp_ms != 0.) p = 1. / invp_ms;
1875  double err_p = p * p * std::sqrt(var_invp);
1876  log << "Multiple scattering momentum estimate = " << p << "+-" << err_p << "\n";
1877  }
1878  }
1879  }
1880 
1881  // Get range momentum estimate.
1882 
1883  KETrack tremom_range;
1884  bool ok_range = false;
1885  if (fFitMomRange) {
1886  ok_range = fitMomentumRange(trg, tremom_range);
1887  if (ok_range) {
1888  KGTrack trg_range(trg);
1889  ok_range = updateMomentum(tremom_range, prop, trg_range);
1890  if (ok_range) {
1891  invp_range = trg_range.startTrack().getVector()(4);
1892  double var_invp = trg_range.startTrack().getError()(4, 4);
1893  double p = 0.;
1894  if (invp_range != 0.) p = 1. / invp_range;
1895  double err_p = p * p * std::sqrt(var_invp);
1896  log << "Range momentum estimate = " << p << "+-" << err_p << "\n";
1897  }
1898  }
1899  }
1900 
1901  bool result = false;
1902  if (ok_range) {
1903  tremom = tremom_range;
1904  result = ok_range;
1905  }
1906  else if (ok_ms) {
1907  tremom = tremom_ms;
1908  result = ok_ms;
1909  }
1910  return result;
1911 }
bool fFitMomRange
Fit momentum using range.
bool fFitMomMS
Fit momentum using multiple scattering.
bool fitMomentumMS(const KGTrack &trg, const Propagator &prop, KETrack &tremom) const
Estimate track momentum using multiple scattering.
bool fitMomentumRange(const KGTrack &trg, KETrack &tremom) const
Estimate track momentum using range.
bool updateMomentum(const KETrack &tremom, const Propagator &prop, KGTrack &trg) const
Set track momentum at each track surface.
bool trkf::KalmanFilterAlg::fitMomentumMS ( const KGTrack trg,
const Propagator prop,
KETrack tremom 
) const

Estimate track momentum using multiple scattering.

Estimate track momentum using multiple scattering.

Arguments:

trg - Global track. prop - Propagator. tremom - Track containing momentum estimate.

Returns: True if success.

This method estimates the momentum of the specified track using multiple scattering. As a result of calling this method, the original global track is not updated, but a KETrack is produced near the starting surface that has the estimated momentum.

The global track passed as argument should have been smoothed prior to calling this method so that all or most fits are optimal. If either the first or last fit is not optimal, return false.

This method assumes that track parameter four is 1/p. This sort of momentum estimation only makes sense if the momentum track parameter is uncorrelated with any other track parameter. The error matrix of the first and last fit is checked for this. If it is found that either error matrix has correlated track parameter four with any other track parameter, this method returns without doing anything (return false).

Definition at line 1659 of file KalmanFilterAlg.cxx.

References util::abs(), trkf::Propagator::err_prop(), fMinSampleDist, trkf::KETrack::getError(), trkf::KFitTrack::getStat(), trkf::KTrack::getSurface(), trkf::KGTrack::getTrackMap(), trkf::KTrack::getVector(), trkf::KTrack::Mass(), trkf::Propagator::noise_prop(), trkf::KFitTrack::OPTIMAL, trkf::Propagator::UNKNOWN, and trkf::Propagator::vec_prop().

Referenced by fitMomentum().

1662 {
1663  // Get iterators pointing to the first and last tracks.
1664 
1665  const std::multimap<double, KHitTrack>& trackmap = trg.getTrackMap();
1666  if (trackmap.size() < 2) return false;
1668  itend[0] = trackmap.begin();
1669  itend[1] = trackmap.end();
1670  --itend[1];
1671 
1672  // Check the fit status and error matrix of the first and last
1673  // track.
1674 
1675  bool result = true;
1676 
1677  for (int i = 0; result && i < 2; ++i) {
1678  const KHitTrack& trh = itend[i]->second;
1679  KFitTrack::FitStatus stat = trh.getStat();
1680  if (stat != KFitTrack::OPTIMAL) result = false;
1681  const TrackError& err = trh.getError();
1682  for (int j = 0; j < 4; ++j) {
1683  if (err(4, j) != 0.) result = false;
1684  }
1685  }
1686  if (!result) return result;
1687 
1688  // We will periodically sample the track trajectory. At each sample
1689  // point, collect the following information.
1690  //
1691  // 1. The path distance at the sample point.
1692  // 2. The original momentum of the track at the sample point and its error.
1693  // 3. One copy of the track (KETrack) that will be propagated without noise
1694  // (infinite momentum track).
1695  // 3. A second copy of the track (KETrack) that will be propagated with the
1696  // minimum allowed momentum (range out track).
1697  // 4. A third copy of the track (KETrack) that will propagated with noise
1698  // with some intermediate momentum (noise track).
1699  //
1700  // Collect the first sample from the maximum path distance track.
1701 
1702  double s_sample = itend[1]->first;
1703  const KETrack& tre = itend[1]->second;
1704  KETrack tre_inf(tre);
1705  KTrack trk_range(tre);
1706  KETrack tre_noise(tre);
1707  tre_inf.getVector()(4) = 0.;
1708  tre_inf.getError()(4, 4) = 0.;
1709  trk_range.getVector()(4) = 100.;
1710  tre_noise.getError()(4, 4) = 0.;
1711  tre_noise.getVector()(4) = 1.;
1712  tre_noise.getError()(4, 4) = 10.;
1713  double invp0 = tre_noise.getVector()(4);
1714  double var_invp0 = tre_noise.getError()(4, 4);
1715 
1716  // Loop over fits, starting at the high path distance (low momentum)
1717  // end.
1718 
1719  for (auto it = trackmap.rbegin(); it != trackmap.rend(); ++it) {
1720  double s = it->first;
1721  const KHitTrack& trh = it->second;
1722 
1723  // Ignore non-optimal fits.
1724 
1725  KFitTrack::FitStatus stat = trh.getStat();
1726  if (stat != KFitTrack::OPTIMAL) continue;
1727 
1728  // See if this track is far enough from the previous sample to
1729  // make a new sample point.
1730 
1731  if (std::abs(s - s_sample) > fMinSampleDist) {
1732 
1733  // Propagate tracks to the current track surface.
1734 
1735  std::shared_ptr<const Surface> psurf = trh.getSurface();
1736  std::optional<double> dist_inf = prop.err_prop(tre_inf, psurf, Propagator::UNKNOWN, false);
1737  std::optional<double> dist_range =
1738  prop.vec_prop(trk_range, psurf, Propagator::UNKNOWN, false);
1739  std::optional<double> dist_noise =
1740  prop.noise_prop(tre_noise, psurf, Propagator::UNKNOWN, true);
1741 
1742  // All propagations should normally succeed. If they don't,
1743  // ignore this sample for the purpose of updating the momentum.
1744 
1745  bool momentum_updated = false;
1746  if (dist_inf && dist_range && dist_noise) {
1747 
1748  // Extract the momentum at the new sample point.
1749 
1750  double invp1 = tre_noise.getVector()(4);
1751  double var_invp1 = tre_noise.getError()(4, 4);
1752 
1753  // Get the average momentum and error for this pair of
1754  // sample points, and other data.
1755 
1756  double invp = 0.5 * (invp0 + invp1);
1757  double var_invp = 0.5 * (var_invp0 + var_invp1);
1758  double mass = tre_inf.Mass();
1759  double beta = std::sqrt(1. + mass * mass * invp * invp);
1760  double invbp = invp / beta;
1761 
1762  // Extract slope subvectors and sub-error-matrices.
1763  // We have the following variables.
1764  //
1765  // slope0 - Predicted slope vector (from infinite momentum track).
1766  // slope1 - Measured slope vector (from new sample point).
1767  // defl - Deflection (slope residual = difference between measured
1768  // and predicted slope vector).
1769  // err0 - Slope error matrix of prediction.
1770  // err1 - Slope error matrix of measurement
1771  // errc - Slope residual error matrix = err0 + err1.
1772  // errn - Noise slope error matrix divided by (1/beta*momentum).
1773 
1774  KVector<2>::type slope0 = project(tre_inf.getVector(), ublas::range(2, 4));
1775  KVector<2>::type slope1 = project(trh.getVector(), ublas::range(2, 4));
1776  KVector<2>::type defl = slope1 - slope0;
1777  KSymMatrix<2>::type err0 =
1778  project(tre_inf.getError(), ublas::range(2, 4), ublas::range(2, 4));
1779  KSymMatrix<2>::type err1 = project(trh.getError(), ublas::range(2, 4), ublas::range(2, 4));
1780  KSymMatrix<2>::type errc = err0 + err1;
1781  KSymMatrix<2>::type errn =
1782  project(tre_noise.getError(), ublas::range(2, 4), ublas::range(2, 4));
1783  errn -= err0;
1784  errn /= invbp;
1785 
1786  // Calculate updated average momentum and error.
1787 
1788  double new_invp = invp;
1789  double new_var_invp = var_invp;
1790  update_momentum(defl, errc, errn, mass, new_invp, new_var_invp);
1791 
1792  // Calculate updated momentum and error at the second sample
1793  // point.
1794 
1795  double dp = 1. / new_invp - 1. / invp;
1796  invp0 = 1. / (1. / invp1 + dp);
1797  var_invp0 = new_var_invp;
1798  momentum_updated = true;
1799 
1800  // Make sure that updated momentum is not less than minimum
1801  // allowed momentum.
1802 
1803  double invp_range = trk_range.getVector()(4);
1804  if (invp0 > invp_range) invp0 = invp_range;
1805  }
1806 
1807  // Update sample.
1808 
1809  if (momentum_updated) {
1810  s_sample = s;
1811  tre_inf = trh;
1812  tre_inf.getVector()(4) = 0.;
1813  tre_inf.getError()(4, 4) = 0.;
1814  double invp_range = trk_range.getVector()(4);
1815  trk_range = trh;
1816  trk_range.getVector()(4) = invp_range;
1817  tre_noise = trh;
1818  tre_noise.getVector()(4) = invp0;
1819  tre_noise.getError()(4, 4) = var_invp0;
1820  }
1821  }
1822  }
1823 
1824  // Propagate noise track to starting (high momentum) track surface
1825  // to get final starting momentum. This propagation should normally
1826  // always succeed, but if it doesn't, don't update the track.
1827 
1828  const KHitTrack& trh0 = itend[0]->second;
1829  std::shared_ptr<const Surface> psurf = trh0.getSurface();
1830  std::optional<double> dist_noise = prop.noise_prop(tre_noise, psurf, Propagator::UNKNOWN, true);
1831  result = dist_noise.has_value();
1832 
1833  // Update momentum-estimating track.
1834 
1835  mf::LogInfo log("KalmanFilterAlg");
1836  if (result) tremom = tre_noise;
1837 
1838  // Done.
1839 
1840  return result;
1841 }
KSymMatrix< 5 >::type TrackError
Track error matrix, dimension 5x5.
constexpr auto abs(T v)
Returns the absolute value of the argument.
intermediate_table::const_iterator const_iterator
FitStatus
Fit status enum.
Definition: KFitTrack.h:39
ublas::symmetric_matrix< double, ublas::lower, ublas::row_major, ublas::bounded_array< double, N *(N+1)/2 > > type
ublas::vector< double, ublas::bounded_array< double, N > > type
double fMinSampleDist
Minimum sample distance (for momentum measurement).
bool trkf::KalmanFilterAlg::fitMomentumRange ( const KGTrack trg,
KETrack tremom 
) const

Estimate track momentum using range.

Estimate track momentum using range.

Arguments:

trg - Global track. prop - Propagator. tremom - Track with momentum estimate.

Returns: True if success.

This method generates a momentum-estimating track by extracting the last track from a global track, and setting its momentum to some small value.

Definition at line 1606 of file KalmanFilterAlg.cxx.

References trkf::KGTrack::endTrack(), trkf::KETrack::getError(), trkf::KFitTrack::getStat(), trkf::KTrack::getVector(), trkf::KFitTrack::INVALID, and trkf::KGTrack::isValid().

Referenced by fitMomentum().

1607 {
1608  if (!trg.isValid()) return false;
1609 
1610  // Extract track with lowest momentum.
1611 
1612  const KHitTrack& trh = trg.endTrack();
1613  if (trh.getStat() == KFitTrack::INVALID)
1614  throw cet::exception("KalmanFilterAlg")
1615  << "KalmanFilterAlg::fitMomentumRange(): invalid end track\n";
1616  tremom = trh;
1617 
1618  // Set track momentum to a small value.
1619 
1620  tremom.getVector()(4) = 100.;
1621  tremom.getError()(4, 0) = 0.;
1622  tremom.getError()(4, 1) = 0.;
1623  tremom.getError()(4, 2) = 0.;
1624  tremom.getError()(4, 3) = 0.;
1625  tremom.getError()(4, 4) = 10000.;
1626 
1627  // Done.
1628 
1629  return true;
1630 }
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
int trkf::KalmanFilterAlg::getPlane ( ) const
inline

Preferred view plane.

Definition at line 75 of file KalmanFilterAlg.h.

bool trkf::KalmanFilterAlg::getTrace ( ) const
inline

Trace config parameters.

Definition at line 74 of file KalmanFilterAlg.h.

void trkf::KalmanFilterAlg::setPlane ( int  plane)
inline

Set preferred view plane.

Definition at line 80 of file KalmanFilterAlg.h.

References dir, and hits().

Referenced by trkf::Track3DKalmanHitAlg::makeKalmanTracks(), and trkf::TrackKalmanCheater::produce().

void trkf::KalmanFilterAlg::setTrace ( bool  trace)
inline

Set trace config parameter.

Definition at line 79 of file KalmanFilterAlg.h.

References trkf::trace().

bool trkf::KalmanFilterAlg::smoothTrack ( KGTrack trg,
KGTrack trg1,
const Propagator prop 
) const

Smooth track.

Smooth track.

Arguments:

trg - Track to be smoothed. trg1 - Track to receive result of unidirectional fit. prop - Propagator.

Returns: True if success.

The starting track should be a global track that has been fit in one direction. Fit status should be optimal at (at least) one end. It is an error if the fit status is not optimal at either end. If the fit status is optimal at both ends, do nothing, but return success.

If the second argument is non-null, save the result of the unidirectional track fit produced as a byproduct of the smoothing operation. This track can be smoothed in order to iterate the Kalman fit, etc.

The Kalman smoothing algorithm starts at the optimal end and fits the track in the reverse direction, calculating optimal track parameters at each measurement surface.

All measurements are included in the reverse fit. No incremental chisquare cut is applied.

If any measurement surface can not be reached because of a measurement error, the entire smoothing operation is considered a failure. In that case, false is returned and the track is left in an undefined state.

Definition at line 845 of file KalmanFilterAlg.cxx.

References trkf::KGTrack::addTrack(), trkf::KFitTrack::BACKWARD, trkf::Propagator::BACKWARD, trkf::KFitTrack::BACKWARD_PREDICTED, trkf::KFitTrack::combineFit(), dir, larg4::dist(), trkf::KGTrack::endTrack(), fCanvases, fGTrace, fMarkerMap, fMaxSmoothIncChisq, trkf::KFitTrack::FORWARD, trkf::Propagator::FORWARD, trkf::KFitTrack::FORWARD_PREDICTED, fTrace, trkf::KFitTrack::getChisq(), trkf::KHitBase::getChisq(), trkf::KHitTrack::getHit(), trkf::KHitBase::getID(), trkf::KFitTrack::getPath(), trkf::KFitTrack::getStat(), trkf::KTrack::getSurface(), trkf::KGTrack::getTrackMap(), trkf::KGTrack::isValid(), trkf::KTrack::isValid(), trkf::Propagator::noise_prop(), trkf::KGTrack::numHits(), trkf::KFitTrack::OPTIMAL, trkf::KHitBase::predict(), trkf::KGTrack::recalibrate(), trkf::KFitTrack::setChisq(), trkf::KFitTrack::setPath(), trkf::KFitTrack::setStat(), trkf::KGTrack::startTrack(), trkf::Propagator::UNKNOWN, and trkf::KHitBase::update().

Referenced by trkf::Track3DKalmanHitAlg::extendandsmoothLoop(), trkf::Track3DKalmanHitAlg::smoothandextendTrack(), and smoothTrackIter().

846 {
847  if (not trg.isValid()) {
848  // It is an error if the KGTrack is not valid.
849  return false;
850  }
851 
852  bool result = false;
853 
854  // Examine the track endpoints and figure out which end of the track
855  // to start from. The fit always starts at the optimal end. It is
856  // an error if neither end point is optimal. Do nothing and return
857  // success if both end points are optimal.
858 
859  const KHitTrack& trh0 = trg.startTrack();
860  const KHitTrack& trh1 = trg.endTrack();
861  KFitTrack::FitStatus stat0 = trh0.getStat();
862  KFitTrack::FitStatus stat1 = trh1.getStat();
863  bool dofit = false;
864 
865  // Remember starting direction, track, and distance.
866 
868  const KTrack* trk = 0;
869  double path = 0.;
870 
871  if (stat0 == KFitTrack::OPTIMAL) {
872  if (stat1 == KFitTrack::OPTIMAL) {
873 
874  // Both ends optimal (do nothing, return success).
875 
876  dofit = false;
877  result = true;
878  }
879  else {
880 
881  // Start optimal.
882 
883  dofit = true;
884  dir = Propagator::FORWARD;
885  trk = &trh0;
886  path = 0.;
887  }
888  }
889  else {
890  if (stat1 == KFitTrack::OPTIMAL) {
891 
892  // End optimal.
893 
894  dofit = true;
895  dir = Propagator::BACKWARD;
896  trk = &trh1;
897  path = trh1.getPath();
898  }
899  else {
900 
901  // Neither end optimal (do nothing, return failure).
902 
903  dofit = false;
904  result = false;
905  }
906  }
907  if (dofit) {
908  if (!trk)
909  throw cet::exception("KalmanFilterAlg") << "KalmanFilterAlg::smoothTrack(): no track!\n";
910 
911  // Cumulative chisquare.
912 
913  double tchisq = 0.;
914 
915  // Construct starting KFitTrack with track information and distance
916  // taken from the optimal end, but with "infinite" errors.
917 
918  TrackError err;
919  trk->getSurface()->getStartingError(err);
920  KETrack tre(*trk, err);
921  KFitTrack trf(tre, path, tchisq);
922 
923  // Make initial reference track to be same as initial fit track.
924 
925  KTrack ref(trf);
926 
927  // Loop over KHitTracks contained in KGTrack.
928 
931  switch (dir) {
932  case Propagator::FORWARD:
933  it = trg.getTrackMap().begin();
934  itend = trg.getTrackMap().end();
935  break;
937  it = trg.getTrackMap().end();
938  itend = trg.getTrackMap().begin();
939  break;
940  default:
941  throw cet::exception("KalmanFilterAlg")
942  << "KalmanFilterAlg::smoothTrack(): invalid direction\n";
943  } // switch
944 
945  mf::LogInfo log("KalmanFilterAlg");
946 
947  // Loop starts here.
948 
949  result = true; // Result success unless we find an error.
950  int step = 0; // Step count.
951  while (dofit && it != itend) {
952  ++step;
953  if (fTrace) {
954  log << "Smooth Step " << step << "\n";
955  log << "Reverse fit track:\n";
956  log << trf;
957  }
958 
959  // For backward fit, decrement iterator at start of loop.
960 
961  if (dir == Propagator::BACKWARD) --it;
962 
963  KHitTrack& trh = (*it).second;
964  if (fTrace) {
965  log << "Forward track:\n";
966  log << trh;
967  }
968 
969  // Extract measurement.
970 
971  const KHitBase& hit = *(trh.getHit());
972 
973  // Propagate KFitTrack to the next track surface.
974 
975  std::shared_ptr<const Surface> psurf = trh.getSurface();
976  std::optional<double> dist = prop.noise_prop(trf, psurf, Propagator::UNKNOWN, true, &ref);
977 
978  // Check if propagation succeeded. If propagation fails, this
979  // measurement will be dropped from the unidirectional fit
980  // track. This measurement will still be in the original
981  // track, but with a status other than optimal.
982 
983  if (dist) {
984 
985  // Propagation succeeded.
986  // Update cumulative path distance and track status.
987 
988  double ds = *dist;
989  path += ds;
990  trf.setPath(path);
991  if (dir == Propagator::FORWARD)
992  trf.setStat(KFitTrack::FORWARD_PREDICTED);
993  else {
994  trf.setStat(KFitTrack::BACKWARD_PREDICTED);
995  }
996  if (fTrace) {
997  log << "Reverse fit track after propagation:\n";
998  log << " Propagation distance = " << ds << "\n";
999  log << trf;
1000  }
1001 
1002  // See if we have the proper information to calculate an optimal track
1003  // at this surface (should normally be possible).
1004 
1005  KFitTrack::FitStatus stat = trh.getStat();
1006  KFitTrack::FitStatus newstat = trf.getStat();
1007 
1008  if ((newstat == KFitTrack::FORWARD_PREDICTED && stat == KFitTrack::BACKWARD) ||
1009  (newstat == KFitTrack::BACKWARD_PREDICTED && stat == KFitTrack::FORWARD)) {
1010 
1011  // Update stored KHitTrack to be optimal.
1012 
1013  bool ok = trh.combineFit(trf);
1014 
1015  // Update the stored path distance to be from the currently fitting track.
1016 
1017  trh.setPath(trf.getPath());
1018 
1019  // Update reference track.
1020 
1021  ref = trh;
1022 
1023  // If combination failed, abandon the fit and return failure.
1024 
1025  if (!ok) {
1026  dofit = false;
1027  result = false;
1028  break;
1029  }
1030  if (fTrace) {
1031  log << "Combined track:\n";
1032  log << trh;
1033  }
1034  }
1035 
1036  // Update measurement predction using current track hypothesis.
1037 
1038  if (not hit.predict(trf, prop, &ref)) {
1039 
1040  // Abandon the fit and return failure.
1041 
1042  dofit = false;
1043  result = false;
1044  break;
1045  }
1046 
1047  // Prediction succeeded. Get incremental chisquare. If this
1048  // hit fails incremental chisquare cut, this hit will be
1049  // dropped from the unidirecitonal Kalman fit track, but may
1050  // still be in the smoothed track.
1051 
1052  double chisq = hit.getChisq();
1053  if (chisq < fMaxSmoothIncChisq) {
1054 
1055  // Update the reverse fitting track using the current measurement
1056  // (both track parameters and status).
1057 
1058  KFitTrack trf0(trf);
1059  hit.update(trf);
1060  bool update_ok = trf.isValid();
1061  if (!update_ok)
1062  trf = trf0;
1063  else {
1064  tchisq += chisq;
1065  trf.setChisq(tchisq);
1066 
1067  if (dir == Propagator::FORWARD)
1068  trf.setStat(KFitTrack::FORWARD);
1069  else {
1070  trf.setStat(KFitTrack::BACKWARD);
1071  }
1072  if (fTrace) {
1073  log << "Reverse fit track after update:\n";
1074  log << trf;
1075  }
1076  if (fGTrace && fCanvases.size() > 0) {
1077  auto marker_it = fMarkerMap.find(hit.getID());
1078  if (marker_it != fMarkerMap.end()) {
1079  TMarker* marker = marker_it->second;
1080  marker->SetMarkerColor(kOrange);
1081  }
1082  fCanvases.back()->Update();
1083  }
1084 
1085  // If unidirectional track pointer is not null, make a
1086  // KHitTrack and save it in the unidirectional track.
1087 
1088  if (trg1 != 0) {
1089  KHitTrack trh1(trf, trh.getHit());
1090  trg1->addTrack(trh1);
1091  }
1092  }
1093  }
1094  }
1095 
1096  // For forward fit, increment iterator at end of loop.
1097 
1098  if (dir == Propagator::FORWARD) ++it;
1099 
1100  } // Loop over KHitTracks.
1101 
1102  // If fit was successful and the unidirectional track pointer
1103  // is not null and the track is valid, set the fit status of
1104  // the last added KHitTrack to optimal.
1105 
1106  if (result && trg1 != 0 && trg1->isValid()) {
1107  if (dir == Propagator::FORWARD)
1108  trg1->endTrack().setStat(KFitTrack::OPTIMAL);
1109  else {
1110  trg1->startTrack().setStat(KFitTrack::OPTIMAL);
1111  }
1112  }
1113 
1114  // Recalibrate track map.
1115 
1116  trg.recalibrate();
1117 
1118  } // Do fit.
1119 
1120  // Get the final chisquare.
1121 
1122  double fchisq = 0.5 * (trg.startTrack().getChisq() + trg.endTrack().getChisq());
1123 
1124  // Summary.
1125 
1126  mf::LogInfo log("KalmanFilterAlg");
1127  log << "KalmanFilterAlg smooth track summary.\n"
1128  << "Smooth direction = " << (dir == Propagator::FORWARD ? "FORWARD" : "BACKWARD") << "\n"
1129  << "Track has " << trg.numHits() << " hits.\n"
1130  << "Track length = " << trg.endTrack().getPath() - trg.startTrack().getPath() << "\n"
1131  << "Track chisquare = " << fchisq << "\n";
1132 
1133  return result;
1134 }
bool fGTrace
Graphical trace flag.
intermediate_table::iterator iterator
std::map< int, TMarker * > fMarkerMap
Markers in current canvas.
KSymMatrix< 5 >::type TrackError
Track error matrix, dimension 5x5.
bool fTrace
Trace flag.
FitStatus
Fit status enum.
Definition: KFitTrack.h:39
std::vector< std::unique_ptr< TCanvas > > fCanvases
Graphical trace canvases.
double fMaxSmoothIncChisq
Maximum incremental chisquare to accept a hit in smooth phase.
Detector simulation of raw signals on wires.
TDirectory * dir
Definition: macro.C:5
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
PropDirection
Propagation direction enum.
Definition: Propagator.h:94
bool trkf::KalmanFilterAlg::smoothTrackIter ( int  nsmooth,
KGTrack trg,
const Propagator prop 
) const

Iteratively smooth a track.

Iteratively smooth a track.

Arguments:

nsmooth - Number of iterations. trg - Track to be smoothed. prop - Propagator.

Returns: True if success.

The initial track should have been unidirectionally fit.

Definition at line 2069 of file KalmanFilterAlg.cxx.

References trkf::KGTrack::getPrefPlane(), and smoothTrack().

Referenced by trkf::TrackKalmanCheater::produce().

2070 {
2071  bool ok = true;
2072 
2073  // Call smoothTrack method in a loop.
2074 
2075  for (int ismooth = 0; ok && ismooth < nsmooth - 1; ++ismooth) {
2076  KGTrack trg1{trg.getPrefPlane()};
2077  ok = smoothTrack(trg, &trg1, prop);
2078  if (ok) trg = trg1;
2079  }
2080 
2081  // Do one final smooth without generating a new unidirectional
2082  // track.
2083 
2084  if (ok) ok = smoothTrack(trg, 0, prop);
2085 
2086  // Done.
2087 
2088  return ok;
2089 }
bool smoothTrack(KGTrack &trg, KGTrack *trg1, const Propagator &prop) const
Smooth track.
bool trkf::KalmanFilterAlg::updateMomentum ( const KETrack tremom,
const Propagator prop,
KGTrack trg 
) const

Set track momentum at each track surface.

Set track momentum at each track surface.

Arguments:

tremom - Track containing momentum estimate. prop - Propagator. trg - Global track to be updated.

The track containing the momentum estimate is propagated to the first or last track fit of the global track (whichever is closer), then the momentum estimate is transfered to that track fit. In similar fashion, the momentum estimate is successively transfered from that track fit to each track fit of the global track.

Only momentum track parameters of the global track fits are updated. Other track parameters and their errors are unmodified. Unreachable track fits are deleted from the global track. Overall failure will occur if the momentum-estimating track can't be propagated to the initial track fit, or if the final global track has no valid track fits.

Definition at line 1934 of file KalmanFilterAlg.cxx.

References util::abs(), larg4::dist(), util::empty(), trkf::KETrack::getError(), trkf::KFitTrack::getStat(), trkf::KTrack::getSurface(), trkf::KGTrack::getTrackMap(), trkf::KTrack::getVector(), trkf::KFitTrack::INVALID, trkf::Propagator::noise_prop(), trkf::Propagator::UNKNOWN, and trkf::Propagator::vec_prop().

Referenced by fitMomentum(), trkf::Track3DKalmanHitAlg::fitnupdateMomentum(), and trkf::TrackKalmanCheater::produce().

1937 {
1938  // Get modifiable track map.
1939 
1940  std::multimap<double, KHitTrack>& trackmap = trg.getTrackMap();
1941 
1942  // If track map is empty, immediately return failure.
1943 
1944  if (trackmap.size() == 0) return false;
1945 
1946  // Make trial propagations to the first and last track fit to figure
1947  // out which track fit is closer to the momentum estimating track.
1948 
1949  // Find distance to first track fit.
1950 
1951  KETrack tre0(tremom);
1952  std::optional<double> dist0 =
1953  prop.vec_prop(tre0, trackmap.begin()->second.getSurface(), Propagator::UNKNOWN, false, 0, 0);
1954  // Find distance to last track fit.
1955 
1956  KETrack tre1(tremom);
1957  std::optional<double> dist1 =
1958  prop.vec_prop(tre1, trackmap.rbegin()->second.getSurface(), Propagator::UNKNOWN, false, 0, 0);
1959 
1960  // Based on distances, make starting iterator and direction flag.
1961 
1962  bool forward = true;
1963  std::multimap<double, KHitTrack>::iterator it = trackmap.begin();
1964  if (dist0) {
1965 
1966  // Propagation to first track succeeded.
1967 
1968  if (dist1) {
1969 
1970  // Propagation to both ends succeeded. If the last track is
1971  // closer, initialize iterator and direction flag for reverse
1972  // direction.
1973 
1974  if (std::abs(*dist0) > std::abs(*dist1)) {
1975  it = trackmap.end();
1976  --it;
1977  forward = false;
1978  }
1979  }
1980  }
1981  else {
1982 
1983  // Propagation to first track failed. Initialize iterator and
1984  // direction flag for reverse direction, provided that the
1985  // propagation to the last track succeeded.
1986 
1987  if (dist1) {
1988  it = trackmap.end();
1989  --it;
1990  forward = false;
1991  }
1992  else {
1993 
1994  // Propagation to both ends failed. Return failure.
1995 
1996  return false;
1997  }
1998  }
1999 
2000  // Loop over track fits in global track.
2001 
2002  KETrack tre(tremom);
2003  bool done = false;
2004  while (!done) {
2005  KHitTrack& trh = it->second;
2006  if (trh.getStat() == KFitTrack::INVALID)
2007  throw cet::exception("KalmanFilterAlg")
2008  << "KalmanFilterAlg::updateMomentum(): invalid track\n";
2009 
2010  // Propagate momentum-estimating track to current track surface
2011  // and update momentum.
2012 
2013  std::optional<double> dist = prop.noise_prop(tre, trh.getSurface(), Propagator::UNKNOWN, true);
2014 
2015  // Copy momentum to global track.
2016 
2017  std::multimap<double, KHitTrack>::iterator erase_it = trackmap.end();
2018  if (dist) {
2019  trh.getVector()(4) = tre.getVector()(4);
2020  trh.getError()(4, 0) = 0.;
2021  trh.getError()(4, 1) = 0.;
2022  trh.getError()(4, 2) = 0.;
2023  trh.getError()(4, 3) = 0.;
2024  trh.getError()(4, 4) = tre.getError()(4, 4);
2025  }
2026  else {
2027 
2028  // If the propagation failed, remember that we are supposed to
2029  // erase this track from the global track.
2030 
2031  erase_it = it;
2032  }
2033 
2034  // Advance the iterator and set the done flag.
2035 
2036  if (forward) {
2037  ++it;
2038  done = (it == trackmap.end());
2039  }
2040  else {
2041  done = (it == trackmap.begin());
2042  if (!done) --it;
2043  }
2044 
2045  // Update momentum-estimating track from just-updated global track
2046  // fit, or erase global track fit.
2047 
2048  if (erase_it == trackmap.end())
2049  tre = trh;
2050  else
2051  trackmap.erase(erase_it);
2052  }
2053 
2054  return not empty(trackmap);
2055 }
intermediate_table::iterator iterator
constexpr auto abs(T v)
Returns the absolute value of the argument.
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
decltype(auto) constexpr empty(T &&obj)
ADL-aware version of std::empty.
Definition: StdUtils.h:109

Member Data Documentation

std::vector<std::unique_ptr<TCanvas> > trkf::KalmanFilterAlg::fCanvases
mutableprivate

Graphical trace canvases.

Definition at line 163 of file KalmanFilterAlg.h.

Referenced by buildTrack(), cleanTrack(), extendTrack(), and smoothTrack().

bool trkf::KalmanFilterAlg::fFitMomMS
private

Fit momentum using multiple scattering.

Definition at line 151 of file KalmanFilterAlg.h.

Referenced by fitMomentum(), and KalmanFilterAlg().

bool trkf::KalmanFilterAlg::fFitMomRange
private

Fit momentum using range.

Definition at line 150 of file KalmanFilterAlg.h.

Referenced by fitMomentum(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fGapDist
private

Minimum gap distance.

Definition at line 147 of file KalmanFilterAlg.h.

Referenced by cleanTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fGoodPErr
private

Pointing error threshold for switching to free propagation.

Definition at line 134 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and KalmanFilterAlg().

bool trkf::KalmanFilterAlg::fGTrace
private

Graphical trace flag.

Definition at line 152 of file KalmanFilterAlg.h.

Referenced by buildTrack(), cleanTrack(), extendTrack(), KalmanFilterAlg(), and smoothTrack().

double trkf::KalmanFilterAlg::fGTraceWH
private

Window height.

Definition at line 154 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fGTraceWW
private

Window width.

Definition at line 153 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fGTraceXMax
private

Graphical trace maximum x.

Definition at line 156 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fGTraceXMin
private

Graphical trace minimum x.

Definition at line 155 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and KalmanFilterAlg().

std::vector<double> trkf::KalmanFilterAlg::fGTraceZMax
private

Graphical trace maximum z for each view.

Definition at line 158 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and KalmanFilterAlg().

std::vector<double> trkf::KalmanFilterAlg::fGTraceZMin
private

Graphical trace minimum z for each view.

Definition at line 157 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and KalmanFilterAlg().

TVirtualPad* trkf::KalmanFilterAlg::fInfoPad
mutableprivate

Information pad.

Definition at line 165 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and extendTrack().

std::map<int, TMarker*> trkf::KalmanFilterAlg::fMarkerMap
mutableprivate

Markers in current canvas.

Definition at line 167 of file KalmanFilterAlg.h.

Referenced by buildTrack(), cleanTrack(), extendTrack(), and smoothTrack().

double trkf::KalmanFilterAlg::fMaxEndChisq
private

Maximum incremental chisquare for endpoint hit.

Definition at line 138 of file KalmanFilterAlg.h.

Referenced by cleanTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fMaxIncChisq
private

Maximum incremental chisquare to accept a hit.

Definition at line 135 of file KalmanFilterAlg.h.

Referenced by extendTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fMaxLDist
private

Maximum distance for linearized propagation.

Definition at line 140 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and KalmanFilterAlg().

int trkf::KalmanFilterAlg::fMaxNoiseHits
private

Maximum number of hits in noise cluster.

Definition at line 148 of file KalmanFilterAlg.h.

Referenced by cleanTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fMaxPErr
private

Maximum pointing error for free propagation.

Definition at line 133 of file KalmanFilterAlg.h.

Referenced by buildTrack(), extendTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fMaxPredDist
private

Maximum prediciton distance to accept a hit.

Definition at line 141 of file KalmanFilterAlg.h.

Referenced by extendTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fMaxPropDist
private

Maximum propagation distance to candidate surface.

Definition at line 143 of file KalmanFilterAlg.h.

Referenced by buildTrack(), extendTrack(), and KalmanFilterAlg().

int trkf::KalmanFilterAlg::fMaxSamePlane
private

Maximum consecutive hits in same plane.

Definition at line 146 of file KalmanFilterAlg.h.

Referenced by buildTrack(), extendTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fMaxSeedIncChisq
private

Maximum incremental chisquare to accept a hit in seed phase.

Definition at line 136 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fMaxSeedPredDist
private

Maximum prediciton distance to accept a hit in seed phase.

Definition at line 142 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fMaxSmoothIncChisq
private

Maximum incremental chisquare to accept a hit in smooth phase.

Definition at line 137 of file KalmanFilterAlg.h.

Referenced by KalmanFilterAlg(), and smoothTrack().

double trkf::KalmanFilterAlg::fMaxSortDist
private

Sort high distance threshold.

Definition at line 145 of file KalmanFilterAlg.h.

Referenced by buildTrack(), extendTrack(), and KalmanFilterAlg().

TPaveText* trkf::KalmanFilterAlg::fMessages
mutableprivate

Message box.

Definition at line 166 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and extendTrack().

int trkf::KalmanFilterAlg::fMinLHits
private

Minimum number of hits to turn off linearized propagation.

Definition at line 139 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fMinSampleDist
private

Minimum sample distance (for momentum measurement).

Definition at line 149 of file KalmanFilterAlg.h.

Referenced by fitMomentumMS(), and KalmanFilterAlg().

double trkf::KalmanFilterAlg::fMinSortDist
private

Sort low distance threshold.

Definition at line 144 of file KalmanFilterAlg.h.

Referenced by KalmanFilterAlg().

std::vector<TVirtualPad*> trkf::KalmanFilterAlg::fPads
mutableprivate

View pads in current canvas.

Definition at line 164 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and extendTrack().

int trkf::KalmanFilterAlg::fPlane
private

Preferred view plane.

Definition at line 162 of file KalmanFilterAlg.h.

Referenced by buildTrack(), and extendTrack().

bool trkf::KalmanFilterAlg::fTrace
private

Trace flag.

Definition at line 132 of file KalmanFilterAlg.h.

Referenced by buildTrack(), extendTrack(), KalmanFilterAlg(), and smoothTrack().


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