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