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