Belle II Software  release-05-01-25
SVDRecoHit2D.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Peter Kvasnicka, Martin Ritter, Moritz Nadler *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <framework/logging/Logger.h>
12 #include <svd/reconstruction/SVDRecoHit2D.h>
13 #include <svd/geometry/SensorInfo.h>
14 #include <vxd/geometry/SensorPlane.h>
15 #include <vxd/geometry/GeoCache.h>
16 
17 #include <genfit/DetPlane.h>
18 #include <TVector3.h>
19 #include <TRandom.h>
20 #include <cmath>
21 
22 using namespace std;
23 using namespace Belle2;
24 
25 SVDRecoHit2D::SVDRecoHit2D():
26  genfit::PlanarMeasurement(HIT_DIMENSIONS), m_sensorID(0), m_trueHit(0), m_uCluster(0), m_vCluster(0),
27  m_energyDep(0)//, m_energyDepError(0)
28 {}
29 
30 SVDRecoHit2D::SVDRecoHit2D(const SVDTrueHit* hit, const genfit::TrackCandHit*, float sigmaU, float sigmaV):
31  genfit::PlanarMeasurement(HIT_DIMENSIONS), m_sensorID(0), m_trueHit(hit), m_uCluster(0), m_vCluster(0),
32  m_energyDep(0)//, m_energyDepError(0)
33 {
34  if (!gRandom) B2FATAL("gRandom not initialized, please set up gRandom first");
35 
36  // Set the sensor UID
37  m_sensorID = hit->getSensorID();
38 
39  //If no error is given, estimate the error by dividing the pixel size by sqrt(12)
40  if (sigmaU < 0 || sigmaV < 0) {
41  const SVD::SensorInfo& geometry = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::get(m_sensorID));
42  sigmaU = geometry.getUPitch(hit->getV()) / sqrt(12);
43  sigmaV = geometry.getVPitch(hit->getV()) / sqrt(12);
44  }
45 
46  // Set positions
47  rawHitCoords_(0) = gRandom->Gaus(hit->getU(), sigmaU);
48  rawHitCoords_(1) = gRandom->Gaus(hit->getV(), sigmaV);
49  // Set the error covariance matrix
50  rawHitCov_(0, 0) = sigmaU * sigmaU;
51  rawHitCov_(0, 1) = 0;
52  rawHitCov_(1, 0) = 0;
53  rawHitCov_(1, 1) = sigmaV * sigmaV;
54  // Set physical parameters
55  m_energyDep = hit->getEnergyDep();
56  // Setup geometry information
57  setDetectorPlane();
58 }
59 
60 SVDRecoHit2D::SVDRecoHit2D(VxdID::baseType vxdid, const double u, const double v, double sigmaU, double sigmaV):
61  genfit::PlanarMeasurement(HIT_DIMENSIONS), m_sensorID(vxdid), m_trueHit(NULL), m_uCluster(0), m_vCluster(0),
62  m_energyDep(0)//, m_energyDepError(0)
63 {
64  //If no error is given, estimate the error by dividing the pixel size by sqrt(12)
65  if (sigmaU < 0 || sigmaV < 0) {
66  const SVD::SensorInfo& geometry = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::get(m_sensorID));
67  sigmaU = geometry.getUPitch(v) / sqrt(12);
68  sigmaV = geometry.getVPitch(v) / sqrt(12);
69  }
70  // Set positions
71  rawHitCoords_(0) = u;
72  rawHitCoords_(1) = v;
73  // Set the error covariance matrix
74  rawHitCov_(0, 0) = sigmaU * sigmaU;
75  rawHitCov_(0, 1) = 0;
76  rawHitCov_(1, 0) = 0;
77  rawHitCov_(1, 1) = sigmaV * sigmaV;
78  // Setup geometry information
79  setDetectorPlane();
80 }
81 
82 SVDRecoHit2D::SVDRecoHit2D(const SVDCluster& uHit, const SVDCluster& vHit):
83  genfit::PlanarMeasurement(HIT_DIMENSIONS), m_trueHit(NULL), m_uCluster(&uHit), m_vCluster(&vHit),
84  m_energyDep(0)
85 {
86  if ((uHit.getRawSensorID() != vHit.getRawSensorID()) || !uHit.isUCluster() || vHit.isUCluster())
87  B2FATAL("Error in SVDRecoHit2D: Incorrect SVDCluster instances on input!");
88 
89  m_sensorID = uHit.getRawSensorID();
90 
91  // Now that we have a v coordinate, we can rescale u.
92  const SVD::SensorInfo& info =
93  dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::get(m_sensorID));
94 
95  double DeltaU =
96  (info.getForwardWidth() - info.getBackwardWidth()) / info.getLength() / info.getWidth(0);
97  double scaleFactorU = 1 + DeltaU * vHit.getPosition();
98  double tan_phi = DeltaU * uHit.getPosition(); // need u at v=0!
99  double one_over_cos_phi_sqr = 1 + tan_phi * tan_phi;
100 
101  rawHitCoords_[0] = uHit.getPosition() * scaleFactorU;
102  rawHitCoords_[1] = vHit.getPosition();
103 
104  double sigmaU = uHit.getPositionSigma() * scaleFactorU;
105  double sigmaU_sq = sigmaU * sigmaU;
106  double sigmaV = vHit.getPositionSigma();
107  double sigmaV_sq = sigmaV * sigmaV;
108 
109  m_energyDep = 0.5 * (uHit.getCharge() + vHit.getCharge());
110 
111  rawHitCov_(0, 0) = sigmaV_sq * tan_phi * tan_phi + sigmaU_sq * one_over_cos_phi_sqr;
112  rawHitCov_(0, 1) = sigmaV_sq * tan_phi;
113  rawHitCov_(1, 0) = sigmaV_sq * tan_phi;
114  rawHitCov_(1, 1) = sigmaV_sq;
115  // Setup geometry information
116  setDetectorPlane();
117 }
118 
119 SVDRecoHit2D::SVDRecoHit2D(const SVDRecoHit& uRecoHit, const SVDRecoHit& vRecoHit):
120  genfit::PlanarMeasurement(HIT_DIMENSIONS), m_trueHit(NULL), m_uCluster(0), m_vCluster(0), m_energyDep(0)
121 {
122  const SVDCluster& uHit = *(uRecoHit.getCluster());
123  const SVDCluster& vHit = *(vRecoHit.getCluster());
124  if ((uHit.getRawSensorID() != vHit.getRawSensorID()) || !uHit.isUCluster() || vHit.isUCluster())
125  B2FATAL("Error in SVDRecoHit2D: Incorrect SVDCluster instances on input!");
126 
127  m_sensorID = uHit.getRawSensorID();
128  m_uCluster = &uHit;
129  m_vCluster = &vHit;
130 
131  // Now that we have a v coordinate, we can rescale u.
132  const SVD::SensorInfo& info =
133  dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::get(m_sensorID));
134 
135  double DeltaU =
136  (info.getForwardWidth() - info.getBackwardWidth()) / info.getLength() / info.getWidth(0);
137  double scaleFactorU = 1 + DeltaU * vHit.getPosition();
138  double tan_phi = DeltaU * uHit.getPosition(); // need u at v=0!
139  double one_over_cos_phi_sqr = 1 + tan_phi * tan_phi;
140 
141  rawHitCoords_[0] = uHit.getPosition() * scaleFactorU;
142  rawHitCoords_[1] = vHit.getPosition();
143 
144  double sigmaU = uHit.getPositionSigma() * scaleFactorU;
145  double sigmaU_sq = sigmaU * sigmaU;
146  double sigmaV = vHit.getPositionSigma();
147  double sigmaV_sq = sigmaV * sigmaV;
148 
149  m_energyDep = 0.5 * (uHit.getCharge() + vHit.getCharge());
150 
151  rawHitCov_(0, 0) = sigmaV_sq * tan_phi * tan_phi + sigmaU_sq * one_over_cos_phi_sqr;
152  rawHitCov_(0, 1) = sigmaV_sq * tan_phi;
153  rawHitCov_(1, 0) = sigmaV_sq * tan_phi;
154  rawHitCov_(1, 1) = sigmaV_sq;
155  // Setup geometry information
156  setDetectorPlane();
157 }
158 
159 void SVDRecoHit2D::setDetectorPlane()
160 {
161  // Construct a finite detector plane and set it.
162  const SVD::SensorInfo& geometry = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::get(m_sensorID));
163 
164  // Construct vectors o, u, v
165  TVector3 origin = geometry.pointToGlobal(TVector3(0, 0, 0), true);
166  TVector3 uGlobal = geometry.vectorToGlobal(TVector3(1, 0, 0), true);
167  TVector3 vGlobal = geometry.vectorToGlobal(TVector3(0, 1, 0), true);
168 
169  //Construct the detector plane
170  genfit::SharedPlanePtr detPlane(new genfit::DetPlane(origin, uGlobal, vGlobal, new VXD::SensorPlane(m_sensorID, 20, 20)));
171  setPlane(detPlane, m_sensorID);
172 }
173 
174 TVectorD SVDRecoHit2D::applyPlanarDeformation(TVectorD rawHit, std::vector<double> planarParameters,
175  const genfit::StateOnPlane& state) const
176 {
177  // Legendre parametrization of deformation
178  auto L1 = [](double x) {return x;};
179  auto L2 = [](double x) {return (3 * x * x - 1) / 2;};
180  auto L3 = [](double x) {return (5 * x * x * x - 3 * x) / 2;};
181  auto L4 = [](double x) {return (35 * x * x * x * x - 30 * x * x + 3) / 8;};
182 
183  const SVD::SensorInfo& geometry = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::get(m_sensorID));
184 
185  double u = rawHit[0]; // U coordinate of hit
186  double v = rawHit[1]; // V coordinate of hit
187  double width = geometry.getWidth(v); // Width of sensor (U side)
188  double length = geometry.getLength(); // Length of sensor (V side)
189  u = u * 2 / width; // Legendre parametrization required U in (-1, 1)
190  v = v * 2 / length; // Legendre parametrization required V in (-1, 1)
191 
192  /* Planar deformation using Legendre parametrization
193  w(u, v) = L_{31} * L2(u) + L_{32} * L1(u) * L1(v) + L_{33} * L2(v) +
194  L_{41} * L3(u) + L_{42} * L2(u) * L1(v) + L_{43} * L1(u) * L2(v) + L_{44} * L3(v) +
195  L_{51} * L4(u) + L_{52} * L3(u) * L1(v) + L_{53} * L2(u) * L2(v) + L_{54} * L1(u) * L3(v) + L_{55} * L4(v); */
196  double dw =
197  planarParameters[0] * L2(u) + planarParameters[1] * L1(u) * L1(v) + planarParameters[2] * L2(v) +
198  planarParameters[3] * L3(u) + planarParameters[4] * L2(u) * L1(v) + planarParameters[5] * L1(u) * L2(v) + planarParameters[6] * L3(
199  v) +
200  planarParameters[7] * L4(u) + planarParameters[8] * L3(u) * L1(v) + planarParameters[9] * L2(u) * L2(v) + planarParameters[10] * L1(
201  u) * L3(v) + planarParameters[11] * L4(v);
202 
203  double du_dw = state.getState()[1]; // slope in U direction
204  double dv_dw = state.getState()[2]; // slope in V direction
205 
206  u = u * width / 2; // from Legendre to Local parametrization
207  v = v * length / 2; // from Legendre to Local parametrization
208 
209  TVectorD pos(2);
210 
211  pos[0] = u + dw * du_dw;
212  pos[1] = v + dw * dv_dw;
213 
214  return pos;
215 }
216 
217 std::vector<genfit::MeasurementOnPlane*> SVDRecoHit2D::constructMeasurementsOnPlane(const genfit::StateOnPlane& state) const
218 {
219  // Apply planar deformation
220  TVectorD pos = applyPlanarDeformation(rawHitCoords_, VXD::GeoCache::get(m_sensorID).getSurfaceParameters(), state);
221 
222  return std::vector<genfit::MeasurementOnPlane*>(1, new genfit::MeasurementOnPlane(pos, rawHitCov_, state.getPlane(),
223  state.getRep(), this->constructHMatrix(state.getRep())));
224 }
225 
226 genfit::AbsMeasurement* SVDRecoHit2D::clone() const
227 {
228  return new SVDRecoHit2D(*this);
229 }
Belle2::SVDRecoHit2D::constructMeasurementsOnPlane
virtual std::vector< genfit::MeasurementOnPlane * > constructMeasurementsOnPlane(const genfit::StateOnPlane &state) const override
Get deposited energy error.
Definition: SVDRecoHit2D.cc:217
Belle2::SVDCluster::getCharge
float getCharge() const
Get collected charge.
Definition: SVDCluster.h:152
Belle2::SVDRecoHit2D::m_vCluster
const SVDCluster * m_vCluster
Pointer to mother vCluster.
Definition: SVDRecoHit2D.h:147
Belle2::SVDRecoHit
SVDRecoHit - an extended form of SVDHit containing geometry information.
Definition: SVDRecoHit.h:57
genfit::SharedPlanePtr
std::shared_ptr< genfit::DetPlane > SharedPlanePtr
Shared Pointer to a DetPlane.
Definition: SharedPlanePtr.h:40
Belle2::SVDRecoHit2D::setDetectorPlane
void setDetectorPlane()
Set up Detector plane information.
Definition: SVDRecoHit2D.cc:159
Belle2::VXD::GeoCache::get
static const SensorInfoBase & get(Belle2::VxdID id)
Return a reference to the SensorInfo of a given SensorID.
Definition: GeoCache.h:141
genfit::StateOnPlane
A state with arbitrary dimension defined in a DetPlane.
Definition: StateOnPlane.h:47
genfit
Defines for I/O streams used for error and debug printing.
Definition: AlignablePXDRecoHit.h:19
Belle2::SVDTrueHit
Class SVDTrueHit - Records of tracks that either enter or leave the sensitive volume.
Definition: SVDTrueHit.h:35
genfit::TrackCandHit
Hit object for use in TrackCand.
Definition: TrackCandHit.h:34
Belle2::SVD::SensorInfo
Specific implementation of SensorInfo for SVD Sensors which provides additional sensor specific infor...
Definition: SensorInfo.h:35
Belle2::SVDRecoHit2D::m_uCluster
const SVDCluster * m_uCluster
Pointer to mother uCluster.
Definition: SVDRecoHit2D.h:146
Belle2::SVDCluster::getPositionSigma
float getPositionSigma() const
Get the error of the reconstructed hit coordinate.
Definition: SVDCluster.h:137
genfit::AbsMeasurement
Contains the measurement and covariance in raw detector coordinates.
Definition: AbsMeasurement.h:42
Belle2::VxdID::baseType
unsigned short baseType
The base integer type for VxdID.
Definition: VxdID.h:46
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::SVDRecoHit2D::m_sensorID
unsigned short m_sensorID
Unique sensor identifier.
Definition: SVDRecoHit2D.h:144
Belle2::SVDRecoHit2D::m_energyDep
float m_energyDep
deposited energy.
Definition: SVDRecoHit2D.h:149
Belle2::SVDRecoHit2D::applyPlanarDeformation
TVectorD applyPlanarDeformation(TVectorD rawHit, std::vector< double > planarParameters, const genfit::StateOnPlane &state) const
Apply planar deformation of sensors.
Definition: SVDRecoHit2D.cc:174
Belle2::SVDCluster::isUCluster
bool isUCluster() const
Get the direction of strips.
Definition: SVDCluster.h:118
Belle2::SVDCluster::getPosition
float getPosition(double v=0) const
Get the coordinate of reconstructed hit.
Definition: SVDCluster.h:125
Belle2::SVDCluster::getRawSensorID
unsigned short getRawSensorID() const
Get raw sensor ID.
Definition: SVDCluster.h:113
Belle2::SVDCluster
The SVD Cluster class This class stores all information about reconstructed SVD clusters.
Definition: SVDCluster.h:38
Belle2::SVDRecoHit2D::clone
genfit::AbsMeasurement * clone() const override
Creating a deep copy of this hit.
Definition: SVDRecoHit2D.cc:226
genfit::DetPlane
Detector plane.
Definition: DetPlane.h:59
Belle2::SVDRecoHit2D::SVDRecoHit2D
SVDRecoHit2D()
Default constructor for ROOT IO.
Definition: SVDRecoHit2D.cc:25
Belle2::SVDRecoHit::getCluster
const SVDCluster * getCluster() const
Get pointer to the Cluster used when creating this RecoHit, can be NULL if created from something els...
Definition: SVDRecoHit.h:97
genfit::MeasurementOnPlane
Measured coordinates on a plane.
Definition: MeasurementOnPlane.h:46
Belle2::VXD::SensorPlane
A Finite plane of one VXD Sensor.
Definition: SensorPlane.h:44