Belle II Software development
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
QualityEstimatorCircleFit Class Reference
Packages » tracking

Class containing the algorithm to perform the simple circle fit. More...

#include <QualityEstimatorCircleFit.h>

Inheritance diagram for QualityEstimatorCircleFit:
QualityEstimatorBase

Public Member Functions

virtual double estimateQuality (std::vector< SpacePoint const * > const &measurements) final
 Minimal implementation of the quality estimation Calculates quality indicator in range [0,1].
 
void setMagneticFieldStrength (double magneticFieldZ=1.5)
 Setter for z component of magnetic field.
 
virtual QualityEstimationResults estimateQualityAndProperties (std::vector< SpacePoint const * > const &measurements)
 Quality estimation providing additional quantities Calculates quality indicator in range [0,1] Optionally returns chi2 and additional extra information.
 

Protected Member Functions

double calcPt (double const radius) const
 Returns a value for the transverse momentum in GeV calculated from a provided radius.
 

Protected Attributes

double m_magneticFieldZ = 1.5
 Member storing the z component of the magnetic field.
 
QualityEstimationResults m_results
 Result of the quality estimation This is stored as a member variable, because some values may be calculated by 'estimateQuality' anyways.
 

Detailed Description

Class containing the algorithm to perform the simple circle fit.

Definition at line 19 of file QualityEstimatorCircleFit.h.

Member Function Documentation

◆ calcPt()

double calcPt ( double const  radius) const
inlineprotectedinherited

Returns a value for the transverse momentum in GeV calculated from a provided radius.

Utilizing m_magneticFieldZ and hard coded speed of light

Definition at line 80 of file QualityEstimatorBase.h.

80{ return m_magneticFieldZ * radius * 0.00299792458; }
Belle2::QualityEstimatorBase::m_magneticFieldZ
double m_magneticFieldZ
Member storing the z component of the magnetic field.
Definition: QualityEstimatorBase.h:84

◆ estimateQuality()

double estimateQuality ( std::vector< SpacePoint const * > const &  measurements)
finalvirtual

Minimal implementation of the quality estimation Calculates quality indicator in range [0,1].

measurements - std::vector<SpacePoint const*> ordered from innermost to outermost measurement

WARNING hardcoded values!

Implements QualityEstimatorBase.

Definition at line 19 of file QualityEstimatorCircleFit.cc.

