CMS 3D CMS Logo

List of all members | Public Types | Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | Private Member Functions
MuonSegFit Class Reference

#include <MuonSegFit.h>

Public Types

typedef std::vector< MuonRecHitPtrMuonRecHitContainer
 
typedef std::shared_ptr< TrackingRecHitMuonRecHitPtr
 
typedef ROOT::Math::SMatrix< double, MaxHits2, 4 > SMatrix12by4
 
typedef ROOT::Math::SMatrix< double, 4 > SMatrix4
 
typedef ROOT::Math::SMatrix< double, MaxHits2, MaxHits2, ROOT::Math::MatRepSym< double, MaxHits2 > > SMatrixSym12
 
typedef ROOT::Math::SMatrix< double, 2, 2, ROOT::Math::MatRepSym< double, 2 > > SMatrixSym2
 
typedef ROOT::Math::SMatrix< double, 4, 4, ROOT::Math::MatRepSym< double, 4 > > SMatrixSym4
 
typedef ROOT::Math::SVector< double, 4 > SVector4
 

Public Member Functions

double chi2 (void) const
 
AlgebraicSymMatrix covarianceMatrix (void)
 
bool fit (void)
 
bool fitdone () const
 
MuonRecHitContainer hits (void) const
 
LocalPoint intercept () const
 
LocalVector localdir () const
 
 MuonSegFit (MuonRecHitContainer hits)
 
int ndof (void) const
 
size_t nhits (void) const
 
float Rdev (float x, float y, float z) const
 
double scaleXError (void) const
 
void setScaleXError (double factor)
 
float xdev (float x, float z) const
 
float xfit (float z) const
 
float ydev (float y, float z) const
 
float yfit (float z) const
 
virtual ~MuonSegFit ()
 

Static Public Attributes

static const int MaxHits2 = 22
 

Protected Member Functions

SMatrix12by4 derivativeMatrix (void)
 
AlgebraicSymMatrix flipErrors (const SMatrixSym4 &)
 
void setOutFromIP (void)
 
SMatrixSym12 weightMatrix (void)
 

Protected Attributes

double chi2_
 
bool fitdone_
 
MuonRecHitContainer hits_
 
LocalPoint intercept_
 
LocalVector localdir_
 
int ndof_
 
double scaleXError_
 
float uslope_
 
float vslope_
 

Private Member Functions

void fit2 (void)
 
void fitlsq (void)
 
void setChi2 (void)
 

Detailed Description

Definition at line 34 of file MuonSegFit.h.

Member Typedef Documentation

◆ MuonRecHitContainer

Definition at line 39 of file MuonSegFit.h.

◆ MuonRecHitPtr

typedef std::shared_ptr<TrackingRecHit> MuonSegFit::MuonRecHitPtr

Definition at line 38 of file MuonSegFit.h.

◆ SMatrix12by4

typedef ROOT::Math::SMatrix<double, MaxHits2, 4> MuonSegFit::SMatrix12by4

Definition at line 46 of file MuonSegFit.h.

◆ SMatrix4

typedef ROOT::Math::SMatrix<double, 4> MuonSegFit::SMatrix4

Definition at line 49 of file MuonSegFit.h.

◆ SMatrixSym12

typedef ROOT::Math::SMatrix<double, MaxHits2, MaxHits2, ROOT::Math::MatRepSym<double, MaxHits2> > MuonSegFit::SMatrixSym12

Definition at line 43 of file MuonSegFit.h.

◆ SMatrixSym2

typedef ROOT::Math::SMatrix<double, 2, 2, ROOT::Math::MatRepSym<double, 2> > MuonSegFit::SMatrixSym2

Definition at line 53 of file MuonSegFit.h.

◆ SMatrixSym4

typedef ROOT::Math::SMatrix<double, 4, 4, ROOT::Math::MatRepSym<double, 4> > MuonSegFit::SMatrixSym4

Definition at line 50 of file MuonSegFit.h.

◆ SVector4

typedef ROOT::Math::SVector<double, 4> MuonSegFit::SVector4

Definition at line 56 of file MuonSegFit.h.

Constructor & Destructor Documentation

◆ MuonSegFit()

MuonSegFit::MuonSegFit ( MuonRecHitContainer  hits)
inline

Definition at line 62 of file MuonSegFit.h.

63  : hits_(hits), uslope_(.0), vslope_(.0), chi2_(.0), ndof_(0), scaleXError_(1.0), fitdone_(false) {}
double chi2_
Definition: MuonSegFit.h:119
float uslope_
Definition: MuonSegFit.h:115
float vslope_
Definition: MuonSegFit.h:116
MuonRecHitContainer hits(void) const
Definition: MuonSegFit.h:86
double scaleXError_
Definition: MuonSegFit.h:121
bool fitdone_
Definition: MuonSegFit.h:122
MuonRecHitContainer hits_
Definition: MuonSegFit.h:114

◆ ~MuonSegFit()

virtual MuonSegFit::~MuonSegFit ( )
inlinevirtual

