CMS 3D CMS Logo

List of all members | Public Member Functions | Protected Member Functions | Protected Attributes
LocalTrackFitter Class Reference

Performs straight-line fit and outlier rejection. More...

#include <LocalTrackFitter.h>

Public Member Functions

bool fit (HitCollection &, const AlignmentGeometry &, LocalTrackFit &) const
 
 LocalTrackFitter ()
 dummy constructor (not to be used) More...
 
 LocalTrackFitter (const edm::ParameterSet &)
 normal constructor More...
 
virtual ~LocalTrackFitter ()
 

Protected Member Functions

void fitAndRemoveOutliers (HitCollection &, const AlignmentGeometry &, LocalTrackFit &, bool &failed, bool &selectionChanged) const
 
void removeInsufficientPots (HitCollection &, bool &selectionChanged) const
 removes the hits of pots with too few planes active More...
 

Protected Attributes

double maxResidualToSigma
 hits with higher ratio residual/sigma will be dropped More...
 
unsigned int minimumHitsPerProjectionPerRP
 minimum of hits to accept data from a RP More...
 
unsigned int verbosity
 verbosity level More...
 

Detailed Description

Performs straight-line fit and outlier rejection.

Definition at line 20 of file LocalTrackFitter.h.

Constructor & Destructor Documentation

◆ LocalTrackFitter() [1/2]

LocalTrackFitter::LocalTrackFitter ( )
inline

dummy constructor (not to be used)

Definition at line 23 of file LocalTrackFitter.h.

23 {}

◆ LocalTrackFitter() [2/2]

LocalTrackFitter::LocalTrackFitter ( const edm::ParameterSet ps)

normal constructor

Definition at line 24 of file LocalTrackFitter.cc.

25  : verbosity(ps.getUntrackedParameter<unsigned int>("verbosity", 0)),
26  minimumHitsPerProjectionPerRP(ps.getParameter<unsigned int>("minimumHitsPerProjectionPerRP")),
27  maxResidualToSigma(ps.getParameter<double>("maxResidualToSigma")) {}
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
unsigned int minimumHitsPerProjectionPerRP
minimum of hits to accept data from a RP
T getUntrackedParameter(std::string const &, T const &) const
unsigned int verbosity
verbosity level
double maxResidualToSigma
hits with higher ratio residual/sigma will be dropped

◆ ~LocalTrackFitter()

virtual LocalTrackFitter::~LocalTrackFitter ( )
inlinevirtual

Definition at line 28 of file LocalTrackFitter.h.

28 {}

Member Function Documentation

◆ fit()

bool LocalTrackFitter::fit ( HitCollection selection,
const AlignmentGeometry geometry,
LocalTrackFit trackFit 
) const

runs the fit and outlier-removal loop returns true in case of success

Definition at line 31 of file LocalTrackFitter.cc.

References fitAndRemoveOutliers(), removeInsufficientPots(), corrVsCorr::selection, and verbosity.

Referenced by trackingPlots.Iteration::modules(), and StraightTrackAlignment::processEvent().

31  {
32  if (verbosity > 5)
33  printf(">> LocalTrackFitter::Fit\n");
34 
35  bool selectionChanged = true;
36  unsigned int loopCounter = 0;
37  while (selectionChanged) {
38  // fit/outlier-removal loop
39  while (selectionChanged) {
40  if (verbosity > 5)
41  printf("* fit loop %u\n", loopCounter++);
42 
43  bool fitFailed = false;
44  fitAndRemoveOutliers(selection, geometry, trackFit, fitFailed, selectionChanged);
45 
46  if (fitFailed) {
47  if (verbosity > 5)
48  printf("\tFIT FAILED\n");
49  return false;
50  }
51  }
52 
53  // remove pots with too few active planes
54  removeInsufficientPots(selection, selectionChanged);
55  }
56 
57  return true;
58 }
selection
main part
Definition: corrVsCorr.py:100
unsigned int verbosity
verbosity level
void removeInsufficientPots(HitCollection &, bool &selectionChanged) const
removes the hits of pots with too few planes active
void fitAndRemoveOutliers(HitCollection &, const AlignmentGeometry &, LocalTrackFit &, bool &failed, bool &selectionChanged) const

◆ fitAndRemoveOutliers()

void LocalTrackFitter::fitAndRemoveOutliers ( HitCollection selection,
const AlignmentGeometry geometry,
LocalTrackFit trackFit,
bool &  failed,
bool &  selectionChanged 
) const
protected

fits the collection of hits and removes hits with too high residual/sigma ratio

