Belle II Software  release-05-01-25
AlignableEKLMRecoHit.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2016 Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Kirill Chilikin *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <alignment/reconstruction/AlignableEKLMRecoHit.h>
12 
13 #include <alignment/GlobalDerivatives.h>
14 #include <alignment/GlobalLabel.h>
15 #include <alignment/reconstruction/AlignableEKLMRecoHit.h>
16 #include <klm/dataobjects/eklm/EKLMElementNumbers.h>
17 #include <klm/dataobjects/eklm/EKLMHit2d.h>
18 #include <klm/dataobjects/KLMDigit.h>
19 #include <klm/dataobjects/KLMElementNumbers.h>
20 #include <klm/dbobjects/eklm/EKLMAlignment.h>
21 #include <klm/eklm/geometry/GeometryData.h>
22 #include <klm/eklm/geometry/TransformDataGlobalAligned.h>
23 
24 using namespace Belle2;
25 using namespace alignment;
26 
28 {
29 }
30 
32  const EKLMAlignmentHit* hit, const genfit::TrackCandHit* trackCandHit) :
33  genfit::PlanarMeasurement(1)
34 {
35  (void)trackCandHit;
36  int digit, plane, segment, strip;
37  const HepGeom::Transform3D* t;
38  CLHEP::HepRotation r;
39  CLHEP::Hep3Vector origin;
40  CLHEP::Hep3Vector u(1, 0, 0);
41  CLHEP::Hep3Vector v(0, 1, 0);
42  B2Vector3D origin2, u2, v2;
44  const KLMElementNumbers* elementNumbers = &(KLMElementNumbers::Instance());
45  const EKLMElementNumbers* eklmElementNumbers =
47  const EKLM::TransformDataGlobalAligned* transformData =
49  RelationVector<EKLMHit2d> hit2ds = hit->getRelationsTo<EKLMHit2d>();
50  if (hit2ds.size() != 1)
51  B2FATAL("Incorrect number of related EKLMHit2ds.");
52  RelationVector<KLMDigit> eklmDigits = hit2ds[0]->getRelationsTo<KLMDigit>();
53  if (eklmDigits.size() != 2)
54  B2FATAL("Incorrect number of related KLMDigits.");
55  digit = hit->getDigitIdentifier();
56  m_Section = eklmDigits[digit]->getSection();
57  m_Layer = eklmDigits[digit]->getLayer();
58  m_Sector = eklmDigits[digit]->getSector();
59  plane = eklmDigits[digit]->getPlane();
60  segment = (eklmDigits[digit]->getStrip() - 1) /
61  eklmElementNumbers->getNStripsSegment() + 1;
62  strip = (segment - 1) * eklmElementNumbers->getNStripsSegment() + 1;
64  m_Segment = eklmElementNumbers->segmentNumber(
65  m_Section, m_Layer, m_Sector, plane, segment);
66  t = transformData->getStripTransform(
67  m_Section, m_Layer, m_Sector, plane, strip);
68  origin = t->getTranslation();
69  origin2.SetX(origin.x() / CLHEP::cm * Unit::cm);
70  origin2.SetY(origin.y() / CLHEP::cm * Unit::cm);
71  origin2.SetZ(origin.z() / CLHEP::cm * Unit::cm);
72  r = t->getRotation();
73  u = r * u;
74  v = r * v;
75  u2.SetX(u.x());
76  u2.SetY(u.y());
77  u2.SetZ(u.z());
78  v2.SetX(v.x());
79  v2.SetY(v.y());
80  v2.SetZ(v.z());
81  t = transformData->getSectorTransform(m_Section, m_Layer, m_Sector);
82  r = t->getRotation().inverse();
83  v = r * v;
84  m_StripV.SetX(v.unit().x());
85  m_StripV.SetY(v.unit().y());
86  genfit::SharedPlanePtr detPlane(new genfit::DetPlane(origin2, u2, v2, 0));
87  setPlane(detPlane, m_Segment);
88  rawHitCoords_[0] = geoDat->getStripGeometry()->getWidth() *
89  ((eklmDigits[digit]->getStrip() - 1) %
90  eklmElementNumbers->getNStripsSegment()) / CLHEP::cm * Unit::cm;
91  rawHitCov_[0][0] = pow(geoDat->getStripGeometry()->getWidth() /
92  CLHEP::cm * Unit::cm, 2) / 12;
93  setStripV();
94 }
95 
97 {
98 }
99 
100 std::pair<std::vector<int>, TMatrixD> AlignableEKLMRecoHit::globalDerivatives(const genfit::StateOnPlane* sop)
101 {
102  std::vector<int> labGlobal;
103  labGlobal.push_back(GlobalLabel::construct<EKLMAlignment>(m_KLMModule, KLMAlignmentData::c_DeltaU)); // dx
104  labGlobal.push_back(GlobalLabel::construct<EKLMAlignment>(m_KLMModule, KLMAlignmentData::c_DeltaV)); // dy
105  labGlobal.push_back(GlobalLabel::construct<EKLMAlignment>(m_KLMModule, KLMAlignmentData::c_DeltaGamma)); // drot
106 
107  /* Local parameters. */
108  const double dalpha = 0;
109  const double dxs = 0;
110  const double dys = 0;
111  const double sinda = sin(dalpha);
112  const double cosda = cos(dalpha);
113  /* Local position in sector coordinates. */
114  HepGeom::Point3D<double> globalPos;
115  HepGeom::Transform3D t;
116  const EKLM::TransformDataGlobalAligned* transformData =
118  t = (*transformData->getSectorTransform(m_Section, m_Layer, m_Sector)).inverse();
119  globalPos.setX(sop->getPos().X() / Unit::cm * CLHEP::cm);
120  globalPos.setY(sop->getPos().Y() / Unit::cm * CLHEP::cm);
121  globalPos.setZ(sop->getPos().Z() / Unit::cm * CLHEP::cm);
122  globalPos = t * globalPos;
123  double x = globalPos.x() / CLHEP::cm * Unit::cm;
124  double y = globalPos.y() / CLHEP::cm * Unit::cm;
125  /*
126  * Matrix of global derivatives of local coordinate (second dimension is
127  * added because of resizing in GblFitterInfo::constructGblPoint()) by
128  * variations of the detector element position.
129  * The sign is inverted to match the convention for other detectors.
130  */
131  TMatrixD derGlobal(2, 3);
132  derGlobal(0, 0) = 0;
133  derGlobal(0, 1) = 0;
134  derGlobal(0, 2) = 0;
135  derGlobal(1, 0) = cosda * m_StripV.X() - sinda * m_StripV.Y();
136  derGlobal(1, 1) = sinda * m_StripV.X() + cosda * m_StripV.Y();
137  derGlobal(1, 2) =
138  -((x - dxs) * (-sinda * m_StripV.X() - cosda * m_StripV.Y()) +
139  (y - dys) * (cosda * m_StripV.X() - sinda * m_StripV.Y()));
140  return alignment::GlobalDerivatives(labGlobal, derGlobal);
141 
142 }
143 
145 {
146  return new AlignableEKLMRecoHit(*this);
147 }
148 
Belle2::RelationVector::size
size_t size() const
Get number of relations.
Definition: RelationVector.h:98
Belle2::AlignableEKLMRecoHit::m_Sector
int m_Sector
Sector number.
Definition: AlignableEKLMRecoHit.h:104
Belle2::Unit::cm
static const double cm
Standard units with the value = 1.
Definition: Unit.h:57
genfit::SharedPlanePtr
std::shared_ptr< genfit::DetPlane > SharedPlanePtr
Shared Pointer to a DetPlane.
Definition: SharedPlanePtr.h:40
Belle2::KLMAlignmentData::c_DeltaV
@ c_DeltaV
Shift in V (EKLM: local Y).
Definition: KLMAlignmentData.h:46
Belle2::EKLMElementNumbers
EKLM element numbers.
Definition: EKLMElementNumbers.h:34
Belle2::EKLM::TransformData::getSectorTransform
const HepGeom::Transform3D * getSectorTransform(int section, int layer, int sector) const
Get sector transformation.
Definition: TransformData.cc:298
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
genfit::TrackCandHit
Hit object for use in TrackCand.
Definition: TrackCandHit.h:34
Belle2::KLMElementNumbers::Instance
static const KLMElementNumbers & Instance()
Instantiation.
Definition: KLMElementNumbers.cc:31
Belle2::AlignableEKLMRecoHit::m_Section
int m_Section
Section number.
Definition: AlignableEKLMRecoHit.h:101
Belle2::AlignableEKLMRecoHit::m_Layer
int m_Layer
Layer number.
Definition: AlignableEKLMRecoHit.h:107
Belle2::EKLMAlignmentHit
This dataobject is used only for EKLM alignment.
Definition: EKLMAlignmentHit.h:39
Belle2::AlignableEKLMRecoHit::m_Segment
uint16_t m_Segment
Segment number.
Definition: AlignableEKLMRecoHit.h:113
Belle2::B2Vector3::SetY
void SetY(DataType y)
set Y/2nd-coordinate
Definition: B2Vector3.h:454
Belle2::KLMDigit
KLM digit (class representing a digitized hit in RPCs or scintillators).
Definition: KLMDigit.h:40
Belle2::B2Vector3< double >
Belle2::EKLM::TransformDataGlobalAligned::Instance
static const TransformDataGlobalAligned & Instance()
Instantiation.
Definition: TransformDataGlobalAligned.cc:17
genfit::AbsMeasurement
Contains the measurement and covariance in raw detector coordinates.
Definition: AbsMeasurement.h:42
Belle2::RelationVector
Class for type safe access to objects that are referred to in relations.
Definition: DataStore.h:38
Belle2::AlignableEKLMRecoHit::m_KLMModule
uint16_t m_KLMModule
KLM module number.
Definition: AlignableEKLMRecoHit.h:110
genfit::PlanarMeasurement::setStripV
void setStripV(bool v=true)
Use if the coordinate for 1D hits measured in V direction.
Definition: PlanarMeasurement.h:70
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::EKLM::TransformData::getStripTransform
const HepGeom::Transform3D * getStripTransform(int section, int layer, int sector, int plane, int strip) const
Get strip transformation.
Definition: TransformData.cc:323
Belle2::EKLMElementNumbers::segmentNumber
int segmentNumber(int section, int layer, int sector, int plane, int segment) const
Get segment number.
Definition: EKLMElementNumbers.cc:191
Belle2::EKLMGeometry::getStripGeometry
const StripGeometry * getStripGeometry() const
Get strip geometry data.
Definition: EKLMGeometry.h:1873
Belle2::AlignableEKLMRecoHit::AlignableEKLMRecoHit
AlignableEKLMRecoHit()
Constructor.
Definition: AlignableEKLMRecoHit.cc:27
Belle2::EKLMElementNumbers::Instance
static const EKLMElementNumbers & Instance()
Instantiation.
Definition: EKLMElementNumbers.cc:20
Belle2::KLMElementNumbers::moduleNumberEKLM
uint16_t moduleNumberEKLM(int section, int sector, int layer) const
Get module number for EKLM.
Definition: KLMElementNumbers.cc:166
Belle2::EKLM::GeometryData
EKLM geometry data.
Definition: GeometryData.h:40
Belle2::AlignableEKLMRecoHit::~AlignableEKLMRecoHit
~AlignableEKLMRecoHit()
Destructor.
Definition: AlignableEKLMRecoHit.cc:96
Belle2::alignment::GlobalDerivatives
Class for easier manipulation with global derivatives (and their labels)
Definition: GlobalDerivatives.h:37
Belle2::AlignableEKLMRecoHit::clone
genfit::AbsMeasurement * clone() const override
Clone.
Definition: AlignableEKLMRecoHit.cc:144
Belle2::EKLMGeometry::StripGeometry::getWidth
double getWidth() const
Get width.
Definition: EKLMGeometry.h:1199
Belle2::EKLMElementNumbers::getNStripsSegment
static constexpr int getNStripsSegment()
Get number of strips in a segment.
Definition: EKLMElementNumbers.h:401
genfit::DetPlane
Detector plane.
Definition: DetPlane.h:59
HepGeom::Point3D< double >
Belle2::KLMElementNumbers
KLM element numbers.
Definition: KLMElementNumbers.h:37
Belle2::B2Vector3::X
DataType X() const
access variable X (= .at(0) without boundary check)
Definition: B2Vector3.h:430
Belle2::KLMAlignmentData::c_DeltaGamma
@ c_DeltaGamma
Rotation in gamma (EKLM: rotation in local plane).
Definition: KLMAlignmentData.h:58
Belle2::AlignableEKLMRecoHit::m_StripV
B2Vector3D m_StripV
V direction.
Definition: AlignableEKLMRecoHit.h:116
Belle2::AlignableEKLMRecoHit::globalDerivatives
virtual std::pair< std::vector< int >, TMatrixD > globalDerivatives(const genfit::StateOnPlane *sop) override
Labels and derivatives of residuals (local measurement coordinates) w.r.t.
Definition: AlignableEKLMRecoHit.cc:100
Belle2::B2Vector3::SetZ
void SetZ(DataType z)
set Z/3rd-coordinate
Definition: B2Vector3.h:456
Belle2::B2Vector3::Y
DataType Y() const
access variable Y (= .at(1) without boundary check)
Definition: B2Vector3.h:432
Belle2::EKLM::GeometryData::Instance
static const GeometryData & Instance(enum DataSource dataSource=c_Database, const GearDir *gearDir=nullptr)
Instantiation.
Definition: GeometryData.cc:35
Belle2::KLMAlignmentData::c_DeltaU
@ c_DeltaU
Shift in U (EKLM: local X).
Definition: KLMAlignmentData.h:43
Belle2::EKLMHit2d
Class for 2d hits handling.
Definition: EKLMHit2d.h:39
Belle2::B2Vector3::SetX
void SetX(DataType x)
set X/1st-coordinate
Definition: B2Vector3.h:452
Belle2::EKLM::TransformDataGlobalAligned
Transformation data (global, aligned): singleton version.
Definition: TransformDataGlobalAligned.h:35