Belle II Software  release-05-02-19
TransformData.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2012 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 /* Own header. */
12 #include <klm/eklm/geometry/TransformData.h>
13 
14 /* KLM headers. */
15 #include <klm/dbobjects/eklm/EKLMAlignment.h>
16 #include <klm/dbobjects/eklm/EKLMSegmentAlignment.h>
17 #include <klm/eklm/geometry/AlignmentChecker.h>
18 #include <klm/eklm/geometry/GeometryData.h>
19 
20 /* Belle 2 headers. */
21 #include <framework/database/DBObjPtr.h>
22 #include <framework/logging/Logger.h>
23 #include <framework/gearbox/Unit.h>
24 
25 using namespace Belle2;
26 
27 EKLM::TransformData::TransformData(bool global, Displacement displacementType)
28 {
29  /* cppcheck-suppress variableScope */
30  int iSection, iLayer, iSector, iPlane, iSegment, iStrip, sector, segment;
31  int nSections, nLayers, nSectors, nPlanes, nStrips, nSegments, nStripsSegment;
32  int nDetectorLayers;
33  AlignmentChecker alignmentChecker(true);
34  m_ElementNumbers = &(EKLMElementNumbers::Instance());
35  m_GeoDat = &(GeometryData::Instance());
36  nSections = m_GeoDat->getNSections();
37  nSectors = m_GeoDat->getNSectors();
38  nLayers = m_GeoDat->getNLayers();
39  nPlanes = m_GeoDat->getNPlanes();
40  nStrips = m_GeoDat->getNStrips();
41  nSegments = m_GeoDat->getNSegments();
42  nStripsSegment = m_ElementNumbers->getNStripsSegment();
43  m_Section = new HepGeom::Transform3D[nSections];
44  m_Layer = new HepGeom::Transform3D*[nSections];
45  m_Sector = new HepGeom::Transform3D** [nSections];
46  m_Plane = new HepGeom::Transform3D** *[nSections];
47  m_PlaneDisplacement = new HepGeom::Transform3D** *[nSections];
48  m_Segment = new HepGeom::Transform3D**** [nSections];
49  m_Strip = new HepGeom::Transform3D**** [nSections];
50  m_StripInverse = new HepGeom::Transform3D**** [nSections];
51  for (iSection = 0; iSection < nSections; iSection++) {
52  m_GeoDat->getSectionTransform(&m_Section[iSection], iSection);
53  nDetectorLayers = m_GeoDat->getNDetectorLayers(iSection + 1);
54  m_Layer[iSection] = new HepGeom::Transform3D[nLayers];
55  m_Sector[iSection] = new HepGeom::Transform3D*[nLayers];
56  m_Plane[iSection] = new HepGeom::Transform3D** [nLayers];
57  m_PlaneDisplacement[iSection] = new HepGeom::Transform3D** [nLayers];
58  m_Segment[iSection] = new HepGeom::Transform3D** *[nLayers];
59  m_Strip[iSection] = new HepGeom::Transform3D** *[nLayers];
60  m_StripInverse[iSection] = new HepGeom::Transform3D** *[nLayers];
61  for (iLayer = 0; iLayer < nLayers; iLayer++) {
62  m_GeoDat->getLayerTransform(&m_Layer[iSection][iLayer], iLayer);
63  m_Sector[iSection][iLayer] = new HepGeom::Transform3D[nSectors];
64  if (iLayer < nDetectorLayers) {
65  m_Plane[iSection][iLayer] = new HepGeom::Transform3D*[nSectors];
66  m_PlaneDisplacement[iSection][iLayer] =
67  new HepGeom::Transform3D*[nSectors];
68  m_Segment[iSection][iLayer] = new HepGeom::Transform3D** [nSectors];
69  m_Strip[iSection][iLayer] = new HepGeom::Transform3D** [nSectors];
70  m_StripInverse[iSection][iLayer] = new HepGeom::Transform3D** [nSectors];
71  }
72  for (iSector = 0; iSector < nSectors; iSector++) {
73  m_GeoDat->getSectorTransform(&m_Sector[iSection][iLayer][iSector],
74  iSector);
75  if (iLayer >= nDetectorLayers)
76  continue;
77  m_Plane[iSection][iLayer][iSector] = new HepGeom::Transform3D[nPlanes];
78  m_PlaneDisplacement[iSection][iLayer][iSector] =
79  new HepGeom::Transform3D[nPlanes];
80  m_Segment[iSection][iLayer][iSector] =
81  new HepGeom::Transform3D*[nPlanes];
82  m_Strip[iSection][iLayer][iSector] = new HepGeom::Transform3D*[nPlanes];
83  m_StripInverse[iSection][iLayer][iSector] =
84  new HepGeom::Transform3D*[nPlanes];
85  for (iPlane = 0; iPlane < nPlanes; iPlane++) {
86  m_GeoDat->getPlaneTransform(
87  &m_Plane[iSection][iLayer][iSector][iPlane], iPlane);
88  m_PlaneDisplacement[iSection][iLayer][iSector][iPlane] =
89  HepGeom::Translate3D(0, 0, 0);
90  m_Segment[iSection][iLayer][iSector][iPlane] =
91  new HepGeom::Transform3D[nSegments];
92  for (iSegment = 0; iSegment < nSegments; iSegment++) {
93  m_Segment[iSection][iLayer][iSector][iPlane][iSegment] =
94  HepGeom::Translate3D(0, 0, 0);
95  }
96  m_Strip[iSection][iLayer][iSector][iPlane] =
97  new HepGeom::Transform3D[nStrips];
98  m_StripInverse[iSection][iLayer][iSector][iPlane] =
99  new HepGeom::Transform3D[nStrips];
100  for (iStrip = 0; iStrip < nStrips; iStrip++) {
101  m_GeoDat->getStripTransform(
102  &m_Strip[iSection][iLayer][iSector][iPlane][iStrip], iStrip);
103  }
104  }
105  }
106  }
107  }
108  /* Read alignment data from the database and modify transformations. */
109  if (displacementType != c_None) {
110  std::string payload, segmentPayload;
111  if (displacementType == c_Displacement) {
112  payload = "EKLMDisplacement";
113  segmentPayload = "EKLMSegmentDisplacement";
114  } else {
115  payload = "EKLMAlignment";
116  segmentPayload = "EKLMSegmentAlignment";
117  }
118  DBObjPtr<EKLMAlignment> alignment(payload);
119  DBObjPtr<EKLMSegmentAlignment> segmentAlignment(segmentPayload);
120  if (!alignment.isValid())
121  B2FATAL("No EKLM displacement (alignment) data.");
122  if (displacementType == c_Displacement) {
123  if (!alignmentChecker.checkAlignment(&(*alignment), &(*segmentAlignment)))
124  B2FATAL("EKLM displacement data are incorrect, overlaps exist.");
125  }
126  for (iSection = 1; iSection <= nSections; iSection++) {
127  nDetectorLayers = m_GeoDat->getNDetectorLayers(iSection);
128  for (iLayer = 1; iLayer <= nDetectorLayers; iLayer++) {
129  for (iSector = 1; iSector <= nSectors; iSector++) {
130  sector = m_ElementNumbers->sectorNumber(iSection, iLayer, iSector);
131  const KLMAlignmentData* sectorAlignment =
132  alignment->getModuleAlignment(sector);
133  if (sectorAlignment == nullptr)
134  B2FATAL("Incomplete EKLM displacement (alignment) data.");
135  for (iPlane = 1; iPlane <= nPlanes; iPlane++) {
136  /* First plane is rotated. */
137  if (iPlane == 1) {
138  m_PlaneDisplacement[iSection - 1][iLayer - 1][iSector - 1][iPlane - 1] =
139  HepGeom::Translate3D(
140  sectorAlignment->getDeltaV() * CLHEP::cm / Unit::cm,
141  sectorAlignment->getDeltaU() * CLHEP::cm / Unit::cm, 0) *
142  HepGeom::RotateZ3D(-sectorAlignment->getDeltaGamma() *
143  CLHEP::rad / Unit::rad);
144  } else {
145  m_PlaneDisplacement[iSection - 1][iLayer - 1][iSector - 1][iPlane - 1] =
146  HepGeom::Translate3D(
147  sectorAlignment->getDeltaU() * CLHEP::cm / Unit::cm,
148  sectorAlignment->getDeltaV() * CLHEP::cm / Unit::cm, 0) *
149  HepGeom::RotateZ3D(sectorAlignment->getDeltaGamma() *
150  CLHEP::rad / Unit::rad);
151  }
152  for (iSegment = 1; iSegment <= nSegments; iSegment++) {
153  segment = m_ElementNumbers->segmentNumber(
154  iSection, iLayer, iSector, iPlane, iSegment);
155  const KLMAlignmentData* segmentAlignmentData =
156  segmentAlignment->getSegmentAlignment(segment);
157  if (segmentAlignmentData == nullptr)
158  B2FATAL("Incomplete EKLM displacement (alignment) data.");
159  m_Segment[iSection - 1][iLayer - 1][iSector - 1][iPlane - 1]
160  [iSegment - 1] =
161  HepGeom::Translate3D(
162  segmentAlignmentData->getDeltaU() * CLHEP::cm / Unit::cm,
163  segmentAlignmentData->getDeltaV() * CLHEP::cm / Unit::cm, 0) *
164  m_Segment[iSection - 1][iLayer - 1][iSector - 1][iPlane - 1]
165  [iSegment - 1] *
166  HepGeom::RotateZ3D(segmentAlignmentData->getDeltaGamma() *
167  CLHEP::rad / Unit::rad);
168  for (iStrip = 1; iStrip <= nStripsSegment; iStrip++) {
169  m_Strip[iSection - 1][iLayer - 1][iSector - 1][iPlane - 1]
170  [nStripsSegment * (iSegment - 1) + iStrip - 1] =
171  HepGeom::Translate3D(
172  segmentAlignmentData->getDeltaU() * CLHEP::cm / Unit::cm,
173  segmentAlignmentData->getDeltaV() * CLHEP::cm / Unit::cm,
174  0) *
175  m_Strip[iSection - 1][iLayer - 1][iSector - 1][iPlane - 1]
176  [nStripsSegment * (iSegment - 1) + iStrip - 1] *
177  HepGeom::RotateZ3D(segmentAlignmentData->getDeltaGamma() *
178  CLHEP::rad / Unit::rad);
179  }
180  }
181  }
182  }
183  }
184  }
185  }
186  if (global)
187  transformsToGlobal();
188 }
189 
190 EKLM::TransformData::~TransformData()
191 {
192  int iSection, iLayer, iSector, iPlane;
193  /* cppcheck-suppress variableScope */
194  int nSections, nLayers, nDetectorLayers, nSectors, nPlanes;
195  nSections = m_GeoDat->getNSections();
196  nLayers = m_GeoDat->getNLayers();
197  nSectors = m_GeoDat->getNSectors();
198  nPlanes = m_GeoDat->getNPlanes();
199  for (iSection = 0; iSection < nSections; iSection++) {
200  nDetectorLayers = m_GeoDat->getNDetectorLayers(iSection + 1);
201  for (iLayer = 0; iLayer < nLayers; iLayer++) {
202  delete[] m_Sector[iSection][iLayer];
203  if (iLayer >= nDetectorLayers)
204  continue;
205  for (iSector = 0; iSector < nSectors; iSector++) {
206  for (iPlane = 0; iPlane < nPlanes; iPlane++) {
207  delete[] m_Segment[iSection][iLayer][iSector][iPlane];
208  delete[] m_Strip[iSection][iLayer][iSector][iPlane];
209  delete[] m_StripInverse[iSection][iLayer][iSector][iPlane];
210  }
211  delete[] m_Plane[iSection][iLayer][iSector];
212  delete[] m_PlaneDisplacement[iSection][iLayer][iSector];
213  delete[] m_Segment[iSection][iLayer][iSector];
214  delete[] m_Strip[iSection][iLayer][iSector];
215  delete[] m_StripInverse[iSection][iLayer][iSector];
216  }
217  delete[] m_Plane[iSection][iLayer];
218  delete[] m_PlaneDisplacement[iSection][iLayer];
219  delete[] m_Segment[iSection][iLayer];
220  delete[] m_Strip[iSection][iLayer];
221  delete[] m_StripInverse[iSection][iLayer];
222  }
223  delete[] m_Layer[iSection];
224  delete[] m_Sector[iSection];
225  delete[] m_PlaneDisplacement[iSection];
226  delete[] m_Plane[iSection];
227  delete[] m_Segment[iSection];
228  delete[] m_Strip[iSection];
229  delete[] m_StripInverse[iSection];
230  }
231  delete[] m_Section;
232  delete[] m_Layer;
233  delete[] m_Sector;
234  delete[] m_PlaneDisplacement;
235  delete[] m_Plane;
236  delete[] m_Segment;
237  delete[] m_Strip;
238  delete[] m_StripInverse;
239 }
240 
241 void EKLM::TransformData::transformsToGlobal()
242 {
243  int iSection, iLayer, iSector, iPlane, iSegment, iStrip;
244  /* cppcheck-suppress variableScope */
245  int nSections, nLayers, nDetectorLayers, nSectors, nPlanes, nSegments, nStrips;
246  nSections = m_GeoDat->getNSections();
247  nLayers = m_GeoDat->getNLayers();
248  nSectors = m_GeoDat->getNSectors();
249  nPlanes = m_GeoDat->getNPlanes();
250  nSegments = m_GeoDat->getNSegments();
251  nStrips = m_GeoDat->getNStrips();
252  for (iSection = 0; iSection < nSections; iSection++) {
253  nDetectorLayers = m_GeoDat->getNDetectorLayers(iSection + 1);
254  for (iLayer = 0; iLayer < nLayers; iLayer++) {
255  m_Layer[iSection][iLayer] = m_Section[iSection] * m_Layer[iSection][iLayer];
256  for (iSector = 0; iSector < nSectors; iSector++) {
257  m_Sector[iSection][iLayer][iSector] =
258  m_Layer[iSection][iLayer] * m_Sector[iSection][iLayer][iSector];
259  if (iLayer >= nDetectorLayers)
260  continue;
261  for (iPlane = 0; iPlane < nPlanes; iPlane++) {
262  m_Plane[iSection][iLayer][iSector][iPlane] =
263  m_Sector[iSection][iLayer][iSector] *
264  m_Plane[iSection][iLayer][iSector][iPlane];
265  for (iSegment = 0; iSegment < nSegments; iSegment++) {
266  m_Segment[iSection][iLayer][iSector][iPlane][iSegment] =
267  m_Plane[iSection][iLayer][iSector][iPlane] *
268  m_PlaneDisplacement[iSection][iLayer][iSector][iPlane] *
269  m_Segment[iSection][iLayer][iSector][iPlane][iSegment];
270  }
271  for (iStrip = 0; iStrip < nStrips; iStrip++) {
272  m_Strip[iSection][iLayer][iSector][iPlane][iStrip] =
273  m_Plane[iSection][iLayer][iSector][iPlane] *
274  m_PlaneDisplacement[iSection][iLayer][iSector][iPlane] *
275  m_Strip[iSection][iLayer][iSector][iPlane][iStrip];
276  m_StripInverse[iSection][iLayer][iSector][iPlane][iStrip] =
277  m_Strip[iSection][iLayer][iSector][iPlane][iStrip].inverse();
278  }
279  }
280  }
281  }
282  }
283 }
284 
285 const HepGeom::Transform3D*
286 EKLM::TransformData::getSectionTransform(int section) const
287 {
288  return &m_Section[section - 1];
289 }
290 
291 const HepGeom::Transform3D*
292 EKLM::TransformData::getLayerTransform(int section, int layer) const
293 {
294  return &m_Layer[section - 1][layer - 1];
295 }
296 
297 const HepGeom::Transform3D* EKLM::TransformData::
298 getSectorTransform(int section, int layer, int sector) const
299 {
300  return &m_Sector[section - 1][layer - 1][sector - 1];
301 }
302 
303 const HepGeom::Transform3D* EKLM::TransformData::
304 getPlaneTransform(int section, int layer, int sector, int plane) const
305 {
306  return &m_Plane[section - 1][layer - 1][sector - 1][plane - 1];
307 }
308 
309 const HepGeom::Transform3D* EKLM::TransformData::
310 getPlaneDisplacement(int section, int layer, int sector, int plane) const
311 {
312  return &m_PlaneDisplacement[section - 1][layer - 1][sector - 1][plane - 1];
313 }
314 
315 const HepGeom::Transform3D* EKLM::TransformData::
316 getSegmentTransform(int section, int layer, int sector, int plane,
317  int segment) const
318 {
319  return &m_Segment[section - 1][layer - 1][sector - 1][plane - 1][segment - 1];
320 }
321 
322 const HepGeom::Transform3D* EKLM::TransformData::
323 getStripTransform(int section, int layer, int sector, int plane, int strip) const
324 {
325  return &m_Strip[section - 1][layer - 1][sector - 1][plane - 1][strip - 1];
326 }
327 
328 const HepGeom::Transform3D*
329 EKLM::TransformData::getStripLocalToGlobal(KLMDigit* hit) const
330 {
331  return &(m_Strip[hit->getSection() - 1][hit->getLayer() - 1]
332  [hit->getSector() - 1][hit->getPlane() - 1][hit->getStrip() - 1]);
333 }
334 
335 const HepGeom::Transform3D*
336 EKLM::TransformData::getStripGlobalToLocal(KLMDigit* hit) const
337 {
338  return &(m_StripInverse[hit->getSection() - 1][hit->getLayer() - 1]
339  [hit->getSector() - 1][hit->getPlane() - 1][hit->getStrip() - 1]);
340 }
341 
342 const HepGeom::Transform3D*
343 EKLM::TransformData::getStripGlobalToLocal(int section, int layer, int sector,
344  int plane, int strip) const
345 {
346  return &(m_StripInverse[section - 1][layer - 1][sector - 1][plane - 1]
347  [strip - 1]);
348 }
349 
350 bool EKLM::TransformData::intersection(KLMDigit* hit1, KLMDigit* hit2,
351  HepGeom::Point3D<double>* cross,
352  double* d1, double* d2, double* sd,
353  bool segments) const
354 {
355  /* Hits must be from the same sector, */
356  if (hit1->getSection() != hit2->getSection())
357  return false;
358  if (hit1->getLayer() != hit2->getLayer())
359  return false;
360  if (hit1->getSector() != hit2->getSector())
361  return false;
362  /* but different planes. */
363  if (hit1->getPlane() == hit2->getPlane())
364  return false;
365  /* Coordinates of strip 1 ends. */
366  double l1 = m_GeoDat->getStripLength(hit1->getStrip());
367  HepGeom::Point3D<double> s1_1(-0.5 * l1, 0.0, 0.0);
368  HepGeom::Point3D<double> s1_2(0.5 * l1, 0.0, 0.0);
369  const HepGeom::Transform3D* tr1 = getStripLocalToGlobal(hit1);
370  HepGeom::Point3D<double> s1_1g = (*tr1) * s1_1;
371  HepGeom::Point3D<double> s1_2g = (*tr1) * s1_2;
372  /* Coordinates of strip 2 ends. */
373  double l2 = m_GeoDat->getStripLength(hit2->getStrip());
374  HepGeom::Point3D<double> s2_1(-0.5 * l2, 0.0, 0.0);
375  HepGeom::Point3D<double> s2_2(0.5 * l2, 0.0, 0.0);
376  const HepGeom::Transform3D* tr2 = getStripLocalToGlobal(hit2);
377  HepGeom::Point3D<double> s2_1g = (*tr2) * s2_1;
378  HepGeom::Point3D<double> s2_2g = (*tr2) * s2_2;
388  HepGeom::Vector3D<double> v1 = s1_2g - s1_1g;
389  HepGeom::Vector3D<double> v2 = s2_2g - s2_1g;
390  HepGeom::Vector3D<double> d = s1_1g - s2_1g;
391  double v1sq = v1.mag2();
392  double v2sq = v2.mag2();
393  double v1dv2 = v1.dot(v2);
394  double ddv1 = d.dot(v1);
395  double ddv2 = d.dot(v2);
396  double den = v1sq * v2sq - v1dv2 * v1dv2;
397  double t1 = (v1dv2 * ddv2 - v2sq * ddv1) / den;
398  double t2 = (- v1dv2 * ddv1 + v1sq * ddv2) / den;
399  /* Segments do not intersect. */
400  if (segments) {
401  if (t1 < 0.0 || t1 > 1.0)
402  return false;
403  if (t2 < 0.0 || t2 > 1.0)
404  return false;
405  }
406  /* Segments intersect, set return values. */
407  HepGeom::Point3D<double> s1_cg = s1_1g + v1 * t1;
408  HepGeom::Point3D<double> s2_cg = s2_1g + v2 * t2;
409  *d1 = s1_2g.distance(s1_cg) / CLHEP::mm * Unit::mm;
410  *d2 = s2_2g.distance(s2_cg) / CLHEP::mm * Unit::mm;
411  *cross = 0.5 * (s1_cg + s2_cg) / CLHEP::mm * Unit::mm;
412  *sd = s1_cg.distance(s2_cg) / CLHEP::mm * Unit::mm;
413  if (s2_cg.mag2() < s1_cg.mag2())
414  *sd = - *sd;
415  return true;
416 }
417 
418 int EKLM::TransformData::getSectorByPosition(
419  int section, const HepGeom::Point3D<double>& position) const
420 {
421  int sector;
422  double x;
423  if (section == 1) {
424  x = position.x();
425  } else {
426  x = -position.x();
427  }
428  if (position.y() > 0) {
429  if (x > 0)
430  sector = 1;
431  else
432  sector = 2;
433  } else {
434  if (x > 0)
435  sector = 4;
436  else
437  sector = 3;
438  }
439  return sector;
440 }
441 
442 int EKLM::TransformData::getStripsByIntersection(
443  const HepGeom::Point3D<double>& intersection, int* strip1, int* strip2) const
444 {
445  /* cppcheck-suppress variableScope */
446  int section, layer, sector, plane, segment, strip, stripSegment, stripGlobal;
447  /* cppcheck-suppress variableScope */
448  int nLayers, nPlanes, nSegments, nStripsSegment, minDistanceSegment;
449  double solenoidCenter, firstLayerCenter, layerShift;
450  /* cppcheck-suppress variableScope */
451  double x, y, z, l, minY, maxY;
452  double minDistance = 0, minDistanceNew, stripWidth;
453  HepGeom::Point3D<double> intersectionClhep, intersectionLocal;
454  intersectionClhep = intersection * CLHEP::cm / Unit::cm;
455  solenoidCenter = m_GeoDat->getSolenoidZ() / CLHEP::cm * Unit::cm;
456  if (intersection.z() < solenoidCenter)
457  section = 1;
458  else
459  section = 2;
460  firstLayerCenter =
461  (m_GeoDat->getSectionPosition()->getZ()
462  - 0.5 * m_GeoDat->getSectionPosition()->getLength()
463  + m_GeoDat->getLayerShiftZ()
464  - 0.5 * m_GeoDat->getLayerPosition()->getLength()) /
465  CLHEP::cm * Unit::cm;
466  layerShift = m_GeoDat->getLayerShiftZ() / CLHEP::cm * Unit::cm;
467  z = fabs(intersection.z() - solenoidCenter);
468  layer = round((z - firstLayerCenter) / layerShift) + 1;
469  if (layer <= 0)
470  layer = 1;
471  nLayers = m_GeoDat->getNDetectorLayers(section);
472  if (layer > nLayers)
473  layer = nLayers;
474  sector = getSectorByPosition(section, intersection);
475  nPlanes = m_GeoDat->getNPlanes();
476  nSegments = m_GeoDat->getNSegments();
477  nStripsSegment = m_ElementNumbers->getNStripsSegment();
478  stripWidth = m_GeoDat->getStripGeometry()->getWidth() / CLHEP::cm * Unit::cm;
479  minY = -stripWidth / 2;
480  maxY = (double(nStripsSegment) - 0.5) * stripWidth;
481  for (plane = 1; plane <= nPlanes; plane++) {
482  minDistanceSegment = 1;
483  for (segment = 1; segment <= nSegments; segment++) {
484  strip = (segment - 1) * nStripsSegment;
485  intersectionLocal = m_StripInverse[section - 1][layer - 1]
486  [sector - 1][plane - 1][strip] * intersectionClhep;
487  y = intersectionLocal.y() / CLHEP::cm * Unit::cm;
488  if (y < minY) {
489  minDistanceNew = minY - y;
490  } else if (y > maxY) {
491  minDistanceNew = y - maxY;
492  } else {
493  minDistance = 0;
494  minDistanceSegment = segment;
495  break;
496  }
497  if (segment == 1) {
498  minDistance = minDistanceNew;
499  } else if (minDistanceNew < minDistance) {
500  minDistance = minDistanceNew;
501  minDistanceSegment = segment;
502  }
503  }
504  /*
505  * The intersection is required to be strictly within a segment,
506  * this condition might be adjusted later.
507  */
508  if (minDistance > 0)
509  return -1;
510  strip = (minDistanceSegment - 1) * nStripsSegment;
511  intersectionLocal = m_StripInverse[section - 1][layer - 1]
512  [sector - 1][plane - 1][strip] * intersectionClhep;
513  y = intersectionLocal.y() / CLHEP::cm * Unit::cm;
514  stripSegment = ceil((y - 0.5 * stripWidth) / stripWidth);
515  if (stripSegment <= 0)
516  stripSegment = 1;
517  else if (stripSegment > nStripsSegment)
518  stripSegment = nStripsSegment;
519  strip = stripSegment + (minDistanceSegment - 1) * nStripsSegment;
520  intersectionLocal = m_StripInverse[section - 1][layer - 1]
521  [sector - 1][plane - 1][strip - 1] * intersectionClhep;
522  x = intersectionLocal.x();
523  l = m_GeoDat->getStripLength(strip);
524  /*
525  * The intersection is required to be strictly within the strip length,
526  * this condition might be adjusted later.
527  */
528  if (fabs(x) > 0.5 * l)
529  return -1;
530  stripGlobal = m_ElementNumbers->stripNumber(
531  section, layer, sector, plane, strip);
532  if (plane == 1)
533  *strip1 = stripGlobal;
534  else
535  *strip2 = stripGlobal;
536  }
537  return 0;
538 }
Belle2::Unit::cm
static const double cm
Standard units with the value = 1.
Definition: Unit.h:57
Belle2::KLMDigit::getPlane
int getPlane() const
Get plane number.
Definition: KLMDigit.h:161
Belle2::EKLMGeometry::getNSegments
int getNSegments() const
Get number of segments.
Definition: EKLMGeometry.h:1725
Belle2::EKLMGeometry::getNLayers
int getNLayers() const
Get number of layers.
Definition: EKLMGeometry.h:1695
Belle2::EKLM::TransformData::getSectorByPosition
int getSectorByPosition(int section, const HepGeom::Point3D< double > &position) const
Get sector by position.
Definition: TransformData.cc:418
Belle2::EKLM::TransformData::m_Segment
HepGeom::Transform3D ***** m_Segment
Segment transformations.
Definition: TransformData.h:235
Belle2::EKLM::TransformData::getSectorTransform
const HepGeom::Transform3D * getSectorTransform(int section, int layer, int sector) const
Get sector transformation.
Definition: TransformData.cc:298
Belle2::EKLM::TransformData::getStripGlobalToLocal
const HepGeom::Transform3D * getStripGlobalToLocal(KLMDigit *hit) const
Get strip global to local transformation by hit.
Definition: TransformData.cc:336
Belle2::EKLMGeometry::getNPlanes
int getNPlanes() const
Get number of planes.
Definition: EKLMGeometry.h:1717
Belle2::EKLM::TransformData::getPlaneDisplacement
const HepGeom::Transform3D * getPlaneDisplacement(int section, int layer, int sector, int plane) const
Get additional displacement for plane internal volumes.
Definition: TransformData.cc:310
Belle2::EKLM::GeometryData::getLayerTransform
void getLayerTransform(HepGeom::Transform3D *t, int n) const
Get layer transformation.
Definition: GeometryData.cc:662
Belle2::KLMAlignmentData::getDeltaV
float getDeltaV() const
Get shift in V.
Definition: KLMAlignmentData.h:104
Belle2::EKLM::TransformData::getSegmentTransform
const HepGeom::Transform3D * getSegmentTransform(int section, int layer, int sector, int plane, int segment) const
Get segment transformation.
Definition: TransformData.cc:316
Belle2::KLMDigit::getSection
int getSection() const
Get section number.
Definition: KLMDigit.h:107
Belle2::EKLM::TransformData::getStripLocalToGlobal
const HepGeom::Transform3D * getStripLocalToGlobal(KLMDigit *hit) const
Get strip local to global transformation by hit.
Definition: TransformData.cc:329
Belle2::EKLM::TransformData::m_PlaneDisplacement
HepGeom::Transform3D **** m_PlaneDisplacement
Plane internal volumes displacements.
Definition: TransformData.h:232
Belle2::KLMDigit::getSector
int getSector() const
Get sector number.
Definition: KLMDigit.h:125
Belle2::KLMDigit::getLayer
int getLayer() const
Get layer number.
Definition: KLMDigit.h:143
Belle2::EKLM::TransformData::getPlaneTransform
const HepGeom::Transform3D * getPlaneTransform(int section, int layer, int sector, int plane) const
Get plane transformation.
Definition: TransformData.cc:304
Belle2::KLMDigit
KLM digit (class representing a digitized hit in RPCs or scintillators).
Definition: KLMDigit.h:40
Belle2::EKLMElementNumbers::sectorNumber
int sectorNumber(int section, int layer, int sector) const
Get sector number.
Definition: EKLMElementNumbers.cc:145
Belle2::DBObjPtr
Class for accessing objects in the database.
Definition: DBObjPtr.h:31
Belle2::EKLM::AlignmentChecker
Class for EKLM alignment checking.
Definition: AlignmentChecker.h:38
Belle2::EKLM::TransformData::m_Layer
HepGeom::Transform3D ** m_Layer
Layer transformations.
Definition: TransformData.h:223
Belle2::Unit::rad
static const double rad
Standard of [angle].
Definition: Unit.h:60
Belle2::KLMAlignmentData::getDeltaU
float getDeltaU() const
Get shift in U.
Definition: KLMAlignmentData.h:87
Belle2::EKLM::TransformData::getSectionTransform
const HepGeom::Transform3D * getSectionTransform(int section) const
Get section transformation.
Definition: TransformData.cc:286
Belle2::EKLMGeometry::getNStrips
int getNStrips() const
Get number of strips.
Definition: EKLMGeometry.h:1741
Belle2::EKLM::GeometryData::getPlaneTransform
void getPlaneTransform(HepGeom::Transform3D *t, int n) const
Get plane transformation.
Definition: GeometryData.cc:690
Belle2::EKLM::TransformData::m_StripInverse
HepGeom::Transform3D ***** m_StripInverse
Inverse strip transformations.
Definition: TransformData.h:241
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::EKLM::TransformData::m_Section
HepGeom::Transform3D * m_Section
Section transformations.
Definition: TransformData.h:220
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::EKLM::TransformData::getLayerTransform
const HepGeom::Transform3D * getLayerTransform(int section, int layer) const
Get layer transformation.
Definition: TransformData.cc:292
Belle2::EKLMElementNumbers::Instance
static const EKLMElementNumbers & Instance()
Instantiation.
Definition: EKLMElementNumbers.cc:20
Belle2::EKLMGeometry::getNDetectorLayers
int getNDetectorLayers(int section) const
Get number of detector layers.
Definition: EKLMGeometry.cc:295
Belle2::KLMAlignmentData::getDeltaGamma
float getDeltaGamma() const
Get rotation in alpha.
Definition: KLMAlignmentData.h:172
Belle2::EKLM::TransformData::m_Strip
HepGeom::Transform3D ***** m_Strip
Strip transformations.
Definition: TransformData.h:238
Belle2::EKLM::TransformData::intersection
bool intersection(KLMDigit *hit1, KLMDigit *hit2, HepGeom::Point3D< double > *cross, double *d1, double *d2, double *sd, bool segments=true) const
Check if strips intersect, and find intersection point if yes.
Definition: TransformData.cc:350
Belle2::EKLM::TransformData::~TransformData
~TransformData()
Destructor.
Definition: TransformData.cc:190
Belle2::EKLM::TransformData::TransformData
TransformData(bool global, Displacement displacementType)
Constructor.
Definition: TransformData.cc:27
Belle2::EKLM::TransformData::m_GeoDat
const GeometryData * m_GeoDat
Geometry data.
Definition: TransformData.h:217
Belle2::EKLM::TransformData::m_Plane
HepGeom::Transform3D **** m_Plane
Plane transformations.
Definition: TransformData.h:229
Belle2::EKLMElementNumbers::getNStripsSegment
static constexpr int getNStripsSegment()
Get number of strips in a segment.
Definition: EKLMElementNumbers.h:401
Belle2::Unit::mm
static const double mm
[millimeters]
Definition: Unit.h:80
Belle2::EKLM::GeometryData::getStripTransform
void getStripTransform(HepGeom::Transform3D *t, int n) const
Get strip transformation.
Definition: GeometryData.cc:703
Belle2::EKLM::TransformData::getStripsByIntersection
int getStripsByIntersection(const HepGeom::Point3D< double > &intersection, int *strip1, int *strip2) const
Find strips by intersection.
Definition: TransformData.cc:442
Belle2::EKLM::AlignmentChecker::checkAlignment
bool checkAlignment(const EKLMAlignment *alignment, const EKLMSegmentAlignment *segmentAlignment) const
Check alignment.
Definition: AlignmentChecker.cc:310
Belle2::EKLM::TransformData::m_ElementNumbers
const EKLMElementNumbers * m_ElementNumbers
Element numbers.
Definition: TransformData.h:214
Belle2::EKLM::TransformData::c_Displacement
@ c_Displacement
Use displacement data (for geometry).
Definition: TransformData.h:46
Belle2::KLMDigit::getStrip
int getStrip() const
Get strip number.
Definition: KLMDigit.h:179
HepGeom::Point3D< double >
Belle2::KLMAlignmentData
KLM Alignment data.
Definition: KLMAlignmentData.h:33
Belle2::EKLMGeometry::getNSections
int getNSections() const
Get number of sections.
Definition: EKLMGeometry.h:1687
Belle2::EKLM::GeometryData::getSectorTransform
void getSectorTransform(HepGeom::Transform3D *t, int n) const
Get sector transformation.
Definition: GeometryData.cc:670
Belle2::EKLMGeometry::getNSectors
int getNSectors() const
Get number of sectors.
Definition: EKLMGeometry.h:1709
Belle2::EKLM::TransformData::m_Sector
HepGeom::Transform3D *** m_Sector
Sector transformations.
Definition: TransformData.h:226
Belle2::EKLM::GeometryData::getSectionTransform
void getSectionTransform(HepGeom::Transform3D *t, int n) const
Get section transformation.
Definition: GeometryData.cc:650
Belle2::EKLM::GeometryData::Instance
static const GeometryData & Instance(enum DataSource dataSource=c_Database, const GearDir *gearDir=nullptr)
Instantiation.
Definition: GeometryData.cc:35
Belle2::EKLM::TransformData::transformsToGlobal
void transformsToGlobal()
Make transformations global from local.
Definition: TransformData.cc:241
Belle2::EKLM::TransformData::c_None
@ c_None
Displacement is not used.
Definition: TransformData.h:45
Belle2::EKLM::TransformData::Displacement
Displacement
Source of displacement (alignment) data.
Definition: TransformData.h:44