Parameters
failedwhether the fit has failed
selectionChangedwhether some hits have been removed

Definition at line 62 of file LocalTrackFitter.cc.

References A, LocalTrackFit::ax, LocalTrackFit::ay, LocalTrackFit::bx, LocalTrackFit::by, LocalTrackFit::chi_sq, ztail::d, hcalRecHitTable_cff::detId, MillePedeFileConverter_cfg::e, runEdmFileComparison::failed, mps_fire::i, dqmiolumiharvest::j, maxResidualToSigma, LocalTrackFit::ndf, printId(), dttmaxenums::R, corrVsCorr::selection, theta(), verbosity, and LocalTrackFit::z0.

Referenced by fit().

66  {
67  if (verbosity > 5)
68  printf(" - LocalTrackFitter::FitAndRemoveOutliers\n");
69 
70  if (selection.empty()) {
71  failed = true;
72  return;
73  }
74 
75  // check if input size is sufficient
76  if (selection.size() < 4) {
77  failed = true;
78  return;
79  }
80 
81  // build matrices and vectors
82  TMatrixD A(selection.size(), 4);
83  TMatrixD Vi(selection.size(), selection.size());
84  TVectorD measVec(selection.size());
85  unsigned int j = 0;
86  for (auto it = selection.begin(); it != selection.end(); ++it, ++j) {
87  const unsigned int &detId = it->id;
88 
89  const DetGeometry &d = geometry.get(detId);
90  const auto &dirData = d.getDirectionData(it->dirIdx);
91 
92  A(j, 0) = it->z * dirData.dx;
93  A(j, 1) = dirData.dx;
94  A(j, 2) = it->z * dirData.dy;
95  A(j, 3) = dirData.dy;
96 
97  measVec(j) = it->position + dirData.s - (it->z - d.z) * dirData.dz; // in mm
98 
99  Vi(j, j) = 1. / it->sigma / it->sigma;
100  }
101 
102  // evaluate local track parameter estimates (h stands for hat)
103  TMatrixD AT(4, selection.size());
104  AT.Transpose(A);
105  TMatrixD ATViA(4, 4);
106  ATViA = AT * Vi * A;
107  TMatrixD ATViAI(ATViA);
108 
109  try {
110  ATViAI = ATViA.Invert();
111  } catch (cms::Exception &e) {
112  failed = true;
113  return;
114  }
115 
116  TVectorD theta(4);
117  theta = ATViAI * AT * Vi * measVec;
118 
119  // residuals
120  TVectorD R(measVec);
121  R -= A * theta;
122 
123  // save results to trackFit
124  trackFit.ax = theta(0);
125  trackFit.bx = theta(1);
126  trackFit.ay = theta(2);
127  trackFit.by = theta(3);
128  trackFit.z0 = geometry.z0;
129  trackFit.ndf = selection.size() - 4;
130  trackFit.chi_sq = 0;
131  for (int i = 0; i < R.GetNrows(); i++)
132  trackFit.chi_sq += R(i) * R(i) * Vi(i, i);
133 
134  if (verbosity > 5) {
135  printf(" ax = %.3f mrad, bx = %.4f mm, ay = %.3f mrad, by = %.4f mm, z0 = %.3f mm\n",
136  trackFit.ax * 1E3,
137  trackFit.bx,
138  trackFit.ay * 1E3,
139  trackFit.by,
140  trackFit.z0);
141  printf(" ndof = %i, chi^2/ndof/si^2 = %.3f\n", trackFit.ndf, trackFit.chi_sq / trackFit.ndf);
142  }
143 
144  // check residuals
145  selectionChanged = false;
146  TVectorD interpolation(A * theta);
147  j = 0;
148  for (auto it = selection.begin(); it != selection.end(); ++j) {
149  if (verbosity > 5) {
150  printf(" %2u, ", j);
151  printId(it->id);
152  printf(", dirIdx=%u: interpol = %+8.1f um, residual = %+6.1f um, residual / sigma = %+6.2f\n",
153  it->dirIdx,
154  interpolation[j] * 1E3,
155  R[j] * 1E3,
156  R[j] / it->sigma);
157  }
158 
159  double resToSigma = R[j] / it->sigma;
160  if (fabs(resToSigma) > maxResidualToSigma) {
161  selection.erase(it);
162  selectionChanged = true;
163  if (verbosity > 5)
164  printf(" Removed\n");
165  } else
166  ++it;
167  }
168 }
double z0
the point where intercepts are measured, in mm
Definition: LocalTrackFit.h:17
A structure to hold relevant geometrical information about one detector/sensor.
selection
main part
Definition: corrVsCorr.py:100
unsigned int verbosity
verbosity level
signed int ndf
the number of degrees of freedom
Definition: LocalTrackFit.h:26
double maxResidualToSigma
hits with higher ratio residual/sigma will be dropped
double chi_sq
the residual sum of squares
Definition: LocalTrackFit.h:29
d
Definition: ztail.py:151
void printId(unsigned int id)
Definition: Utilities.cc:26
Definition: APVGainStruct.h:7
Geom::Theta< T > theta() const
double bx
intercepts in mm
Definition: LocalTrackFit.h:23
double ax
slopes in rad
Definition: LocalTrackFit.h:20

