Belle II Software development
CDCSZFitter Class Reference

Class implementing the z coordinate over travel distance line fit. More...

#include <CDCSZFitter.h>

Public Member Functions

CDCTrajectorySZ fitWithStereoHits (const CDCTrack &track) const
 Returns the fitted sz trajectory of the track with the z-information of all stereo hits of the number of stereo hits is big enough.
 
CDCTrajectorySZ fit (const CDCSegment2D &stereoSegment, const CDCTrajectory2D &axialTrajectory2D) const
 Returns a fitted trajectory.
 
CDCTrajectorySZ fit (const CDCSegment3D &segment3D) const
 Fits a linear sz trajectory to the z and s coordinates in the stereo segment.
 
CDCTrajectorySZ fit (CDCSZObservations observationsSZ) const
 Fits a linear sz trajectory to the s and z coordinates given in the observations.
 
CDCTrajectorySZ fit (const CDCObservations2D &observations2D) const
 Legacy - Fits a linear sz trajectory to the x and y coordinates interpreted as sz space.
 
void update (const CDCSegmentPair &segmentPair) const
 Updates the trajectory of the axial stereo segment pair inplace.
 
void update (CDCTrajectorySZ &trajectorySZ, const CDCSegment2D &stereoSegment, const CDCTrajectory2D &axialTrajectory2D) const
 Update the given sz trajectory reconstructing the stereo segment with a near by axial segment.
 
void update (CDCTrajectorySZ &trajectorySZ, CDCSZObservations &observationsSZ) const
 Update the trajectory with a fit to the observations.
 

Static Public Member Functions

static const CDCSZFittergetFitter ()
 Getter for a standard sz line fitter instance.
 

Detailed Description

Class implementing the z coordinate over travel distance line fit.

Definition at line 27 of file CDCSZFitter.h.

Member Function Documentation

◆ fit() [1/4]

CDCTrajectorySZ fit ( CDCSZObservations  observationsSZ) const

Fits a linear sz trajectory to the s and z coordinates given in the observations.

Definition at line 156 of file CDCSZFitter.cc.

157{
158 CDCTrajectorySZ trajectorySZ;
159 update(trajectorySZ, observationsSZ);
160 return trajectorySZ;
161}
void update(const CDCSegmentPair &segmentPair) const
Updates the trajectory of the axial stereo segment pair inplace.
Definition: CDCSZFitter.cc:163
Linear trajectory in sz space.

◆ fit() [2/4]

CDCTrajectorySZ fit ( const CDCObservations2D observations2D) const

Legacy - Fits a linear sz trajectory to the x and y coordinates interpreted as sz space.

Definition at line 182 of file CDCSZFitter.cc.

183{
184 CDCSZObservations szObservations;
185 for (size_t i = 0; i < observations2D.size(); ++i) {
186 const double s = observations2D.getX(i);
187 const double z = observations2D.getY(i);
188 szObservations.fill(s, z);
189 }
190 return fit(std::move(szObservations));
191}
double getX(int iObservation) const
Getter for the x value of the observation at the given index.
double getY(int iObservation) const
Getter for the y value of the observation at the given index.
std::size_t size() const
Returns the number of observations stored.
CDCTrajectorySZ fit(const CDCSegment2D &stereoSegment, const CDCTrajectory2D &axialTrajectory2D) const
Returns a fitted trajectory.
Definition: CDCSZFitter.cc:139
Class serving as a storage of observed sz positions to present to the sz line fitters.
std::size_t fill(double s, double z, double weight=1.0)
Appends the observed position.

◆ fit() [3/4]

CDCTrajectorySZ fit ( const CDCSegment2D stereoSegment,
const CDCTrajectory2D axialTrajectory2D 
) const

Returns a fitted trajectory.

Definition at line 139 of file CDCSZFitter.cc.

141{
142 B2ASSERT("Expected stereo segment", not stereoSegment.isAxial());
143
144 CDCTrajectorySZ trajectorySZ;
145 update(trajectorySZ, stereoSegment, axialTrajectory2D);
146 return trajectorySZ;
147}
bool isAxial() const
Indicator if the underlying wires are axial.
Definition: CDCSegment.h:45

◆ fit() [4/4]

CDCTrajectorySZ fit ( const CDCSegment3D segment3D) const

Fits a linear sz trajectory to the z and s coordinates in the stereo segment.

Definition at line 149 of file CDCSZFitter.cc.

150{
151 CDCSZObservations observationsSZ;
152 observationsSZ.appendRange(segment3D);
153 return fit(std::move(observationsSZ));
154}
std::size_t appendRange(const std::vector< CDCRecoHit3D > &recoHit3Ds)
Appends all reconstructed hits from the three dimensional track.

◆ fitWithStereoHits()

CDCTrajectorySZ fitWithStereoHits ( const CDCTrack track) const

Returns the fitted sz trajectory of the track with the z-information of all stereo hits of the number of stereo hits is big enough.

Else return the basic assumption.

Definition at line 125 of file CDCSZFitter.cc.

