Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TGeoBBox.cxx
Go to the documentation of this file.
1// @(#)root/geom:$Id$// Author: Andrei Gheata 24/10/01
2
3// Contains() and DistFromOutside/Out() implemented by Mihaela Gheata
4
5/*************************************************************************
6 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
7 * All rights reserved. *
8 * *
9 * For the licensing terms see $ROOTSYS/LICENSE. *
10 * For the list of contributors see $ROOTSYS/README/CREDITS. *
11 *************************************************************************/
12
13/** \class TGeoBBox
14\ingroup Shapes_classes
15\brief Box class.
16
17 - [Building boxes](\ref GEOB00)
18 - [Creation of boxes](\ref GEOB01)
19 - [Divisions of boxes](\ref GEOB02)
20
21All shape primitives inherit from this, their
22 constructor filling automatically the parameters of the box that bounds
23 the given shape. Defined by 6 parameters :
24
25```
26TGeoBBox(Double_t dx,Double_t dy,Double_t dz,Double_t *origin=0);
27```
28
29 - `fDX`, `fDY`, `fDZ` : half lengths on X, Y and Z axis
30 - `fOrigin[3]` : position of box origin
31
32\anchor GEOB00
33### Building boxes
34
35Normally a box has to be built only with 3 parameters: `DX,DY,DZ`
36representing the half-lengths on X, Y and Z-axes. In this case, the
37origin of the box will match the one of its reference frame and the box
38will range from: `-DX` to `DX` on X-axis, from `-DY` to `DY` on Y and
39from `-DZ` to `DZ` on Z. On the other hand, any other shape needs to
40compute and store the parameters of their minimal bounding box. The
41bounding boxes are essential to optimize navigation algorithms.
42Therefore all other primitives derive from **`TGeoBBox`**. Since the
43minimal bounding box is not necessary centered in the origin, any box
44allows an origin translation `(Ox`,`Oy`,`Oz)`. All primitive
45constructors automatically compute the bounding box parameters. Users
46should be aware that building a translated box that will represent a
47primitive shape by itself would affect any further positioning of other
48shapes inside. Therefore it is highly recommendable to build
49non-translated boxes as primitives and translate/rotate their
50corresponding volumes only during positioning stage.
51
52\anchor GEOB01
53#### Creation of boxes
54
55```
56 TGeoBBox *box = new TGeoBBox("BOX", 20, 30, 40);
57```
58
59Begin_Macro
60{
61 TCanvas *c = new TCanvas("c", "c",0,0,600,600);
62 new TGeoManager("box", "poza1");
63 TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
64 TGeoMedium *med = new TGeoMedium("MED",1,mat);
65 TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
66 gGeoManager->SetTopVolume(top);
67 TGeoVolume *vol = gGeoManager->MakeBox("BOX",med, 20,30,40);
68 vol->SetLineWidth(2);
69 top->AddNode(vol,1);
70 gGeoManager->CloseGeometry();
71 gGeoManager->SetNsegments(80);
72 top->Draw();
73 TView *view = gPad->GetView();
74 if (view) view->ShowAxis();
75}
76End_Macro
77
78A volume having a box shape can be built in one step:
79
80```
81 TGeoVolume *vbox = gGeoManager->MakeBox("vbox", ptrMed, 20,30,40);
82```
83
84\anchor GEOB02
85#### Divisions of boxes
86
87 Volumes having box shape can be divided with equal-length slices on
88X, Y or Z axis. The following options are supported:
89
90 - Dividing the full range of one axis in N slices
91```
92 TGeoVolume *divx = vbox->Divide("SLICEX", 1, N);
93```
94here `1` stands for the division axis (1-X, 2-Y, 3-Z)
95
96Begin_Macro
97{
98 TCanvas *c = new TCanvas("c", "c",0,0,600,600);
99 new TGeoManager("box", "poza1");
100 TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
101 TGeoMedium *med = new TGeoMedium("MED",1,mat);
102 TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
103 gGeoManager->SetTopVolume(top);
104 TGeoVolume *vol = gGeoManager->MakeBox("BOX",med, 20,30,40);
105 vol->SetLineWidth(2);
106 top->AddNode(vol,1);
107 TGeoVolume *divx = vol->Divide("SLICE",1,8,0,0);
108 gGeoManager->CloseGeometry();
109 gGeoManager->SetNsegments(80);
110 top->Draw();
111 TView *view = gPad->GetView();
112 if (view) view->ShowAxis();
113}
114End_Macro
115
116 - Dividing in a limited range - general case.
117```
118 TGeoVolume *divy = vbox->Divide("SLICEY",2,N,start,step);
119```
120 - start = starting offset within (-fDY, fDY)
121 - step = slicing step
122
123Begin_Macro
124{
125 TCanvas *c = new TCanvas("c", "c",0,0,600,600);
126 new TGeoManager("box", "poza1");
127 TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
128 TGeoMedium *med = new TGeoMedium("MED",1,mat);
129 TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
130 gGeoManager->SetTopVolume(top);
131 TGeoVolume *vol = gGeoManager->MakeBox("BOX",med, 20,30,40);
132 vol->SetLineWidth(2);
133 top->AddNode(vol,1);
134 TGeoVolume *divx = vol->Divide("SLICE",2,8,2,3);
135 gGeoManager->CloseGeometry();
136 gGeoManager->SetNsegments(80);
137 top->Draw();
138 TView *view = gPad->GetView();
139 if (view) view->ShowAxis();
140}
141End_Macro
142
143Both cases are supported by all shapes.
144See also class TGeoShape for utility methods provided by any particular shape.
145*/
146
147#include <iostream>
148
149#include "TGeoManager.h"
150#include "TGeoMatrix.h"
151#include "TGeoVolume.h"
152#include "TVirtualGeoPainter.h"
153#include "TGeoBBox.h"
154#include "TBuffer3D.h"
155#include "TBuffer3DTypes.h"
156#include "TMath.h"
157#include "TRandom.h"
158
159
160////////////////////////////////////////////////////////////////////////////////
161/// Default constructor
162
164{
166 fDX = fDY = fDZ = 0;
167 fOrigin[0] = fOrigin[1] = fOrigin[2] = 0.0;
168}
169
170////////////////////////////////////////////////////////////////////////////////
171/// Constructor where half-lengths are provided.
172
179
180////////////////////////////////////////////////////////////////////////////////
181/// Constructor with shape name.
182
189
190////////////////////////////////////////////////////////////////////////////////
191/// Constructor based on the array of parameters.
192/// - param[0] - half-length in x
193/// - param[1] - half-length in y
194/// - param[2] - half-length in z
195
197{
199 fOrigin[0] = fOrigin[1] = fOrigin[2] = 0.0;
200 SetDimensions(param);
201}
202
203////////////////////////////////////////////////////////////////////////////////
204/// Destructor
205
207
208////////////////////////////////////////////////////////////////////////////////
209/// Computes capacity of the shape in [length^3].
210
212{
213 return (8. * fDX * fDY * fDZ);
214}
215
216////////////////////////////////////////////////////////////////////////////////
217/// Computes normal to closest surface from POINT.
218
219void TGeoBBox::ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm) const
220{
221 memset(norm, 0, 3 * sizeof(Double_t));
222 Double_t saf[3];
223 Int_t i;
224 saf[0] = TMath::Abs(TMath::Abs(point[0] - fOrigin[0]) - fDX);
225 saf[1] = TMath::Abs(TMath::Abs(point[1] - fOrigin[1]) - fDY);
226 saf[2] = TMath::Abs(TMath::Abs(point[2] - fOrigin[2]) - fDZ);
227 i = TMath::LocMin(3, saf);
228 norm[i] = (dir[i] > 0) ? 1 : (-1);
229}
230
231////////////////////////////////////////////////////////////////////////////////
232/// Decides fast if the bounding box could be crossed by a vector.
233
234Bool_t TGeoBBox::CouldBeCrossed(const Double_t *point, const Double_t *dir) const
235{
236 Double_t mind = fDX;
237 if (fDY < mind)
238 mind = fDY;
239 if (fDZ < mind)
240 mind = fDZ;
241 Double_t dx = fOrigin[0] - point[0];
242 Double_t dy = fOrigin[1] - point[1];
243 Double_t dz = fOrigin[2] - point[2];
244 Double_t do2 = dx * dx + dy * dy + dz * dz;
245 if (do2 <= (mind * mind))
246 return kTRUE;
247 Double_t rmax2 = fDX * fDX + fDY * fDY + fDZ * fDZ;
248 if (do2 <= rmax2)
249 return kTRUE;
250 // inside bounding sphere
251 Double_t doct = dx * dir[0] + dy * dir[1] + dz * dir[2];
252 // leaving ray
253 if (doct <= 0)
254 return kFALSE;
255 Double_t dirnorm = dir[0] * dir[0] + dir[1] * dir[1] + dir[2] * dir[2];
256 if ((doct * doct) >= (do2 - rmax2) * dirnorm)
257 return kTRUE;
258 return kFALSE;
259}
260
261////////////////////////////////////////////////////////////////////////////////
262/// Compute closest distance from point px,py to each corner.
263
265{
266 const Int_t numPoints = 8;
267 return ShapeDistancetoPrimitive(numPoints, px, py);
268}
269
270////////////////////////////////////////////////////////////////////////////////
271/// Divide this box shape belonging to volume "voldiv" into ndiv equal volumes
272/// called divname, from start position with the given step. Returns pointer
273/// to created division cell volume. In case a wrong division axis is supplied,
274/// returns pointer to volume to be divided.
275
278{
279 TGeoShape *shape; //--- shape to be created
280 TGeoVolume *vol; //--- division volume to be created
281 TGeoVolumeMulti *vmulti; //--- generic divided volume
282 TGeoPatternFinder *finder; //--- finder to be attached
283 TString opt = ""; //--- option to be attached
284 Double_t end = start + ndiv * step;
285 switch (iaxis) {
286 case 1: //--- divide on X
287 shape = new TGeoBBox(step / 2., fDY, fDZ);
288 finder = new TGeoPatternX(voldiv, ndiv, start, end);
289 opt = "X";
290 break;
291 case 2: //--- divide on Y
292 shape = new TGeoBBox(fDX, step / 2., fDZ);
293 finder = new TGeoPatternY(voldiv, ndiv, start, end);
294 opt = "Y";
295 break;
296 case 3: //--- divide on Z
297 shape = new TGeoBBox(fDX, fDY, step / 2.);
298 finder = new TGeoPatternZ(voldiv, ndiv, start, end);
299 opt = "Z";
300 break;
301 default: Error("Divide", "Wrong axis type for division"); return nullptr;
302 }
303 vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
305 vmulti->AddVolume(vol);
306 voldiv->SetFinder(finder);
307 finder->SetDivIndex(voldiv->GetNdaughters());
308 for (Int_t ic = 0; ic < ndiv; ic++) {
309 voldiv->AddNodeOffset(vol, ic, start + step / 2. + ic * step, opt.Data());
310 ((TGeoNodeOffset *)voldiv->GetNodes()->At(voldiv->GetNdaughters() - 1))->SetFinder(finder);
311 }
312 return vmulti;
313}
314
315////////////////////////////////////////////////////////////////////////////////
316/// Compute bounding box - nothing to do in this case.
317
319
320////////////////////////////////////////////////////////////////////////////////
321/// Test if point is inside this shape.
322
324{
325 if (TMath::Abs(point[2] - fOrigin[2]) > fDZ)
326 return kFALSE;
327 if (TMath::Abs(point[0] - fOrigin[0]) > fDX)
328 return kFALSE;
329 if (TMath::Abs(point[1] - fOrigin[1]) > fDY)
330 return kFALSE;
331 return kTRUE;
332}
333
334////////////////////////////////////////////////////////////////////////////////
335/// Static method to check if point[3] is located inside a box of having dx, dy, dz
336/// as half-lengths.
337
339{
340 if (TMath::Abs(point[2] - origin[2]) > dz)
341 return kFALSE;
342 if (TMath::Abs(point[0] - origin[0]) > dx)
343 return kFALSE;
344 if (TMath::Abs(point[1] - origin[1]) > dy)
345 return kFALSE;
346 return kTRUE;
347}
348
349////////////////////////////////////////////////////////////////////////////////
350/// Compute distance from inside point to surface of the box.
351/// Boundary safe algorithm.
352
355{
356 Double_t s, smin, saf[6];
357 Double_t newpt[3];
358 Int_t i;
359 for (i = 0; i < 3; i++)
360 newpt[i] = point[i] - fOrigin[i];
361 saf[0] = fDX + newpt[0];
362 saf[1] = fDX - newpt[0];
363 saf[2] = fDY + newpt[1];
364 saf[3] = fDY - newpt[1];
365 saf[4] = fDZ + newpt[2];
366 saf[5] = fDZ - newpt[2];
367 if (iact < 3 && safe) {
368 smin = saf[0];
369 // compute safe distance
370 for (i = 1; i < 6; i++)
371 if (saf[i] < smin)
372 smin = saf[i];
373 *safe = smin;
374 if (smin < 0)
375 *safe = 0.0;
376 if (iact == 0)
377 return TGeoShape::Big();
378 if (iact == 1 && step < *safe)
379 return TGeoShape::Big();
380 }
381 // compute distance to surface
383 for (i = 0; i < 3; i++) {
384 if (dir[i] != 0) {
385 s = (dir[i] > 0) ? (saf[(i << 1) + 1] / dir[i]) : (-saf[i << 1] / dir[i]);
386 if (s < 0)
387 return 0.0;
388 if (s < smin)
389 smin = s;
390 }
391 }
392 return smin;
393}
394
395////////////////////////////////////////////////////////////////////////////////
396/// Compute distance from inside point to surface of the box.
397/// Boundary safe algorithm.
398
400 const Double_t *origin, Double_t /*stepmax*/)
401{
402 Double_t s, smin, saf[6];
403 Double_t newpt[3];
404 Int_t i;
405 for (i = 0; i < 3; i++)
406 newpt[i] = point[i] - origin[i];
407 saf[0] = dx + newpt[0];
408 saf[1] = dx - newpt[0];
409 saf[2] = dy + newpt[1];
410 saf[3] = dy - newpt[1];
411 saf[4] = dz + newpt[2];
412 saf[5] = dz - newpt[2];
413 // compute distance to surface
415 for (i = 0; i < 3; i++) {
416 if (dir[i] != 0) {
417 s = (dir[i] > 0) ? (saf[(i << 1) + 1] / dir[i]) : (-saf[i << 1] / dir[i]);
418 if (s < 0)
419 return 0.0;
420 if (s < smin)
421 smin = s;
422 }
423 }
424 return smin;
425}
426
427////////////////////////////////////////////////////////////////////////////////
428/// Compute distance from outside point to surface of the box.
429/// Boundary safe algorithm.
430
433{
434 Bool_t in = kTRUE;
435 Double_t saf[3];
436 Double_t par[3];
437 Double_t newpt[3];
438 Int_t i, j;
439 for (i = 0; i < 3; i++)
440 newpt[i] = point[i] - fOrigin[i];
441 par[0] = fDX;
442 par[1] = fDY;
443 par[2] = fDZ;
444 for (i = 0; i < 3; i++) {
445 saf[i] = TMath::Abs(newpt[i]) - par[i];
446 if (saf[i] >= step)
447 return TGeoShape::Big();
448 if (in && saf[i] > 0)
449 in = kFALSE;
450 }
451 if (iact < 3 && safe) {
452 // compute safe distance
453 if (in) {
454 *safe = 0.0;
455 } else {
456 *safe = saf[0];
457 if (saf[1] > *safe)
458 *safe = saf[1];
459 if (saf[2] > *safe)
460 *safe = saf[2];
461 }
462 if (iact == 0)
463 return TGeoShape::Big();
464 if (iact == 1 && step < *safe)
465 return TGeoShape::Big();
466 }
467 // compute distance from point to box
468 // protection in case point is actually inside box
469 if (in) {
470 j = 0;
471 Double_t ss = saf[0];
472 if (saf[1] > ss) {
473 ss = saf[1];
474 j = 1;
475 }
476 if (saf[2] > ss)
477 j = 2;
478 if (newpt[j] * dir[j] > 0)
479 return TGeoShape::Big(); // in fact exiting
480 return 0.0;
481 }
482 for (i = 0; i < 3; i++) {
483 if (saf[i] < 0)
484 continue;
485 if (newpt[i] * dir[i] >= 0)
486 continue;
487 Double_t snxt = saf[i] / TMath::Abs(dir[i]);
488 Int_t ibreak = 0;
489 for (j = 0; j < 3; j++) {
490 if (j == i)
491 continue;
492 Double_t coord = newpt[j] + snxt * dir[j];
493 if (TMath::Abs(coord) > par[j]) {
494 ibreak = 1;
495 break;
496 }
497 }
498 if (!ibreak)
499 return snxt;
500 }
501 return TGeoShape::Big();
502}
503
504////////////////////////////////////////////////////////////////////////////////
505/// Compute distance from outside point to surface of the box.
506/// Boundary safe algorithm.
507
510{
511 Bool_t in = kTRUE;
512 Double_t saf[3];
513 Double_t par[3];
514 Double_t newpt[3];
515 Int_t i, j;
516 for (i = 0; i < 3; i++)
517 newpt[i] = point[i] - origin[i];
518 par[0] = dx;
519 par[1] = dy;
520 par[2] = dz;
521 for (i = 0; i < 3; i++) {
522 saf[i] = TMath::Abs(newpt[i]) - par[i];
523 if (saf[i] >= stepmax)
524 return TGeoShape::Big();
525 if (in && saf[i] > 0)
526 in = kFALSE;
527 }
528 // In case point is inside return ZERO
529 if (in)
530 return 0.0;
532 Int_t ibreak = 0;
533 for (i = 0; i < 3; i++) {
534 if (saf[i] < 0)
535 continue;
536 if (newpt[i] * dir[i] >= 0)
537 continue;
538 snxt = saf[i] / TMath::Abs(dir[i]);
539 ibreak = 0;
540 for (j = 0; j < 3; j++) {
541 if (j == i)
542 continue;
543 coord = newpt[j] + snxt * dir[j];
544 if (TMath::Abs(coord) > par[j]) {
545 ibreak = 1;
546 break;
547 }
548 }
549 if (!ibreak)
550 return snxt;
551 }
552 return TGeoShape::Big();
553}
554
555////////////////////////////////////////////////////////////////////////////////
556/// Returns name of axis IAXIS.
557
559{
560 switch (iaxis) {
561 case 1: return "X";
562 case 2: return "Y";
563 case 3: return "Z";
564 default: return "UNDEFINED";
565 }
566}
567
568////////////////////////////////////////////////////////////////////////////////
569/// Get range of shape for a given axis.
570
572{
573 xlo = 0;
574 xhi = 0;
575 Double_t dx = 0;
576 switch (iaxis) {
577 case 1:
578 xlo = fOrigin[0] - fDX;
579 xhi = fOrigin[0] + fDX;
580 dx = 2 * fDX;
581 return dx;
582 case 2:
583 xlo = fOrigin[1] - fDY;
584 xhi = fOrigin[1] + fDY;
585 dx = 2 * fDY;
586 return dx;
587 case 3:
588 xlo = fOrigin[2] - fDZ;
589 xhi = fOrigin[2] + fDZ;
590 dx = 2 * fDZ;
591 return dx;
592 }
593 return dx;
594}
595
596////////////////////////////////////////////////////////////////////////////////
597/// Fill vector param[4] with the bounding cylinder parameters. The order
598/// is the following : Rmin, Rmax, Phi1, Phi2
599
601{
602 param[0] = 0.; // Rmin
603 param[1] = fDX * fDX + fDY * fDY; // Rmax
604 param[2] = 0.; // Phi1
605 param[3] = 360.; // Phi2
606}
607
608////////////////////////////////////////////////////////////////////////////////
609/// Get area in internal units of the facet with a given index.
610/// Possible index values:
611/// - 0 - all facets together
612/// - 1 to 6 - facet index from bottom to top Z
613
615{
616 Double_t area = 0.;
617 switch (index) {
618 case 0: area = 8. * (fDX * fDY + fDX * fDZ + fDY * fDZ); return area;
619 case 1:
620 case 6: area = 4. * fDX * fDY; return area;
621 case 2:
622 case 4: area = 4. * fDX * fDZ; return area;
623 case 3:
624 case 5: area = 4. * fDY * fDZ; return area;
625 }
626 return area;
627}
628
629////////////////////////////////////////////////////////////////////////////////
630/// Fills array with n random points located on the surface of indexed facet.
631/// The output array must be provided with a length of minimum 3*npoints. Returns
632/// true if operation succeeded.
633/// Possible index values:
634/// - 0 - all facets together
635/// - 1 to 6 - facet index from bottom to top Z
636
638{
640 return kFALSE;
641 Double_t surf[6];
642 Double_t area = 0.;
643 if (index == 0) {
644 for (Int_t isurf = 0; isurf < 6; isurf++) {
646 if (isurf > 0)
647 surf[isurf] += surf[isurf - 1];
648 }
649 area = surf[5];
650 }
651
652 for (Int_t i = 0; i < npoints; i++) {
653 // Generate randomly a surface index if needed.
654 Double_t *point = &array[3 * i];
656 if (surfindex == 0) {
657 Double_t val = area * gRandom->Rndm();
658 surfindex = 2 + TMath::BinarySearch(6, surf, val);
659 if (surfindex > 6)
660 surfindex = 6;
661 }
662 switch (surfindex) {
663 case 1:
664 point[0] = -fDX + 2 * fDX * gRandom->Rndm();
665 point[1] = -fDY + 2 * fDY * gRandom->Rndm();
666 point[2] = -fDZ;
667 break;
668 case 2:
669 point[0] = -fDX + 2 * fDX * gRandom->Rndm();
670 point[1] = -fDY;
671 point[2] = -fDZ + 2 * fDZ * gRandom->Rndm();
672 break;
673 case 3:
674 point[0] = -fDX;
675 point[1] = -fDY + 2 * fDY * gRandom->Rndm();
676 point[2] = -fDZ + 2 * fDZ * gRandom->Rndm();
677 break;
678 case 4:
679 point[0] = -fDX + 2 * fDX * gRandom->Rndm();
680 point[1] = fDY;
681 point[2] = -fDZ + 2 * fDZ * gRandom->Rndm();
682 break;
683 case 5:
684 point[0] = fDX;
685 point[1] = -fDY + 2 * fDY * gRandom->Rndm();
686 point[2] = -fDZ + 2 * fDZ * gRandom->Rndm();
687 break;
688 case 6:
689 point[0] = -fDX + 2 * fDX * gRandom->Rndm();
690 point[1] = -fDY + 2 * fDY * gRandom->Rndm();
691 point[2] = fDZ;
692 break;
693 }
694 }
695 return kTRUE;
696}
697
698////////////////////////////////////////////////////////////////////////////////
699/// Fills array with n random points located on the line segments of the shape mesh.
700/// The output array must be provided with a length of minimum 3*npoints. Returns
701/// true if operation is implemented.
702
704{
705 if (npoints < GetNmeshVertices()) {
706 Error("GetPointsOnSegments", "You should require at least %d points", GetNmeshVertices());
707 return kFALSE;
708 }
710 Int_t npnts = buff.NbPnts();
711 Int_t nsegs = buff.NbSegs();
712 // Copy buffered points in the array
713 memcpy(array, buff.fPnts, 3 * npnts * sizeof(Double_t));
715 Int_t icrt = 3 * npnts;
717 Double_t *p0, *p1;
718 Double_t x, y, z, dx, dy, dz;
719 for (Int_t i = 0; i < nsegs; i++) {
720 p0 = &array[3 * buff.fSegs[3 * i + 1]];
721 p1 = &array[3 * buff.fSegs[3 * i + 2]];
722 if (i == (nsegs - 1))
724 dx = (p1[0] - p0[0]) / (nperseg + 1);
725 dy = (p1[1] - p0[1]) / (nperseg + 1);
726 dz = (p1[2] - p0[2]) / (nperseg + 1);
727 for (Int_t j = 0; j < nperseg; j++) {
728 x = p0[0] + (j + 1) * dx;
729 y = p0[1] + (j + 1) * dy;
730 z = p0[2] + (j + 1) * dz;
731 array[icrt++] = x;
732 array[icrt++] = y;
733 array[icrt++] = z;
734 ipoints--;
735 }
736 }
737 return kTRUE;
738}
739
740////////////////////////////////////////////////////////////////////////////////
741/// Fills real parameters of a positioned box inside this one. Returns 0 if successful.
742
744{
745 dx = dy = dz = 0;
746 if (mat->IsRotation()) {
747 Error("GetFittingBox", "cannot handle parametrized rotated volumes");
748 return 1; // ### rotation not accepted ###
749 }
750 //--> translate the origin of the parametrized box to the frame of this box.
751 Double_t origin[3];
752 mat->LocalToMaster(parambox->GetOrigin(), origin);
753 if (!Contains(origin)) {
754 Error("GetFittingBox", "wrong matrix - parametrized box is outside this");
755 return 1; // ### wrong matrix ###
756 }
757 //--> now we have to get the valid range for all parametrized axis
758 Double_t xlo = 0, xhi = 0;
759 Double_t dd[3];
760 dd[0] = parambox->GetDX();
761 dd[1] = parambox->GetDY();
762 dd[2] = parambox->GetDZ();
763 for (Int_t iaxis = 0; iaxis < 3; iaxis++) {
764 if (dd[iaxis] >= 0)
765 continue;
766 TGeoBBox::GetAxisRange(iaxis + 1, xlo, xhi);
767 //-> compute best fitting parameter
768 dd[iaxis] = TMath::Min(origin[iaxis] - xlo, xhi - origin[iaxis]);
769 if (dd[iaxis] < 0) {
770 Error("GetFittingBox", "wrong matrix");
771 return 1;
772 }
773 }
774 dx = dd[0];
775 dy = dd[1];
776 dz = dd[2];
777 return 0;
778}
779
780////////////////////////////////////////////////////////////////////////////////
781/// In case shape has some negative parameters, these has to be computed
782/// in order to fit the mother
783
785{
787 return nullptr;
788 Double_t dx, dy, dz;
789 Int_t ierr = mother->GetFittingBox(this, mat, dx, dy, dz);
790 if (ierr) {
791 Error("GetMakeRuntimeShape", "cannot fit this to mother");
792 return nullptr;
793 }
794 return (new TGeoBBox(dx, dy, dz));
795}
796
797////////////////////////////////////////////////////////////////////////////////
798/// Returns numbers of vertices, segments and polygons composing the shape mesh.
799
801{
802 nvert = 8;
803 nsegs = 12;
804 npols = 6;
805}
806
807////////////////////////////////////////////////////////////////////////////////
808/// Prints shape parameters
809
811{
812 printf("*** Shape %s: TGeoBBox ***\n", GetName());
813 printf(" dX = %11.5f\n", fDX);
814 printf(" dY = %11.5f\n", fDY);
815 printf(" dZ = %11.5f\n", fDZ);
816 printf(" origin: x=%11.5f y=%11.5f z=%11.5f\n", fOrigin[0], fOrigin[1], fOrigin[2]);
817}
818
819////////////////////////////////////////////////////////////////////////////////
820/// Creates a TBuffer3D describing *this* shape.
821/// Coordinates are in local reference frame.
822
824{
825 TBuffer3D *buff = new TBuffer3D(TBuffer3DTypes::kGeneric, 8, 24, 12, 36, 6, 36);
826 if (buff) {
827 SetPoints(buff->fPnts);
829 }
830
831 return buff;
832}
833
834////////////////////////////////////////////////////////////////////////////////
835/// Fills TBuffer3D structure for segments and polygons.
836
838{
840
841 buff.fSegs[0] = c;
842 buff.fSegs[1] = 0;
843 buff.fSegs[2] = 1;
844 buff.fSegs[3] = c + 1;
845 buff.fSegs[4] = 1;
846 buff.fSegs[5] = 2;
847 buff.fSegs[6] = c + 1;
848 buff.fSegs[7] = 2;
849 buff.fSegs[8] = 3;
850 buff.fSegs[9] = c;
851 buff.fSegs[10] = 3;
852 buff.fSegs[11] = 0;
853 buff.fSegs[12] = c + 2;
854 buff.fSegs[13] = 4;
855 buff.fSegs[14] = 5;
856 buff.fSegs[15] = c + 2;
857 buff.fSegs[16] = 5;
858 buff.fSegs[17] = 6;
859 buff.fSegs[18] = c + 3;
860 buff.fSegs[19] = 6;
861 buff.fSegs[20] = 7;
862 buff.fSegs[21] = c + 3;
863 buff.fSegs[22] = 7;
864 buff.fSegs[23] = 4;
865 buff.fSegs[24] = c;
866 buff.fSegs[25] = 0;
867 buff.fSegs[26] = 4;
868 buff.fSegs[27] = c + 2;
869 buff.fSegs[28] = 1;
870 buff.fSegs[29] = 5;
871 buff.fSegs[30] = c + 1;
872 buff.fSegs[31] = 2;
873 buff.fSegs[32] = 6;
874 buff.fSegs[33] = c + 3;
875 buff.fSegs[34] = 3;
876 buff.fSegs[35] = 7;
877
878 buff.fPols[0] = c;
879 buff.fPols[1] = 4;
880 buff.fPols[2] = 0;
881 buff.fPols[3] = 9;
882 buff.fPols[4] = 4;
883 buff.fPols[5] = 8;
884 buff.fPols[6] = c + 1;
885 buff.fPols[7] = 4;
886 buff.fPols[8] = 1;
887 buff.fPols[9] = 10;
888 buff.fPols[10] = 5;
889 buff.fPols[11] = 9;
890 buff.fPols[12] = c;
891 buff.fPols[13] = 4;
892 buff.fPols[14] = 2;
893 buff.fPols[15] = 11;
894 buff.fPols[16] = 6;
895 buff.fPols[17] = 10;
896 buff.fPols[18] = c + 1;
897 buff.fPols[19] = 4;
898 buff.fPols[20] = 3;
899 buff.fPols[21] = 8;
900 buff.fPols[22] = 7;
901 buff.fPols[23] = 11;
902 buff.fPols[24] = c + 2;
903 buff.fPols[25] = 4;
904 buff.fPols[26] = 0;
905 buff.fPols[27] = 3;
906 buff.fPols[28] = 2;
907 buff.fPols[29] = 1;
908 buff.fPols[30] = c + 3;
909 buff.fPols[31] = 4;
910 buff.fPols[32] = 4;
911 buff.fPols[33] = 5;
912 buff.fPols[34] = 6;
913 buff.fPols[35] = 7;
914}
915
916////////////////////////////////////////////////////////////////////////////////
917/// Computes the closest distance from given point to this shape.
918
920{
922 if (in) {
923 safe = fDX - TMath::Abs(point[0] - fOrigin[0]);
924 safy = fDY - TMath::Abs(point[1] - fOrigin[1]);
925 safz = fDZ - TMath::Abs(point[2] - fOrigin[2]);
926 if (safy < safe)
927 safe = safy;
928 if (safz < safe)
929 safe = safz;
930 } else {
931 safe = -fDX + TMath::Abs(point[0] - fOrigin[0]);
932 safy = -fDY + TMath::Abs(point[1] - fOrigin[1]);
933 safz = -fDZ + TMath::Abs(point[2] - fOrigin[2]);
934 if (safy > safe)
935 safe = safy;
936 if (safz > safe)
937 safe = safz;
938 }
939 return safe;
940}
941
942////////////////////////////////////////////////////////////////////////////////
943/// Save a primitive as a C++ statement(s) on output stream "out".
944
945void TGeoBBox::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
946{
948 return;
949 out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
950 out << " dx = " << fDX << ";" << std::endl;
951 out << " dy = " << fDY << ";" << std::endl;
952 out << " dz = " << fDZ << ";" << std::endl;
955 out << " origin[0] = " << fOrigin[0] << ";" << std::endl;
956 out << " origin[1] = " << fOrigin[1] << ";" << std::endl;
957 out << " origin[2] = " << fOrigin[2] << ";" << std::endl;
958 out << " TGeoShape *" << GetPointerName() << " = new TGeoBBox(\"" << GetName() << "\", dx,dy,dz,origin);"
959 << std::endl;
960 } else {
961 out << " TGeoShape *" << GetPointerName() << " = new TGeoBBox(\"" << GetName() << "\", dx,dy,dz);" << std::endl;
962 }
964}
965
966////////////////////////////////////////////////////////////////////////////////
967/// Set parameters of the box.
968
970{
971 fDX = dx;
972 fDY = dy;
973 fDZ = dz;
974 if (origin) {
975 fOrigin[0] = origin[0];
976 fOrigin[1] = origin[1];
977 fOrigin[2] = origin[2];
978 }
981 return;
982 if ((fDX < 0) || (fDY < 0) || (fDZ < 0))
984}
985
986////////////////////////////////////////////////////////////////////////////////
987/// Set dimensions based on the array of parameters
988/// param[0] - half-length in x
989/// param[1] - half-length in y
990/// param[2] - half-length in z
991
993{
994 if (!param) {
995 Error("SetDimensions", "null parameters");
996 return;
997 }
998 fDX = param[0];
999 fDY = param[1];
1000 fDZ = param[2];
1003 return;
1004 if ((fDX < 0) || (fDY < 0) || (fDZ < 0))
1006}
1007
1008////////////////////////////////////////////////////////////////////////////////
1009/// Fill box vertices to an array.
1010
1015
1016////////////////////////////////////////////////////////////////////////////////
1017/// Fill box points.
1018
1020{
1021 if (!points)
1022 return;
1023 Double_t xmin, xmax, ymin, ymax, zmin, zmax;
1024 xmin = -fDX + fOrigin[0];
1025 xmax = fDX + fOrigin[0];
1026 ymin = -fDY + fOrigin[1];
1027 ymax = fDY + fOrigin[1];
1028 zmin = -fDZ + fOrigin[2];
1029 zmax = fDZ + fOrigin[2];
1030 points[0] = xmin;
1031 points[1] = ymin;
1032 points[2] = zmin;
1033 points[3] = xmin;
1034 points[4] = ymax;
1035 points[5] = zmin;
1036 points[6] = xmax;
1037 points[7] = ymax;
1038 points[8] = zmin;
1039 points[9] = xmax;
1040 points[10] = ymin;
1041 points[11] = zmin;
1042 points[12] = xmin;
1043 points[13] = ymin;
1044 points[14] = zmax;
1045 points[15] = xmin;
1046 points[16] = ymax;
1047 points[17] = zmax;
1048 points[18] = xmax;
1049 points[19] = ymax;
1050 points[20] = zmax;
1051 points[21] = xmax;
1052 points[22] = ymin;
1053 points[23] = zmax;
1054}
1055
1056////////////////////////////////////////////////////////////////////////////////
1057/// Fill box points.
1058
1060{
1061 if (!points)
1062 return;
1063 Double_t xmin, xmax, ymin, ymax, zmin, zmax;
1064 xmin = -fDX + fOrigin[0];
1065 xmax = fDX + fOrigin[0];
1066 ymin = -fDY + fOrigin[1];
1067 ymax = fDY + fOrigin[1];
1068 zmin = -fDZ + fOrigin[2];
1069 zmax = fDZ + fOrigin[2];
1070 points[0] = xmin;
1071 points[1] = ymin;
1072 points[2] = zmin;
1073 points[3] = xmin;
1074 points[4] = ymax;
1075 points[5] = zmin;
1076 points[6] = xmax;
1077 points[7] = ymax;
1078 points[8] = zmin;
1079 points[9] = xmax;
1080 points[10] = ymin;
1081 points[11] = zmin;
1082 points[12] = xmin;
1083 points[13] = ymin;
1084 points[14] = zmax;
1085 points[15] = xmin;
1086 points[16] = ymax;
1087 points[17] = zmax;
1088 points[18] = xmax;
1089 points[19] = ymax;
1090 points[20] = zmax;
1091 points[21] = xmax;
1092 points[22] = ymin;
1093 points[23] = zmax;
1094}
1095
1096////////////////////////////////////////////////////////////////////////////////
1097////// fill size of this 3-D object
1098//// TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter();
1099//// if (painter) painter->AddSize3D(8, 12, 6);
1100
1101void TGeoBBox::Sizeof3D() const {}
1102
1103////////////////////////////////////////////////////////////////////////////////
1104/// Fills a static 3D buffer and returns a reference.
1105
1107{
1108 static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
1109
1111
1112 // TODO: A box itself has has nothing more as already described
1113 // by bounding box. How will viewer interpret?
1115 if (buffer.SetRawSizes(8, 3 * 8, 12, 3 * 12, 6, 6 * 6)) {
1117 }
1118 }
1120 SetPoints(buffer.fPnts);
1121 if (!buffer.fLocalFrame) {
1122 TransformPoints(buffer.fPnts, buffer.NbPnts());
1123 }
1124
1125 SetSegsAndPols(buffer);
1127 }
1128
1129 return buffer;
1130}
1131
1132////////////////////////////////////////////////////////////////////////////////
1133/// Fills the supplied buffer, with sections in desired frame
1134/// See TBuffer3D.h for explanation of sections, frame etc.
1135
1137{
1139
1141 Double_t halfLengths[3] = {fDX, fDY, fDZ};
1143
1144 if (!buffer.fLocalFrame) {
1145 TransformPoints(buffer.fBBVertex[0], 8);
1146 }
1148 }
1149}
1150
1151////////////////////////////////////////////////////////////////////////////////
1152/// Check the inside status for each of the points in the array.
1153/// Input: Array of point coordinates + vector size
1154/// Output: Array of Booleans for the inside of each point
1155
1157{
1158 for (Int_t i = 0; i < vecsize; i++)
1159 inside[i] = Contains(&points[3 * i]);
1160}
1161
1162////////////////////////////////////////////////////////////////////////////////
1163/// Compute the normal for an array o points so that norm.dot.dir is positive
1164/// Input: Arrays of point coordinates and directions + vector size
1165/// Output: Array of normal directions
1166
1168{
1169 for (Int_t i = 0; i < vecsize; i++)
1170 ComputeNormal(&points[3 * i], &dirs[3 * i], &norms[3 * i]);
1171}
1172
1173////////////////////////////////////////////////////////////////////////////////
1174/// Compute distance from array of input points having directions specified by dirs. Store output in dists
1175
1177 Double_t *step) const
1178{
1179 for (Int_t i = 0; i < vecsize; i++)
1180 dists[i] = DistFromInside(&points[3 * i], &dirs[3 * i], 3, step[i]);
1181}
1182
1183////////////////////////////////////////////////////////////////////////////////
1184/// Compute distance from array of input points having directions specified by dirs. Store output in dists
1185
1187 Double_t *step) const
1188{
1189 for (Int_t i = 0; i < vecsize; i++)
1190 dists[i] = DistFromOutside(&points[3 * i], &dirs[3 * i], 3, step[i]);
1191}
1192
1193////////////////////////////////////////////////////////////////////////////////
1194/// Compute safe distance from each of the points in the input array.
1195/// Input: Array of point coordinates, array of statuses for these points, size of the arrays
1196/// Output: Safety values
1197
1199{
1200 for (Int_t i = 0; i < vecsize; i++)
1201 safe[i] = Safety(&points[3 * i], inside[i]);
1202}
#define c(i)
Definition RSha256.hxx:101
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 TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t index
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t points
char name[80]
Definition TGX11.cxx:110
R__EXTERN TGeoManager * gGeoManager
float xmin
float ymin
float xmax
float ymax
R__EXTERN TRandom * gRandom
Definition TRandom.h:62
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
@ kBoundingBox
Definition TBuffer3D.h:51
void SetSectionsValid(UInt_t mask)
Definition TBuffer3D.h:65
Bool_t fLocalFrame
Definition TBuffer3D.h:90
void SetAABoundingBox(const Double_t origin[3], const Double_t halfLengths[3])
Set fBBVertex in kBoundingBox section to a axis aligned (local) BB using supplied origin and box half...
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
Double_t fBBVertex[8][3]
Definition TBuffer3D.h:108
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...
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save a primitive as a C++ statement(s) on output stream "out".
Definition TGeoBBox.cxx:945
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 box.
Definition TGeoBBox.cxx:354
Double_t fDX
Definition TGeoBBox.h:20
const char * GetAxisName(Int_t iaxis) const override
Returns name of axis IAXIS.
Definition TGeoBBox.cxx:558
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...
void GetBoundingCylinder(Double_t *param) const override
Fill vector param[4] with the bounding cylinder parameters.
Definition TGeoBBox.cxx:600
void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm) const override
Computes normal to closest surface from POINT.
Definition TGeoBBox.cxx:219
virtual Bool_t GetPointsOnFacet(Int_t index, Int_t npoints, Double_t *array) const
Fills array with n random points located on the surface of indexed facet.
Definition TGeoBBox.cxx:637
void SetPoints(Double_t *points) const override
Fill box points.
void SetBoxDimensions(Double_t dx, Double_t dy, Double_t dz, Double_t *origin=nullptr)
Set parameters of the box.
Definition TGeoBBox.cxx:969
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 Capacity() const override
Computes capacity of the shape in [length^3].
Definition TGeoBBox.cxx:211
virtual Double_t GetFacetArea(Int_t index=0) const
Get area in internal units of the facet with a given index.
Definition TGeoBBox.cxx:614
Int_t DistancetoPrimitive(Int_t px, Int_t py) override
Compute closest distance from point px,py to each corner.
Definition TGeoBBox.cxx:264
Bool_t CouldBeCrossed(const Double_t *point, const Double_t *dir) const override
Decides fast if the bounding box could be crossed by a vector.
Definition TGeoBBox.cxx:234
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.
void SetBoxPoints(Double_t *points) const
Fill box vertices to an array.
Bool_t GetPointsOnSegments(Int_t npoints, Double_t *array) const override
Fills array with n random points located on the line segments of the shape mesh.
Definition TGeoBBox.cxx:703
Int_t GetFittingBox(const TGeoBBox *parambox, TGeoMatrix *mat, Double_t &dx, Double_t &dy, Double_t &dz) const override
Fills real parameters of a positioned box inside this one. Returns 0 if successful.
Definition TGeoBBox.cxx:743
Double_t fOrigin[3]
Definition TGeoBBox.h:23
Double_t GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const override
Get range of shape for a given axis.
Definition TGeoBBox.cxx:571
TBuffer3D * MakeBuffer3D() const override
Creates a TBuffer3D describing this shape.
Definition TGeoBBox.cxx:823
void InspectShape() const override
Prints shape parameters.
Definition TGeoBBox.cxx:810
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...
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
void GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const override
Returns numbers of vertices, segments and polygons composing the shape mesh.
Definition TGeoBBox.cxx:800
void Sizeof3D() const override
TGeoShape * GetMakeRuntimeShape(TGeoShape *mother, TGeoMatrix *mat) const override
In case shape has some negative parameters, these has to be computed in order to fit the mother.
Definition TGeoBBox.cxx:784
Double_t fDY
Definition TGeoBBox.h:21
void SetSegsAndPols(TBuffer3D &buffer) const override
Fills TBuffer3D structure for segments and polygons.
Definition TGeoBBox.cxx:837
void SetDimensions(Double_t *param) override
Set dimensions based on the array of parameters param[0] - half-length in x param[1] - half-length in...
Definition TGeoBBox.cxx:992
const TBuffer3D & GetBuffer3D(Int_t reqSections, Bool_t localFrame) const override
Fills a static 3D buffer and returns a reference.
Int_t GetNmeshVertices() const override
Definition TGeoBBox.h:75
Double_t fDZ
Definition TGeoBBox.h:22
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...
~TGeoBBox() override
Destructor.
Definition TGeoBBox.cxx:206
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 ComputeBBox() override
Compute bounding box - nothing to do in this case.
Definition TGeoBBox.cxx:318
TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step) override
Divide this box shape belonging to volume "voldiv" into ndiv equal volumes called divname,...
Definition TGeoBBox.cxx:277
TGeoBBox()
Default constructor.
Definition TGeoBBox.cxx:163
TGeoVolumeMulti * MakeVolumeMulti(const char *name, TGeoMedium *medium)
Make a TGeoVolumeMulti handling a list of volumes.
Geometrical transformation package.
Definition TGeoMatrix.h:38
Node containing an offset.
Definition TGeoNode.h:184
base finder class for patterns. A pattern is specifying a division type
a Y axis divison pattern
a Z axis divison pattern
Base abstract class for all shapes.
Definition TGeoShape.h:25
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.
virtual void FillBuffer3D(TBuffer3D &buffer, Int_t reqSections, Bool_t localFrame) const
Fill the supplied buffer, with sections in desired frame See TBuffer3D.h for explanation of sections,...
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
@ kGeoRunTimeShape
Definition TGeoShape.h:40
static Double_t Tolerance()
Definition TGeoShape.h:97
Bool_t TestShapeBit(UInt_t f) const
Definition TGeoShape.h:175
Volume families.
Definition TGeoVolume.h:266
TGeoVolume, TGeoVolumeMulti, TGeoVolumeAssembly are the volume classes.
Definition TGeoVolume.h:43
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
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
Double_t Rndm() override
Machine independent random number generator.
Definition TRandom.cxx:558
Basic string class.
Definition TString.h:138
const char * Data() const
Definition TString.h:384
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
Long64_t LocMin(Long64_t n, const T *a)
Returns index of array with the minimum element.
Definition TMath.h:993
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