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