126{
127 const bool onlyStereo = true;
128 CDCSZObservations observationsSZ(EFitVariance::c_Proper, onlyStereo);
129 observationsSZ.appendRange(track);
130 if (observationsSZ.size() > 3) {
131 CDCTrajectorySZ szTrajectory;
132 update(szTrajectory, observationsSZ);
133 return szTrajectory;
134 } else {
136 }
137}
static CDCTrajectorySZ basicAssumption()
Constructs a basic assumption, what the z0 start position and the sz slope are, including some broad ...

◆ getFitter()

const CDCSZFitter & getFitter ( )
static

Getter for a standard sz line fitter instance.

Definition at line 36 of file CDCSZFitter.cc.

37{
38 static CDCSZFitter szFitter;
39 return szFitter;
40}
Class implementing the z coordinate over travel distance line fit.
Definition: CDCSZFitter.h:27

◆ update() [1/3]

void update ( CDCTrajectorySZ trajectorySZ,
CDCSZObservations observationsSZ 
) const

Update the trajectory with a fit to the observations.

Definition at line 209 of file CDCSZFitter.cc.

210{
211 trajectorySZ.clear();
212 if (observationsSZ.size() < 3) {
213 return;
214 }
215
216 // Determine NDF : Line fit eats up 2 degrees of freedom.
217 size_t ndf = observationsSZ.size() - 2;
218
219 // Matrix of weighted sums
220 Eigen::Matrix<double, 3, 3> sumMatrixWSZ = getWSZSumMatrix(observationsSZ);
221 UncertainSZLine uncertainSZLine = fitZ(sumMatrixWSZ);
222
223 uncertainSZLine.setNDF(ndf);
224 trajectorySZ.setSZLine(uncertainSZLine);
225}
std::size_t size() const
Returns the number of observations stored.
void setSZLine(const UncertainSZLine &szLine)
Setter for the line in sz space.
void clear()
Clears all information from this trajectory line.
A line in sz where s is the transverse travel distance as seen in the xy projection with uncertaintie...
void setNDF(std::size_t ndf)
Setter for the number of degrees of freediom used in the line fit.

◆ update() [2/3]

void update ( CDCTrajectorySZ trajectorySZ,
const CDCSegment2D stereoSegment,
const CDCTrajectory2D axialTrajectory2D 
) const

Update the given sz trajectory reconstructing the stereo segment with a near by axial segment.

Definition at line 193 of file CDCSZFitter.cc.

196{
197 B2ASSERT("Expected stereo segment", not stereoSegment.isAxial());
198
199 // recostruct the stereo segment
200 CDCSZObservations observationsSZ;
201 for (const CDCRecoHit2D& recoHit2D : stereoSegment) {
202 CDCRecoHit3D recoHit3D = CDCRecoHit3D::reconstruct(recoHit2D, axialTrajectory2D);
203 observationsSZ.append(recoHit3D);
204 }
205
206 update(trajectorySZ, observationsSZ);
207}
Class representing a two dimensional reconstructed hit in the central drift chamber.
Definition: CDCRecoHit2D.h:47
Class representing a three dimensional reconstructed hit.
Definition: CDCRecoHit3D.h:52
static CDCRecoHit3D reconstruct(const CDCRecoHit2D &recoHit2D, const CDCTrajectory2D &trajectory2D)
Reconstructs the three dimensional hit from the two dimensional and the two dimensional trajectory.
Definition: CDCRecoHit3D.cc:56
std::size_t append(const CDCRecoHit3D &recoHit3D)
Appends the observed position.

◆ update() [3/3]

void update ( const CDCSegmentPair segmentPair) const

Updates the trajectory of the axial stereo segment pair inplace.

Definition at line 163 of file CDCSZFitter.cc.

164{
165 const CDCSegment2D* ptrStereoSegment = segmentPair.getStereoSegment();
166 const CDCSegment2D* ptrAxialSegment = segmentPair.getAxialSegment();
167
168 assert(ptrStereoSegment);
169 assert(ptrAxialSegment);
170
171 const CDCSegment2D& stereoSegment = *ptrStereoSegment;
172 const CDCSegment2D& axialSegment = *ptrAxialSegment;
173 const CDCTrajectory2D& axialTrajectory2D = axialSegment.getTrajectory2D();
174
175 CDCTrajectorySZ trajectorySZ;
176 update(trajectorySZ, stereoSegment, axialTrajectory2D);
177
178 CDCTrajectory3D trajectory3D(axialTrajectory2D, trajectorySZ);
179 segmentPair.setTrajectory3D(trajectory3D);
180}
A reconstructed sequence of two dimensional hits in one super layer.
Definition: CDCSegment2D.h:39
const CDCSegment2D * getAxialSegment() const
Getter for the axial segment.
void setTrajectory3D(const CDCTrajectory3D &trajectory3D) const
Setter for the three dimensional trajectory.
const CDCSegment2D * getStereoSegment() const
Getter for the stereo segment.
CDCTrajectory2D & getTrajectory2D() const
Getter for the two dimensional trajectory fitted to the segment.
Definition: CDCSegment.h:69
Particle trajectory as it is seen in xy projection represented as a circle.
Particle full three dimensional trajectory.

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