Definition at line 65 of file MuonSegFit.h.

65 {}

Member Function Documentation

◆ chi2()

double MuonSegFit::chi2 ( void  ) const
inline

Definition at line 89 of file MuonSegFit.h.

References chi2_.

89 { return chi2_; }
double chi2_
Definition: MuonSegFit.h:119

◆ covarianceMatrix()

AlgebraicSymMatrix MuonSegFit::covarianceMatrix ( void  )

Definition at line 359 of file MuonSegFit.cc.

References A, derivativeMatrix(), flipErrors(), convertSQLiteXML::ok, mps_fire::result, weightMatrix(), and hltDeepSecondaryVertexTagInfosPFPuppi_cfi::weights.

359  {
362 #ifdef EDM_ML_DEBUG
363  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::covarianceMatrix] weights matrix W: \n" << weights;
364  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::covarianceMatrix] derivatives matrix A: \n" << A;
365 #endif
366 
367  // (AT W A)^-1
368  // e.g. See http://www.phys.ufl.edu/~avery/fitting.html, part I
369 
370  bool ok;
371  SMatrixSym4 result = ROOT::Math::SimilarityT(A, weights);
372 #ifdef EDM_ML_DEBUG
373  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::covarianceMatrix] (AT W A): \n" << result;
374 #endif
375 
376  ok = result.Invert(); // inverts in place
377  if (!ok) {
378  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::calculateError] Failed to invert matrix: \n" << result;
379  // return ok; //@@ SHOULD PASS THIS BACK TO CALLER?
380  }
381 #ifdef EDM_ML_DEBUG
382  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::covarianceMatrix] (AT W A)^-1: \n" << result;
383 #endif
384 
385  // reorder components to match TrackingRecHit interface (GEMSegment isa TrackingRecHit)
386  // i.e. slopes first, then positions
388 
389  return flipped;
390 }
Log< level::Info, true > LogVerbatim
SMatrixSym12 weightMatrix(void)
Definition: MuonSegFit.cc:295
SMatrix12by4 derivativeMatrix(void)
Definition: MuonSegFit.cc:320
ROOT::Math::SMatrix< double, MaxHits2, 4 > SMatrix12by4
Definition: MuonSegFit.h:46
ROOT::Math::SMatrix< double, 4, 4, ROOT::Math::MatRepSym< double, 4 > > SMatrixSym4
Definition: MuonSegFit.h:50
AlgebraicSymMatrix flipErrors(const SMatrixSym4 &)
Definition: MuonSegFit.cc:392
ROOT::Math::SMatrix< double, MaxHits2, MaxHits2, ROOT::Math::MatRepSym< double, MaxHits2 > > SMatrixSym12
Definition: MuonSegFit.h:43
CLHEP::HepSymMatrix AlgebraicSymMatrix
Definition: APVGainStruct.h:7

◆ derivativeMatrix()

MuonSegFit::SMatrix12by4 MuonSegFit::derivativeMatrix ( void  )
protected

Definition at line 320 of file MuonSegFit.cc.

References hits_, ALPAKA_ACCELERATOR_NAMESPACE::vertexFinder::it, makeMuonMisalignmentScenario::matrix, z, and PV3DBase< T, PVType, FrameType >::z().

Referenced by covarianceMatrix().

320  {
321  SMatrix12by4 matrix; // 12x4, init to 0
322  int row = 0;
323 
324  for (MuonRecHitContainer::const_iterator it = hits_.begin(); it != hits_.end(); ++it) {
325  LocalPoint lp = (*it)->localPosition();
326 
327  float z = lp.z();
328 
329  matrix(row, 0) = 1.;
330  matrix(row, 2) = z;
331  ++row;
332  matrix(row, 1) = 1.;
333  matrix(row, 3) = z;
334  ++row;
335  }
336  return matrix;
337 }
T z() const
Definition: PV3DBase.h:61
ROOT::Math::SMatrix< double, MaxHits2, 4 > SMatrix12by4
Definition: MuonSegFit.h:46
MuonRecHitContainer hits_
Definition: MuonSegFit.h:114

◆ fit()

bool MuonSegFit::fit ( void  )

Definition at line 13 of file MuonSegFit.cc.

References fit2(), fitdone(), fitdone_, fitlsq(), MaxHits2, dqmiodumpmetadata::n, and nhits().

Referenced by trackingPlots.Iteration::modules().

13  {
14  if (fitdone())
15  return fitdone_; // don't redo fit unnecessarily
16  short n = nhits();
17  if (n < 2) {
18  edm::LogVerbatim("MuonSegFit") << "[MuonSegFit::fit] - cannot fit just 1 hit!!";
19  } else if (n == 2) {
20  fit2();
21  } else if (2 * n <= MaxHits2) {
22  fitlsq();
23  } else {
24  edm::LogVerbatim("MuonSegFit") << "[MuonSegFit::fit] - cannot fit more than " << MaxHits2 / 2 << " hits!!";
25  }
26  return fitdone_;
27 }
static const int MaxHits2
Definition: MuonSegFit.h:41
Log< level::Info, true > LogVerbatim
size_t nhits(void) const
Definition: MuonSegFit.h:88
void fitlsq(void)
Definition: MuonSegFit.cc:76
bool fitdone() const
Definition: MuonSegFit.h:93
void fit2(void)
Definition: MuonSegFit.cc:29
bool fitdone_
Definition: MuonSegFit.h:122

