Belle II Software  release-05-01-25
SegmentTrackVarSet.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2015 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Nils Braun *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 #include <tracking/trackFindingCDC/filters/segmentTrack/SegmentTrackVarSet.h>
11 
12 #include <tracking/trackFindingCDC/fitting/CDCObservations2D.h>
13 #include <tracking/trackFindingCDC/fitting/CDCRiemannFitter.h>
14 #include <tracking/trackFindingCDC/fitting/CDCSZFitter.h>
15 
16 #include <tracking/trackFindingCDC/eventdata/tracks/CDCTrack.h>
17 #include <tracking/trackFindingCDC/eventdata/segments/CDCSegment2D.h>
18 #include <tracking/trackFindingCDC/eventdata/hits/CDCWireHit.h>
19 #include <tracking/trackFindingCDC/eventdata/trajectories/CDCTrajectorySZ.h>
20 
21 #include <tracking/trackFindingCDC/topology/CDCWire.h>
22 
23 #include <tracking/trackFindingCDC/numerics/ToFinite.h>
24 
25 using namespace Belle2;
26 using namespace TrackFindingCDC;
27 
29 {
30  const CDCTrack* track = testPair->getFrom();
31  const CDCSegment2D* segment = testPair->getTo();
32 
33  assert(segment);
34  assert(track);
35 
36  double maxmimumTrajectoryDistanceFront = 0;
37  double maxmimumTrajectoryDistanceBack = 0;
38  double maxmimumHitDistanceFront = 0;
39  double maxmimumHitDistanceBack = 0;
40  double outOfCDC = 0; // 0 means no, 1 means yes
41  double hitsInSameRegion = 0;
42  double hitsInCommon = 0;
43 
44  const CDCRecoHit2D& front = segment->front();
45  const CDCRecoHit2D& back = segment->back();
46 
47  // Calculate distances
48  const CDCTrajectory3D& trajectoryTrack3D = track->getStartTrajectory3D();
49  const CDCTrajectory2D& trajectoryTrack2D = trajectoryTrack3D.getTrajectory2D();
50  const CDCTrajectorySZ& szTrajectoryTrack = trajectoryTrack3D.getTrajectorySZ();
51  double radius = trajectoryTrack2D.getGlobalCircle().absRadius();
52 
53  maxmimumTrajectoryDistanceFront = trajectoryTrack2D.getDist2D(front.getWireHit().getRefPos2D());
54  maxmimumTrajectoryDistanceBack = trajectoryTrack2D.getDist2D(back.getWireHit().getRefPos2D());
55 
56  var<named("z_distance")>() = 0;
57  var<named("theta_distance")>() = 0;
58 
59  if (segment->getStereoKind() == EStereoKind::c_Axial) {
60  CDCTrajectory2D& trajectorySegment = segment->getTrajectory2D();
61  if (not trajectoryTrack2D.isFitted()) {
63  fitter.update(trajectorySegment, *segment);
64  }
65  } else {
66  CDCObservations2D observations;
67  for (const CDCRecoHit2D& recoHit : *segment) {
68  const CDCRLWireHit& rlWireHit = recoHit.getRLWireHit();
69  CDCRecoHit3D recoHit3D = CDCRecoHit3D::reconstruct(rlWireHit, trajectoryTrack2D);
70  double s = recoHit3D.getArcLength2D();
71  double z = recoHit3D.getRecoZ();
72  observations.fill(s, z);
73  }
74 
75  if (observations.size() > 3) {
76  const CDCSZFitter& fitter = CDCSZFitter::getFitter();
77  const CDCTrajectorySZ& szTrajectorySegments = fitter.fit(observations);
78 
79  double startZTrack = szTrajectoryTrack.getZ0();
80  double startZSegments = szTrajectorySegments.getZ0();
81 
82  var<named("z_distance")>() = startZTrack - startZSegments;
83  var<named("theta_distance")>() = szTrajectoryTrack.getTanLambda() - szTrajectorySegments.getTanLambda();
84  }
85  }
86 
87  // Calculate if it is out of the CDC
88  Vector3D frontRecoPos3D = front.reconstruct3D(trajectoryTrack2D);
89  Vector3D backRecoPos3D = back.reconstruct3D(trajectoryTrack2D);
90 
91  if (segment->getStereoKind() != EStereoKind::c_Axial) {
92  double forwardZ = front.getWire().getWireLine().forwardZ();
93  double backwardZ = front.getWire().getWireLine().backwardZ();
94 
95  if (frontRecoPos3D.z() > forwardZ or frontRecoPos3D.z() < backwardZ or backRecoPos3D.z() > forwardZ
96  or backRecoPos3D.z() < backwardZ) {
97  outOfCDC = 1.0;
98  }
99  }
100 
101  // Get perpS of track in the beginning and the end
102  double perpSOfFront = trajectoryTrack2D.calcArcLength2D(segment->front().getRecoPos2D());
103  double perpSOfBack = trajectoryTrack2D.calcArcLength2D(segment->back().getRecoPos2D());
104 
105  double perpSMinimum = std::min(perpSOfFront, perpSOfBack);
106  double perpSMaximum = std::max(perpSOfFront, perpSOfBack);
107 
108  // Count number of hits in the same region
109  for (const CDCRecoHit3D& recoHit : *track) {
110  if (recoHit.getArcLength2D() < 0.8 * perpSMinimum or
111  recoHit.getArcLength2D() > 1.2 * perpSMaximum) {
112  continue;
113  }
114  if (recoHit.getISuperLayer() == segment->getISuperLayer()) {
115  hitsInSameRegion++;
116  } else if (abs(recoHit.getISuperLayer() - segment->getISuperLayer()) == 1) {
117  double distanceFront = (front.getWireHit().getRefPos2D() - recoHit.getRecoPos2D()).norm();
118  if (distanceFront > maxmimumHitDistanceFront) {
119  maxmimumHitDistanceFront = distanceFront;
120  }
121  double distanceBack = (back.getWireHit().getRefPos2D() - recoHit.getRecoPos2D()).norm();
122  if (distanceBack > maxmimumHitDistanceBack) {
123  maxmimumHitDistanceBack = distanceBack;
124  }
125  }
126  }
127 
128  // Count number of common hits
129  for (const CDCRecoHit3D& trackHit : *track) {
130  if (std::find_if(segment->begin(), segment->end(), [&trackHit](const CDCRecoHit2D & segmentHit) {
131  return segmentHit.getWireHit().getHit() == trackHit.getWireHit().getHit();
132  }) != segment->end()) {
133  hitsInCommon += 1;
134  }
135  }
136 
137  // Make a fit with all the hits and one with only the hits in the near range
138  CDCObservations2D observationsFull;
139  CDCObservations2D observationsNeigh;
140 
141  // Collect the observations
142  bool isAxialSegment = segment->getStereoKind() != EStereoKind::c_Axial;
143 
144  for (const CDCRecoHit3D& recoHit : *track) {
145  if (isAxialSegment and recoHit.getStereoKind() == EStereoKind::c_Axial) {
146  observationsFull.fill(recoHit.getWireHit().getRefPos2D());
147  if (abs(recoHit.getISuperLayer() - segment->getISuperLayer()) < 3) {
148  observationsNeigh.fill(recoHit.getWireHit().getRefPos2D());
149  }
150  } else if (not isAxialSegment and recoHit.getStereoKind() != EStereoKind::c_Axial) {
151  double s = recoHit.getArcLength2D();
152  double z = recoHit.getRecoZ();
153  observationsFull.fill(s, z);
154  if (abs(recoHit.getISuperLayer() - segment->getISuperLayer()) < 3) {
155  observationsNeigh.fill(s, z);
156  }
157  }
158  }
159 
160  const CDCTrajectorySZ& trajectorySZ = track->getStartTrajectory3D().getTrajectorySZ();
161  double tanLambda = trajectorySZ.getTanLambda();
162 
163  bool hasZInformation = tanLambda != 0;
164  double max_hit_z_distance = -1;
165  double sum_hit_z_distance = 0;
166  double stereo_quad_tree_distance = 0;
167 
168  if (hasZInformation) {
169  double thetaFirstSegmentHit = -10;
170 
171  for (const CDCRecoHit2D& recoHit2D : *segment) {
172  Vector3D reconstructedPosition = recoHit2D.reconstruct3D(trajectoryTrack2D);
173  const Vector2D& recoPos2D = recoHit2D.getRecoPos2D();
174  double perpS = trajectoryTrack2D.calcArcLength2D(recoPos2D);
175 
176 
177  double current_z_distance = std::abs(trajectorySZ.getZDist(perpS, reconstructedPosition.z()));
178  if (std::isnan(current_z_distance)) {
179  continue;
180  }
181 
182  if (thetaFirstSegmentHit == -10) {
183  thetaFirstSegmentHit = reconstructedPosition.theta();
184  }
185  sum_hit_z_distance += current_z_distance;
186  if (current_z_distance > max_hit_z_distance) {
187  max_hit_z_distance = current_z_distance;
188  }
189  }
190 
191  double thetaTrack = trajectoryTrack3D.getFlightDirection3DAtSupport().theta();
192  stereo_quad_tree_distance = thetaTrack - thetaFirstSegmentHit;
193  }
194 
195 
196  for (const CDCRecoHit2D& recoHit : *segment) {
197  if (isAxialSegment) {
198  observationsFull.fill(recoHit.getRecoPos2D());
199  observationsNeigh.fill(recoHit.getRecoPos2D());
200  } else {
201  const CDCRLWireHit& rlWireHit = recoHit.getRLWireHit();
202  CDCRecoHit3D recoHit3D = CDCRecoHit3D::reconstruct(rlWireHit, trajectoryTrack2D);
203  double s = recoHit3D.getArcLength2D();
204  double z = recoHit3D.getRecoZ();
205  observationsFull.fill(s, z);
206  observationsNeigh.fill(s, z);
207  }
208  }
209 
210  // Do the fit
211  var<named("fit_neigh")>() = 0;
212  var<named("fit_full")>() = 0;
213  if (segment->getStereoKind() == EStereoKind::c_Axial) {
215  var<named("fit_full")>() = fitter.fit(observationsFull).getPValue();
216  } else {
217  const CDCSZFitter& fitter = CDCSZFitter::getFitter();
218  var<named("fit_full")>() = toFinite(fitter.fit(observationsFull).getPValue(), 0);
219 
220  if (observationsNeigh.size() > 3) {
221  var<named("fit_neigh")>() = toFinite(fitter.fit(observationsNeigh).getPValue(), 0);
222  } else {
223  var<named("fit_neigh")>() = 0;
224  }
225  }
226 
227  if (observationsFull.size() == observationsNeigh.size()) {
228  var<named("fit_neigh")>() = -1;
229  }
230 
231  var<named("is_stereo")>() = segment->getStereoKind() != EStereoKind::c_Axial;
232  var<named("segment_size")>() = segment->size();
233  var<named("track_size")>() = track->size();
234  var<named("mean_hit_z_distance")>() = sum_hit_z_distance;
235  var<named("max_hit_z_distance")>() = max_hit_z_distance;
236  var<named("stereo_quad_tree_distance")>() = toFinite(stereo_quad_tree_distance, 0);
237 
238  var<named("pt_of_track")>() = toFinite(std::isnan(trajectoryTrack2D.getAbsMom2D()) ? 0.0 : trajectoryTrack2D.getAbsMom2D(), 0);
239  var<named("track_is_curler")>() = trajectoryTrack2D.getExit().hasNAN();
240 
241  var<named("superlayer_already_full")>() = not trajectoryTrack2D.getOuterExit().hasNAN() and hitsInSameRegion > 5;
242 
243  var<named("maxmimum_trajectory_distance_front")>() = toFinite(maxmimumTrajectoryDistanceFront, 999);
244  var<named("maxmimum_trajectory_distance_back")>() = toFinite(maxmimumTrajectoryDistanceBack, 999);
245 
246  var<named("maxmimum_hit_distance_front")>() = maxmimumHitDistanceFront;
247  var<named("maxmimum_hit_distance_back")>() = maxmimumHitDistanceBack;
248 
249  var<named("out_of_CDC")>() = outOfCDC;
250  var<named("hits_in_same_region")>() = hitsInSameRegion;
251 
252  var<named("number_of_hits_in_common")>() = hitsInCommon;
253 
254  var<named("segment_super_layer")>() = segment->getISuperLayer();
255 
256  double phiBetweenTrackAndSegment = trajectoryTrack2D.getMom2DAtSupport().angleWith(segment->front().getRecoPos2D());
257 
258  var<named("phi_between_track_and_segment")>() = toFinite(phiBetweenTrackAndSegment, 0);
259  var<named("perp_s_of_front")>() = toFinite(perpSOfFront / radius, 0);
260  var<named("perp_s_of_back")>() = toFinite(perpSOfBack / radius, 0);
261 
262  return true;
263 }
Belle2::TrackFindingCDC::SegmentTrackVarSet::extract
bool extract(const BaseSegmentTrackFilter::Object *testPair) final
Generate and assign the contained variables.
Definition: SegmentTrackVarSet.cc:28
Belle2::TrackFindingCDC::CDCObservations2D::size
std::size_t size() const
Returns the number of observations stored.
Definition: CDCObservations2D.h:88
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::CDCTrack
Class representing a sequence of three dimensional reconstructed hits.
Definition: CDCTrack.h:51
Belle2::TrackFindingCDC::CDCTrajectory2D::getOuterExit
Vector2D getOuterExit(double factor=1) const
Calculates the point where the trajectory meets the outer wall of the CDC.
Definition: CDCTrajectory2D.cc:316
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::CDCWireHit::getRefPos2D
const Vector2D & getRefPos2D() const
The two dimensional reference position (z=0) of the underlying wire.
Definition: CDCWireHit.cc:214
Belle2::TrackFindingCDC::CDCRecoHit3D::getRecoZ
double getRecoZ() const
Getter for the z coordinate of the reconstructed position.
Definition: CDCRecoHit3D.h:313
Belle2::TrackFindingCDC::CDCTrajectorySZ
Linear trajectory in sz space.
Definition: CDCTrajectorySZ.h:41
Belle2::TrackFindingCDC::CDCRecoHit3D::getArcLength2D
double getArcLength2D() const
Getter for the travel distance in the xy projection.
Definition: CDCRecoHit3D.h:380
Belle2::TrackFindingCDC::CDCTrajectory2D
Particle trajectory as it is seen in xy projection represented as a circle.
Definition: CDCTrajectory2D.h:46
Belle2::TrackFindingCDC::Vector2D::hasNAN
bool hasNAN() const
Checks if one of the coordinates is NAN.
Definition: Vector2D.h:163
Belle2::TrackFindingCDC::CDCRecoHit2D::getWire
const CDCWire & getWire() const
Getter for the wire the reconstructed hit assoziated to.
Definition: CDCRecoHit2D.h:185
Belle2::TrackFindingCDC::CDCTrajectory3D::getTrajectory2D
CDCTrajectory2D getTrajectory2D() const
Getter for the two dimensional trajectory.
Definition: CDCTrajectory3D.cc:336
Belle2::TrackFindingCDC::CDCObservations2D
Class serving as a storage of observed drift circles to present to the Riemann fitter.
Definition: CDCObservations2D.h:53
Belle2::TrackFindingCDC::CDCTrajectory3D::getTrajectorySZ
CDCTrajectorySZ getTrajectorySZ() const
Getter for the sz trajectory.
Definition: CDCTrajectory3D.cc:343
Belle2::TrackFindingCDC::CDCRecoHit2D
Class representing a two dimensional reconstructed hit in the central drift chamber.
Definition: CDCRecoHit2D.h:57
Belle2::TrackFindingCDC::CDCTrajectorySZ::getZ0
double getZ0() const
Getter for the z coordinate at zero travel distance.
Definition: CDCTrajectorySZ.h:118
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TrackFindingCDC::CDCRLWireHit
Class representing an oriented hit wire including a hypotheses whether the causing track passes left ...
Definition: CDCRLWireHit.h:51
Belle2::TrackFindingCDC::Vector3D::theta
double theta() const
Getter for the polar angle.
Definition: Vector3D.h:538
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::CDCWire::getWireLine
const WireLine & getWireLine() const
Getter for the wire line represenation of the wire.
Definition: CDCWire.h:190
Belle2::TrackFindingCDC::WireLine::forwardZ
double forwardZ() const
Gives the forward z coodinate.
Definition: WireLine.h:140
Belle2::TrackFindingCDC::CDCTrajectory2D::isFitted
bool isFitted() const
Checks if the circle is already set to a valid value.
Definition: CDCTrajectory2D.cc:85
Belle2::TrackFindingCDC::CDCTrajectory3D::getFlightDirection3DAtSupport
Vector3D getFlightDirection3DAtSupport() const
Get the unit momentum at the start point of the trajectory.
Definition: CDCTrajectory3D.h:222
Belle2::TrackFindingCDC::CDCRiemannFitter::getFitter
static const CDCRiemannFitter & getFitter()
Static getter for a general Riemann fitter.
Definition: CDCRiemannFitter.cc:22
Belle2::TrackFindingCDC::WireLine::backwardZ
double backwardZ() const
Gives the backward z coodinate.
Definition: WireLine.h:144
Belle2::TrackFindingCDC::CDCTrajectorySZ::getZDist
double getZDist(const double s, const double z) const
Calculates the distance along between the given point an the sz trajectory.
Definition: CDCTrajectorySZ.h:78
Belle2::TrackFindingCDC::CDCSegment2D
A reconstructed sequence of two dimensional hits in one super layer.
Definition: CDCSegment2D.h:40
Belle2::TrackFindingCDC::VarSet< SegmentTrackVarNames >::named
constexpr static int named(const char *name)
Getter for the index from the name.
Definition: VarSet.h:88
Belle2::TrackFindingCDC::CDCTrajectory2D::getDist2D
double getDist2D(const Vector2D &point) const
Calculates the distance from the point to the trajectory as seen from the xy projection.
Definition: CDCTrajectory2D.h:419
Belle2::TrackFindingCDC::Vector2D::angleWith
double angleWith(const Vector2D &rhs) const
The angle between this and rhs.
Definition: Vector2D.h:211
Belle2::TrackFindingCDC::Vector3D::z
double z() const
Getter for the z coordinate.
Definition: Vector3D.h:488
Belle2::TrackFindingCDC::VarSet< SegmentTrackVarNames >::var
Float_t & var()
Reference getter for the value of the ith variable. Static version.
Definition: VarSet.h:103
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::CDCTrajectory2D::getExit
Vector2D getExit() const
Calculates the point where the trajectory leaves the CDC.
Definition: CDCTrajectory2D.cc:336
Belle2::TrackFindingCDC::CDCTrajectory2D::getMom2DAtSupport
Vector2D getMom2DAtSupport(const double bZ) const
Get the momentum at the support point of the trajectory.
Definition: CDCTrajectory2D.h:308
Belle2::TrackFindingCDC::CDCTrajectorySZ::getTanLambda
double getTanLambda() const
Getter for the slope over the travel distance coordinate.
Definition: CDCTrajectorySZ.h:114
Belle2::TrackFindingCDC::CDCObservations2D::fill
std::size_t fill(double x, double y, double signedRadius=0.0, double weight=1.0)
Appends the observed position.
Belle2::TrackFindingCDC::CDCTrajectory2D::calcArcLength2D
double calcArcLength2D(const Vector2D &point) const
Calculate the travel distance from the start position of the trajectory.
Definition: CDCTrajectory2D.h:270
Belle2::TrackFindingCDC::CDCRecoHit2D::reconstruct3D
Vector3D reconstruct3D(const CDCTrajectory2D &trajectory2D, const double z=0) const
Reconstruct the three dimensional position (especially of stereo hits) by determinating the z coordin...
Definition: CDCRecoHit2D.cc:133
Belle2::TrackFindingCDC::CDCTrajectory3D
Particle full three dimensional trajectory.
Definition: CDCTrajectory3D.h:47
Belle2::TrackFindingCDC::PerigeeCircle::absRadius
double absRadius() const
Gives the signed radius of the circle. If it was a line this will be infinity.
Definition: PerigeeCircle.h:350
Belle2::TrackFindingCDC::CDCRecoHit2D::getWireHit
const CDCWireHit & getWireHit() const
Getter for the wire hit assoziated with the reconstructed hit.
Definition: CDCRecoHit2D.h:203
Belle2::TrackFindingCDC::Filter::Object
AObject Object
Type of the object to be analysed.
Definition: Filter.dcl.h:43
Belle2::TrackFindingCDC::CDCTrajectory2D::getGlobalCircle
PerigeeCircle getGlobalCircle() const
Getter for the circle in global coordinates.
Definition: CDCTrajectory2D.h:451
Belle2::TrackFindingCDC::CDCTrajectory2D::getAbsMom2D
double getAbsMom2D(double bZ) const
Get the estimation for the absolute value of the transvers momentum.