Belle II Software  release-05-01-25
TrackQualityTools.cc
1 #include <tracking/trackFindingCDC/processing/TrackQualityTools.h>
2 #include <tracking/trackFindingCDC/eventdata/tracks/CDCTrack.h>
3 #include <tracking/trackFindingCDC/eventdata/hits/CDCWireHit.h>
4 #include <tracking/trackFindingCDC/eventdata/trajectories/CDCTrajectory3D.h>
5 #include <tracking/trackFindingCDC/eventdata/trajectories/CDCTrajectory2D.h>
6 
7 #include <tracking/trackFindingCDC/topology/CDCWire.h>
8 #include <tracking/trackFindingCDC/utilities/ReversedRange.h>
9 
10 #include <TVector2.h>
11 
12 using namespace Belle2;
13 using namespace TrackFindingCDC;
14 
15 
17 {
18  for (const CDCRecoHit3D& recoHit : track) {
19  if (recoHit.getStereoKind() == EStereoKind::c_Axial) {
20  break;
21  } else {
22  recoHit.getWireHit().getAutomatonCell().setAssignedFlag();
23  }
24  }
25 }
26 
27 void TrackQualityTools::splitSecondHalfOfTrack(CDCTrack& track, std::vector<CDCTrack>& tracks)
28 {
29  const CDCTrajectory3D& trajectory3D = track.getStartTrajectory3D();
30  const CDCTrajectory2D& trajectory2D = trajectory3D.getTrajectory2D();
31  const double radius = trajectory2D.getLocalCircle()->absRadius();
32  const Vector2D& apogee = trajectory2D.getGlobalCircle().apogee();
33  double arcLength2DOfApogee = trajectory2D.calcArcLength2D(apogee);
34  if (arcLength2DOfApogee < 0) {
35  arcLength2DOfApogee += 2 * TMath::Pi() * radius;
36  }
37 
38  CDCTrack splittedCDCTrack;
39  splittedCDCTrack.setStartTrajectory3D(trajectory3D);
40 
41  for (CDCRecoHit3D& recoHit : track) {
42  double currentArcLength2D = recoHit.getArcLength2D();
43  if (currentArcLength2D < 0) B2INFO("Below 0");
44  if (currentArcLength2D > arcLength2DOfApogee) {
45  splittedCDCTrack.push_back(recoHit);
46  recoHit.getWireHit().getAutomatonCell().setAssignedFlag();
47  }
48  }
49 
50  if (splittedCDCTrack.size() > 0) {
51  tracks.push_back(splittedCDCTrack);
52  }
53 }
54 
56 {
57  if (track.empty()) return;
58 
59  CDCTrajectory3D trajectory3D = track.getStartTrajectory3D();
60 
61  // We reset the trajectory here to start (later) at the newStartPosition of the first hit
62  const Vector3D startPosition(0, 0, 0);
63  trajectory3D.setLocalOrigin(startPosition);
64  trajectory3D.setFlightTime(0);
65 
66  CDCTrajectory2D trajectory2D = trajectory3D.getTrajectory2D();
67 
68  // Check if we have to reverse the trajectory. This is done by counting the number of hits
69  // with positive and with negative arc length
70  unsigned int numberOfPositiveHits = 0;
71  for (const CDCRecoHit3D& recoHit : track) {
72  const double currentArcLength = trajectory2D.calcArcLength2D(recoHit.getRecoPos2D());
73  if (currentArcLength > 0) {
74  numberOfPositiveHits++;
75  }
76  }
77  const bool reverseTrajectory = 2 * numberOfPositiveHits < track.size();
78 
79  if (reverseTrajectory) {
80  trajectory3D.reverse();
81  trajectory2D.reverse();
82  }
83 
84  double arcLength2DPeriod = trajectory2D.getArcLength2DPeriod();
85  for (CDCRecoHit3D& recoHit : track) {
86  Vector2D recoPos2D = recoHit.getRLWireHit().reconstruct3D(trajectory2D).xy();
87  double arcLength2D = trajectory2D.calcArcLength2D(recoPos2D);
88 
89  if (arcLength2D < 0) {
90  arcLength2D += arcLength2DPeriod;
91  }
92  recoHit.setArcLength2D(arcLength2D);
93  recoHit.getWireHit().getAutomatonCell().unsetAssignedFlag();
94  recoHit.getWireHit().getAutomatonCell().setTakenFlag();
95  }
96 
97  // We can now sort by arc length
98  track.sortByArcLength2D();
99 
100  // Set the position to the first hit and let the hits start at arc length of 0
101  Vector3D frontPosition = track.front().getRLWireHit().reconstruct3D(trajectory2D);
102  double arcLengthOffset = trajectory3D.setLocalOrigin(frontPosition);
103  track.setStartTrajectory3D(trajectory3D);
104  for (CDCRecoHit3D& recoHit : track) {
105  recoHit.shiftArcLength2D(-arcLengthOffset);
106  }
107 
108  // Set the back trajectory to start at the last hit position (and shift if necessary)
109  Vector3D backPosition = track.back().getRLWireHit().reconstruct3D(trajectory2D);
110  double backArcLength2D = trajectory3D.setLocalOrigin(backPosition);
111  if (backArcLength2D < 0) {
112  trajectory3D.shiftPeriod(1);
113  }
114  track.setEndTrajectory3D(trajectory3D);
115 
116 }
117 
118 void TrackQualityTools::removeHitsAfterCDCWall(CDCTrack& track, double outerCylindricalRFactor)
119 {
120  const CDCTrajectory2D& trajectory2D = track.getStartTrajectory3D().getTrajectory2D();
121  const double radius = trajectory2D.getLocalCircle()->absRadius();
122 
123  // Curler are allowed to have hits on both arms
124  if (trajectory2D.isCurler(outerCylindricalRFactor)) {
125  return;
126  }
127 
128  const Vector2D& outerExitWithFactor = trajectory2D.getOuterExit(outerCylindricalRFactor);
129 
130  double arcLength2DOfExitWithFactor = trajectory2D.calcArcLength2D(outerExitWithFactor);
131  if (arcLength2DOfExitWithFactor < 0) {
132  arcLength2DOfExitWithFactor += 2 * TMath::Pi() * radius;
133  }
134  bool removeAfterThis = false;
135 
136  for (const CDCRecoHit3D& recoHit : track) {
137  if (removeAfterThis) {
138  recoHit.getWireHit().getAutomatonCell().setAssignedFlag();
139  continue;
140  }
141 
142  double currentArcLength2D = recoHit.getArcLength2D();
143  if (currentArcLength2D < 0) {
144  currentArcLength2D += 2 * TMath::Pi() * radius;
145  }
146 
147  if (currentArcLength2D > arcLength2DOfExitWithFactor) {
148  recoHit.getWireHit().getAutomatonCell().setAssignedFlag();
149  removeAfterThis = true;
150  }
151  }
152 }
153 
155 {
156  ILayer lastLayer = -1;
157  Vector2D lastWirePosition;
158 
159  std::vector<std::vector<const CDCRecoHit3D*>> trackletList;
160  trackletList.reserve(3);
161 
162  std::vector<const CDCRecoHit3D*>* currentTracklet = nullptr;
163 
164  for (const CDCRecoHit3D& recoHit : track) {
165  if (currentTracklet == nullptr) {
166  trackletList.emplace_back();
167  currentTracklet = &(trackletList.back());
168  }
169 
170  const ILayer currentLayer = recoHit.getWire().getICLayer();
171  const Vector2D& currentPosition = recoHit.getRecoPos2D();
172  if (lastLayer != -1) {
173  const ILayer delta = currentLayer - lastLayer;
174  const double distance = (currentPosition - lastWirePosition).norm();
175  if (abs(delta) > 4 or distance > 50) {
176  trackletList.emplace_back();
177  currentTracklet = &(trackletList.back());
178  }
179  }
180 
181  lastWirePosition = currentPosition;
182  lastLayer = currentLayer;
183 
184  currentTracklet->push_back(&recoHit);
185  }
186 
187  if (trackletList.size() > 1) {
188  for (const std::vector<const CDCRecoHit3D*>& tracklet : trackletList) {
189  if (tracklet.size() < 5) {
190  for (const CDCRecoHit3D* recoHit : tracklet) {
191  recoHit->getWireHit().getAutomatonCell().setAssignedFlag();
192  }
193  }
194  }
195  }
196 }
197 
198 void TrackQualityTools::removeHitsAfterLayerBreak(CDCTrack& track, double m_maximumArcLength2DDistance)
199 {
200  const CDCTrajectory3D& trajectory3D = track.getStartTrajectory3D();
201  const CDCTrajectory2D& trajectory2D = trajectory3D.getTrajectory2D();
202  const double radius = trajectory2D.getLocalCircle()->absRadius();
203 
204  if (std::isnan(radius)) {
205  return;
206  }
207 
208  std::vector<std::vector<const CDCRecoHit3D*>> trackletList;
209 
210  double lastArcLength2D = NAN;
211 
212  std::vector<const CDCRecoHit3D*>* currentTracklet = nullptr;
213 
214  for (const CDCRecoHit3D& recoHit : track) {
215  if (currentTracklet == nullptr) {
216  trackletList.emplace_back();
217  currentTracklet = &(trackletList.back());
218  }
219 
220  const double currentArcLength2D = recoHit.getArcLength2D();
221  if (not std::isnan(lastArcLength2D)) {
222  const double delta = (currentArcLength2D - lastArcLength2D);
223  if (std::fabs(delta) > m_maximumArcLength2DDistance) {
224  trackletList.emplace_back();
225  currentTracklet = &(trackletList.back());
226  }
227  }
228 
229  lastArcLength2D = currentArcLength2D;
230 
231  currentTracklet->push_back(&recoHit);
232  }
233 
234  if (trackletList.size() > 1) {
235  // Throw away the ends if they are too small
236  while (trackletList.size() > 0) {
237  if (trackletList.back().size() > 4) {
238  break;
239  }
240 
241  for (const CDCRecoHit3D* recoHit : trackletList.back()) {
242  recoHit->getWireHit().getAutomatonCell().setAssignedFlag();
243  }
244 
245  trackletList.pop_back();
246  }
247 
248  std::reverse(trackletList.begin(), trackletList.end());
249 
250  while (trackletList.size() > 0) {
251  if (trackletList.back().size() > 4) {
252  break;
253  }
254 
255  for (const CDCRecoHit3D* recoHit : trackletList.back()) {
256  recoHit->getWireHit().getAutomatonCell().setAssignedFlag();
257  }
258 
259  trackletList.pop_back();
260  }
261  }
262 }
263 
264 void TrackQualityTools::removeHitsIfSmall(CDCTrack& track, unsigned int minimalHits)
265 {
266  const bool deleteTrack = track.size() < minimalHits;
267 
268  if (deleteTrack) {
269  for (const CDCRecoHit3D& recoHit : track) {
270  recoHit.getWireHit().getAutomatonCell().setAssignedFlag();
271  }
272  }
273 }
274 
276 {
277  double lastAngle = NAN;
278  bool removeAfterThis = false;
279 
280  for (const CDCRecoHit3D& recoHit : reversedRange(track)) {
281  if (removeAfterThis) {
282  recoHit.getWireHit().getAutomatonCell().setAssignedFlag();
283  continue;
284  }
285 
286  const double currentAngle = recoHit.getRecoPos2D().phi();
287  if (not std::isnan(lastAngle)) {
288  const double delta = currentAngle - lastAngle;
289  const double normalizedDelta = std::min(TVector2::Phi_0_2pi(delta), TVector2::Phi_0_2pi(-delta));
290  if (fabs(normalizedDelta) > maximalAngle) {
291  removeAfterThis = true;
292  recoHit.getWireHit().getAutomatonCell().setAssignedFlag();
293  }
294  }
295 
296  lastAngle = currentAngle;
297  }
298 }
299 
300 
302 {
303  ISuperLayer lastLayer = -1;
304  bool deleteTrack = true;
305 
306  for (const CDCRecoHit3D& recoHit : track) {
307  const ISuperLayer currentLayer = recoHit.getISuperLayer();
308  if (lastLayer != -1 and lastLayer != currentLayer) {
309  deleteTrack = false;
310  break;
311  }
312 
313  lastLayer = currentLayer;
314  }
315 
316  if (deleteTrack) {
317  for (const CDCRecoHit3D& recoHit : track) {
318  recoHit.getWireHit().getAutomatonCell().setAssignedFlag();
319  }
320  }
321 }
322 
324 {
325  const CDCTrajectory3D& trajectory3D = track.getStartTrajectory3D();
326  const CDCTrajectory2D& trajectory2D = trajectory3D.getTrajectory2D();
327 
328  // Curler are allowed to have negative arc lengths
329  if (trajectory2D.isCurler()) {
330  return;
331  }
332  for (const CDCRecoHit3D& recoHit : track) {
333  if (trajectory2D.calcArcLength2D(recoHit.getRecoPos2D()) < 0) {
334  recoHit.getWireHit().getAutomatonCell().setAssignedFlag();
335  }
336  }
337 }
338 
339 void TrackQualityTools::removeArcLength2DHoles(CDCTrack& track, double m_maximumArcLength2DDistance)
340 {
341  const CDCTrajectory3D& trajectory3D = track.getStartTrajectory3D();
342  const CDCTrajectory2D& trajectory2D = trajectory3D.getTrajectory2D();
343  const double radius = trajectory2D.getLocalCircle()->absRadius();
344 
345  if (std::isnan(radius)) {
346  return;
347  }
348 
349  double lastArcLength2D = NAN;
350  bool removeAfterThis = false;
351 
352  for (const CDCRecoHit3D& recoHit : track) {
353  if (removeAfterThis) {
354  recoHit.getWireHit().getAutomatonCell().setAssignedFlag();
355  continue;
356  }
357 
358  const double currentArcLength2D = recoHit.getArcLength2D();
359  if (not std::isnan(lastArcLength2D)) {
360  const double delta = (currentArcLength2D - lastArcLength2D) / radius;
361  if (fabs(delta) > m_maximumArcLength2DDistance) {
362  removeAfterThis = true;
363  recoHit.getWireHit().getAutomatonCell().setAssignedFlag();
364  continue;
365  }
366  }
367 
368  lastArcLength2D = currentArcLength2D;
369  }
370 }
Belle2::TrackFindingCDC::TrackQualityTools::removeArcLength2DHoles
static void removeArcLength2DHoles(CDCTrack &track, double m_maximumArcLength2DDistance=10)
Remove all hits that come after a large hole in the two dimensional arc length.
Definition: TrackQualityTools.cc:339
Belle2::TrackFindingCDC::TrackQualityTools::removeHitsOnTheWrongSide
static void removeHitsOnTheWrongSide(CDCTrack &track)
Remove all hits that are on the wrong side of the detector (so to say: "beyond the IP").
Definition: TrackQualityTools.cc:323
Belle2::TrackFindingCDC::CDCRecoHit3D
Class representing a three dimensional reconstructed hit.
Definition: CDCRecoHit3D.h:62
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::CDCTrajectory2D::getLocalCircle
const UncertainPerigeeCircle & getLocalCircle() const
Getter for the cirlce in local coordinates.
Definition: CDCTrajectory2D.h:466
Belle2::TrackFindingCDC::CDCTrajectory2D::reverse
void reverse()
Reverses the trajectory in place.
Definition: CDCTrajectory2D.cc:97
Belle2::TrackFindingCDC::CDCTrajectory3D::setLocalOrigin
double setLocalOrigin(const Vector3D &localOrigin)
Setter for the origin of the local coordinate system.
Definition: CDCTrajectory3D.cc:369
Belle2::TrackFindingCDC::TrackQualityTools::splitSecondHalfOfTrack
static void splitSecondHalfOfTrack(CDCTrack &track, std::vector< CDCTrack > &tracks)
Trasan did output curlers in split two halves - this method can be used to mimic this.
Definition: TrackQualityTools.cc:27
Belle2::TrackFindingCDC::TrackQualityTools::removeHitsAfterCDCWall
static void removeHitsAfterCDCWall(CDCTrack &track, double outerCylindricalRFactor=1.1)
Remove all hits which can not belong to the track, as the particle can not exit and enter the CDC aga...
Definition: TrackQualityTools.cc:118
Belle2::TrackFindingCDC::CDCTrajectory3D::shiftPeriod
double shiftPeriod(int nPeriods)
Adjusts the z0 to the one that lies n periods forward.
Definition: CDCTrajectory3D.cc:329
Belle2::TrackFindingCDC::CDCTrajectory2D::isCurler
bool isCurler(double factor=1) const
Checks if the trajectory leaves the outer radius of the CDC times the given tolerance factor.
Definition: CDCTrajectory2D.cc:267
Belle2::TrackFindingCDC::TrackQualityTools::normalizeHitsAndResetTrajectory
static void normalizeHitsAndResetTrajectory(CDCTrack &track)
Update all hits to have a positive perpS, a taken flag and no background flag Also set the trajectory...
Definition: TrackQualityTools.cc:55
Belle2::TrackFindingCDC::TrackQualityTools::removeHitsAfterLayerBreak
static void removeHitsAfterLayerBreak(CDCTrack &track, double m_maximumArcLength2DDistance=10)
Delete all hits after a large layer break.
Definition: TrackQualityTools.cc:198
Belle2::TrackFindingCDC::CDCTrajectory2D::getArcLength2DPeriod
double getArcLength2DPeriod() const
Getter for the arc length for one round trip around the trajectory.
Definition: CDCTrajectory2D.h:289
Belle2::TrackFindingCDC::CDCTrajectory2D
Particle trajectory as it is seen in xy projection represented as a circle.
Definition: CDCTrajectory2D.h:46
Belle2::TrackFindingCDC::TrackQualityTools::moveToNextAxialLayer
static void moveToNextAxialLayer(CDCTrack &track)
Delete hits of the first superlayer if it is a stereo one (fitting does not work very well when start...
Definition: TrackQualityTools.cc:16
Belle2::TrackFindingCDC::AutomatonCell::setAssignedFlag
void setAssignedFlag(bool setTo=true)
Sets the already assigned marker flag to the given value. Default value true.
Definition: AutomatonCell.h:138
Belle2::TrackFindingCDC::CDCTrajectory3D::getTrajectory2D
CDCTrajectory2D getTrajectory2D() const
Getter for the two dimensional trajectory.
Definition: CDCTrajectory3D.cc:336
Belle2::TrackFindingCDC::PerigeeCircle::apogee
Vector2D apogee() const
Getter for the apogee of the circle. If it was a line both components will be infinity.
Definition: PerigeeCircle.h:314
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TrackFindingCDC::TrackQualityTools::removeHitsAfterLayerBreak2
static void removeHitsAfterLayerBreak2(CDCTrack &track)
Delete all hits after a large layer break.
Definition: TrackQualityTools.cc:154
Belle2::TrackFindingCDC::Vector3D
A three dimensional vector.
Definition: Vector3D.h:34
Belle2::TrackFindingCDC::TrackQualityTools::removeHitsIfSmall
static void removeHitsIfSmall(CDCTrack &track, unsigned int minimalHits=7)
Delete a track fully of the number of hits is below minimalHits.
Definition: TrackQualityTools.cc:264
Belle2::TrackFindingCDC::CDCTrajectory3D::setFlightTime
void setFlightTime(double flightTime)
Setter for the time when the particle reached the support point position.
Definition: CDCTrajectory3D.h:372
Belle2::TrackFindingCDC::TrackQualityTools::removeHitsIfOnlyOneSuperLayer
static void removeHitsIfOnlyOneSuperLayer(CDCTrack &track)
Remove the whole track if it only consists of one superlayer.
Definition: TrackQualityTools.cc:301
Belle2::TrackFindingCDC::CDCTrajectory3D::reverse
void reverse()
Reverses the trajectory in place.
Definition: CDCTrajectory3D.h:151
Belle2::TrackFindingCDC::CDCTrack::getAutomatonCell
AutomatonCell & getAutomatonCell() const
Mutable getter for the automaton cell.
Definition: CDCTrack.h:135
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::CDCTrajectory3D
Particle full three dimensional trajectory.
Definition: CDCTrajectory3D.h:47
Belle2::TrackFindingCDC::CDCTrack::setStartTrajectory3D
void setStartTrajectory3D(const CDCTrajectory3D &startTrajectory3D)
Setter for the two dimensional trajectory.
Definition: CDCTrack.h:108
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::CDCTrajectory2D::getGlobalCircle
PerigeeCircle getGlobalCircle() const
Getter for the circle in global coordinates.
Definition: CDCTrajectory2D.h:451
Belle2::TrackFindingCDC::TrackQualityTools::removeHitsInTheBeginningIfAngleLarge
static void removeHitsInTheBeginningIfAngleLarge(CDCTrack &track, double maximalAngle=0.7)
If the angle between two following hits is larger than maximalAngle, delete all hits before (!...
Definition: TrackQualityTools.cc:275