◆ removeInsufficientPots()

void LocalTrackFitter::removeInsufficientPots ( HitCollection selection,
bool &  selectionChanged 
) const
protected

removes the hits of pots with too few planes active

Definition at line 172 of file LocalTrackFitter.cc.

References minimumHitsPerProjectionPerRP, TotemRPDetId::plane(), CTPPSDetId::rpId(), CTPPSDetId::sdTrackingPixel, CTPPSDetId::sdTrackingStrip, corrVsCorr::selection, DetId::subdetId(), and verbosity.

Referenced by fit().

172  {
173  if (verbosity > 5)
174  printf(" - RemoveInsufficientPots\n");
175 
176  selectionChanged = false;
177 
178  // map: RP id -> (active planes in projection 1, active planes in projection 2)
179  map<unsigned int, pair<set<unsigned int>, set<unsigned int> > > planeMap;
180  for (auto it = selection.begin(); it != selection.end(); ++it) {
181  CTPPSDetId senId(it->id);
182  const unsigned int rpId = senId.rpId();
183 
184  if (senId.subdetId() == CTPPSDetId::sdTrackingStrip) {
185  const unsigned int plane = TotemRPDetId(it->id).plane();
186  if ((plane % 2) == 0)
187  planeMap[rpId].first.insert(senId);
188  else
189  planeMap[rpId].second.insert(senId);
190  }
191 
192  if (senId.subdetId() == CTPPSDetId::sdTrackingPixel) {
193  planeMap[rpId].first.insert(senId);
194  planeMap[rpId].second.insert(senId);
195  }
196  }
197 
198  // remove RPs with insufficient information
199  selectionChanged = false;
200 
201  for (auto it = planeMap.begin(); it != planeMap.end(); ++it) {
202  if (it->second.first.size() < minimumHitsPerProjectionPerRP ||
203  it->second.second.size() < minimumHitsPerProjectionPerRP) {
204  if (verbosity > 5)
205  printf("\tRP %u: projection1 = %lu, projection 2 = %lu\n",
206  it->first,
207  it->second.first.size(),
208  it->second.second.size());
209 
210  // remove all hits from that RP
211  for (auto dit = selection.begin(); dit != selection.end();) {
212  if (it->first == CTPPSDetId(dit->id).rpId()) {
213  if (verbosity > 5)
214  printf("\t\tremoving %u\n", dit->id);
215  selection.erase(dit);
216  selectionChanged = true;
217  } else {
218  dit++;
219  }
220  }
221  }
222  }
223 }
Detector ID class for TOTEM Si strip detectors.
Definition: TotemRPDetId.h:30
unsigned int minimumHitsPerProjectionPerRP
minimum of hits to accept data from a RP
selection
main part
Definition: corrVsCorr.py:100
uint32_t plane() const
Definition: TotemRPDetId.h:46
CTPPSDetId rpId() const
Definition: CTPPSDetId.h:84
unsigned int verbosity
verbosity level
Base class for CTPPS detector IDs.
Definition: CTPPSDetId.h:32

Member Data Documentation

◆ maxResidualToSigma

double LocalTrackFitter::maxResidualToSigma
protected

hits with higher ratio residual/sigma will be dropped

Definition at line 42 of file LocalTrackFitter.h.

Referenced by fitAndRemoveOutliers().

◆ minimumHitsPerProjectionPerRP

unsigned int LocalTrackFitter::minimumHitsPerProjectionPerRP
protected

minimum of hits to accept data from a RP

Definition at line 39 of file LocalTrackFitter.h.

Referenced by removeInsufficientPots().

◆ verbosity

unsigned int LocalTrackFitter::verbosity
protected

verbosity level

Definition at line 36 of file LocalTrackFitter.h.

Referenced by fit(), fitAndRemoveOutliers(), and removeInsufficientPots().