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