Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TGeoXtru.cxx
Go to the documentation of this file.
1// @(#)root/geom:$Id$
2// Author: Mihaela Gheata 24/01/04
3
4/*************************************************************************
5 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12/** \class TGeoXtru
13\ingroup Shapes_classes
14A TGeoXtru shape is represented by the extrusion of an arbitrary
15polygon with fixed outline between several Z sections. Each Z section is
16a scaled version of the same "blueprint" polygon. Different global XY
17translations are allowed from section to section. Corresponding polygon
18vertices from consecutive sections are connected.
19
20An extruded polygon can be created using the constructor:
21
22~~~ {.cpp}
23TGeoXtru::TGeoXtru(Int_t nplanes);
24~~~
25
26 - `nplanes: `number of Z sections (minimum 2)
27
28Begin_Macro
29{
30 new TGeoManager("xtru", "poza12");
31 TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
32 TGeoMedium *med = new TGeoMedium("MED",1,mat);
33 TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
34 gGeoManager->SetTopVolume(top);
35 TGeoVolume *vol = gGeoManager->MakeXtru("XTRU",med,4);
36 TGeoXtru *xtru = (TGeoXtru*)vol->GetShape();
37 Double_t x[8] = {-30,-30,30,30,15,15,-15,-15};
38 Double_t y[8] = {-30,30,30,-30,-30,15,15,-30};
39 xtru->DefinePolygon(8,x,y);
40 xtru->DefineSection(0,-40, -20., 10., 1.5);
41 xtru->DefineSection(1, 10, 0., 0., 0.5);
42 xtru->DefineSection(2, 10, 0., 0., 0.7);
43 xtru->DefineSection(3, 40, 10., 20., 0.9);
44 top->AddNode(vol,1);
45 gGeoManager->CloseGeometry();
46 gGeoManager->SetNsegments(80);
47 top->Draw();
48 if (gPad) {
49 TView *view = gPad->GetView();
50 if (view) view->ShowAxis();
51 }
52}
53End_Macro
54
55The lists of X and Y positions for all vertices have to be provided for
56the "blueprint" polygon:
57
58~~~{.cpp}
59TGeoXtru::DefinePolygon (Int_t nvertices, Double_t *xv,
60Double_t *yv);
61~~~
62
63 - `nvertices: `number of vertices of the polygon
64 - `xv,yv: `arrays of X and Y coordinates for polygon vertices
65
66The method creates an object of the class **`TGeoPolygon`** for which
67the convexity is automatically determined . The polygon is decomposed
68into convex polygons if needed.
69
70Next step is to define the Z positions for each section plane as well as
71the XY offset and scaling for the corresponding polygons.
72
73~~~{.cpp}
74TGeoXtru::DefineSection(Int_t snum,Double_t zsection,Double_t x0,
75Double_t y0, Double_t scale);
76~~~
77
78 - `snum: `Z section index (0, nplanes-1). The section with
79 `snum = nplanes-1` must be defined last and triggers the computation
80 of the bounding box for the whole shape
81 - `zsection: `Z position of section `snum`. Sections must be defined
82 in increasing order of Z (e.g. `snum=0` correspond to the minimum Z
83 and `snum=nplanes-1` to the maximum one).
84 - `x0,y0: `offset of section `snum` with respect to the local shape
85 reference frame `T`
86 - `scale: `factor that multiplies the X/Y coordinates of each vertex
87 of the polygon at section `snum`:
88 - `x[ivert] = x0 + scale*xv[ivert]`
89 - `y[ivert] = y0 + scale*yv[ivert]`
90*/
91
92#include "TGeoXtru.h"
93
94#include <iostream>
95
96#include "TBuffer3D.h"
97#include "TBuffer3DTypes.h"
98#include "TMath.h"
99
100#include "TVirtualGeoPainter.h"
101#include "TGeoManager.h"
102#include "TGeoVolume.h"
103#include "TGeoPolygon.h"
104
105
106////////////////////////////////////////////////////////////////////////////////
107/// Constructor.
108
109TGeoXtru::ThreadData_t::ThreadData_t() : fSeg(0), fIz(0), fXc(nullptr), fYc(nullptr), fPoly(nullptr) {}
110
111////////////////////////////////////////////////////////////////////////////////
112/// Destructor.
113
115{
116 delete[] fXc;
117 delete[] fYc;
118 delete fPoly;
119}
120
121////////////////////////////////////////////////////////////////////////////////
122
124{
125 if (!fThreadSize)
126 ((TGeoXtru *)this)->CreateThreadData(1);
128 return *fThreadData[tid];
129}
130
131////////////////////////////////////////////////////////////////////////////////
132
134{
135 std::lock_guard<std::mutex> guard(fMutex);
136 std::vector<ThreadData_t *>::iterator i = fThreadData.begin();
137 while (i != fThreadData.end()) {
138 delete *i;
139 ++i;
140 }
141 fThreadData.clear();
142 fThreadSize = 0;
143}
144
145////////////////////////////////////////////////////////////////////////////////
146/// Create thread data for n threads max.
147
149{
150 std::lock_guard<std::mutex> guard(fMutex);
151 fThreadData.resize(nthreads);
153 for (Int_t tid = 0; tid < nthreads; tid++) {
154 if (fThreadData[tid] == nullptr) {
157 td.fXc = new Double_t[fNvert];
158 td.fYc = new Double_t[fNvert];
159 memcpy(td.fXc, fX, fNvert * sizeof(Double_t));
160 memcpy(td.fYc, fY, fNvert * sizeof(Double_t));
161 td.fPoly = new TGeoPolygon(fNvert);
162 td.fPoly->SetXY(td.fXc, td.fYc); // initialize with current coordinates
163 td.fPoly->FinishPolygon();
164 if (tid == 0 && td.fPoly->IsIllegalCheck()) {
165 Error("DefinePolygon", "Shape %s of type XTRU has an illegal polygon.", GetName());
166 }
167 }
168 }
169}
170
171////////////////////////////////////////////////////////////////////////////////
172/// Set current z-plane.
173
175{
176 GetThreadData().fIz = iz;
177}
178////////////////////////////////////////////////////////////////////////////////
179/// Set current segment.
180
185
186////////////////////////////////////////////////////////////////////////////////
187/// dummy ctor
188
190 : TGeoBBox(),
191 fNvert(0),
192 fNz(0),
193 fZcurrent(0.),
194 fX(nullptr),
195 fY(nullptr),
196 fZ(nullptr),
197 fScale(nullptr),
198 fX0(nullptr),
199 fY0(nullptr),
200 fThreadData(0),
201 fThreadSize(0)
202{
204}
205
206////////////////////////////////////////////////////////////////////////////////
207/// Default constructor
208
210 : TGeoBBox(0, 0, 0),
211 fNvert(0),
212 fNz(nz),
213 fZcurrent(0.),
214 fX(nullptr),
215 fY(nullptr),
216 fZ(new Double_t[nz]),
217 fScale(new Double_t[nz]),
218 fX0(new Double_t[nz]),
219 fY0(new Double_t[nz]),
220 fThreadData(0),
221 fThreadSize(0)
222{
224 if (nz < 2) {
225 Error("ctor", "Cannot create TGeoXtru %s with less than 2 Z planes", GetName());
227 return;
228 }
229}
230
231////////////////////////////////////////////////////////////////////////////////
232/// Default constructor in GEANT3 style
233/// - param[0] = nz // number of z planes
234///
235/// - param[1] = z1 // Z position of first plane
236/// - param[2] = x1 // X position of first plane
237/// - param[3] = y1 // Y position of first plane
238/// - param[4] = scale1 // scale factor for first plane
239/// ...
240/// - param[4*(nz-1]+1] = zn
241/// - param[4*(nz-1)+2] = xn
242/// - param[4*(nz-1)+3] = yn
243/// - param[4*(nz-1)+4] = scalen
244
246 : TGeoBBox(0, 0, 0),
247 fNvert(0),
248 fNz(0),
249 fZcurrent(0.),
250 fX(nullptr),
251 fY(nullptr),
252 fZ(nullptr),
253 fScale(nullptr),
254 fX0(nullptr),
255 fY0(nullptr),
256 fThreadData(0),
257 fThreadSize(0)
258{
260 SetDimensions(param);
261}
262
263////////////////////////////////////////////////////////////////////////////////
264/// destructor
265
267{
268 if (fX) {
269 delete[] fX;
270 fX = nullptr;
271 }
272 if (fY) {
273 delete[] fY;
274 fY = nullptr;
275 }
276 if (fZ) {
277 delete[] fZ;
278 fZ = nullptr;
279 }
280 if (fScale) {
281 delete[] fScale;
282 fScale = nullptr;
283 }
284 if (fX0) {
285 delete[] fX0;
286 fX0 = nullptr;
287 }
288 if (fY0) {
289 delete[] fY0;
290 fY0 = nullptr;
291 }
293}
294
295////////////////////////////////////////////////////////////////////////////////
296/// Compute capacity [length^3] of this shape.
297
299{
301 Int_t iz;
302 Double_t capacity = 0;
304 TGeoXtru *xtru = (TGeoXtru *)this;
305 xtru->SetCurrentVertices(0., 0., 1.);
306 area = td.fPoly->Area();
307 for (iz = 0; iz < fNz - 1; iz++) {
308 dz = fZ[iz + 1] - fZ[iz];
310 continue;
311 sc1 = fScale[iz];
312 sc2 = fScale[iz + 1];
313 capacity += (area * dz / 3.) * (sc1 * sc1 + sc1 * sc2 + sc2 * sc2);
314 }
315 return capacity;
316}
317
318////////////////////////////////////////////////////////////////////////////////
319/// compute bounding box of the pcon
320
322{
324 if (!fX || !fZ || !fNvert) {
325 Error("ComputeBBox", "In shape %s polygon not defined", GetName());
327 return;
328 }
329 Double_t zmin = fZ[0];
330 Double_t zmax = fZ[fNz - 1];
335 for (Int_t i = 0; i < fNz; i++) {
336 SetCurrentVertices(fX0[i], fY0[i], fScale[i]);
337 for (Int_t j = 0; j < fNvert; j++) {
338 if (td.fXc[j] < xmin)
339 xmin = td.fXc[j];
340 if (td.fXc[j] > xmax)
341 xmax = td.fXc[j];
342 if (td.fYc[j] < ymin)
343 ymin = td.fYc[j];
344 if (td.fYc[j] > ymax)
345 ymax = td.fYc[j];
346 }
347 }
348 fOrigin[0] = 0.5 * (xmin + xmax);
349 fOrigin[1] = 0.5 * (ymin + ymax);
350 fOrigin[2] = 0.5 * (zmin + zmax);
351 fDX = 0.5 * (xmax - xmin);
352 fDY = 0.5 * (ymax - ymin);
353 fDZ = 0.5 * (zmax - zmin);
354}
355
356////////////////////////////////////////////////////////////////////////////////
357/// Compute normal to closest surface from POINT.
358
359void TGeoXtru::ComputeNormal(const Double_t * /*point*/, const Double_t *dir, Double_t *norm) const
360{
362 if (td.fIz < 0) {
363 memset(norm, 0, 3 * sizeof(Double_t));
364 norm[2] = (dir[2] > 0) ? 1 : -1;
365 return;
366 }
367 Double_t vert[12];
368 GetPlaneVertices(td.fIz, td.fSeg, vert);
370 Double_t ndotd = norm[0] * dir[0] + norm[1] * dir[1] + norm[2] * dir[2];
371 if (ndotd < 0) {
372 norm[0] = -norm[0];
373 norm[1] = -norm[1];
374 norm[2] = -norm[2];
375 }
376}
377
378////////////////////////////////////////////////////////////////////////////////
379/// test if point is inside this shape
380
382{
384 // Check Z range
385 TGeoXtru *xtru = (TGeoXtru *)this;
386 if (point[2] < fZ[0])
387 return kFALSE;
388 if (point[2] > fZ[fNz - 1])
389 return kFALSE;
390 Int_t iz = TMath::BinarySearch(fNz, fZ, point[2]);
391 if (iz < 0 || iz == fNz - 1)
392 return kFALSE;
393 if (TGeoShape::IsSameWithinTolerance(point[2], fZ[iz])) {
394 xtru->SetIz(-1);
395 xtru->SetCurrentVertices(fX0[iz], fY0[iz], fScale[iz]);
396 if (td.fPoly->Contains(point))
397 return kTRUE;
398 if (iz > 1 && TGeoShape::IsSameWithinTolerance(fZ[iz], fZ[iz - 1])) {
399 xtru->SetCurrentVertices(fX0[iz - 1], fY0[iz - 1], fScale[iz - 1]);
400 return td.fPoly->Contains(point);
401 } else if (iz < fNz - 2 && TGeoShape::IsSameWithinTolerance(fZ[iz], fZ[iz + 1])) {
402 xtru->SetCurrentVertices(fX0[iz + 1], fY0[iz + 1], fScale[iz + 1]);
403 return td.fPoly->Contains(point);
404 }
405 }
406 xtru->SetCurrentZ(point[2], iz);
407 if (TMath::Abs(point[2] - fZ[iz]) < TGeoShape::Tolerance() ||
408 TMath::Abs(fZ[iz + 1] - point[2]) < TGeoShape::Tolerance())
409 xtru->SetIz(-1);
410 // Now td.fXc,fYc represent the vertices of the section at point[2]
411 return td.fPoly->Contains(point);
412}
413
414////////////////////////////////////////////////////////////////////////////////
415/// compute closest distance from point px,py to each corner
416
418{
419 const Int_t numPoints = fNvert * fNz;
420 return ShapeDistancetoPrimitive(numPoints, px, py);
421}
422
423////////////////////////////////////////////////////////////////////////////////
424/// Draw the section polygon.
425
427{
429 if (td.fPoly)
430 td.fPoly->Draw(option);
431}
432
433////////////////////////////////////////////////////////////////////////////////
434/// Compute distance to a Xtru lateral surface.
435
437 Bool_t in) const
438{
441 Double_t vert[12];
442 Double_t norm[3];
444 Double_t pt[3];
446 if (TGeoShape::IsSameWithinTolerance(fZ[iz], fZ[iz + 1]) && !in) {
447 TGeoXtru *xtru = (TGeoXtru *)this;
448 snext = (fZ[iz] - point[2]) / dir[2];
449 if (snext < 0)
450 return TGeoShape::Big();
451 pt[0] = point[0] + snext * dir[0];
452 pt[1] = point[1] + snext * dir[1];
453 pt[2] = point[2] + snext * dir[2];
454 if (dir[2] < 0.)
455 xtru->SetCurrentVertices(fX0[iz], fY0[iz], fScale[iz]);
456 else
457 xtru->SetCurrentVertices(fX0[iz + 1], fY0[iz + 1], fScale[iz + 1]);
458 if (!td.fPoly->Contains(pt))
459 return TGeoShape::Big();
460 return snext;
461 }
464 Double_t ndotd = norm[0] * dir[0] + norm[1] * dir[1] + norm[2] * dir[2];
465 if (in) {
466 if (ndotd <= 0)
467 return TGeoShape::Big();
468 safe = (vert[0] - point[0]) * norm[0] + (vert[1] - point[1]) * norm[1] + (vert[2] - point[2]) * norm[2];
469 if (safe < -1.E-8)
470 return TGeoShape::Big(); // direction outwards plane
471 } else {
472 ndotd = -ndotd;
473 if (ndotd <= 0)
474 return TGeoShape::Big();
475 safe = (point[0] - vert[0]) * norm[0] + (point[1] - vert[1]) * norm[1] + (point[2] - vert[2]) * norm[2];
476 if (safe < -1.E-8)
477 return TGeoShape::Big(); // direction outwards plane
478 }
479 snext = safe / ndotd;
480 if (snext > stepmax)
481 return TGeoShape::Big();
482 if (fZ[iz] < fZ[iz + 1]) {
483 znew = point[2] + snext * dir[2];
484 if (znew < fZ[iz])
485 return TGeoShape::Big();
486 if (znew > fZ[iz + 1])
487 return TGeoShape::Big();
488 }
489 pt[0] = point[0] + snext * dir[0];
490 pt[1] = point[1] + snext * dir[1];
491 pt[2] = point[2] + snext * dir[2];
493 return TGeoShape::Big();
494 return TMath::Max(snext, 0.);
495}
496
497////////////////////////////////////////////////////////////////////////////////
498/// compute distance from inside point to surface of the polycone
499/// locate Z segment
500
503{
505 if (iact < 3 && safe) {
506 *safe = Safety(point, kTRUE);
507 if (iact == 0)
508 return TGeoShape::Big();
509 if (iact == 1 && step < *safe)
510 return TGeoShape::Big();
511 }
512 TGeoXtru *xtru = (TGeoXtru *)this;
513 Int_t iz = TMath::BinarySearch(fNz, fZ, point[2]);
514 if (iz < 0) {
515 if (dir[2] <= 0) {
516 xtru->SetIz(-1);
517 return 0.;
518 }
519 iz = 0;
520 }
521 if (iz == fNz - 1) {
522 if (dir[2] >= 0) {
523 xtru->SetIz(-1);
524 return 0.;
525 }
526 iz--;
527 } else {
528 if (iz > 0) {
529 if (TGeoShape::IsSameWithinTolerance(point[2], fZ[iz])) {
530 if (TGeoShape::IsSameWithinTolerance(fZ[iz], fZ[iz + 1]) && dir[2] < 0)
531 iz++;
532 else if (TGeoShape::IsSameWithinTolerance(fZ[iz], fZ[iz - 1]) && dir[2] > 0)
533 iz--;
534 }
535 }
536 }
537 Bool_t convex = td.fPoly->IsConvex();
538 // Double_t stepmax = step;
539 // if (stepmax>TGeoShape::Big()) stepmax = TGeoShape::Big();
541 Double_t dist, sz;
542 Double_t pt[3];
543 Int_t iv, ipl, inext;
544 // we treat the special case when dir[2]=0
545 if (TGeoShape::IsSameWithinTolerance(dir[2], 0)) {
546 for (iv = 0; iv < fNvert; iv++) {
547 xtru->SetIz(-1);
548 dist = DistToPlane(point, dir, iz, iv, TGeoShape::Big(), kTRUE);
549 if (dist < snext) {
550 snext = dist;
551 xtru->SetSeg(iv);
552 if (convex)
553 return snext;
554 }
555 }
556 if (snext < 1.E10)
557 return snext;
558 return TGeoShape::Tolerance();
559 }
560
561 // normal case
562 Int_t incseg = (dir[2] > 0) ? 1 : -1;
563 Int_t iznext = iz;
565 while (iz >= 0 && iz < fNz - 1) {
566 // find the distance to current segment end Z surface
567 ipl = iz + ((incseg + 1) >> 1); // next plane
568 inext = ipl + incseg; // next next plane
569 sz = (fZ[ipl] - point[2]) / dir[2];
570 if (sz < snext) {
571 iznext += incseg;
572 // we cross the next Z section before stepmax
573 pt[0] = point[0] + sz * dir[0];
574 pt[1] = point[1] + sz * dir[1];
575 xtru->SetCurrentVertices(fX0[ipl], fY0[ipl], fScale[ipl]);
576 if (td.fPoly->Contains(pt)) {
577 // ray gets through next polygon - is it the last one?
578 if (ipl == 0 || ipl == fNz - 1) {
579 xtru->SetIz(-1);
580 if (convex)
581 return sz;
582 zexit = kTRUE;
583 snext = sz;
584 }
585 // maybe a Z discontinuity - check this
587 xtru->SetCurrentVertices(fX0[inext], fY0[inext], fScale[inext]);
588 // if we do not cross the next polygone, we are out
589 if (!td.fPoly->Contains(pt)) {
590 xtru->SetIz(-1);
591 if (convex)
592 return sz;
593 zexit = kTRUE;
594 snext = sz;
595 } else {
596 iznext = inext;
597 }
598 }
599 }
600 } else {
601 iznext = fNz - 1; // stop
602 }
603 // ray may cross the lateral surfaces of section iz
604 for (iv = 0; iv < fNvert; iv++) {
605 dist = DistToPlane(point, dir, iz, iv, TGeoShape::Big(), kTRUE);
606 if (dist < snext) {
607 xtru->SetIz(iz);
608 xtru->SetSeg(iv);
609 snext = dist;
610 if (convex)
611 return snext;
612 zexit = kTRUE;
613 }
614 }
615 if (zexit)
616 return snext;
617 iz = iznext;
618 }
619 return TGeoShape::Tolerance();
620}
621
622////////////////////////////////////////////////////////////////////////////////
623/// compute distance from outside point to surface of the tube
624/// Warning("DistFromOutside", "not implemented");
625
628{
630 if (iact < 3 && safe) {
631 *safe = Safety(point, kTRUE);
632 if (iact == 0)
633 return TGeoShape::Big();
634 if (iact == 1 && step < *safe)
635 return TGeoShape::Big();
636 }
637 // Check if the bounding box is crossed within the requested distance
638 Double_t sdist = TGeoBBox::DistFromOutside(point, dir, fDX, fDY, fDZ, fOrigin, step);
639 if (sdist >= step)
640 return TGeoShape::Big();
641 Double_t stepmax = step;
642 if (stepmax > TGeoShape::Big())
644 Double_t snext = 0.;
645 Int_t i, iv;
646 Double_t pt[3];
647 memcpy(pt, point, 3 * sizeof(Double_t));
648 TGeoXtru *xtru = (TGeoXtru *)this;
649 // We might get out easy with Z checks
650 Int_t iz = TMath::BinarySearch(fNz, fZ, point[2]);
651 if (iz < 0) {
652 if (dir[2] <= 0)
653 return TGeoShape::Big();
654 // propagate to first Z plane
655 snext = (fZ[0] - point[2]) / dir[2];
656 if (snext > stepmax)
657 return TGeoShape::Big();
658 for (i = 0; i < 3; i++)
659 pt[i] = point[i] + snext * dir[i];
660 xtru->SetCurrentVertices(fX0[0], fY0[0], fScale[0]);
661 if (td.fPoly->Contains(pt)) {
662 xtru->SetIz(-1);
663 return snext;
664 }
665 iz = 0; // valid starting value = first segment
666 stepmax -= snext;
667 } else {
668 if (iz == fNz - 1) {
669 if (dir[2] >= 0)
670 return TGeoShape::Big();
671 // propagate to last Z plane
672 snext = (fZ[fNz - 1] - point[2]) / dir[2];
673 if (snext > stepmax)
674 return TGeoShape::Big();
675 for (i = 0; i < 3; i++)
676 pt[i] = point[i] + snext * dir[i];
677 xtru->SetCurrentVertices(fX0[fNz - 1], fY0[fNz - 1], fScale[fNz - 1]);
678 if (td.fPoly->Contains(pt)) {
679 xtru->SetIz(-1);
680 return snext;
681 }
682 iz = fNz - 2; // valid value = last segment
683 stepmax -= snext;
684 }
685 }
686 // Check if the bounding box is missed by the track
687 if (!TGeoBBox::Contains(pt)) {
688 Double_t dist = TGeoBBox::DistFromOutside(pt, dir, 3);
689 if (dist > stepmax)
690 return TGeoShape::Big();
691 if (dist > 1E-6)
692 dist -= 1E-6; // decrease snext to make sure we do not cross the xtru
693 else
694 dist = 0;
695 for (i = 0; i < 3; i++)
696 pt[i] += dist * dir[i]; // we are now closer
697 iz = TMath::BinarySearch(fNz, fZ, pt[2]);
698 if (iz < 0)
699 iz = 0;
700 else if (iz == fNz - 1)
701 iz = fNz - 2;
702 snext += dist;
703 stepmax -= dist;
704 }
705 // not the case - we have to do some work...
706 // Start tracking from current iz
707 // - first solve particular case dir[2]=0
708 Bool_t convex = td.fPoly->IsConvex();
709 Bool_t hit = kFALSE;
710 if (TGeoShape::IsSameWithinTolerance(dir[2], 0)) {
711 // loop lateral planes to see if we cross something
712 xtru->SetIz(iz);
713 for (iv = 0; iv < fNvert; iv++) {
714 Double_t dist = DistToPlane(pt, dir, iz, iv, stepmax, kFALSE);
715 if (dist < stepmax) {
716 xtru->SetSeg(iv);
717 if (convex)
718 return (snext + dist);
719 stepmax = dist;
720 hit = kTRUE;
721 }
722 }
723 if (hit)
724 return (snext + stepmax);
725 return TGeoShape::Big();
726 }
727 // general case
728 Int_t incseg = (dir[2] > 0) ? 1 : -1;
729 while (iz >= 0 && iz < fNz - 1) {
730 // compute distance to lateral planes
731 xtru->SetIz(iz);
732 if (TGeoShape::IsSameWithinTolerance(fZ[iz], fZ[iz + 1]))
733 xtru->SetIz(-1);
734 for (iv = 0; iv < fNvert; iv++) {
735 Double_t dist = DistToPlane(pt, dir, iz, iv, stepmax, kFALSE);
736 if (dist < stepmax) {
737 // HIT
738 xtru->SetSeg(iv);
739 if (convex)
740 return (snext + dist);
741 stepmax = dist;
742 hit = kTRUE;
743 }
744 }
745 if (hit)
746 return (snext + stepmax);
747 iz += incseg;
748 }
749 return TGeoShape::Big();
750}
751
752////////////////////////////////////////////////////////////////////////////////
753/// Creates the polygon representing the blueprint of any Xtru section.
754/// - nvert = number of vertices >2
755/// - xv[nvert] = array of X vertex positions
756/// - yv[nvert] = array of Y vertex positions
757///
758/// *NOTE* should be called before DefineSection or ctor with 'param'
759
761{
762 if (nvert < 3) {
763 Error("DefinePolygon", "In shape %s cannot create polygon with less than 3 vertices", GetName());
765 return kFALSE;
766 }
767 for (Int_t i = 0; i < nvert - 1; i++) {
768 for (Int_t j = i + 1; j < nvert; j++) {
770 Error("DefinePolygon", "In shape %s 2 vertices cannot be identical", GetName());
772 // return kFALSE;
773 }
774 }
775 }
776 fNvert = nvert;
777 if (fX)
778 delete[] fX;
779 fX = new Double_t[nvert];
780 if (fY)
781 delete[] fY;
782 fY = new Double_t[nvert];
783 memcpy(fX, xv, nvert * sizeof(Double_t));
784 memcpy(fY, yv, nvert * sizeof(Double_t));
785
787
788 return kTRUE;
789}
790
791////////////////////////////////////////////////////////////////////////////////
792/// defines z position of a section plane, rmin and rmax at this z.
793
795{
796 if ((snum < 0) || (snum >= fNz))
797 return;
798 fZ[snum] = z;
799 fX0[snum] = x0;
800 fY0[snum] = y0;
801 fScale[snum] = scale;
802 if (snum) {
803 if (fZ[snum] < fZ[snum - 1]) {
804 Warning("DefineSection",
805 "In shape: %s, Z position of section "
806 "%i, z=%e, not in increasing order, %i, z=%e",
807 GetName(), snum, fZ[snum], snum - 1, fZ[snum - 1]);
808 return;
809 }
810 }
811 if (snum == (fNz - 1)) {
812 ComputeBBox();
814 InspectShape();
815 }
816}
817
818////////////////////////////////////////////////////////////////////////////////
819/// Return the Z coordinate for segment ipl.
820
822{
823 if (ipl < 0 || ipl > (fNz - 1)) {
824 Error("GetZ", "In shape %s, ipl=%i out of range (0,%i)", GetName(), ipl, fNz - 1);
825 return 0.;
826 }
827 return fZ[ipl];
828}
829////////////////////////////////////////////////////////////////////////////////
830/// Returns normal vector to the planar quadrilateral defined by vector VERT.
831/// The normal points outwards the xtru.
832
834{
835 Double_t cross = 0.;
836 Double_t v1[3], v2[3];
837 v1[0] = vert[9] - vert[0];
838 v1[1] = vert[10] - vert[1];
839 v1[2] = vert[11] - vert[2];
840 v2[0] = vert[3] - vert[0];
841 v2[1] = vert[4] - vert[1];
842 v2[2] = vert[5] - vert[2];
843 norm[0] = v1[1] * v2[2] - v1[2] * v2[1];
844 cross += norm[0] * norm[0];
845 norm[1] = v1[2] * v2[0] - v1[0] * v2[2];
846 cross += norm[1] * norm[1];
847 norm[2] = v1[0] * v2[1] - v1[1] * v2[0];
848 cross += norm[2] * norm[2];
849 if (cross < TGeoShape::Tolerance())
850 return;
851 cross = 1. / TMath::Sqrt(cross);
852 for (Int_t i = 0; i < 3; i++)
853 norm[i] *= cross;
854}
855
856////////////////////////////////////////////////////////////////////////////////
857/// Returns (x,y,z) of 3 vertices of the surface defined by Z sections (iz, iz+1)
858/// and polygon vertices (ivert, ivert+1). No range check.
859
861{
863 Double_t x, y, z1, z2;
864 Int_t iv1 = (ivert + 1) % fNvert;
865 Int_t icrt = 0;
866 z1 = fZ[iz];
867 z2 = fZ[iz + 1];
868 if (td.fPoly->IsClockwise()) {
869 x = fX[ivert] * fScale[iz] + fX0[iz];
870 y = fY[ivert] * fScale[iz] + fY0[iz];
871 vert[icrt++] = x;
872 vert[icrt++] = y;
873 vert[icrt++] = z1;
874 x = fX[iv1] * fScale[iz] + fX0[iz];
875 y = fY[iv1] * fScale[iz] + fY0[iz];
876 vert[icrt++] = x;
877 vert[icrt++] = y;
878 vert[icrt++] = z1;
879 x = fX[iv1] * fScale[iz + 1] + fX0[iz + 1];
880 y = fY[iv1] * fScale[iz + 1] + fY0[iz + 1];
881 vert[icrt++] = x;
882 vert[icrt++] = y;
883 vert[icrt++] = z2;
884 x = fX[ivert] * fScale[iz + 1] + fX0[iz + 1];
885 y = fY[ivert] * fScale[iz + 1] + fY0[iz + 1];
886 vert[icrt++] = x;
887 vert[icrt++] = y;
888 vert[icrt++] = z2;
889 } else {
890 x = fX[iv1] * fScale[iz] + fX0[iz];
891 y = fY[iv1] * fScale[iz] + fY0[iz];
892 vert[icrt++] = x;
893 vert[icrt++] = y;
894 vert[icrt++] = z1;
895 x = fX[ivert] * fScale[iz] + fX0[iz];
896 y = fY[ivert] * fScale[iz] + fY0[iz];
897 vert[icrt++] = x;
898 vert[icrt++] = y;
899 vert[icrt++] = z1;
900 x = fX[ivert] * fScale[iz + 1] + fX0[iz + 1];
901 y = fY[ivert] * fScale[iz + 1] + fY0[iz + 1];
902 vert[icrt++] = x;
903 vert[icrt++] = y;
904 vert[icrt++] = z2;
905 x = fX[iv1] * fScale[iz + 1] + fX0[iz + 1];
906 y = fY[iv1] * fScale[iz + 1] + fY0[iz + 1];
907 vert[icrt++] = x;
908 vert[icrt++] = y;
909 vert[icrt++] = z2;
910 }
911}
912////////////////////////////////////////////////////////////////////////////////
913/// Check if the quadrilateral defined by VERT contains a coplanar POINT.
914
916{
917 Double_t v1[3], v2[3];
918 Double_t cross;
919 Int_t j, k;
920 for (Int_t i = 0; i < 4; i++) { // loop vertices
921 j = 3 * i;
922 k = 3 * ((i + 1) % 4);
923 v1[0] = point[0] - vert[j];
924 v1[1] = point[1] - vert[j + 1];
925 v1[2] = point[2] - vert[j + 2];
926 v2[0] = vert[k] - vert[j];
927 v2[1] = vert[k + 1] - vert[j + 1];
928 v2[2] = vert[k + 2] - vert[j + 2];
929 cross = (v1[1] * v2[2] - v1[2] * v2[1]) * norm[0] + (v1[2] * v2[0] - v1[0] * v2[2]) * norm[1] +
930 (v1[0] * v2[1] - v1[1] * v2[0]) * norm[2];
931 if (cross < 0)
932 return kFALSE;
933 }
934 return kTRUE;
935}
936
937////////////////////////////////////////////////////////////////////////////////
938/// Print actual Xtru parameters.
939
941{
942 printf("*** Shape %s: TGeoXtru ***\n", GetName());
943 printf(" Nz = %i\n", fNz);
944 printf(" List of (x,y) of polygon vertices:\n");
945 for (Int_t ivert = 0; ivert < fNvert; ivert++)
946 printf(" x = %11.5f y = %11.5f\n", fX[ivert], fY[ivert]);
947 for (Int_t ipl = 0; ipl < fNz; ipl++)
948 printf(" plane %i: z=%11.5f x0=%11.5f y0=%11.5f scale=%11.5f\n", ipl, fZ[ipl], fX0[ipl], fY0[ipl],
949 fScale[ipl]);
950 printf(" Bounding box:\n");
952}
953
954////////////////////////////////////////////////////////////////////////////////
955/// Creates a TBuffer3D describing *this* shape.
956/// Coordinates are in local reference frame.
957
959{
960 Int_t nz = GetNz();
961 Int_t nvert = GetNvert();
962 Int_t nbPnts = nz * nvert;
963 Int_t nbSegs = nvert * (2 * nz - 1);
964 Int_t nbPols = nvert * (nz - 1) + 2;
965
967 6 * (nbPols - 2) + 2 * (2 + nvert));
968 if (buff) {
969 SetPoints(buff->fPnts);
971 }
972
973 return buff;
974}
975
976////////////////////////////////////////////////////////////////////////////////
977/// Fill TBuffer3D structure for segments and polygons.
978
980{
981 Int_t nz = GetNz();
982 Int_t nvert = GetNvert();
984
985 Int_t i, j;
986 Int_t indx = 0, indx2, k;
987 for (i = 0; i < nz; i++) {
988 // loop Z planes
989 indx2 = i * nvert;
990 // loop polygon segments
991 for (j = 0; j < nvert; j++) {
992 k = (j + 1) % nvert;
993 buff.fSegs[indx++] = c;
994 buff.fSegs[indx++] = indx2 + j;
995 buff.fSegs[indx++] = indx2 + k;
996 }
997 } // total: nz*nvert polygon segments
998 for (i = 0; i < nz - 1; i++) {
999 // loop Z planes
1000 indx2 = i * nvert;
1001 // loop polygon segments
1002 for (j = 0; j < nvert; j++) {
1003 k = j + nvert;
1004 buff.fSegs[indx++] = c;
1005 buff.fSegs[indx++] = indx2 + j;
1006 buff.fSegs[indx++] = indx2 + k;
1007 }
1008 } // total (nz-1)*nvert lateral segments
1009
1010 indx = 0;
1011
1012 // fill lateral polygons
1013 for (i = 0; i < nz - 1; i++) {
1014 indx2 = i * nvert;
1015 for (j = 0; j < nvert; j++) {
1016 k = (j + 1) % nvert;
1017 buff.fPols[indx++] = c + j % 3;
1018 buff.fPols[indx++] = 4;
1019 buff.fPols[indx++] = indx2 + j;
1020 buff.fPols[indx++] = nz * nvert + indx2 + k;
1021 buff.fPols[indx++] = indx2 + nvert + j;
1022 buff.fPols[indx++] = nz * nvert + indx2 + j;
1023 }
1024 } // total (nz-1)*nvert polys
1025 buff.fPols[indx++] = c + 2;
1026 buff.fPols[indx++] = nvert;
1027 indx2 = 0;
1028 for (j = nvert - 1; j >= 0; --j) {
1029 buff.fPols[indx++] = indx2 + j;
1030 }
1031
1032 buff.fPols[indx++] = c;
1033 buff.fPols[indx++] = nvert;
1034 indx2 = (nz - 1) * nvert;
1035
1036 for (j = 0; j < nvert; j++) {
1037 buff.fPols[indx++] = indx2 + j;
1038 }
1039}
1040
1041////////////////////////////////////////////////////////////////////////////////
1042/// Compute safety to sector iz, returning also the closest segment index.
1043
1045{
1048 Bool_t in1, in2;
1049 Int_t iseg;
1050 // segment-break case
1051 if (TGeoShape::IsSameWithinTolerance(fZ[iz], fZ[iz + 1])) {
1052 safz = TMath::Abs(point[2] - fZ[iz]);
1053 if (safz > safmin)
1054 return TGeoShape::Big();
1055 SetCurrentVertices(fX0[iz], fY0[iz], fScale[iz]);
1056 saf1 = td.fPoly->Safety(point, iseg);
1057 in1 = td.fPoly->Contains(point);
1058 // if (!in1 && saf1>safmin) return TGeoShape::Big();
1059 SetCurrentVertices(fX0[iz + 1], fY0[iz + 1], fScale[iz + 1]);
1060 saf2 = td.fPoly->Safety(point, iseg);
1061 in2 = td.fPoly->Contains(point);
1062 if ((in1 & !in2) | (in2 & !in1)) {
1063 safe = safz;
1064 } else {
1067 }
1068 if (safe > safmin)
1069 return TGeoShape::Big();
1070 return safe;
1071 }
1072 // normal case
1073 safz = fZ[iz] - point[2];
1074 if (safz > safmin)
1075 return TGeoShape::Big();
1076 if (safz < 0) {
1077 saf1 = point[2] - fZ[iz + 1];
1078 if (saf1 > safmin)
1079 return TGeoShape::Big();
1080 if (saf1 < 0) {
1081 safz = TMath::Max(safz, saf1); // we are in between the 2 Z segments - we ignore safz
1082 } else {
1083 safz = saf1;
1084 }
1085 }
1086
1087 // loop segments
1088 Bool_t found = kFALSE;
1089 Double_t vert[12];
1090 Double_t norm[3];
1091 // printf("plane %d: safz=%f in=%d\n", iz, safz, in);
1092 for (iseg = 0; iseg < fNvert; iseg++) {
1095 saf1 = (point[0] - vert[0]) * norm[0] + (point[1] - vert[1]) * norm[1] + (point[2] - vert[2]) * norm[2];
1096 if (in)
1097 saf1 = -saf1;
1098 // printf("segment %d: (%f,%f)-(%f,%f) norm=(%f,%f,%f): saf1=%f\n", iseg,
1099 // vert[0],vert[1],vert[3],vert[4],norm[0],norm[1],norm[2],saf1);
1100 if (saf1 < -1.E-8)
1101 continue;
1103 safe = TMath::Abs(safe);
1104 if (safe > safmin)
1105 continue;
1106 safmin = safe;
1107 found = kTRUE;
1108 }
1109 if (found)
1110 return safmin;
1111 return TGeoShape::Big();
1112}
1113
1114////////////////////////////////////////////////////////////////////////////////
1115/// computes the closest distance from given point to this shape, according
1116/// to option. The matching point on the shape is stored in spoint.
1117///> localize the Z segment
1118
1120{
1122 Double_t safe;
1124 TGeoXtru *xtru = (TGeoXtru *)this;
1125 Int_t iz;
1126 if (in) {
1127 safmin = TMath::Min(point[2] - fZ[0], fZ[fNz - 1] - point[2]);
1128 for (iz = 0; iz < fNz - 1; iz++) {
1129 safe = xtru->SafetyToSector(point, iz, safmin, in);
1130 if (safe < safmin)
1131 safmin = safe;
1132 }
1133 return safmin;
1134 }
1135 // Accurate safety is expensive, use the bounding box
1136 if (!TGeoBBox::Contains(point))
1137 return TGeoBBox::Safety(point, in);
1138 iz = TMath::BinarySearch(fNz, fZ, point[2]);
1139 if (iz < 0) {
1140 iz = 0;
1141 safz = fZ[0] - point[2];
1142 } else {
1143 if (iz == fNz - 1) {
1144 iz = fNz - 2;
1145 safz = point[2] - fZ[fNz - 1];
1146 }
1147 }
1148 // loop segments from iz up
1149 Int_t i;
1150 for (i = iz; i < fNz - 1; i++) {
1151 safe = xtru->SafetyToSector(point, i, safmin, in);
1152 if (safe < safmin)
1153 safmin = safe;
1154 }
1155 // loop segments from iz-1 down
1156 for (i = iz - 1; i >= 0; i--) {
1157 safe = xtru->SafetyToSector(point, i, safmin, in);
1158 if (safe < safmin)
1159 safmin = safe;
1160 }
1162 return safe;
1163}
1164
1165////////////////////////////////////////////////////////////////////////////////
1166/// Save a primitive as a C++ statement(s) on output stream "out".
1167
1168void TGeoXtru::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
1169{
1171 return;
1172 out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
1173 out << " auto " << GetPointerName() << " = new TGeoXtru(" << fNz << ");" << std::endl;
1174 out << " " << GetPointerName() << "->SetName(\"" << GetName() << "\");" << std::endl;
1175 for (Int_t i = 0; i < fNvert; i++) {
1176 out << " xvert[" << i << "] = " << fX[i] << "; yvert[" << i << "] = " << fY[i] << ";" << std::endl;
1177 }
1178 out << " " << GetPointerName() << "->DefinePolygon(" << fNvert << ", xvert, yvert);" << std::endl;
1179 for (Int_t i = 0; i < fNz; i++)
1180 out << " " << GetPointerName() << "->DefineSection(" << i << ", " << fZ[i] << ", " << fX0[i] << ", " << fY0[i]
1181 << ", " << fScale[i] << ");" << std::endl;
1183}
1184
1185////////////////////////////////////////////////////////////////////////////////
1186/// Recompute current section vertices for a given Z position within range of section iz.
1187
1189{
1190 Double_t x0, y0, scale, a, b;
1191 Int_t ind1, ind2;
1192 ind1 = iz;
1193 ind2 = iz + 1;
1194 Double_t invdz = 1. / (fZ[ind2] - fZ[ind1]);
1195 a = (fX0[ind1] * fZ[ind2] - fX0[ind2] * fZ[ind1]) * invdz;
1196 b = (fX0[ind2] - fX0[ind1]) * invdz;
1197 x0 = a + b * z;
1198 a = (fY0[ind1] * fZ[ind2] - fY0[ind2] * fZ[ind1]) * invdz;
1199 b = (fY0[ind2] - fY0[ind1]) * invdz;
1200 y0 = a + b * z;
1201 a = (fScale[ind1] * fZ[ind2] - fScale[ind2] * fZ[ind1]) * invdz;
1202 b = (fScale[ind2] - fScale[ind1]) * invdz;
1203 scale = a + b * z;
1205}
1206
1207////////////////////////////////////////////////////////////////////////////////
1208/// Set current vertex coordinates according X0, Y0 and SCALE.
1209
1211{
1213 for (Int_t i = 0; i < fNvert; i++) {
1214 td.fXc[i] = scale * fX[i] + x0;
1215 td.fYc[i] = scale * fY[i] + y0;
1216 }
1217}
1218
1219////////////////////////////////////////////////////////////////////////////////
1220/// - param[0] = nz // number of z planes
1221///
1222/// - param[1] = z1 // Z position of first plane
1223/// - param[2] = x1 // X position of first plane
1224/// - param[3] = y1 // Y position of first plane
1225/// - param[4] = scale1 // scale factor for first plane
1226/// ...
1227/// - param[4*(nz-1]+1] = zn
1228/// - param[4*(nz-1)+2] = xn
1229/// - param[4*(nz-1)+3] = yn
1230/// - param[4*(nz-1)+4] = scalen
1231
1233{
1234 fNz = (Int_t)param[0];
1235 if (fNz < 2) {
1236 Error("SetDimensions", "Cannot create TGeoXtru %s with less than 2 Z planes", GetName());
1238 return;
1239 }
1240 if (fZ)
1241 delete[] fZ;
1242 if (fScale)
1243 delete[] fScale;
1244 if (fX0)
1245 delete[] fX0;
1246 if (fY0)
1247 delete[] fY0;
1248 fZ = new Double_t[fNz];
1249 fScale = new Double_t[fNz];
1250 fX0 = new Double_t[fNz];
1251 fY0 = new Double_t[fNz];
1252
1253 for (Int_t i = 0; i < fNz; i++)
1254 DefineSection(i, param[1 + 4 * i], param[2 + 4 * i], param[3 + 4 * i], param[4 + 4 * i]);
1255}
1256
1257////////////////////////////////////////////////////////////////////////////////
1258/// create polycone mesh points
1259
1261{
1263 Int_t i, j;
1264 Int_t indx = 0;
1265 TGeoXtru *xtru = (TGeoXtru *)this;
1266 if (points) {
1267 for (i = 0; i < fNz; i++) {
1268 xtru->SetCurrentVertices(fX0[i], fY0[i], fScale[i]);
1269 if (td.fPoly->IsClockwise()) {
1270 for (j = 0; j < fNvert; j++) {
1271 points[indx++] = td.fXc[j];
1272 points[indx++] = td.fYc[j];
1273 points[indx++] = fZ[i];
1274 }
1275 } else {
1276 for (j = 0; j < fNvert; j++) {
1277 points[indx++] = td.fXc[fNvert - 1 - j];
1278 points[indx++] = td.fYc[fNvert - 1 - j];
1279 points[indx++] = fZ[i];
1280 }
1281 }
1282 }
1283 }
1284}
1285
1286////////////////////////////////////////////////////////////////////////////////
1287/// create polycone mesh points
1288
1290{
1292 Int_t i, j;
1293 Int_t indx = 0;
1294 TGeoXtru *xtru = (TGeoXtru *)this;
1295 if (points) {
1296 for (i = 0; i < fNz; i++) {
1297 xtru->SetCurrentVertices(fX0[i], fY0[i], fScale[i]);
1298 if (td.fPoly->IsClockwise()) {
1299 for (j = 0; j < fNvert; j++) {
1300 points[indx++] = td.fXc[j];
1301 points[indx++] = td.fYc[j];
1302 points[indx++] = fZ[i];
1303 }
1304 } else {
1305 for (j = 0; j < fNvert; j++) {
1306 points[indx++] = td.fXc[fNvert - 1 - j];
1307 points[indx++] = td.fYc[fNvert - 1 - j];
1308 points[indx++] = fZ[i];
1309 }
1310 }
1311 }
1312 }
1313}
1314
1315////////////////////////////////////////////////////////////////////////////////
1316/// Returns numbers of vertices, segments and polygons composing the shape mesh.
1317
1319{
1320 Int_t nz = GetNz();
1321 Int_t nv = GetNvert();
1322 nvert = nz * nv;
1323 nsegs = nv * (2 * nz - 1);
1324 npols = nv * (nz - 1) + 2;
1325}
1326
1327////////////////////////////////////////////////////////////////////////////////
1328/// Return number of vertices of the mesh representation
1329
1331{
1332 Int_t numPoints = fNz * fNvert;
1333 return numPoints;
1334}
1335
1336////////////////////////////////////////////////////////////////////////////////
1337/// fill size of this 3-D object
1338
1339void TGeoXtru::Sizeof3D() const {}
1340
1341////////////////////////////////////////////////////////////////////////////////
1342/// Fills a static 3D buffer and returns a reference.
1343
1345{
1346 static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
1347
1349
1351 Int_t nz = GetNz();
1352 Int_t nvert = GetNvert();
1353 Int_t nbPnts = nz * nvert;
1354 Int_t nbSegs = nvert * (2 * nz - 1);
1355 Int_t nbPols = nvert * (nz - 1) + 2;
1356 if (buffer.SetRawSizes(nbPnts, 3 * nbPnts, nbSegs, 3 * nbSegs, nbPols, 6 * (nbPols - 2) + 2 * (2 + nvert))) {
1358 }
1359 }
1360 // TODO: Push down to TGeoShape?
1362 SetPoints(buffer.fPnts);
1363 if (!buffer.fLocalFrame) {
1364 TransformPoints(buffer.fPnts, buffer.NbPnts());
1365 }
1366
1367 SetSegsAndPols(buffer);
1369 }
1370
1371 return buffer;
1372}
1373
1374////////////////////////////////////////////////////////////////////////////////
1375/// Check the inside status for each of the points in the array.
1376/// Input: Array of point coordinates + vector size
1377/// Output: Array of Booleans for the inside of each point
1378
1380{
1381 for (Int_t i = 0; i < vecsize; i++)
1382 inside[i] = Contains(&points[3 * i]);
1383}
1384
1385////////////////////////////////////////////////////////////////////////////////
1386/// Compute the normal for an array o points so that norm.dot.dir is positive
1387/// Input: Arrays of point coordinates and directions + vector size
1388/// Output: Array of normal directions
1389
1391{
1392 for (Int_t i = 0; i < vecsize; i++)
1393 ComputeNormal(&points[3 * i], &dirs[3 * i], &norms[3 * i]);
1394}
1395
1396////////////////////////////////////////////////////////////////////////////////
1397/// Compute distance from array of input points having directions specified by dirs. Store output in dists
1398
1400 Double_t *step) const
1401{
1402 for (Int_t i = 0; i < vecsize; i++)
1403 dists[i] = DistFromInside(&points[3 * i], &dirs[3 * i], 3, step[i]);
1404}
1405
1406////////////////////////////////////////////////////////////////////////////////
1407/// Compute distance from array of input points having directions specified by dirs. Store output in dists
1408
1410 Double_t *step) const
1411{
1412 for (Int_t i = 0; i < vecsize; i++)
1413 dists[i] = DistFromOutside(&points[3 * i], &dirs[3 * i], 3, step[i]);
1414}
1415
1416////////////////////////////////////////////////////////////////////////////////
1417/// Compute safe distance from each of the points in the input array.
1418/// Input: Array of point coordinates, array of statuses for these points, size of the arrays
1419/// Output: Safety values
1420
1422{
1423 for (Int_t i = 0; i < vecsize; i++)
1424 safe[i] = Safety(&points[3 * i], inside[i]);
1425}
#define b(i)
Definition RSha256.hxx:100
#define c(i)
Definition RSha256.hxx:101
#define a(i)
Definition RSha256.hxx:99
int Int_t
Signed integer 4 bytes (int)
Definition RtypesCore.h:59
float Float_t
Float 4 bytes (float)
Definition RtypesCore.h:71
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
double Double_t
Double 8 bytes.
Definition RtypesCore.h:73
constexpr Bool_t kTRUE
Definition RtypesCore.h:107
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t option
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t points
float xmin
float ymin
float xmax
float ymax
Generic 3D primitive description class.
Definition TBuffer3D.h:18
UInt_t NbPnts() const
Definition TBuffer3D.h:80
Bool_t SectionsValid(UInt_t mask) const
Definition TBuffer3D.h:67
void SetSectionsValid(UInt_t mask)
Definition TBuffer3D.h:65
Bool_t fLocalFrame
Definition TBuffer3D.h:90
Bool_t SetRawSizes(UInt_t reqPnts, UInt_t reqPntsCapacity, UInt_t reqSegs, UInt_t reqSegsCapacity, UInt_t reqPols, UInt_t reqPolsCapacity)
Set kRaw tessellation section of buffer with supplied sizes.
Double_t * fPnts
Definition TBuffer3D.h:113
Box class.
Definition TGeoBBox.h:17
void FillBuffer3D(TBuffer3D &buffer, Int_t reqSections, Bool_t localFrame) const override
Fills the supplied buffer, with sections in desired frame See TBuffer3D.h for explanation of sections...
Double_t fDX
Definition TGeoBBox.h:20
Double_t DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=nullptr) const override
Compute distance from outside point to surface of the box.
Definition TGeoBBox.cxx:432
Double_t fOrigin[3]
Definition TGeoBBox.h:23
void InspectShape() const override
Prints shape parameters.
Definition TGeoBBox.cxx:810
Bool_t Contains(const Double_t *point) const override
Test if point is inside this shape.
Definition TGeoBBox.cxx:323
Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const override
Computes the closest distance from given point to this shape.
Definition TGeoBBox.cxx:919
Double_t fDY
Definition TGeoBBox.h:21
Double_t fDZ
Definition TGeoBBox.h:22
static Int_t ThreadId()
Translates the current thread id to an ordinal number.
An arbitrary polygon defined by vertices.
Definition TGeoPolygon.h:19
static Double_t Big()
Definition TGeoShape.h:94
Int_t GetBasicColor() const
Get the basic color (0-7).
void TransformPoints(Double_t *points, UInt_t NbPoints) const
Tranform a set of points (LocalToMaster)
void SetShapeBit(UInt_t f, Bool_t set)
Equivalent of TObject::SetBit.
static Bool_t IsSameWithinTolerance(Double_t a, Double_t b)
Check if two numbers differ with less than a tolerance.
const char * GetPointerName() const
Provide a pointer name containing uid.
Int_t ShapeDistancetoPrimitive(Int_t numpoints, Int_t px, Int_t py) const
Returns distance to shape primitive mesh.
const char * GetName() const override
Get the shape name.
@ kGeoSavePrimitive
Definition TGeoShape.h:64
static Double_t Tolerance()
Definition TGeoShape.h:97
Bool_t TestShapeBit(UInt_t f) const
Definition TGeoShape.h:175
A TGeoXtru shape is represented by the extrusion of an arbitrary polygon with fixed outline between s...
Definition TGeoXtru.h:22
Double_t * fScale
Definition TGeoXtru.h:46
Int_t GetNmeshVertices() const override
Return number of vertices of the mesh representation.
Double_t * fX0
Definition TGeoXtru.h:47
void DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t *step) const override
Compute distance from array of input points having directions specified by dirs. Store output in dist...
TGeoXtru()
dummy ctor
Definition TGeoXtru.cxx:189
ThreadData_t & GetThreadData() const
Definition TGeoXtru.cxx:123
void ClearThreadData() const override
Definition TGeoXtru.cxx:133
Double_t Capacity() const override
Compute capacity [length^3] of this shape.
Definition TGeoXtru.cxx:298
Int_t fNz
Definition TGeoXtru.h:41
Double_t * fZ
Definition TGeoXtru.h:45
void ComputeBBox() override
compute bounding box of the pcon
Definition TGeoXtru.cxx:321
void Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const override
Check the inside status for each of the points in the array.
Bool_t Contains(const Double_t *point) const override
test if point is inside this shape
Definition TGeoXtru.cxx:381
Int_t fThreadSize
Navigation data per thread.
Definition TGeoXtru.h:51
Int_t fNvert
Definition TGeoXtru.h:40
void SetCurrentVertices(Double_t x0, Double_t y0, Double_t scale)
Set current vertex coordinates according X0, Y0 and SCALE.
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save a primitive as a C++ statement(s) on output stream "out".
void GetPlaneVertices(Int_t iz, Int_t ivert, Double_t *vert) const
Returns (x,y,z) of 3 vertices of the surface defined by Z sections (iz, iz+1) and polygon vertices (i...
Definition TGeoXtru.cxx:860
void InspectShape() const override
Print actual Xtru parameters.
Definition TGeoXtru.cxx:940
Int_t DistancetoPrimitive(Int_t px, Int_t py) override
compute closest distance from point px,py to each corner
Definition TGeoXtru.cxx:417
Double_t * GetZ() const
Definition TGeoXtru.h:102
Bool_t IsPointInsidePlane(const Double_t *point, Double_t *vert, Double_t *norm) const
Check if the quadrilateral defined by VERT contains a coplanar POINT.
Definition TGeoXtru.cxx:915
void Sizeof3D() const override
fill size of this 3-D object
~TGeoXtru() override
destructor
Definition TGeoXtru.cxx:266
const TBuffer3D & GetBuffer3D(Int_t reqSections, Bool_t localFrame) const override
Fills a static 3D buffer and returns a reference.
Int_t GetNvert() const
Definition TGeoXtru.h:96
Double_t * fY0
Definition TGeoXtru.h:48
Bool_t DefinePolygon(Int_t nvert, const Double_t *xv, const Double_t *yv)
Creates the polygon representing the blueprint of any Xtru section.
Definition TGeoXtru.cxx:760
void SetCurrentZ(Double_t z, Int_t iz)
Recompute current section vertices for a given Z position within range of section iz.
virtual void DefineSection(Int_t snum, Double_t z, Double_t x0=0., Double_t y0=0., Double_t scale=1.)
defines z position of a section plane, rmin and rmax at this z.
Definition TGeoXtru.cxx:794
Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const override
computes the closest distance from given point to this shape, according to option.
void ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize) override
Compute the normal for an array o points so that norm.dot.dir is positive Input: Arrays of point coor...
void SetDimensions(Double_t *param) override
void DrawPolygon(Option_t *option="")
Draw the section polygon.
Definition TGeoXtru.cxx:426
void SetPoints(Double_t *points) const override
create polycone mesh points
Double_t SafetyToSector(const Double_t *point, Int_t iz, Double_t safmin, Bool_t in)
Compute safety to sector iz, returning also the closest segment index.
Double_t * fY
Definition TGeoXtru.h:44
Double_t DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=nullptr) const override
compute distance from inside point to surface of the polycone locate Z segment
Definition TGeoXtru.cxx:502
Double_t * fX
Definition TGeoXtru.h:43
void DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t *step) const override
Compute distance from array of input points having directions specified by dirs. Store output in dist...
void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm) const override
Compute normal to closest surface from POINT.
Definition TGeoXtru.cxx:359
TBuffer3D * MakeBuffer3D() const override
Creates a TBuffer3D describing this shape.
Definition TGeoXtru.cxx:958
std::vector< ThreadData_t * > fThreadData
Definition TGeoXtru.h:50
Double_t DistToPlane(const Double_t *point, const Double_t *dir, Int_t iz, Int_t ivert, Double_t stepmax, Bool_t in) const
Compute distance to a Xtru lateral surface.
Definition TGeoXtru.cxx:436
void SetSegsAndPols(TBuffer3D &buff) const override
Fill TBuffer3D structure for segments and polygons.
Definition TGeoXtru.cxx:979
void SetSeg(Int_t iseg)
Set current segment.
Definition TGeoXtru.cxx:181
Double_t DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=nullptr) const override
compute distance from outside point to surface of the tube Warning("DistFromOutside",...
Definition TGeoXtru.cxx:627
std::mutex fMutex
size of thread-specific array
Definition TGeoXtru.h:52
void GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const override
Returns numbers of vertices, segments and polygons composing the shape mesh.
void GetPlaneNormal(const Double_t *vert, Double_t *norm) const
Returns normal vector to the planar quadrilateral defined by vector VERT.
Definition TGeoXtru.cxx:833
void SetIz(Int_t iz)
Set current z-plane.
Definition TGeoXtru.cxx:174
Int_t GetNz() const
Definition TGeoXtru.h:95
void Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const override
Compute safe distance from each of the points in the input array.
void CreateThreadData(Int_t nthreads) override
Create thread data for n threads max.
Definition TGeoXtru.cxx:148
Double_t fZcurrent
Definition TGeoXtru.h:42
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition TObject.h:202
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:226
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:1057
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition TObject.cxx:864
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1071
TPaveText * pt
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
Short_t Max(Short_t a, Short_t b)
Returns the largest of a and b.
Definition TMathBase.h:251
Double_t Sqrt(Double_t x)
Returns the square root of x.
Definition TMath.h:673
Short_t Min(Short_t a, Short_t b)
Returns the smallest of a and b.
Definition TMathBase.h:199
Long64_t BinarySearch(Long64_t n, const T *array, T value)
Binary search in an array of n values to locate value.
Definition TMathBase.h:348
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:124
~ThreadData_t()
Destructor.
Definition TGeoXtru.cxx:114
ThreadData_t()
Constructor.
Definition TGeoXtru.cxx:109