Belle II Software  release-05-02-19
CDCSZObservations.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 *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 #include <tracking/trackFindingCDC/fitting/CDCSZObservations.h>
11 
12 #include <tracking/trackFindingCDC/fitting/EigenObservationMatrix.h>
13 
14 #include <tracking/trackFindingCDC/eventdata/tracks/CDCTrack.h>
15 #include <tracking/trackFindingCDC/eventdata/segments/CDCSegment3D.h>
16 #include <tracking/trackFindingCDC/eventdata/hits/CDCRecoHit3D.h>
17 
18 #include <tracking/trackFindingCDC/topology/CDCWire.h>
19 
20 #include <framework/logging/Logger.h>
21 
22 using namespace Belle2;
23 using namespace TrackFindingCDC;
24 
25 std::size_t CDCSZObservations::fill(double s, double z, double weight)
26 {
27  if (std::isnan(s)) return 0;
28  if (std::isnan(z)) return 0;
29  if (std::isnan(weight)) {
30  B2WARNING("Weight is nan. Skipping observation");
31  return 0;
32  }
33 
34  m_szObservations.push_back(s);
35  m_szObservations.push_back(z);
36  m_szObservations.push_back(weight);
37  return 1;
38 }
39 
40 std::size_t CDCSZObservations::append(const CDCRecoHit3D& recoHit3D)
41 {
42  if (m_onlyStereo and recoHit3D.isAxial()) return 0;
43 
44  const double s = recoHit3D.getArcLength2D();
45  const double z = recoHit3D.getRecoPos3D().z();
46 
47  double weight = 1.0;
48  if (m_fitVariance == EFitVariance::c_Unit) {
49  weight = 1;
50  } else {
51  // Translate the drift length uncertainty to a uncertainty in z
52  // by the taking the projected wire vector part parallel to the displacement
53  // as a proportionality factor to the z direction.
54  const CDCWire& wire = recoHit3D.getWire();
55  const Vector3D& wireVector = wire.getWireVector();
56  const Vector2D disp2D = recoHit3D.getRecoDisp2D();
57  const double driftlengthVariance = recoHit3D.getRecoDriftLengthVariance();
58 
59  double dispNorm = disp2D.norm();
60 
61  double zeta = 1.0;
62  if (dispNorm == 0.0) {
63  zeta = wireVector.xy().norm() / wireVector.z();
64  } else {
65  zeta = wireVector.xy().dot(disp2D) / wireVector.z() / dispNorm;
66  }
67 
68  weight = zeta * zeta / driftlengthVariance;
69  }
70 
71  size_t appendedHit = fill(s, z, weight);
72  return appendedHit;
73 }
74 
75 std::size_t CDCSZObservations::appendRange(const std::vector<CDCRecoHit3D>& recoHit3Ds)
76 {
77  std::size_t nAppendedHits = 0;
78  for (const CDCRecoHit3D& recoHit3D : recoHit3Ds) {
79  nAppendedHits += append(recoHit3D);
80  }
81  return nAppendedHits;
82 }
83 
84 std::size_t CDCSZObservations::appendRange(const CDCSegment3D& segment3D)
85 {
86  const std::vector<CDCRecoHit3D> recoHit3Ds = segment3D;
87  return this->appendRange(recoHit3Ds);
88 }
89 
90 std::size_t CDCSZObservations::appendRange(const CDCTrack& track)
91 {
92  const std::vector<CDCRecoHit3D> recoHit3Ds = track;
93  return this->appendRange(recoHit3Ds);
94 }
95 
97 {
98  std::size_t n = size();
99  if (n == 0) return Vector2D(NAN, NAN);
100  std::size_t i = n / 2;
101 
102  if (isEven(n)) {
103  // For even number of observations use the middle one with the bigger distance from IP
104  Vector2D center1(getS(i), getZ(i));
105  Vector2D center2(getS(i - 1), getZ(i - 1));
106  return center1.normSquared() > center2.normSquared() ? center1 : center2;
107  } else {
108  Vector2D center1(getS(i), getZ(i));
109  return center1;
110  }
111 }
112 
114 {
115  Eigen::Matrix<double, 1, 2> eigenOrigin(origin.x(), origin.y());
116  EigenObservationMatrix eigenObservations = getEigenObservationMatrix(this);
117  eigenObservations.leftCols<2>().rowwise() -= eigenOrigin;
118 }
119 
121 {
122  // Pick an observation at the center
123  Vector2D centralPoint = getCentralPoint();
124  passiveMoveBy(centralPoint);
125  return centralPoint;
126 }
Belle2::TrackFindingCDC::CDCRecoHit3D::getRecoDisp2D
Vector2D getRecoDisp2D() const
Gets the displacement from the wire position in the xy plain at the reconstructed position.
Definition: CDCRecoHit3D.cc:154
Belle2::TrackFindingCDC::CDCRecoHit3D::getWire
const CDCWire & getWire() const
Getter for the wire.
Definition: CDCRecoHit3D.h:236
Belle2::TrackFindingCDC::CDCSZObservations::centralize
Vector2D centralize()
Picks one observation as a reference point and transform all observations to that new origin.
Definition: CDCSZObservations.cc:120
Belle2::TrackFindingCDC::CDCRecoHit3D
Class representing a three dimensional reconstructed hit.
Definition: CDCRecoHit3D.h:62
Belle2::TrackFindingCDC::CDCSZObservations::getCentralPoint
Vector2D getCentralPoint() const
Extracts the observation center that is at the index in the middle.
Definition: CDCSZObservations.cc:96
Belle2::TrackFindingCDC::CDCTrack
Class representing a sequence of three dimensional reconstructed hits.
Definition: CDCTrack.h:51
Belle2::TrackFindingCDC::Vector2D
A two dimensional vector which is equipped with functions for correct handeling of orientation relat...
Definition: Vector2D.h:37
Belle2::TrackFindingCDC::Vector2D::normSquared
double normSquared() const
Calculates .
Definition: Vector2D.h:183
Belle2::TrackFindingCDC::CDCSZObservations::m_fitVariance
EFitVariance m_fitVariance
Indicator which variance information should preferably be extracted from hits in calls to append.
Definition: CDCSZObservations.h:164
Belle2::TrackFindingCDC::Vector2D::y
double y() const
Getter for the y coordinate.
Definition: Vector2D.h:619
Belle2::TrackFindingCDC::CDCRecoHit3D::isAxial
bool isAxial() const
Indicator if the underlying wire is axial.
Definition: CDCRecoHit3D.h:224
Belle2::TrackFindingCDC::CDCSZObservations::append
std::size_t append(const CDCRecoHit3D &recoHit3D)
Appends the observed position.
Definition: CDCSZObservations.cc:40
Belle2::TrackFindingCDC::CDCSZObservations::size
std::size_t size() const
Returns the number of observations stored.
Definition: CDCSZObservations.h:55
Belle2::TrackFindingCDC::CDCSegment3D
A segment consisting of three dimensional reconstructed hits.
Definition: CDCSegment3D.h:36
Belle2::TrackFindingCDC::CDCWire::getWireVector
Vector3D getWireVector() const
Getter for the vector pointing from the back end ofthe wire to the front end of the wire.
Definition: CDCWire.h:250
Belle2::TrackFindingCDC::CDCRecoHit3D::getArcLength2D
double getArcLength2D() const
Getter for the travel distance in the xy projection.
Definition: CDCRecoHit3D.h:380
Belle2::TrackFindingCDC::Vector2D::dot
double dot(const Vector2D &rhs) const
Calculates the two dimensional dot product.
Definition: Vector2D.h:172
Belle2::TrackFindingCDC::CDCSZObservations::passiveMoveBy
void passiveMoveBy(const Vector2D &origin)
Moves all observations passively such that the given vector becomes to origin of the new coordinate s...
Definition: CDCSZObservations.cc:113
Belle2::TrackFindingCDC::CDCRecoHit3D::getRecoDriftLengthVariance
double getRecoDriftLengthVariance() const
Returns the drift length variance next to the reconstructed position.
Definition: CDCRecoHit3D.h:368
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
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::Vector3D::xy
const Vector2D & xy() const
Getter for the xy projected vector ( reference ! )
Definition: Vector3D.h:500
Belle2::TrackFindingCDC::Vector2D::x
double x() const
Getter for the x coordinate.
Definition: Vector2D.h:609
Belle2::TrackFindingCDC::CDCSZObservations::appendRange
std::size_t appendRange(const std::vector< CDCRecoHit3D > &recoHit3Ds)
Appends all reconstructed hits from the three dimensional track.
Definition: CDCSZObservations.cc:75
Belle2::TrackFindingCDC::Vector2D::norm
double norm() const
Calculates the length of the vector.
Definition: Vector2D.h:189
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::CDCWire
Class representing a sense wire in the central drift chamber.
Definition: CDCWire.h:60
Belle2::TrackFindingCDC::CDCSZObservations::m_szObservations
std::vector< double > m_szObservations
Memory for the individual observations.
Definition: CDCSZObservations.h:157
Belle2::TrackFindingCDC::Vector3D::z
double z() const
Getter for the z coordinate.
Definition: Vector3D.h:488
Belle2::TrackFindingCDC::CDCRecoHit3D::getRecoPos3D
const Vector3D & getRecoPos3D() const
Getter for the 3d position of the hit.
Definition: CDCRecoHit3D.h:295
Belle2::TrackFindingCDC::CDCSZObservations::fill
std::size_t fill(double s, double z, double weight=1.0)
Appends the observed position.
Definition: CDCSZObservations.cc:25
Belle2::TrackFindingCDC::CDCSZObservations::m_onlyStereo
bool m_onlyStereo
Switch to only use information from stereo hits.
Definition: CDCSZObservations.h:167