Belle II Software  release-05-01-25
CDCRobustSZFitter.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2016 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Oliver Frost, Nils Braun *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 #include <tracking/trackFindingCDC/fitting/CDCRobustSZFitter.h>
11 
12 #include <tracking/trackFindingCDC/fitting/CDCSZFitter.h>
13 
14 #include <tracking/trackFindingCDC/fitting/CDCSZObservations.h>
15 
16 #include <tracking/trackFindingCDC/eventdata/trajectories/CDCTrajectorySZ.h>
17 
18 #include <tracking/trackFindingCDC/numerics/Median.h>
19 #include <tracking/trackFindingCDC/numerics/WithWeight.h>
20 
21 using namespace Belle2;
22 using namespace TrackFindingCDC;
23 
25 {
26  // This seems to be some other algorithm
27 
28 
29  CDCTrajectorySZ trajectorySZ;
30  CDCSZFitter szFitter;
31 
32  if (observationsSZ.size() > 4) {
33  CDCSZObservations observationsSZFiltered;
34 
35  double meanTanLambda = 0;
36  double meanStartZ = 0;
37 
38  for (unsigned int i = 0; i < observationsSZ.size(); i++) {
39  for (unsigned int j = 0; j < observationsSZ.size(); j++) {
40  if (i != j) {
41  observationsSZFiltered.fill(observationsSZ.getS(j),
42  observationsSZ.getZ(j),
43  observationsSZ.getWeight(j));
44  }
45  } // for j
46 
47  szFitter.update(trajectorySZ, observationsSZFiltered);
48  meanTanLambda += trajectorySZ.getTanLambda();
49  meanStartZ += trajectorySZ.getZ0();
50  } // for i
51 
52  return CDCTrajectorySZ(meanTanLambda / observationsSZ.size(),
53  meanStartZ / observationsSZ.size());
54  } else {
56  }
57 }
58 
60 {
61  std::vector<double> tanLambdas;
62  tanLambdas.reserve(szObservations.size() * (szObservations.size() - 1) / 2);
63  for (unsigned int i = 0; i < szObservations.size(); i++) {
64  for (unsigned int j = i + 1; j < szObservations.size(); j++) {
65  double s1 = szObservations.getS(i);
66  double s2 = szObservations.getS(j);
67  if (s1 == s2) continue;
68  double z1 = szObservations.getZ(i);
69  double z2 = szObservations.getZ(j);
70  tanLambdas.push_back((z2 - z1) / (s2 - s1));
71  }
72  }
73 
74  const double tanLambda = median(std::move(tanLambdas));
75  const double z0 = getMedianZ0(szObservations, tanLambda);
76 
77  CDCTrajectorySZ trajectorySZ(tanLambda, z0);
78  trajectorySZ.setNDF(szObservations.size());
79  return trajectorySZ;
80 }
81 
83 {
84  std::vector<WithWeight<double> > weightedTanLambdas;
85  Weight totalWeight = 0;
86  weightedTanLambdas.reserve(szObservations.size() * (szObservations.size() - 1) / 2);
87  for (unsigned int i = 0; i < szObservations.size(); i++) {
88  for (unsigned int j = i + 1; j < szObservations.size(); j++) {
89  double s1 = szObservations.getS(i);
90  double s2 = szObservations.getS(j);
91  if (s1 == s2) continue;
92  double z1 = szObservations.getZ(i);
93  double z2 = szObservations.getZ(j);
94 
95  double w1 = szObservations.getWeight(i);
96  double w2 = szObservations.getWeight(j);
97 
98  // Longer legs receive proportionally longer weight.
99  Weight weight = std::abs(s2 - s1) * hypot2(w1, w2);
100  weightedTanLambdas.emplace_back((z2 - z1) / (s2 - s1), weight);
101  totalWeight += weight;
102  }
103  }
104 
105  for (WithWeight<double>& weightedTanLambda : weightedTanLambdas) {
106  weightedTanLambda.weight() /= totalWeight;
107  }
108 
109  const double tanLambda = weightedMedian(std::move(weightedTanLambdas));
110  const double z0 = getMedianZ0(szObservations, tanLambda);
111 
112  CDCTrajectorySZ trajectorySZ(tanLambda, z0);
113  trajectorySZ.setNDF(szObservations.size());
114  return trajectorySZ;
115 }
116 
117 double CDCRobustSZFitter::getMedianZ0(const CDCSZObservations& szObservations, double tanLambda) const
118 {
119  std::vector<double> z0s;
120  z0s.reserve(szObservations.size());
121  for (unsigned int i = 0; i < szObservations.size(); i++) {
122  double s = szObservations.getS(i);
123  double z = szObservations.getZ(i);
124  z0s.push_back(z - s * tanLambda);
125  }
126  const double z0 = median(std::move(z0s));
127  return z0;
128 }
Belle2::TrackFindingCDC::CDCTrajectorySZ::basicAssumption
static CDCTrajectorySZ basicAssumption()
Constucts a basic assumption, what the z0 start position and the sz slope are, including some broad v...
Definition: CDCTrajectorySZ.cc:29
Belle2::TrackFindingCDC::CDCSZObservations::size
std::size_t size() const
Returns the number of observations stored.
Definition: CDCSZObservations.h:55
Belle2::TrackFindingCDC::CDCRobustSZFitter::fitWeightedTheilSen
CDCTrajectorySZ fitWeightedTheilSen(const CDCSZObservations &szObservations) const
Implements the weighted Theil-Sen line fit algorithm.
Definition: CDCRobustSZFitter.cc:82
Belle2::TrackFindingCDC::CDCTrajectorySZ
Linear trajectory in sz space.
Definition: CDCTrajectorySZ.h:41
Belle2::TrackFindingCDC::CDCTrajectorySZ::setNDF
void setNDF(std::size_t ndf)
Setter for the number of degrees of freedom of the line fit.
Definition: CDCTrajectorySZ.h:138
Belle2::TrackFindingCDC::CDCRobustSZFitter::getMedianZ0
double getMedianZ0(const CDCSZObservations &szObservations, double tanLambda) const
Compute the median z0 intercept from the given observations and an estimated slope.
Definition: CDCRobustSZFitter.cc:117
Belle2::TrackFindingCDC::CDCSZObservations::getZ
double getZ(int iObservation) const
Getter for the z value of the observation at the given index.
Definition: CDCSZObservations.h:91
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::CDCRobustSZFitter::fitTheilSen
CDCTrajectorySZ fitTheilSen(const CDCSZObservations &szObservations) const
Implements the original Theil-Sen line fit algorithm.
Definition: CDCRobustSZFitter.cc:59
Belle2::TrackFindingCDC::CDCSZFitter::update
void update(const CDCSegmentPair &segmentPair) const
Updates the trajectory of the axial stereo segment pair inplace.
Definition: CDCSZFitter.cc:165
Belle2::TrackFindingCDC::CDCSZFitter
Class implementing the z coordinate over travel distance line fit.
Definition: CDCSZFitter.h:37
Belle2::TrackFindingCDC::WithWeight
A mixin class to attach a weight to an object.
Definition: WithWeight.h:34
Belle2::TrackFindingCDC::CDCRobustSZFitter::fitUsingSimplifiedTheilSen
CDCTrajectorySZ fitUsingSimplifiedTheilSen(const CDCSZObservations &szObservations) const
Fit a linear sz trajectory to the reconstructed stereo segment.
Definition: CDCRobustSZFitter.cc:24
Belle2::TrackFindingCDC::CDCSZObservations
Class serving as a storage of observed sz positions to present to the sz line fitters.
Definition: CDCSZObservations.h:36
Belle2::TrackFindingCDC::CDCSZObservations::getS
double getS(int iObservation) const
Getter for the arc length value of the observation at the given index.
Definition: CDCSZObservations.h:85
Belle2::TrackFindingCDC::CDCSZObservations::getWeight
double getWeight(int iObservation) const
Getter for the weight / inverse variance of the observation at the given index.
Definition: CDCSZObservations.h:97
Belle2::TrackFindingCDC::CDCTrajectorySZ::getTanLambda
double getTanLambda() const
Getter for the slope over the travel distance coordinate.
Definition: CDCTrajectorySZ.h:114
Belle2::TrackFindingCDC::CDCSZObservations::fill
std::size_t fill(double s, double z, double weight=1.0)
Appends the observed position.
Definition: CDCSZObservations.cc:25