20{
21 if (measurements.size() < 4) {
22 if (measurements.size() == 3) return 0.2; // Arbitrary value to prevent excluding measurements with 3 hits later on.
23 else return 0;
24 }
25 // Calculates Curvature: True means clockwise, False means counterclockwise.
26 // TODO this is not an optimized approach; just to get things to work.
27 // CalcCurvature could be integrated into the looping over the hits which CircleFit does anyhow.
28 m_results.curvatureSign = calcCurvatureSignum(measurements);
29 bool clockwise = *(m_results.curvatureSign) >= 0;
30
31 double stopper = 0.000000001;
32 double meanX = 0, meanY = 0, meanX2 = 0, meanY2 = 0, meanR2 = 0, meanR4 = 0, meanXR2 = 0, meanYR2 = 0, meanXY = 0; //mean values
33 double sumWeights = 0, divisor/*, weightNormalizer = 0*/; // sumWeights is sum of weights, divisor is 1/sumWeights;
34 double tuningParameter = 1.; //0.02; // this parameter is for internal tuning of the weights, 1 means no influence of parameter.
35
36 // looping over all hits and do the division afterwards
37 for (const SpacePoint* hit : measurements) {
38 double weight = 1. / (sqrt(hit->getPositionError().X() * hit->getPositionError().X() + hit->getPositionError().Y() *
39 hit->getPositionError().Y()) * tuningParameter);
40 sumWeights += weight;
41 if (std::isnan(weight) or std::isinf(weight)) {
42 B2ERROR("QualityEstimators::circleFit, chosen sigma is 'nan': " << weight << ", setting arbitrary error: "
43 << stopper << ")"); weight = stopper;
44 }
45 double x = hit->getPosition().X();
46 double y = hit->getPosition().Y();
47 double x2 = x * x;
48 double y2 = y * y;
49 double r2 = x2 + y2;
50 meanX += x * weight;
51 meanY += y * weight;
52 meanXY += x * y * weight;
53 meanX2 += x2 * weight;
54 meanY2 += y2 * weight;
55 meanXR2 += x * r2 * weight;
56 meanYR2 += y * r2 * weight;
57 meanR2 += r2 * weight;
58 meanR4 += r2 * r2 * weight;
59 }
60 divisor = 1. / sumWeights;
61 meanX *= divisor;
62 meanY *= divisor;
63 meanXY *= divisor;
64 meanY2 *= divisor;
65 meanX2 *= divisor;
66 meanXR2 *= divisor;
67 meanYR2 *= divisor;
68 meanR2 *= divisor;
69 meanR4 *= divisor;
70
71 // covariances:
72 double covXX = meanX2 - meanX * meanX;
73 double covXY = meanXY - meanX * meanY;
74 double covYY = meanY2 - meanY * meanY;
75 double covXR2 = meanXR2 - meanX * meanR2;
76 double covYR2 = meanYR2 - meanY * meanR2;
77 double covR2R2 = meanR4 - meanR2 * meanR2;
78
79 if (covR2R2 == 0) {
80 return 0; // TODO could be problematic if it is pretty near to 0
81 }
82
83 // q1, q2: helping variables, to make the code more readable
84 double q1 = covR2R2 * covXY - covXR2 * covYR2;
85 double q2 = covR2R2 * (covXX - covYY) - covXR2 * covXR2 + covYR2 * covYR2;
86
87 // physical meaning: phi value of the point of closest approach of the fitted circle to the origin
88 double pocaPhi = 0.5 * atan2(2. * q1, q2);
89
90 double sinPhi = sin(pocaPhi);
91 double cosPhi = cos(pocaPhi);
92 double kappa = (sinPhi * covXR2 - cosPhi * covYR2) / covR2R2;
93 double delta = -kappa * meanR2 + sinPhi * meanX - cosPhi * meanY;
94 double rootTerm = sqrt(1. - 4.*delta * kappa);
95 // rho = curvature in X-Y-plane = 1/radius of fitting circle, used for pT-calculation
96 double curvature = 2.*kappa / (rootTerm);
97 double pocaD = 2.*delta / (1. + rootTerm);
98
99 // Checks if the random Curvature of CircleFit corresponds to CalcCurvature and adjust the results accordingly.
100 if ((curvature < 0 && clockwise) || (curvature > 0 && !clockwise)) {
101 // this is according to eq. 23 in the paper of Karimäki
102 curvature = -curvature;
103 pocaD = -pocaD;
104 //TODO ..and swap correlation Terms V_rho_phi and V_rho_d (which are not implemented anyway)
105 }
106 double radius = 1. / curvature;
107 double absRadius = fabs(radius);
108 double pt = calcPt(absRadius);
109 m_results.pt = pt;
110 m_results.pocaD = pocaD;
111
112 // Estimating theta of the track with the theta of the innermost hit.
113 // Otherwise the track is solely a circle without theta information.
114 double innermostHitTheta = 0.;
115 double innermostHitRadiusSquared = 100000000;
116 for (const SpacePoint* hit : measurements) {
117 if (hit->getPosition().Perp2() < innermostHitRadiusSquared) {
118 innermostHitRadiusSquared = hit->getPosition().Perp2();
119 innermostHitTheta = hit->getPosition().Theta();
120 }
121 }
122 // Account for precision when checking equality of innermostHitTheta with 0
123 double pz = (innermostHitTheta <= 1e-6 ? 0 : pt / tan(innermostHitTheta));
124 // Both p and pmag are only rough estimates based and by far not accutate!
125 m_results.p = B2Vector3D(pt * cosPhi, pt * sinPhi, pz);
126 m_results.pmag = sqrt(pt * pt + pz * pz);
127
128 double chi2 = sumWeights * (1. + curvature * pocaD) * (1. + curvature * pocaD)
129 * (sinPhi * sinPhi * covXX - 2.*sinPhi * cosPhi * covXY + cosPhi * cosPhi * covYY - kappa * kappa * covR2R2);
130
131 m_results.chiSquared = chi2;
132
133 return TMath::Prob(chi2, measurements.size() - 3);
134}
Belle2::QualityEstimatorBase::m_results
QualityEstimationResults m_results
Result of the quality estimation This is stored as a member variable, because some values may be calc...
Definition: QualityEstimatorBase.h:90
Belle2::QualityEstimatorBase::calcPt
double calcPt(double const radius) const
Returns a value for the transverse momentum in GeV calculated from a provided radius.
Definition: QualityEstimatorBase.h:80
Belle2::SpacePoint
SpacePoint typically is build from 1 PXDCluster or 1-2 SVDClusters.
Definition: SpacePoint.h:42
Belle2::B2Vector3D
B2Vector3< double > B2Vector3D
typedef for common usage with double
Definition: B2Vector3.h:516
Belle2::sqrt
double sqrt(double a)
sqrt for double
Definition: beamHelpers.h:28
Belle2::calcCurvatureSignum
short calcCurvatureSignum(std::vector< SpacePoint const * > const &measurements)
Calculate curvature based on triplets of measurements.
Definition: CalcCurvatureSignum.h:22
Belle2::QualityEstimationResults::pt
std::optional< double > pt
transverse momentum estimate from the QE
Definition: QualityEstimatorBase.h:29
Belle2::QualityEstimationResults::curvatureSign
std::optional< short > curvatureSign
direction of curvature as obtained by the QE
Definition: QualityEstimatorBase.h:27
Belle2::QualityEstimationResults::pocaD
std::optional< double > pocaD
distance to the z-axis of the POCA
Definition: QualityEstimatorBase.h:28
Belle2::QualityEstimationResults::chiSquared
std::optional< double > chiSquared
chi squared value obtained by the fit of the QE
Definition: QualityEstimatorBase.h:26
Belle2::QualityEstimationResults::p
std::optional< B2Vector3D > p
momentum vector estimate from the QE
Definition: QualityEstimatorBase.h:31
Belle2::QualityEstimationResults::pmag
std::optional< double > pmag
momentum magnitute estimate from the QE
Definition: QualityEstimatorBase.h:30