◆ fit2()

void MuonSegFit::fit2 ( void  )
private

Definition at line 29 of file MuonSegFit.cc.

References chi2_, PVValHelper::dz, fitdone_, hits_, intercept_, ndof_, setOutFromIP(), uslope_, vslope_, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by fit().

29  {
30  // Just join the two points
31  // Equation of straight line between (x1, y1) and (x2, y2) in xy-plane is
32  // y = mx + c
33  // with m = (y2-y1)/(x2-x1)
34  // and c = (y1*x2-x2*y1)/(x2-x1)
35  //
36  // Now we will make two straight lines
37  // one in xz-plane, another in yz-plane
38  // x = uz + c1
39  // y = vz + c2
40  // 1) Check whether hits are on the same layer
41  // should be done before, so removed
42  // -------------------------------------------
43 
44  MuonRecHitContainer::const_iterator ih = hits_.begin();
45  // 2) Global Positions of hit 1 and 2 and
46  // Local Positions of hit 1 and 2 w.r.t. reference GEM Eta Partition
47  // ---------------------------------------------------------------------
48  LocalPoint h1pos = (*ih)->localPosition();
49  ++ih;
50  LocalPoint h2pos = (*ih)->localPosition();
51 
52  // 3) Now make straight line between the two points in local coords
53  // ----------------------------------------------------------------
54  float dz = h2pos.z() - h1pos.z();
55  if (dz != 0.0) {
56  uslope_ = (h2pos.x() - h1pos.x()) / dz;
57  vslope_ = (h2pos.y() - h1pos.y()) / dz;
58  }
59 
60  float uintercept = (h1pos.x() * h2pos.z() - h2pos.x() * h1pos.z()) / dz;
61  float vintercept = (h1pos.y() * h2pos.z() - h2pos.y() * h1pos.z()) / dz;
62 
63  // calculate local position (variable: intercept_)
64  intercept_ = LocalPoint(uintercept, vintercept, 0.);
65 
66  // calculate the local direction (variable: localdir_)
67  setOutFromIP();
68 
69  //@@ NOT SURE WHAT IS SENSIBLE FOR THESE...
70  chi2_ = 0.;
71  ndof_ = 0;
72 
73  fitdone_ = true;
74 }
Point3DBase< Scalar, LocalTag > LocalPoint
Definition: Definitions.h:30
double chi2_
Definition: MuonSegFit.h:119
LocalPoint intercept_
Definition: MuonSegFit.h:117
T z() const
Definition: PV3DBase.h:61
float uslope_
Definition: MuonSegFit.h:115
T x() const
Definition: PV3DBase.h:59
T y() const
Definition: PV3DBase.h:60
float vslope_
Definition: MuonSegFit.h:116
void setOutFromIP(void)
Definition: MuonSegFit.cc:339
bool fitdone_
Definition: MuonSegFit.h:122
MuonRecHitContainer hits_
Definition: MuonSegFit.h:114

◆ fitdone()

bool MuonSegFit::fitdone ( ) const
inline

Definition at line 93 of file MuonSegFit.h.

References fitdone_.

Referenced by fit().

93 { return fitdone_; }
bool fitdone_
Definition: MuonSegFit.h:122

◆ fitlsq()

void MuonSegFit::fitlsq ( void  )
private

Definition at line 76 of file MuonSegFit.cc.

References B, fitdone_, hits_, intercept_, LogTrace, convertSQLiteXML::ok, AlCaHLTBitMon_ParallelJobs::p, setChi2(), setOutFromIP(), uslope_, findQualityFiles::v, vslope_, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), z, and PV3DBase< T, PVType, FrameType >::z().

Referenced by fit().

