Belle II Software  release-05-01-25
CDCAxialStereoFusion.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2014 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Oliver Frost *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 #include <tracking/trackFindingCDC/fitting/CDCAxialStereoFusion.h>
11 
12 #include <tracking/trackFindingCDC/fitting/CDCSZFitter.h>
13 #include <tracking/trackFindingCDC/fitting/CDCRiemannFitter.h>
14 
15 #include <tracking/trackFindingCDC/eventdata/tracks/CDCSegmentPair.h>
16 #include <tracking/trackFindingCDC/eventdata/segments/CDCSegment3D.h>
17 #include <tracking/trackFindingCDC/eventdata/segments/CDCSegment2D.h>
18 #include <tracking/trackFindingCDC/eventdata/trajectories/CDCTrajectory3D.h>
19 #include <tracking/trackFindingCDC/eventdata/trajectories/CDCTrajectory2D.h>
20 #include <tracking/trackFindingCDC/eventdata/trajectories/CDCTrajectorySZ.h>
21 
22 #include <tracking/trackFindingCDC/topology/CDCWire.h>
23 
24 using namespace Belle2;
25 using namespace TrackFindingCDC;
26 
27 namespace {
28  CDCSegment3D reconstruct(const CDCSegment2D& segment2D,
29  const CDCTrajectory3D& trajectory3D)
30  {
31  CDCSegment3D result;
32  CDCTrajectory2D trajectory2D = trajectory3D.getTrajectory2D();
33  CDCTrajectorySZ trajectorySZ = trajectory3D.getTrajectorySZ();
34 
35  result.reserve(segment2D.size());
36  for (const CDCRecoHit2D& recoHit2D : segment2D) {
37  result.push_back(CDCRecoHit3D::reconstruct(recoHit2D, trajectory2D, trajectorySZ));
38  }
39  return result;
40  }
41 }
42 
44 {
45  const CDCSegment2D* ptrFromSegment = segmentPair.getFromSegment();
46  const CDCSegment2D* ptrToSegment = segmentPair.getToSegment();
47 
48  if (not ptrFromSegment) {
49  B2WARNING("From segment unset.");
50  return;
51  }
52 
53  if (not ptrToSegment) {
54  B2WARNING("To segment unset.");
55  return;
56  }
57 
58  const CDCSegment2D& fromSegment = *ptrFromSegment;
59  const CDCSegment2D& toSegment = *ptrToSegment;
60 
61  CDCTrajectory3D trajectory3D = reconstructFuseTrajectories(fromSegment, toSegment);
62  segmentPair.setTrajectory3D(trajectory3D);
63 }
64 
66 {
67  const CDCSegment2D* ptrFromSegment = segmentPair.getFromSegment();
68  const CDCSegment2D* ptrToSegment = segmentPair.getToSegment();
69 
70  if (not ptrFromSegment) {
71  B2WARNING("From segment unset.");
72  return;
73  }
74 
75  if (not ptrToSegment) {
76  B2WARNING("To segment unset.");
77  return;
78  }
79  const CDCSegment2D& fromSegment = *ptrFromSegment;
80  const CDCSegment2D& toSegment = *ptrToSegment;
81 
82  CDCTrajectory3D trajectory3D = fusePreliminary(fromSegment, toSegment);
83  segmentPair.setTrajectory3D(trajectory3D);
84 }
85 
86 
88  const CDCSegment2D& toSegment2D)
89 {
90  CDCTrajectory3D preliminaryTrajectory3D = fusePreliminary(fromSegment2D, toSegment2D);
91  return reconstructFuseTrajectories(fromSegment2D, toSegment2D, preliminaryTrajectory3D);
92 }
93 
95  const CDCSegment2D& toSegment2D)
96 {
97  if (fromSegment2D.empty()) {
98  B2WARNING("From segment is empty.");
99  return CDCTrajectory3D();
100  }
101 
102  if (toSegment2D.empty()) {
103  B2WARNING("To segment is empty.");
104  return CDCTrajectory3D();
105  }
106 
107  bool fromIsAxial = fromSegment2D.isAxial();
108  const CDCSegment2D& axialSegment2D = fromIsAxial ? fromSegment2D : toSegment2D;
109  const CDCSegment2D& stereoSegment2D = not fromIsAxial ? fromSegment2D : toSegment2D;
110 
111  CDCTrajectory2D axialTrajectory2D = axialSegment2D.getTrajectory2D();
112 
113  Vector2D localOrigin2D = (fromIsAxial ? fromSegment2D.back() : toSegment2D.front()).getRecoPos2D();
114  axialTrajectory2D.setLocalOrigin(localOrigin2D);
115 
116  CDCSegment3D stereoSegment3D = CDCSegment3D::reconstruct(stereoSegment2D, axialTrajectory2D);
117 
118  CDCTrajectorySZ trajectorySZ;
119 
120  CDCSZFitter szFitter = CDCSZFitter::getFitter();
121  trajectorySZ = szFitter.fit(stereoSegment3D);
122  if (not trajectorySZ.isFitted()) {
123  CDCTrajectory3D result;
124  result.setChi2(NAN);
125  return result;
126  }
127 
128  CDCTrajectory3D preliminaryTrajectory3D(axialTrajectory2D, trajectorySZ);
129  Vector3D localOrigin3D(localOrigin2D, 0.0);
130  preliminaryTrajectory3D.setLocalOrigin(localOrigin3D);
131  return preliminaryTrajectory3D;
132 }
133 
135  const CDCSegment2D& toSegment2D,
136  const CDCTrajectory3D& preliminaryTrajectory3D)
137 {
138  Vector3D localOrigin3D = preliminaryTrajectory3D.getLocalOrigin();
139  Vector2D localOrigin2D = localOrigin3D.xy();
140 
141  CDCRiemannFitter riemannFitter;
142  //riemannFitter.useOnlyOrientation();
143  riemannFitter.useOnlyPosition();
144 
145  CDCSegment3D fromSegment3D = reconstruct(fromSegment2D, preliminaryTrajectory3D);
146  CDCSegment3D toSegment3D = reconstruct(toSegment2D, preliminaryTrajectory3D);
147 
149  double tanLambda = preliminaryTrajectory3D.getTanLambda();
150  m_driftLengthEstimator.updateDriftLength(fromSegment3D, tanLambda);
151  m_driftLengthEstimator.updateDriftLength(toSegment3D, tanLambda);
152  }
153 
154  CDCTrajectory2D fromTrajectory2D = riemannFitter.fit(fromSegment3D);
155  CDCTrajectory2D toTrajectory2D = riemannFitter.fit(toSegment3D);
156 
157  fromTrajectory2D.setLocalOrigin(localOrigin2D);
158  toTrajectory2D.setLocalOrigin(localOrigin2D);
159 
160  SZParameters refSZ = preliminaryTrajectory3D.getLocalSZLine().szParameters();
161 
162  const UncertainPerigeeCircle& fromCircle = fromTrajectory2D.getLocalCircle();
163  const UncertainPerigeeCircle& toCircle = toTrajectory2D.getLocalCircle();
164 
165  JacobianMatrix<3, 5> fromH = calcAmbiguity(fromSegment3D, fromTrajectory2D);
166  JacobianMatrix<3, 5> toH = calcAmbiguity(toSegment3D, toTrajectory2D);
167 
168  UncertainHelix resultHelix = UncertainHelix::average(fromCircle, fromH, toCircle, toH, refSZ);
169  return CDCTrajectory3D(localOrigin3D, resultHelix);
170 }
171 
173  const CDCTrajectory2D& trajectory2D)
174 {
175  size_t nHits = segment3D.size();
176 
177  const Vector2D& localOrigin2D = trajectory2D.getLocalOrigin();
178  const UncertainPerigeeCircle& localCircle = trajectory2D.getLocalCircle();
179 
180  double zeta = 0;
181  for (const CDCRecoHit3D& recoHit3D : segment3D) {
182  const Vector2D& recoPos2D = recoHit3D.getRecoPos2D();
183  const Vector2D localRecoPos2D = recoPos2D - localOrigin2D;
184  const Vector2D normal = localCircle->normal(localRecoPos2D);
185  const CDCWire& wire = recoHit3D.getWire();
186  zeta += wire.getWireLine().sagMovePerZ(recoHit3D.getRecoZ()).dot(normal);
187  }
188  zeta /= nHits;
189 
191 
192  using namespace NHelixParameterIndices;
193  result(c_Curv, c_Curv) = 1.0;
194  result(c_Phi0, c_Phi0) = 1.0;
195  result(c_I, c_I) = 1.0;
196 
197  result(c_Phi0, c_TanL) = zeta;
198  result(c_I, c_Z0) = - zeta;
199 
200  return result;
201 }
Belle2::TrackFindingCDC::UncertainHelix::average
static UncertainHelix average(const UncertainHelix &fromHelix, const UncertainHelix &toHelix)
Construct the averages of the two given helices by properly considering their covariance matrix.
Definition: UncertainHelix.cc:25
Belle2::TrackFindingCDC::CDCTrajectory3D::getTanLambda
double getTanLambda() const
Getter for the slope of z over the transverse travel distance s.
Definition: CDCTrajectory3D.h:299
Belle2::TrackFindingCDC::CDCRecoHit3D
Class representing a three dimensional reconstructed hit.
Definition: CDCRecoHit3D.h:62
Belle2::TrackFindingCDC::CDCRiemannFitter
Class implementing the Riemann fit for two dimensional trajectory circle.
Definition: CDCRiemannFitter.h:34
Belle2::TrackFindingCDC::CDCSegmentPair
Class representing a pair of one reconstructed axial segement and one stereo segment in adjacent supe...
Definition: CDCSegmentPair.h:44
Belle2::TrackFindingCDC::Vector2D
A two dimensional vector which is equipped with functions for correct handeling of orientation relat...
Definition: Vector2D.h:37
Belle2::TrackFindingCDC::CDCSZFitter::getFitter
static const CDCSZFitter & getFitter()
Getter for a standard sz line fitter instance.
Definition: CDCSZFitter.cc:38
Belle2::TrackFindingCDC::WireLine::sagMovePerZ
Vector2D sagMovePerZ(const double z) const
Gives the two dimensional position with wire sag effect of the line at the given z value.
Definition: WireLine.h:85
Belle2::TrackFindingCDC::CDCTrajectorySZ::isFitted
bool isFitted() const
Indicates if the line has been fitted.
Definition: CDCTrajectorySZ.h:106
Belle2::TrackFindingCDC::CDCTrajectory2D::getLocalCircle
const UncertainPerigeeCircle & getLocalCircle() const
Getter for the cirlce in local coordinates.
Definition: CDCTrajectory2D.h:466
Belle2::TrackFindingCDC::CDCFitter2D::useOnlyPosition
void useOnlyPosition()
Setup the fitter to use only the reconstructed positions of the hits.
Definition: CDCFitter2D.icc.h:209
Belle2::TrackFindingCDC::CDCTrajectory3D::getLocalSZLine
UncertainSZLine getLocalSZLine() const
Getter for the sz line starting from the local origin.
Definition: CDCTrajectory3D.cc:364
Belle2::TrackFindingCDC::CDCAxialStereoFusion::reconstructFuseTrajectories
void reconstructFuseTrajectories(const CDCSegmentPair &segmentPair)
Combine the two trajectories of the segments in the pair and assign the resulting three dimensional t...
Definition: CDCAxialStereoFusion.cc:43
Belle2::TrackFindingCDC::CDCTrajectory3D::setLocalOrigin
double setLocalOrigin(const Vector3D &localOrigin)
Setter for the origin of the local coordinate system.
Definition: CDCTrajectory3D.cc:369
Belle2::TrackFindingCDC::CDCSegment3D::reconstruct
static CDCSegment3D reconstruct(const CDCSegment2D &segment2D, const CDCTrajectory2D &trajectory2D)
Reconstructs a two dimensional stereo segment by shifting each hit onto the given two dimensional tra...
Definition: CDCSegment3D.cc:22
Belle2::TrackFindingCDC::CDCTrajectory2D::setLocalOrigin
double setLocalOrigin(const Vector2D &localOrigin)
Setter for the origin of the local coordinate system.
Definition: CDCTrajectory2D.cc:357
Belle2::TrackFindingCDC::CDCTrajectorySZ
Linear trajectory in sz space.
Definition: CDCTrajectorySZ.h:41
Belle2::TrackFindingCDC::CDCSegment3D
A segment consisting of three dimensional reconstructed hits.
Definition: CDCSegment3D.h:36
Belle2::TrackFindingCDC::Vector2D::dot
double dot(const Vector2D &rhs) const
Calculates the two dimensional dot product.
Definition: Vector2D.h:172
Belle2::TrackFindingCDC::CDCTrajectory2D
Particle trajectory as it is seen in xy projection represented as a circle.
Definition: CDCTrajectory2D.h:46
Belle2::TrackFindingCDC::DriftLengthEstimator::updateDriftLength
double updateDriftLength(CDCRecoHit2D &recoHit2D)
Update the drift length of the reconstructed hit in place.
Definition: DriftLengthEstimator.cc:58
Belle2::TrackFindingCDC::CDCSegmentPair::getToSegment
const CDCSegment2D * getToSegment() const
Getter for the to segment.
Definition: CDCSegmentPair.h:132
Belle2::TrackFindingCDC::CDCAxialStereoFusion::m_reestimateDriftLength
bool m_reestimateDriftLength
Swtich to reestimate the drift length.
Definition: CDCAxialStereoFusion.h:93
Belle2::TrackFindingCDC::CDCTrajectory2D::getLocalOrigin
const Vector2D & getLocalOrigin() const
Getter for the origin of the local coordinate system.
Definition: CDCTrajectory2D.h:508
Belle2::TrackFindingCDC::CDCSegment::isAxial
bool isAxial() const
Indicator if the underlying wires are axial.
Definition: CDCSegment.h:55
Belle2::TrackFindingCDC::UncertainPerigeeCircle
Adds an uncertainty matrix to the circle in perigee parameterisation.
Definition: UncertainPerigeeCircle.h:39
Belle2::TrackFindingCDC::CDCSegment::getTrajectory2D
CDCTrajectory2D & getTrajectory2D() const
Getter for the two dimensional trajectory fitted to the segment.
Definition: CDCSegment.h:79
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TrackFindingCDC::Vector3D
A three dimensional vector.
Definition: Vector3D.h:34
Belle2::TrackFindingCDC::CDCSZFitter
Class implementing the z coordinate over travel distance line fit.
Definition: CDCSZFitter.h:37
Belle2::TrackFindingCDC::CDCFitter2D::fit
CDCTrajectory2D fit(const CDCObservations2D &observations2D) const
Fits a collection of observation drift circles.
Definition: CDCFitter2D.icc.h:48
Belle2::TrackFindingCDC::CDCAxialStereoFusion::calcAmbiguity
PerigeeHelixAmbiguity calcAmbiguity(const CDCSegment3D &segment3D, const CDCTrajectory2D &trajectory2D)
Calculate the ambiguity of the helix parameters relative to the three circle parameters given the hit...
Definition: CDCAxialStereoFusion.cc:172
Belle2::TrackFindingCDC::CDCTrajectory3D::getLocalOrigin
const Vector3D & getLocalOrigin() const
Getter for the origin of the local coordinate system.
Definition: CDCTrajectory3D.h:353
Belle2::TrackFindingCDC::CDCWire::getWireLine
const WireLine & getWireLine() const
Getter for the wire line represenation of the wire.
Definition: CDCWire.h:190
Belle2::TrackFindingCDC::Vector3D::xy
const Vector2D & xy() const
Getter for the xy projected vector ( reference ! )
Definition: Vector3D.h:500
Belle2::TrackFindingCDC::UncertainSZLine::szParameters
SZParameters szParameters() const
Getter for the sz parameters in the order defined by ESZParameter.h.
Definition: UncertainSZLine.h:113
Belle2::TrackFindingCDC::CDCWire
Class representing a sense wire in the central drift chamber.
Definition: CDCWire.h:60
Belle2::TrackFindingCDC::CDCSegment2D
A reconstructed sequence of two dimensional hits in one super layer.
Definition: CDCSegment2D.h:40
Belle2::TrackFindingCDC::CDCAxialStereoFusion::fusePreliminary
void fusePreliminary(const CDCSegmentPair &segmentPair)
Fit the given segment pair using the preliminary helix fit without proper covariance matrix.
Definition: CDCAxialStereoFusion.cc:65
Belle2::TrackFindingCDC::PlainMatrix
A matrix implementation to be used as an interface typ through out the track finder.
Definition: PlainMatrix.h:50
Belle2::TrackFindingCDC::CDCRecoHit3D::reconstruct
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:58
Belle2::TrackFindingCDC::UncertainHelix
A general helix class including a covariance matrix.
Definition: UncertainHelix.h:44
Belle2::TrackFindingCDC::CDCSegmentPair::getFromSegment
const CDCSegment2D * getFromSegment() const
Getter for the from segment.
Definition: CDCSegmentPair.h:120
Belle2::TrackFindingCDC::PerigeeCircle::normal
Vector2D normal(const Vector2D &point) const
Unit normal vector from the circle to the given point.
Definition: PerigeeCircle.h:284
Belle2::TrackFindingCDC::HelixUtil::defaultPerigeeAmbiguity
static PerigeeAmbiguity defaultPerigeeAmbiguity()
Initialse a default covariance matrix to zero.
Definition: HelixParameters.cc:67
Belle2::TrackFindingCDC::CDCSegmentPair::setTrajectory3D
void setTrajectory3D(const CDCTrajectory3D &trajectory3D) const
Setter for the three dimensional trajectory.
Definition: CDCSegmentPair.h:207
Belle2::TrackFindingCDC::CDCTrajectory3D
Particle full three dimensional trajectory.
Definition: CDCTrajectory3D.h:47
Belle2::TrackFindingCDC::CDCSZFitter::fit
CDCTrajectorySZ fit(const CDCSegment2D &stereoSegment, const CDCTrajectory2D &axialTrajectory2D) const
Returns a fitted trajectory.
Definition: CDCSZFitter.cc:141
Belle2::TrackFindingCDC::CDCAxialStereoFusion::m_driftLengthEstimator
DriftLengthEstimator m_driftLengthEstimator
Helper object to carry out the drift length estimation.
Definition: CDCAxialStereoFusion.h:96