◆ estimateQualityAndProperties()

virtual QualityEstimationResults estimateQualityAndProperties ( std::vector< SpacePoint const * > const &  measurements)
inlinevirtualinherited

Quality estimation providing additional quantities Calculates quality indicator in range [0,1] Optionally returns chi2 and additional extra information.

Eg. momentum estimation.

measurements - std::vector<SpacePoint const*> ordered from innermost to outermost measurement

Reimplemented in QualityEstimatorTripletFit, and QualityEstimatorMC.

Definition at line 68 of file QualityEstimatorBase.h.

69 {
70 m_results = QualityEstimationResults();
71 m_results.qualityIndicator = estimateQuality(measurements);
72 return m_results;
73 }
Belle2::QualityEstimatorBase::estimateQuality
virtual double estimateQuality(std::vector< SpacePoint const * > const &measurements)=0
Minimal implementation of the quality estimation Calculates quality indicator in range [0,...
Belle2::QualityEstimationResults::qualityIndicator
double qualityIndicator
return value of the quality estimator
Definition: QualityEstimatorBase.h:25

◆ setMagneticFieldStrength()

void setMagneticFieldStrength ( double  magneticFieldZ = 1.5)
inlineinherited

Setter for z component of magnetic field.

Parameters
magneticFieldZ: value to set it to

Definition at line 53 of file QualityEstimatorBase.h.

53{m_magneticFieldZ = magneticFieldZ;}

Member Data Documentation

◆ m_magneticFieldZ

double m_magneticFieldZ = 1.5
protectedinherited

Member storing the z component of the magnetic field.

Definition at line 84 of file QualityEstimatorBase.h.

◆ m_results

QualityEstimationResults m_results
protectedinherited

Result of the quality estimation This is stored as a member variable, because some values may be calculated by 'estimateQuality' anyways.

Therefore they don't need to be calculated explicitly in 'estimateQualityAndProperties'.

Definition at line 90 of file QualityEstimatorBase.h.

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