76  {
77  // Linear least-squares fit to up to 6 GEM rechits, one per layer in a GEM chamber.
78  // Comments adapted from Tim Cox' comments in the original GEMSegAlgoSK algorithm.
79 
80  // Fit to the local x, y rechit coordinates in z projection
81  // The strip measurement controls the precision of x
82  // The wire measurement controls the precision of y.
83  // Typical precision: u (strip, sigma~200um), v (wire, sigma~1cm)
84 
85  // Set up the normal equations for the least-squares fit as a matrix equation
86 
87  // We have a vector of measurements m, which is a 2n x 1 dim matrix
88  // The transpose mT is (u1, v1, u2, v2, ..., un, vn) where
89  // ui is the strip-associated measurement and
90  // vi is the wire-associated measurement
91  // for a given rechit i.
92 
93  // The fit is to
94  // u = u0 + uz * z
95  // v = v0 + vz * z
96  // where u0, uz, v0, vz are the parameters to be obtained from the fit.
97 
98  // These are contained in a vector p which is a 4x1 dim matrix, and
99  // its transpose pT is (u0, v0, uz, vz). Note the ordering!
100 
101  // The covariance matrix for each pair of measurements is 2 x 2 and
102  // the inverse of this is the error matrix E.
103  // The error matrix for the whole set of n measurements is a diagonal
104  // matrix with diagonal elements the individual 2 x 2 error matrices
105  // (because the inverse of a diagonal matrix is a diagonal matrix
106  // with each element the inverse of the original.)
107 
108  // In function 'weightMatrix()', the variable 'matrix' is filled with this
109  // block-diagonal overall covariance matrix. Then 'matrix' is inverted to the
110  // block-diagonal error matrix, and returned.
111 
112  // Define the matrix A as
113  // 1 0 z1 0
114  // 0 1 0 z1
115  // 1 0 z2 0
116  // 0 1 0 z2
117  // .. .. .. ..
118  // 1 0 zn 0
119  // 0 1 0 zn
120 
121  // This matrix A is set up and returned by function 'derivativeMatrix()'.
122 
123  // Then the normal equations are described by the matrix equation
124  //
125  // (AT E A)p = (AT E)m
126  //
127  // where AT is the transpose of A.
128 
129  // Call the combined matrix on the LHS, M, and that on the RHS, B:
130  // M p = B
131 
132  // We solve this for the parameter vector, p.
133  // The elements of M and B then involve sums over the hits
134 
135  // The covariance matrix of the parameters is obtained by
136  // (AT E A)^-1 calculated in 'covarianceMatrix()'.
137 
138  SMatrix4 M; // 4x4, init to 0
139  SVector4 B; // 4x1, init to 0;
140 
141  MuonRecHitContainer::const_iterator ih = hits_.begin();
142 
143  // Loop over the GEMRecHits and make small (2x2) matrices used to fill the blockdiagonal covariance matrix E^-1
144  for (ih = hits_.begin(); ih != hits_.end(); ++ih) {
145  LocalPoint lp = (*ih)->localPosition();
146  // Local position of hit w.r.t. chamber
147  double u = lp.x();
148  double v = lp.y();
149  double z = lp.z();
150 
151  // Covariance matrix of local errors
152  SMatrixSym2 IC; // 2x2, init to 0
153 
154  IC(0, 0) = (*ih)->localPositionError().xx();
155  IC(1, 1) = (*ih)->localPositionError().yy();
156  //@@ NOT SURE WHICH OFF-DIAGONAL ELEMENT MUST BE DEFINED BUT (1,0) WORKS
157  //@@ (and SMatrix enforces symmetry)
158  IC(1, 0) = (*ih)->localPositionError().xy();
159  // IC(0,1) = IC(1,0);
160 
161 #ifdef EDM_ML_DEBUG
162  edm::LogVerbatim("MuonSegFitMatrixDetails")
163  << "[MuonSegFit::fit] 2x2 covariance matrix for this GEMRecHit :: [[" << IC(0, 0) << ", " << IC(0, 1) << "]["
164  << IC(1, 0) << "," << IC(1, 1) << "]]";
165 #endif
166 
167  // Invert covariance matrix (and trap if it fails!)
168  bool ok = IC.Invert();
169  if (!ok) {
170  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::fit] Failed to invert covariance matrix: \n" << IC;
171  // return ok; //@@ SHOULD PASS THIS BACK TO CALLER?
172  }
173 
174  M(0, 0) += IC(0, 0);
175  M(0, 1) += IC(0, 1);
176  M(0, 2) += IC(0, 0) * z;
177  M(0, 3) += IC(0, 1) * z;
178  B(0) += u * IC(0, 0) + v * IC(0, 1);
179 
180  M(1, 0) += IC(1, 0);
181  M(1, 1) += IC(1, 1);
182  M(1, 2) += IC(1, 0) * z;
183  M(1, 3) += IC(1, 1) * z;
184  B(1) += u * IC(1, 0) + v * IC(1, 1);
185 
186  M(2, 0) += IC(0, 0) * z;
187  M(2, 1) += IC(0, 1) * z;
188  M(2, 2) += IC(0, 0) * z * z;
189  M(2, 3) += IC(0, 1) * z * z;
190  B(2) += (u * IC(0, 0) + v * IC(0, 1)) * z;
191 
192  M(3, 0) += IC(1, 0) * z;
193  M(3, 1) += IC(1, 1) * z;
194  M(3, 2) += IC(1, 0) * z * z;
195  M(3, 3) += IC(1, 1) * z * z;
196  B(3) += (u * IC(1, 0) + v * IC(1, 1)) * z;
197  }
198 
199  SVector4 p;
200  bool ok = M.Invert();
201  if (!ok) {
202  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::fit] Failed to invert matrix: \n" << M;
203  // return ok; //@@ SHOULD PASS THIS BACK TO CALLER?
204  } else {
205  p = M * B;
206  }
207 
208 #ifdef EDM_ML_DEBUG
209  LogTrace("MuonSegFitMatrixDetails") << "[MuonSegFit::fit] p = " << p(0) << ", " << p(1) << ", " << p(2) << ", "
210  << p(3);
211 #endif
212 
213  // fill member variables (note origin has local z = 0)
214  // intercept_
215  intercept_ = LocalPoint(p(0), p(1), 0.);
216 
217  // localdir_ - set so segment points outwards from IP
218  uslope_ = p(2);
219  vslope_ = p(3);
220  setOutFromIP();
221 
222  // calculate chi2 of fit
223  setChi2();
224 
225  // flag fit has been done
226  fitdone_ = true;
227 }
Log< level::Info, true > LogVerbatim
Point3DBase< Scalar, LocalTag > LocalPoint
Definition: Definitions.h:30
Definition: APVGainStruct.h:7
LocalPoint intercept_
Definition: MuonSegFit.h:117
T z() const
Definition: PV3DBase.h:61
float uslope_
Definition: MuonSegFit.h:115
ROOT::Math::SMatrix< double, 2, 2, ROOT::Math::MatRepSym< double, 2 > > SMatrixSym2
Definition: MuonSegFit.h:53
#define LogTrace(id)
ROOT::Math::SMatrix< double, 4 > SMatrix4
Definition: MuonSegFit.h:49
T x() const
Definition: PV3DBase.h:59
T y() const
Definition: PV3DBase.h:60
float vslope_
Definition: MuonSegFit.h:116
void setChi2(void)
Definition: MuonSegFit.cc:229
ROOT::Math::SVector< double, 4 > SVector4
Definition: MuonSegFit.h:56
void setOutFromIP(void)
Definition: MuonSegFit.cc:339
bool fitdone_
Definition: MuonSegFit.h:122
MuonRecHitContainer hits_
Definition: MuonSegFit.h:114

