Belle II Software  release-05-01-25
RiemannsMethod Class Reference

Class implementing the Riemann fit for two dimensional trajectory circle. More...

#include <RiemannsMethod.h>

Public Member Functions

 RiemannsMethod ()
 Default constructor.
 
void update (CDCTrajectory2D &trajectory2D, CDCObservations2D &observations2D) const
 Executes the fit and updates the trajectory parameters. This may render the information in the observation object.
 
bool isLineConstrained () const
 Getter for the indictor that lines should be fitted by this fitter.
 
bool isOriginConstrained () const
 Getter for the indictor that curves through the origin should be fitted by this fitter.
 
void setLineConstrained (bool constrained=true)
 Indicator if this fitter is setup to fit lines.
 
void setOriginConstrained (bool constrained=true)
 Indicator if this fitter is setup to fit curves through the origin.
 

Private Member Functions

void updateWithoutDriftLength (CDCTrajectory2D &trajectory2D, CDCObservations2D &observations2D) const
 Executes the fit without using the drift length information. More...
 
void updateWithDriftLength (CDCTrajectory2D &trajectory2D, CDCObservations2D &observations2D) const
 Executes the fit using the drift length information. More...
 

Private Attributes

bool m_lineConstrained
 Memory for the flag indicating that lines should be fitted.
 
bool m_originConstrained
 Memory for the flag indicating that curves through the origin shall be fitter.
 

Detailed Description

Class implementing the Riemann fit for two dimensional trajectory circle.

Definition at line 30 of file RiemannsMethod.h.

Member Function Documentation

◆ updateWithDriftLength()

void updateWithDriftLength ( CDCTrajectory2D trajectory2D,
CDCObservations2D observations2D 
) const
private

Executes the fit using the drift length information.

This method is used if there is no drift length information is available from the observations.

Definition at line 146 of file RiemannsMethod.cc.

147 {
148  EigenObservationMatrix eigenObservation = getEigenObservationMatrix(&observations2D);
149  assert(eigenObservation.cols() == 4);
150  size_t nObservations = observations2D.size();
151 
152  //cout << "updateWithRightLeft : " << endl;
153  //observations always have the structure
154  /*
155  observables[0][0] == x of first point
156  observables[0][1] == y of first point
157  observables[0][2] == desired distance of first point
158  */
159 
160  Eigen::Matrix< double, Eigen::Dynamic, 1 > distances = eigenObservation.col(2);
161 
162  //cout << "distances : " << endl << distances << endl;
163 
164  if ((isLineConstrained()) && (isOriginConstrained())) {
165 
166 
167  Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > projectedPoints(nObservations, 2);
168 
169  //all coordiates
170  //projectedPoints.col(0) = Eigen::Matrix<double,Eigen::Dynamic,1>::Constant(nObservations,1.0);
171  projectedPoints.col(0) = eigenObservation.col(0);
172  projectedPoints.col(1) = eigenObservation.col(1);
173  //projectedPoints.col(2) = eigenObservation.leftCols<2>().rowwise().squaredNorm() - distances.rowwise().squaredNorm();
174 
175  Eigen::Matrix<double, 2, 1> parameters =
176  projectedPoints.jacobiSvd(Eigen::ComputeThinU | Eigen::ComputeThinV).solve(distances);
177 
178  trajectory2D.setGlobalCircle(UncertainPerigeeCircle(PerigeeCircle::fromN(0.0, parameters(0), parameters(1), 0.0)));
179  }
180 
181  else if ((! isLineConstrained()) && (isOriginConstrained())) {
182 
183  Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > projectedPoints(nObservations, 3);
184 
185  //all coordiates
186  //projectedPoints.col(0) = Eigen::Matrix<double,Eigen::Dynamic,1>::Constant(1.0);
187  projectedPoints.col(0) = eigenObservation.col(0);
188  projectedPoints.col(1) = eigenObservation.col(1);
189  projectedPoints.col(2) = eigenObservation.leftCols<2>().rowwise().squaredNorm() - distances.rowwise().squaredNorm();
190 
191  Eigen::Matrix<double, 3, 1> parameters =
192  projectedPoints.jacobiSvd(Eigen::ComputeThinU | Eigen::ComputeThinV).solve(distances);
193 
194  trajectory2D.setGlobalCircle(UncertainPerigeeCircle(PerigeeCircle::fromN(0.0, parameters(0), parameters(1), parameters(2))));
195  }
196 
197  else if ((isLineConstrained()) && (! isOriginConstrained())) {
198 
199  Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > projectedPoints(nObservations, 3);
200 
201  //all coordiates
202  projectedPoints.col(0) = Eigen::Matrix<double, Eigen::Dynamic, 1>::Constant(nObservations, 1.0);
203  projectedPoints.col(1) = eigenObservation.col(0);
204  projectedPoints.col(2) = eigenObservation.col(1);
205  //projectedPoints.col(3) = eigenObservation.leftCols<2>().rowwise().squaredNorm()- distances.rowwise().squaredNorm();
206 
207  Eigen::Matrix<double, 3, 1> parameters =
208  projectedPoints.jacobiSvd(Eigen::ComputeThinU | Eigen::ComputeThinV).solve(distances);
209 
210  trajectory2D.setGlobalCircle(UncertainPerigeeCircle(PerigeeCircle::fromN(parameters(0), parameters(1), parameters(2), 0.0)));
211  //fit.setParameters(parameters(0),parameters(1),parameters(2),0.0);
212 
213  }
214 
215  else if ((! isLineConstrained()) && (! isOriginConstrained())) {
216 
217  Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > projectedPoints(nObservations, 4);
218 
219  //all coordiates
220  projectedPoints.col(0) = Eigen::Matrix<double, Eigen::Dynamic, 1>::Constant(nObservations, 1.0);
221  projectedPoints.col(1) = eigenObservation.col(0);
222  projectedPoints.col(2) = eigenObservation.col(1);
223  projectedPoints.col(3) = eigenObservation.leftCols<2>().rowwise().squaredNorm() - distances.rowwise().squaredNorm();
224 
225  //cout << "points : " << endl << projectedPoints << endl;
226 
227  Eigen::Matrix<double, 4, 1> parameters =
228  projectedPoints.jacobiSvd(Eigen::ComputeThinU | Eigen::ComputeThinV).solve(distances);
229 
230  trajectory2D.setGlobalCircle(UncertainPerigeeCircle(PerigeeCircle::fromN(parameters(0), parameters(1), parameters(2),
231  parameters(3))));
232 
233  }
234 }

◆ updateWithoutDriftLength()

void updateWithoutDriftLength ( CDCTrajectory2D trajectory2D,
CDCObservations2D observations2D 
) const
private

Executes the fit without using the drift length information.

This method is used if there is drift length information is available from the observations.

Definition at line 49 of file RiemannsMethod.cc.


The documentation for this class was generated from the following files:
Belle2::TrackFindingCDC::PerigeeCircle::fromN
static PerigeeCircle fromN(double n0, double n1, double n2, double n3=0)
Constructor with the four parameters of the generalized circle.
Belle2::TrackFindingCDC::RiemannsMethod::isOriginConstrained
bool isOriginConstrained() const
Getter for the indictor that curves through the origin should be fitted by this fitter.
Definition: RiemannsMethod.h:64
Belle2::TrackFindingCDC::RiemannsMethod::isLineConstrained
bool isLineConstrained() const
Getter for the indictor that lines should be fitted by this fitter.
Definition: RiemannsMethod.h:58