◆ flipErrors()

AlgebraicSymMatrix MuonSegFit::flipErrors ( const SMatrixSym4 a)
protected

Definition at line 392 of file MuonSegFit.cc.

References a, mps_fire::i, and dqmiolumiharvest::j.

Referenced by covarianceMatrix().

392  {
393  // The GEMSegment needs the error matrix re-arranged to match
394  // parameters in order (uz, vz, u0, v0)
395  // where uz, vz = slopes, u0, v0 = intercepts
396 
397 #ifdef EDM_ML_DEBUG
398  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::flipErrors] input: \n" << a;
399 #endif
400 
401  AlgebraicSymMatrix hold(4, 0.);
402 
403  for (short j = 0; j != 4; ++j) {
404  for (short i = 0; i != 4; ++i) {
405  hold(i + 1, j + 1) = a(i, j); // SMatrix counts from 0, AlgebraicMatrix from 1
406  }
407  }
408 
409 #ifdef EDM_ML_DEBUG
410  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::flipErrors] after copy:";
411  edm::LogVerbatim("MuonSegFitMatrixDetails")
412  << "(" << hold(1, 1) << " " << hold(1, 2) << " " << hold(1, 3) << " " << hold(1, 4);
413  edm::LogVerbatim("MuonSegFitMatrixDetails")
414  << " " << hold(2, 1) << " " << hold(2, 2) << " " << hold(2, 3) << " " << hold(2, 4);
415  edm::LogVerbatim("MuonSegFitMatrixDetails")
416  << " " << hold(3, 1) << " " << hold(3, 2) << " " << hold(3, 3) << " " << hold(3, 4);
417  edm::LogVerbatim("MuonSegFitMatrixDetails")
418  << " " << hold(4, 1) << " " << hold(4, 2) << " " << hold(4, 3) << " " << hold(4, 4) << ")";
419 #endif
420 
421  // errors on slopes into upper left
422  hold(1, 1) = a(2, 2);
423  hold(1, 2) = a(2, 3);
424  hold(2, 1) = a(3, 2);
425  hold(2, 2) = a(3, 3);
426 
427  // errors on positions into lower right
428  hold(3, 3) = a(0, 0);
429  hold(3, 4) = a(0, 1);
430  hold(4, 3) = a(1, 0);
431  hold(4, 4) = a(1, 1);
432 
433  // must also interchange off-diagonal elements of off-diagonal 2x2 submatrices
434  hold(4, 1) = a(1, 2);
435  hold(3, 2) = a(0, 3);
436  hold(2, 3) = a(3, 0); // = a(0,3)
437  hold(1, 4) = a(2, 1); // = a(1,2)
438 
439 #ifdef EDM_ML_DEBUG
440  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::flipErrors] after flip:";
441  edm::LogVerbatim("MuonSegFitMatrixDetails")
442  << "(" << hold(1, 1) << " " << hold(1, 2) << " " << hold(1, 3) << " " << hold(1, 4);
443  edm::LogVerbatim("MuonSegFitMatrixDetails")
444  << " " << hold(2, 1) << " " << hold(2, 2) << " " << hold(2, 3) << " " << hold(2, 4);
445  edm::LogVerbatim("MuonSegFitMatrixDetails")
446  << " " << hold(3, 1) << " " << hold(3, 2) << " " << hold(3, 3) << " " << hold(3, 4);
447  edm::LogVerbatim("MuonSegFitMatrixDetails")
448  << " " << hold(4, 1) << " " << hold(4, 2) << " " << hold(4, 3) << " " << hold(4, 4) << ")";
449 #endif
450 
451  return hold;
452 }
Log< level::Info, true > LogVerbatim
double a
Definition: hdecay.h:121
CLHEP::HepSymMatrix AlgebraicSymMatrix

◆ hits()

MuonRecHitContainer MuonSegFit::hits ( void  ) const
inline

Definition at line 86 of file MuonSegFit.h.

References hits_.

86 { return hits_; }
MuonRecHitContainer hits_
Definition: MuonSegFit.h:114

◆ intercept()

LocalPoint MuonSegFit::intercept ( ) const
inline

Definition at line 91 of file MuonSegFit.h.

References intercept_.

91 { return intercept_; }
LocalPoint intercept_
Definition: MuonSegFit.h:117

◆ localdir()

LocalVector MuonSegFit::localdir ( ) const
inline

Definition at line 92 of file MuonSegFit.h.

References localdir_.

92 { return localdir_; }
LocalVector localdir_
Definition: MuonSegFit.h:118

◆ ndof()

int MuonSegFit::ndof ( void  ) const
inline

Definition at line 90 of file MuonSegFit.h.

References ndof_.

90 { return ndof_; }

◆ nhits()

size_t MuonSegFit::nhits ( void  ) const
inline

Definition at line 88 of file MuonSegFit.h.

References hits_.

Referenced by fit().

88 { return hits_.size(); }
MuonRecHitContainer hits_
Definition: MuonSegFit.h:114

◆ Rdev()

float MuonSegFit::Rdev ( float  x,
float  y,
float  z 
) const

Definition at line 466 of file MuonSegFit.cc.

References mathSSE::sqrt(), x, xdev(), y, ydev(), and z.

466  {
467  return sqrt(xdev(x, z) * xdev(x, z) + ydev(y, z) * ydev(y, z));
468 }
T sqrt(T t)
Definition: SSEVec.h:23
float xdev(float x, float z) const
Definition: MuonSegFit.cc:462
float ydev(float y, float z) const
Definition: MuonSegFit.cc:464

◆ scaleXError()

double MuonSegFit::scaleXError ( void  ) const
inline

Definition at line 87 of file MuonSegFit.h.

References scaleXError_.

Referenced by weightMatrix().

87 { return scaleXError_; }
double scaleXError_
Definition: MuonSegFit.h:121

◆ setChi2()

void MuonSegFit::setChi2 ( void  )
private

Definition at line 229 of file MuonSegFit.cc.

References chi2_, hits_, intercept_, localdir_, ndof_, convertSQLiteXML::ok, uslope_, findQualityFiles::v, vslope_, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), z, and PV3DBase< T, PVType, FrameType >::z().

Referenced by fitlsq().

229  {
230  double chsq = 0.;
231 
232  MuonRecHitContainer::const_iterator ih;
233 
234  for (ih = hits_.begin(); ih != hits_.end(); ++ih) {
235  LocalPoint lp = (*ih)->localPosition();
236 
237  double u = lp.x();
238  double v = lp.y();
239  double z = lp.z();
240 
241  double du = intercept_.x() + uslope_ * z - u;
242  double dv = intercept_.y() + vslope_ * z - v;
243 
244  SMatrixSym2 IC; // 2x2, init to 0
245 
246  IC(0, 0) = (*ih)->localPositionError().xx();
247  // IC(0,1) = (*ih)->localPositionError().xy();
248  IC(1, 0) = (*ih)->localPositionError().xy();
249  IC(1, 1) = (*ih)->localPositionError().yy();
250  // IC(1,0) = IC(0,1);
251 
252 #ifdef EDM_ML_DEBUG
253  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::setChi2] IC before = \n" << IC;
254 #endif
255 
256  // Invert covariance matrix
257  bool ok = IC.Invert();
258  if (!ok) {
259  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::setChi2] Failed to invert covariance matrix: \n"
260  << IC;
261  // return ok;
262  }
263  chsq += du * du * IC(0, 0) + 2. * du * dv * IC(0, 1) + dv * dv * IC(1, 1);
264 
265 #ifdef EDM_ML_DEBUG
266  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::setChi2] IC after = \n" << IC;
267  edm::LogVerbatim("MuonSegFit")
268  << "[GEM RecHit ] Contribution to Chi^2 of this hit :: du^2*Cov(0,0) + 2*du*dv*Cov(0,1) + dv^2*IC(1,1) = "
269  << du * du << "*" << IC(0, 0) << " + 2.*" << du << "*" << dv << "*" << IC(0, 1) << " + " << dv * dv << "*"
270  << IC(1, 1) << " = " << chsq;
271 #endif
272  }
273 
274  // fill member variables
275  chi2_ = chsq;
276  ndof_ = 2. * hits_.size() - 4;
277 
278 #ifdef EDM_ML_DEBUG
279  edm::LogVerbatim("MuonSegFit") << "[MuonSegFit::setChi2] chi2/ndof = " << chi2_ << "/" << ndof_;
280  edm::LogVerbatim("MuonSegFit") << "-----------------------------------------------------";
281 
282  // check fit quality ... maybe write a separate function later on
283  // that is there only for debugging issues
284 
285  edm::LogVerbatim("MuonSegFit") << "[GEM Segment with Local Direction = " << localdir_ << " and Local Position "
286  << intercept_ << "] can be written as:";
287  edm::LogVerbatim("MuonSegFit") << "[ x ] = " << localdir_.x() << " * t + " << intercept_.x();
288  edm::LogVerbatim("MuonSegFit") << "[ y ] = " << localdir_.y() << " * t + " << intercept_.y();
289  edm::LogVerbatim("MuonSegFit") << "[ z ] = " << localdir_.z() << " * t + " << intercept_.z();
290  edm::LogVerbatim("MuonSegFit")
291  << "Now extrapolate to each of the GEMRecHits XY plane (so constant z = RH LP.z()) to obtain [x1,y1]";
292 #endif
293 }
Log< level::Info, true > LogVerbatim
LocalVector localdir_
Definition: MuonSegFit.h:118
double chi2_
Definition: MuonSegFit.h:119
LocalPoint intercept_
Definition: MuonSegFit.h:117
T z() const
Definition: PV3DBase.h:61
float uslope_
Definition: MuonSegFit.h:115
ROOT::Math::SMatrix< double, 2, 2, ROOT::Math::MatRepSym< double, 2 > > SMatrixSym2
Definition: MuonSegFit.h:53
T x() const
Definition: PV3DBase.h:59
T y() const
Definition: PV3DBase.h:60
float vslope_
Definition: MuonSegFit.h:116
MuonRecHitContainer hits_
Definition: MuonSegFit.h:114

◆ setOutFromIP()

void MuonSegFit::setOutFromIP ( void  )
protected

Definition at line 339 of file MuonSegFit.cc.

References PVValHelper::dx, dxdz, PVValHelper::dy, dydz, PVValHelper::dz, localdir_, PV3DBase< T, PVType, FrameType >::phi(), mathSSE::sqrt(), unit(), uslope_, and vslope_.

Referenced by fit2(), and fitlsq().

339  {
340  // Set direction of segment to point from IP outwards
341  // (Incorrect for particles not coming from IP, of course.)
342 
343  double dxdz = uslope_;
344  double dydz = vslope_;
345  double dz = 1. / sqrt(1. + dxdz * dxdz + dydz * dydz);
346  double dx = dz * dxdz;
347  double dy = dz * dydz;
348  LocalVector localDir(dx, dy, dz);
349 
350  localdir_ = (localDir).unit();
351 #ifdef EDM_ML_DEBUG
352  edm::LogVerbatim("MuonSegFit") << "[MuonSegFit::setOutFromIP] :: dxdz = uslope_ = " << std::setw(9) << uslope_
353  << " dydz = vslope_ = " << std::setw(9) << vslope_ << " local dir = " << localDir;
354  edm::LogVerbatim("MuonSegFit") << "[MuonSegFit::setOutFromIP] :: ==> local dir = " << localdir_
355  << " localdir.phi = " << localdir_.phi();
356 #endif
357 }
Log< level::Info, true > LogVerbatim
float dydz
LocalVector localdir_
Definition: MuonSegFit.h:118
float dxdz
float uslope_
Definition: MuonSegFit.h:115
Geom::Phi< T > phi() const
Definition: PV3DBase.h:66
float vslope_
Definition: MuonSegFit.h:116
T sqrt(T t)
Definition: SSEVec.h:23
Basic3DVector unit() const

◆ setScaleXError()

void MuonSegFit::setScaleXError ( double  factor)
inline

Definition at line 74 of file MuonSegFit.h.

References scaleXError_.

◆ weightMatrix()

MuonSegFit::SMatrixSym12 MuonSegFit::weightMatrix ( void  )
protected

Definition at line 295 of file MuonSegFit.cc.

References hits_, ALPAKA_ACCELERATOR_NAMESPACE::vertexFinder::it, makeMuonMisalignmentScenario::matrix, convertSQLiteXML::ok, and scaleXError().

Referenced by covarianceMatrix().

295  {
296  bool ok = true;
297 
298  SMatrixSym12 matrix = ROOT::Math::SMatrixIdentity(); // 12x12, init to 1's on diag
299 
300  int row = 0;
301 
302  for (MuonRecHitContainer::const_iterator it = hits_.begin(); it != hits_.end(); ++it) {
303  // Note scaleXError allows rescaling the x error if necessary
304  matrix(row, row) = scaleXError() * (*it)->localPositionError().xx();
305  matrix(row, row + 1) = (*it)->localPositionError().xy();
306  ++row;
307  matrix(row, row - 1) = (*it)->localPositionError().xy();
308  matrix(row, row) = (*it)->localPositionError().yy();
309  ++row;
310  }
311 
312  ok = matrix.Invert(); // invert in place
313  if (!ok) {
314  edm::LogVerbatim("MuonSegFitMatrixDetails") << "[MuonSegFit::weightMatrix] Failed to invert matrix: \n" << matrix;
315  // return ok; //@@ SHOULD PASS THIS BACK TO CALLER?
316  }
317  return matrix;
318 }
Log< level::Info, true > LogVerbatim
double scaleXError(void) const
Definition: MuonSegFit.h:87
MuonRecHitContainer hits_
Definition: MuonSegFit.h:114
ROOT::Math::SMatrix< double, MaxHits2, MaxHits2, ROOT::Math::MatRepSym< double, MaxHits2 > > SMatrixSym12
Definition: MuonSegFit.h:43

◆ xdev()

float MuonSegFit::xdev ( float  x,
float  z 
) const

Definition at line 462 of file MuonSegFit.cc.

References intercept_, uslope_, x, PV3DBase< T, PVType, FrameType >::x(), and z.

Referenced by Rdev().

462 { return intercept_.x() + uslope_ * z - x; }
LocalPoint intercept_
Definition: MuonSegFit.h:117
float uslope_
Definition: MuonSegFit.h:115
T x() const
Definition: PV3DBase.h:59

◆ xfit()

float MuonSegFit::xfit ( float  z) const

Definition at line 454 of file MuonSegFit.cc.

References intercept_, uslope_, PV3DBase< T, PVType, FrameType >::x(), and z.

454  {
455  //@@ ADD THIS TO EACH ACCESSOR OF FIT RESULTS?
456  // if ( !fitdone() ) fit();
457  return intercept_.x() + uslope_ * z;
458 }
LocalPoint intercept_
Definition: MuonSegFit.h:117
float uslope_
Definition: MuonSegFit.h:115
T x() const
Definition: PV3DBase.h:59

◆ ydev()

float MuonSegFit::ydev ( float  y,
float  z 
) const

Definition at line 464 of file MuonSegFit.cc.

References intercept_, vslope_, y, PV3DBase< T, PVType, FrameType >::y(), and z.

Referenced by Rdev().

464 { return intercept_.y() + vslope_ * z - y; }
LocalPoint intercept_
Definition: MuonSegFit.h:117
T y() const
Definition: PV3DBase.h:60
float vslope_
Definition: MuonSegFit.h:116

◆ yfit()

float MuonSegFit::yfit ( float  z) const

Definition at line 460 of file MuonSegFit.cc.

References intercept_, vslope_, PV3DBase< T, PVType, FrameType >::y(), and z.

460 { return intercept_.y() + vslope_ * z; }
LocalPoint intercept_
Definition: MuonSegFit.h:117
T y() const
Definition: PV3DBase.h:60
float vslope_
Definition: MuonSegFit.h:116

Member Data Documentation

◆ chi2_

double MuonSegFit::chi2_
protected

Definition at line 119 of file MuonSegFit.h.

Referenced by chi2(), fit2(), and setChi2().

◆ fitdone_

bool MuonSegFit::fitdone_
protected

Definition at line 122 of file MuonSegFit.h.

Referenced by fit(), fit2(), fitdone(), and fitlsq().

◆ hits_

MuonRecHitContainer MuonSegFit::hits_
protected

Definition at line 114 of file MuonSegFit.h.

Referenced by derivativeMatrix(), fit2(), fitlsq(), hits(), nhits(), setChi2(), and weightMatrix().

◆ intercept_

LocalPoint MuonSegFit::intercept_
protected

Definition at line 117 of file MuonSegFit.h.

Referenced by fit2(), fitlsq(), intercept(), setChi2(), xdev(), xfit(), ydev(), and yfit().

◆ localdir_

LocalVector MuonSegFit::localdir_
protected

Definition at line 118 of file MuonSegFit.h.

Referenced by localdir(), setChi2(), and setOutFromIP().

◆ MaxHits2

const int MuonSegFit::MaxHits2 = 22
static

Definition at line 41 of file MuonSegFit.h.

Referenced by fit().

◆ ndof_

int MuonSegFit::ndof_
protected

Definition at line 120 of file MuonSegFit.h.

Referenced by fit2(), ndof(), and setChi2().

◆ scaleXError_

double MuonSegFit::scaleXError_
protected

Definition at line 121 of file MuonSegFit.h.

Referenced by scaleXError(), and setScaleXError().

◆ uslope_

float MuonSegFit::uslope_
protected

Definition at line 115 of file MuonSegFit.h.

Referenced by fit2(), fitlsq(), setChi2(), setOutFromIP(), xdev(), and xfit().

◆ vslope_

float MuonSegFit::vslope_
protected

Definition at line 116 of file MuonSegFit.h.

Referenced by fit2(), fitlsq(), setChi2(), setOutFromIP(), ydev(), and yfit().