Logo ROOT   6.08/07
Reference Guide
TGeoShape.cxx
Go to the documentation of this file.
1 // @(#)root/geom:$Id$
2 // Author: Andrei Gheata 31/01/02
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 /** \class TGeoShape
13 \ingroup Geometry_classes
14 Base abstract class for all shapes.
15 
16  Shapes are geometrical objects that provide the basic modelling
17 functionality. They provide the definition of the LOCAL frame of coordinates,
18 with respect to which they are defined. Any implementation of a shape deriving
19 from the base TGeoShape class has to provide methods for :
20 
21  - finding out if a point defined in their local frame is or not contained
22  inside;
23  - computing the distance from a local point to getting outside/entering the
24  shape, given a known direction;
25  - computing the maximum distance in any direction from a local point that
26  does NOT result in a boundary crossing of the shape (safe distance);
27  - computing the cosines of the normal vector to the crossed shape surface,
28  given a starting local point and an ongoing direction.
29  All the features above are globally managed by the modeller in order to
30  provide navigation functionality. In addition to those, shapes have also to
31  implement additional specific abstract methods :
32  - computation of the minimal box bounding the shape, given that this box have
33  to be aligned with the local coordinates;
34  - algorithms for dividing the shape along a given axis and producing resulting
35  divisions volumes.
36 
37  The modeler currently provides a set of 16 basic shapes, which we will call
38 primitives. It also provides a special class allowing the creation of shapes
39 made as a result of boolean operations between primitives. These are called
40 composite shapes and the composition operation can be recursive (composition
41 of composites). This allows the creation of a quite large number of different
42 shape topologies and combinations.
43 
44  Shapes are named objects and register themselves to the manager class at
45 creation time. This is responsible for their final deletion. Shapes
46 can be created without name if their retrieval by name is no needed. Generally
47 shapes are objects that are useful only at geometry creation stage. The pointer
48 to a shape is in fact needed only when referring to a given volume and it is
49 always accessible at that level. A shape may be referenced by several volumes,
50 therefore its deletion is not possible once volumes were defined based on it.
51 
52 ### Creating shapes
53 
54  Shape objects embed only the minimum set of parameters that are fully
55 describing a valid physical shape. For instance, a tube is represented by
56 its half length, the minimum radius and the maximum radius. Shapes are used
57 together with media in order to create volumes, which in their turn
58 are the main components of the geometrical tree. A specific shape can be created
59 stand-alone :
60 
61 ~~~ {.cpp}
62  TGeoBBox *box = new TGeoBBox("s_box", halfX, halfY, halfZ); // named
63  TGeoTube *tub = new TGeoTube(rmin, rmax, halfZ); // no name
64  ... (see each specific shape constructors)
65 ~~~
66 
67  Sometimes it is much easier to create a volume having a given shape in one
68 step, since shapes are not directly linked in the geometrical tree but volumes
69 are :
70 
71 ~~~ {.cpp}
72  TGeoVolume *vol_box = gGeoManager->MakeBox("BOX_VOL", "mat1", halfX, halfY, halfZ);
73  TGeoVolume *vol_tub = gGeoManager->MakeTube("TUB_VOL", "mat2", rmin, rmax, halfZ);
74  ... (see MakeXXX() utilities in TGeoManager class)
75 ~~~
76 
77 ### Shape queries
78 
79 Note that global queries related to a geometry are handled by the manager class.
80 However, shape-related queries might be sometimes useful.
81 
82 #### `Bool_t TGeoShape::Contains(const Double_t *point[3])`
83 
84 this method returns true if POINT is actually inside the shape. The point
85 has to be defined in the local shape reference. For instance, for a box having
86 DX, DY and DZ half-lengths a point will be considered inside if :
87 
88 ~~~ {.cpp}
89  | -DX <= point[0] <= DX
90  | -DY <= point[1] <= DY
91  | -DZ <= point[2] <= DZ
92 ~~~
93 
94 #### `Double_t TGeoShape::DistFromInside(Double_t *point[3], Double_t *dir[3], Int_t iact, Double_t step, Double_t *safe)`
95 
96 computes the distance to exiting a shape from a given point INSIDE, along
97 a given direction. The direction is given by its director cosines with respect
98 to the local shape coordinate system. This method provides additional
99 information according the value of IACT input parameter :
100 
101  - IACT = 0 => compute only safe distance and fill it at the location
102  given by SAFE
103  - IACT = 1 => a proposed STEP is supplied. The safe distance is computed
104  first. If this is bigger than STEP than the proposed step
105  is approved and returned by the method since it does not
106  cross the shape boundaries. Otherwise, the distance to
107  exiting the shape is computed and returned.
108  - IACT = 2 => compute both safe distance and distance to exiting, ignoring
109  the proposed step.
110  - IACT > 2 => compute only the distance to exiting, ignoring anything else.
111 
112 #### `Double_t TGeoShape::DistFromOutside(Double_t *point[3], Double_t *dir[3], Int_t iact, Double_t step, Double_t *safe)`
113 
114 computes the distance to entering a shape from a given point OUTSIDE. Acts
115 in the same way as B).
116 
117 #### `Double_t Safety(const Double_t *point[3], Bool_t inside)`
118 
119 compute maximum shift of a point in any direction that does not change its
120 INSIDE/OUTSIDE state (does not cross shape boundaries). The state of the point
121 have to be properly supplied.
122 
123 #### `Double_t *Normal(Double_t *point[3], Double_t *dir[3], Bool_t inside)`
124 
125 returns director cosines of normal to the crossed shape surface from a
126 given point towards a direction. One has to specify if the point is inside
127 or outside shape. According to this, the normal will be outwards or inwards
128 shape respectively. Normal components are statically stored by shape class,
129 so it has to be copied after retrieval in a different array.
130 
131 ### Dividing shapes
132 
133 Shapes can generally be divided along a given axis. Supported axis are
134 X, Y, Z, Rxy, Phi, Rxyz. A given shape cannot be divided however on any axis.
135 The general rule is that that divisions are possible on whatever axis that
136 produces still known shapes as slices. The division of shapes should not be
137 performed by TGeoShape::Divide() calls, but rather by TGeoVolume::Divide().
138 The algorithm for dividing a specific shape is known by the shape object, but
139 is always invoked in a generic way from the volume level. Details on how to
140 do that can be found in TGeoVolume class. One can see how all division options
141 are interpreted and which is their result inside specific shape classes.
142 
143 \image html geom_t_shape.png
144 */
145 
146 #include "TObjArray.h"
147 #include "TEnv.h"
148 #include "TError.h"
149 
150 #include "TGeoMatrix.h"
151 #include "TGeoManager.h"
152 #include "TGeoVolume.h"
153 #include "TGeoShape.h"
154 #include "TVirtualGeoPainter.h"
155 #include "TBuffer3D.h"
156 #include "TBuffer3DTypes.h"
157 #include "TMath.h"
158 
160 
162 Double_t TGeoShape::fgEpsMch = 2.220446049250313e-16;
163 
164 ////////////////////////////////////////////////////////////////////////////////
165 /// Default constructor
166 
168 {
169  fShapeBits = 0;
170  fShapeId = 0;
171  if (!gGeoManager) {
172  gGeoManager = new TGeoManager("Geometry", "default geometry");
173  // gROOT->AddGeoManager(gGeoManager);
174  }
175 // fShapeId = gGeoManager->GetListOfShapes()->GetSize();
176 // gGeoManager->AddShape(this);
177 }
178 
179 ////////////////////////////////////////////////////////////////////////////////
180 /// Default constructor
181 
183  :TNamed(name, "")
184 {
185  fShapeBits = 0;
186  fShapeId = 0;
187  if (!gGeoManager) {
188  gGeoManager = new TGeoManager("Geometry", "default geometry");
189  // gROOT->AddGeoManager(gGeoManager);
190  }
192  gGeoManager->AddShape(this);
193 }
194 
195 ////////////////////////////////////////////////////////////////////////////////
196 /// Destructor
197 
199 {
201 }
202 
203 ////////////////////////////////////////////////////////////////////////////////
204 /// Test for shape navigation methods. Summary for test numbers:
205 ///
206 /// 1: DistFromInside/Outside. Sample points inside the shape. Generate
207 /// directions randomly in cos(theta). Compute DistFromInside and move the
208 /// point with bigger distance. Compute DistFromOutside back from new point.
209 /// Plot d-(d1+d2)
210 
211 void TGeoShape::CheckShape(Int_t testNo, Int_t nsamples, Option_t *option)
212 {
213  if (!gGeoManager) {
214  Error("CheckShape","No geometry manager");
215  return;
216  }
217  TGeoShape *shape = (TGeoShape*)this;
218  gGeoManager->CheckShape(shape, testNo, nsamples, option);
219 }
220 
221 ////////////////////////////////////////////////////////////////////////////////
222 /// Compute machine round-off double precision error as the smallest number that
223 /// if added to 1.0 is different than 1.0.
224 
226 {
227  Double_t temp1 = 1.0;
228  Double_t temp2 = 1.0 + temp1;
229  Double_t mchEps = 0.;
230  while (temp2>1.0) {
231  mchEps = temp1;
232  temp1 /= 2;
233  temp2 = 1.0 + temp1;
234  }
235  fgEpsMch = mchEps;
236  return fgEpsMch;
237 }
238 
239 ////////////////////////////////////////////////////////////////////////////////
240 ///static function returning the machine round-off error
241 
243 {
244  return fgEpsMch;
245 }
246 
247 ////////////////////////////////////////////////////////////////////////////////
248 /// Get the shape name.
249 
250 const char *TGeoShape::GetName() const
251 {
252  if (!fName[0]) {
253  return ((TObject *)this)->ClassName();
254  }
255  return TNamed::GetName();
256 }
257 
258 ////////////////////////////////////////////////////////////////////////////////
259 /// Returns distance to shape primitive mesh.
260 
262 {
264  if (!painter) return 9999;
265  return painter->ShapeDistancetoPrimitive(this, numpoints, px, py);
266 }
267 
268 ////////////////////////////////////////////////////////////////////////////////
269 /// True if point is closer than epsil to one of the phi planes defined by c1,s1 or c2,s2
270 
272 {
273  Double_t saf1 = TGeoShape::Big();
274  Double_t saf2 = TGeoShape::Big();
275  if (point[0]*c1+point[1]*s1 >= 0) saf1 = TMath::Abs(-point[0]*s1 + point[1]*c1);
276  if (point[0]*c2+point[1]*s2 >= 0) saf2 = TMath::Abs(point[0]*s2 - point[1]*c2);
277  Double_t saf = TMath::Min(saf1,saf2);
278  if (saf<epsil) return kTRUE;
279  return kFALSE;
280 }
281 
282 ////////////////////////////////////////////////////////////////////////////////
283 /// Static method to check if a point is in the phi range (phi1, phi2) [degrees]
284 
286 {
287  Double_t phi = TMath::ATan2(point[1], point[0]) * TMath::RadToDeg();
288  while (phi<phi1) phi+=360.;
289  Double_t ddp = phi-phi1;
290  if (ddp>phi2-phi1) return kFALSE;
291  return kTRUE;
292 }
293 
294 ////////////////////////////////////////////////////////////////////////////////
295 /// Compute distance from POINT to semiplane defined by PHI angle along DIR. Computes
296 /// also radius at crossing point. This might be negative in case the crossing is
297 /// on the other side of the semiplane.
298 
300 {
301  snext = rxy = TGeoShape::Big();
302  Double_t nx = -sphi;
303  Double_t ny = cphi;
304  Double_t rxy0 = point[0]*cphi+point[1]*sphi;
305  Double_t rdotn = point[0]*nx + point[1]*ny;
306  if (TMath::Abs(rdotn)<TGeoShape::Tolerance()) {
307  snext = 0.0;
308  rxy = rxy0;
309  return kTRUE;
310  }
311  if (rdotn<0) {
312  rdotn = -rdotn;
313  } else {
314  nx = -nx;
315  ny = -ny;
316  }
317  Double_t ddotn = dir[0]*nx + dir[1]*ny;
318  if (ddotn<=0) return kFALSE;
319  snext = rdotn/ddotn;
320  rxy = rxy0+snext*(dir[0]*cphi+dir[1]*sphi);
321  if (rxy<0) return kFALSE;
322  return kTRUE;
323 }
324 
325 ////////////////////////////////////////////////////////////////////////////////
326 /// Check if two numbers differ with less than a tolerance.
327 
329 {
330  if (TMath::Abs(a-b)<1.E-10) return kTRUE;
331  return kFALSE;
332 }
333 
334 ////////////////////////////////////////////////////////////////////////////////
335 /// Check if segments (A,B) and (C,D) are crossing,
336 /// where: A(x1,y1), B(x2,y2), C(x3,y3), D(x4,y4)
337 
339 {
341  Bool_t stand1 = kFALSE;
342  Double_t dx1 = x2-x1;
343  Bool_t stand2 = kFALSE;
344  Double_t dx2 = x4-x3;
345  Double_t xm = 0.;
346  Double_t ym = 0.;
347  Double_t a1 = 0.;
348  Double_t b1 = 0.;
349  Double_t a2 = 0.;
350  Double_t b2 = 0.;
351  if (TMath::Abs(dx1) < eps) stand1 = kTRUE;
352  if (TMath::Abs(dx2) < eps) stand2 = kTRUE;
353  if (!stand1) {
354  a1 = (x2*y1-x1*y2)/dx1;
355  b1 = (y2-y1)/dx1;
356  }
357  if (!stand2) {
358  a2 = (x4*y3-x3*y4)/dx2;
359  b2 = (y4-y3)/dx2;
360  }
361  if (stand1 && stand2) {
362  // Segments parallel and vertical
363  if (TMath::Abs(x1-x3)<eps) {
364  // Check if segments are overlapping
365  if ((y3-y1)*(y3-y2)<-eps || (y4-y1)*(y4-y2)<-eps ||
366  (y1-y3)*(y1-y4)<-eps || (y2-y3)*(y2-y4)<-eps) return kTRUE;
367  return kFALSE;
368  }
369  // Different x values
370  return kFALSE;
371  }
372 
373  if (stand1) {
374  // First segment vertical
375  xm = x1;
376  ym = a2+b2*xm;
377  } else {
378  if (stand2) {
379  // Second segment vertical
380  xm = x3;
381  ym = a1+b1*xm;
382  } else {
383  // Normal crossing
384  if (TMath::Abs(b1-b2)<eps) {
385  // Parallel segments, are they aligned
386  if (TMath::Abs(y3-(a1+b1*x3))>eps) return kFALSE;
387  // Aligned segments, are they overlapping
388  if ((x3-x1)*(x3-x2)<-eps || (x4-x1)*(x4-x2)<-eps ||
389  (x1-x3)*(x1-x4)<-eps || (x2-x3)*(x2-x4)<-eps) return kTRUE;
390  return kFALSE;
391  }
392  xm = (a1-a2)/(b2-b1);
393  ym = (a1*b2-a2*b1)/(b2-b1);
394  }
395  }
396  // Check if crossing point is both between A,B and C,D
397  Double_t check = (xm-x1)*(xm-x2)+(ym-y1)*(ym-y2);
398  if (check > -eps) return kFALSE;
399  check = (xm-x3)*(xm-x4)+(ym-y3)*(ym-y4);
400  if (check > -eps) return kFALSE;
401  return kTRUE;
402 }
403 
404 ////////////////////////////////////////////////////////////////////////////////
405 /// compute distance from point (inside phi) to both phi planes. Return minimum.
406 
408  Double_t s2, Double_t c2, Double_t sm, Double_t cm, Bool_t in)
409 {
410  Double_t sfi1=TGeoShape::Big();
411  Double_t sfi2=TGeoShape::Big();
412  Double_t s=0;
413  Double_t un = dir[0]*s1-dir[1]*c1;
414  if (!in) un=-un;
415  if (un>0) {
416  s=-point[0]*s1+point[1]*c1;
417  if (!in) s=-s;
418  if (s>=0) {
419  s /= un;
420  if (((point[0]+s*dir[0])*sm-(point[1]+s*dir[1])*cm)>=0) sfi1=s;
421  }
422  }
423  un = -dir[0]*s2+dir[1]*c2;
424  if (!in) un=-un;
425  if (un>0) {
426  s=point[0]*s2-point[1]*c2;
427  if (!in) s=-s;
428  if (s>=0) {
429  s /= un;
430  if ((-(point[0]+s*dir[0])*sm+(point[1]+s*dir[1])*cm)>=0) sfi2=s;
431  }
432  }
433  return TMath::Min(sfi1, sfi2);
434 }
435 
436 ////////////////////////////////////////////////////////////////////////////////
437 /// Static method to compute normal to phi planes.
438 
440 {
441  Double_t saf1 = TGeoShape::Big();
442  Double_t saf2 = TGeoShape::Big();
443  if (point[0]*c1+point[1]*s1 >= 0) saf1 = TMath::Abs(-point[0]*s1 + point[1]*c1);
444  if (point[0]*c2+point[1]*s2 >= 0) saf2 = TMath::Abs(point[0]*s2 - point[1]*c2);
445  Double_t c,s;
446  if (saf1<saf2) {
447  c=c1;
448  s=s1;
449  } else {
450  c=c2;
451  s=s2;
452  }
453  norm[2] = 0;
454  norm[0] = -s;
455  norm[1] = c;
456  if (dir[0]*norm[0]+dir[1]*norm[1] < 0) {
457  norm[0] = s;
458  norm[1] = -c;
459  }
460 }
461 
462 ////////////////////////////////////////////////////////////////////////////////
463 /// Static method to compute safety w.r.t a phi corner defined by cosines/sines
464 /// of the angles phi1, phi2.
465 
467 {
468  Bool_t inphi = TGeoShape::IsInPhiRange(point, phi1, phi2);
469  if (inphi && !in) return -TGeoShape::Big();
470  phi1 *= TMath::DegToRad();
471  phi2 *= TMath::DegToRad();
472  Double_t c1 = TMath::Cos(phi1);
473  Double_t s1 = TMath::Sin(phi1);
474  Double_t c2 = TMath::Cos(phi2);
475  Double_t s2 = TMath::Sin(phi2);
476  Double_t rsq = point[0]*point[0]+point[1]*point[1];
477  Double_t rproj = point[0]*c1+point[1]*s1;
478  Double_t safsq = rsq-rproj*rproj;
479  if (safsq<0) return 0.;
480  Double_t saf1 = (rproj<0)?TGeoShape::Big():TMath::Sqrt(safsq);
481  rproj = point[0]*c2+point[1]*s2;
482  safsq = rsq-rproj*rproj;
483  if (safsq<0) return 0.;
484  Double_t saf2 = (rproj<0)?TGeoShape::Big():TMath::Sqrt(safsq);
485  Double_t safe = TMath::Min(saf1, saf2); // >0
486  if (safe>1E10) {
487  if (in) return TGeoShape::Big();
488  return -TGeoShape::Big();
489  }
490  return safe;
491 }
492 
493 ////////////////////////////////////////////////////////////////////////////////
494 /// Compute distance from point of coordinates (r,z) to segment (r1,z1):(r2,z2)
495 
497 {
498  Double_t crossp = (z2-z1)*(r-r1)-(z-z1)*(r2-r1);
499  crossp *= (outer) ? 1. : -1.;
500  // Positive crossp means point on the requested side of the (1,2) segment
501  if (crossp < -TGeoShape::Tolerance()) {
502 // if (((z-z1)*(z2-z)) > -1.E-10) return 0;
503  if (outer) return TGeoShape::Big();
504  else return 0.;
505  }
506  // Compute (1,P) dot (1,2)
507  Double_t c1 = (z-z1)*(z2-z1)+(r-r1)*(r2-r1);
508  // Negative c1 means point (1) is closest
509  if (c1<1.E-10) return TMath::Sqrt((r-r1)*(r-r1)+(z-z1)*(z-z1));
510  // Compute (2,P) dot (1,2)
511  Double_t c2 = (z-z2)*(z2-z1)+(r-r2)*(r2-r1);
512  // Positive c2 means point (2) is closest
513  if (c2>-1.E-10) return TMath::Sqrt((r-r2)*(r-r2)+(z-z2)*(z-z2));
514  // The closest point is between (1) and (2)
515  c2 = (z2-z1)*(z2-z1)+(r2-r1)*(r2-r1);
516  // projected length factor with respect to (1,2) length
517  Double_t alpha = c1/c2;
518  Double_t rp = r1 + alpha*(r2-r1);
519  Double_t zp = z1 + alpha*(z2-z1);
520  return TMath::Sqrt((r-rp)*(r-rp)+(z-zp)*(z-zp));
521 }
522 
523 ////////////////////////////////////////////////////////////////////////////////
524 /// Equivalent of TObject::SetBit.
525 
527 {
528  if (set) {
529  SetShapeBit(f);
530  } else {
531  ResetShapeBit(f);
532  }
533 }
534 
535 ////////////////////////////////////////////////////////////////////////////////
536 /// Returns current transformation matrix that applies to shape.
537 
539 {
540  return fgTransform;
541 }
542 
543 ////////////////////////////////////////////////////////////////////////////////
544 /// Set current transformation matrix that applies to shape.
545 
547 {
548  fgTransform = matrix;
549 }
550 
551 ////////////////////////////////////////////////////////////////////////////////
552 /// Tranform a set of points (LocalToMaster)
553 
555 {
556  UInt_t i,j;
557  Double_t dmaster[3];
558  if (fgTransform) {
559  for (j = 0; j < NbPnts; j++) {
560  i = 3*j;
561  fgTransform->LocalToMaster(&points[i], dmaster);
562  points[i] = dmaster[0];
563  points[i+1] = dmaster[1];
564  points[i+2] = dmaster[2];
565  }
566  return;
567  }
568  if (!gGeoManager) return;
569  Bool_t bomb = (gGeoManager->GetBombMode()==0)?kFALSE:kTRUE;
570 
571  for (j = 0; j < NbPnts; j++) {
572  i = 3*j;
574  TGeoHMatrix *glmat = gGeoManager->GetGLMatrix();
575  if (bomb) glmat->LocalToMasterBomb(&points[i], dmaster);
576  else glmat->LocalToMaster(&points[i], dmaster);
577  } else {
578  if (bomb) gGeoManager->LocalToMasterBomb(&points[i], dmaster);
579  else gGeoManager->LocalToMaster(&points[i],dmaster);
580  }
581  points[i] = dmaster[0];
582  points[i+1] = dmaster[1];
583  points[i+2] = dmaster[2];
584  }
585 }
586 
587 ////////////////////////////////////////////////////////////////////////////////
588 /// Fill the supplied buffer, with sections in desired frame
589 /// See TBuffer3D.h for explanation of sections, frame etc.
590 
591 void TGeoShape::FillBuffer3D(TBuffer3D & buffer, Int_t reqSections, Bool_t localFrame) const
592 {
593  // Catch this common potential error here
594  // We have to set kRawSize (unless already done) to allocate buffer space
595  // before kRaw can be filled
596  if (reqSections & TBuffer3D::kRaw) {
597  if (!(reqSections & TBuffer3D::kRawSizes) && !buffer.SectionsValid(TBuffer3D::kRawSizes)) {
598  R__ASSERT(kFALSE);
599  }
600  }
601 
602  if (reqSections & TBuffer3D::kCore) {
603  // If writing core section all others will be invalid
604  buffer.ClearSectionsValid();
605 
606  // Check/grab some objects we need
607  if (!gGeoManager) {
608  R__ASSERT(kFALSE);
609  return;
610  }
611  const TGeoVolume * paintVolume = gGeoManager->GetPaintVolume();
612  if (!paintVolume) paintVolume = gGeoManager->GetTopVolume();
613  if (!paintVolume) {
614  buffer.fID = const_cast<TGeoShape *>(this);
615  buffer.fColor = 0;
616  buffer.fTransparency = 0;
617 // R__ASSERT(kFALSE);
618 // return;
619  } else {
620  buffer.fID = const_cast<TGeoVolume *>(paintVolume);
621  buffer.fColor = paintVolume->GetLineColor();
622 
623  buffer.fTransparency = paintVolume->GetTransparency();
624  Double_t visdensity = gGeoManager->GetVisDensity();
625  if (visdensity>0 && paintVolume->GetMedium()) {
626  if (paintVolume->GetMaterial()->GetDensity() < visdensity) {
627  buffer.fTransparency = 90;
628  }
629  }
630  }
631 
632  buffer.fLocalFrame = localFrame;
633  Bool_t r1,r2=kFALSE;
635  if (paintVolume && paintVolume->GetShape()) {
636  if (paintVolume->GetShape()->IsReflected()) {
637  // Temporary trick to deal with reflected shapes.
638  // Still lighting gets wrong...
639  if (buffer.Type() < TBuffer3DTypes::kTube) r2 = kTRUE;
640  }
641  }
642  buffer.fReflection = ((r1&(!r2))|(r2&!(r1)));
643 
644  // Set up local -> master translation matrix
645  if (localFrame) {
646  TGeoMatrix * localMasterMat = 0;
647  if (TGeoShape::GetTransform()) {
648  localMasterMat = TGeoShape::GetTransform();
649  } else {
650  localMasterMat = gGeoManager->GetCurrentMatrix();
651 
652  // For overlap drawing the correct matrix needs to obtained in
653  // from GetGLMatrix() - this should not be applied in the case
654  // of composite shapes
656  localMasterMat = gGeoManager->GetGLMatrix();
657  }
658  }
659  if (!localMasterMat) {
660  R__ASSERT(kFALSE);
661  return;
662  }
663  localMasterMat->GetHomogenousMatrix(buffer.fLocalMaster);
664  } else {
665  buffer.SetLocalMasterIdentity();
666  }
667 
668  buffer.SetSectionsValid(TBuffer3D::kCore);
669  }
670 }
671 
672 ////////////////////////////////////////////////////////////////////////////////
673 /// Get the basic color (0-7).
674 
676 {
677  Int_t basicColor = 0; // TODO: Check on sensible fallback
678  if (gGeoManager) {
679  const TGeoVolume * volume = gGeoManager->GetPaintVolume();
680  if (volume) {
681  basicColor = ((volume->GetLineColor() %8) -1) * 4;
682  if (basicColor < 0) basicColor = 0;
683  }
684  }
685  return basicColor;
686 }
687 
688 ////////////////////////////////////////////////////////////////////////////////
689 /// Stub implementation to avoid forcing implementation at this stage
690 
691 const TBuffer3D &TGeoShape::GetBuffer3D(Int_t /*reqSections*/, Bool_t /*localFrame*/) const
692 {
693  static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
694  Warning("GetBuffer3D", "this must be implemented for shapes in a TGeoPainter hierarchy. This will be come a pure virtual fn eventually.");
695  return buffer;
696 }
697 
698 ////////////////////////////////////////////////////////////////////////////////
699 /// Provide a pointer name containing uid.
700 
701 const char *TGeoShape::GetPointerName() const
702 {
703  static TString name;
704  Int_t uid = GetUniqueID();
705  if (uid) name = TString::Format("p%s_%d", GetName(),uid);
706  else name = TString::Format("p%s", GetName());
707  return name.Data();
708 }
709 
710 ////////////////////////////////////////////////////////////////////////////////
711 /// Execute mouse actions on this shape.
712 
714 {
715  if (!gGeoManager) return;
717  painter->ExecuteShapeEvent(this, event, px, py);
718 }
719 
720 ////////////////////////////////////////////////////////////////////////////////
721 /// Draw this shape.
722 
724 {
726  if (option && option[0]) {
727  painter->DrawShape(this, option);
728  } else {
729  painter->DrawShape(this, gEnv->GetValue("Viewer3D.DefaultDrawOption",""));
730  }
731 }
732 
733 ////////////////////////////////////////////////////////////////////////////////
734 /// Paint this shape.
735 
737 {
739  if (option && option[0]) {
740  painter->PaintShape(this, option);
741  } else {
742  painter->PaintShape(this, gEnv->GetValue("Viewer3D.DefaultDrawOption",""));
743  }
744 }
const int nx
Definition: kalman.C:16
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
virtual UInt_t GetUniqueID() const
Return the unique object id.
Definition: TObject.cxx:434
void CheckShape(Int_t testNo, Int_t nsamples=10000, Option_t *option="")
Test for shape navigation methods.
Definition: TGeoShape.cxx:211
TGeoVolume * GetPaintVolume() const
Definition: TGeoManager.h:193
virtual ~TGeoShape()
Destructor.
Definition: TGeoShape.cxx:198
void CheckShape(TGeoShape *shape, Int_t testNo, Int_t nsamples, Option_t *option)
Test for shape navigation methods.
The manager class for any TGeo geometry.
Definition: TGeoManager.h:38
#define snext(osub1, osub2)
Definition: triangle.c:1167
static Double_t ComputeEpsMch()
Compute machine round-off double precision error as the smallest number that if added to 1...
Definition: TGeoShape.cxx:225
virtual Bool_t IsReflected() const
Definition: TGeoShape.h:139
virtual void ExecuteEvent(Int_t event, Int_t px, Int_t py)
Execute mouse actions on this shape.
Definition: TGeoShape.cxx:713
static Double_t EpsMch()
static function returning the machine round-off error
Definition: TGeoShape.cxx:242
virtual Double_t GetDensity() const
Definition: TGeoMaterial.h:93
return c
const char Option_t
Definition: RtypesCore.h:62
Geometrical transformation package.
Definition: TGeoMatrix.h:40
return c1
Definition: legend1.C:41
TVirtualGeoPainter * GetPainter() const
Definition: TGeoManager.h:187
Double_t DegToRad()
Definition: TMath.h:50
TGeoVolume, TGeoVolumeMulti, TGeoVolumeAssembly are the volume classes.
Definition: TGeoVolume.h:61
Double_t fLocalMaster[16]
Definition: TBuffer3D.h:94
void SetLocalMasterIdentity()
Set kRaw tessellation section of buffer with supplied sizes.
Definition: TBuffer3D.cxx:296
#define R__ASSERT(e)
Definition: TError.h:98
virtual TObject * Remove(TObject *obj)
Remove object from array.
Definition: TObjArray.cxx:653
Double_t GetVisDensity() const
Definition: TGeoManager.h:195
void ClearSectionsValid()
Clear any sections marked valid.
Definition: TBuffer3D.cxx:286
Double_t RadToDeg()
Definition: TMath.h:49
Basic string class.
Definition: TString.h:137
Matrix class used for computing global transformations Should NOT be used for node definition...
Definition: TGeoMatrix.h:410
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:170
Bool_t IsMatrixReflection() const
Definition: TGeoManager.h:387
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
TArc * a
Definition: textangle.C:12
const Bool_t kFALSE
Definition: Rtypes.h:92
static Double_t SafetyPhi(const Double_t *point, Bool_t in, Double_t phi1, Double_t phi2)
Static method to compute safety w.r.t a phi corner defined by cosines/sines of the angles phi1...
Definition: TGeoShape.cxx:466
Short_t Abs(Short_t d)
Definition: TMathBase.h:110
static Bool_t IsSameWithinTolerance(Double_t a, Double_t b)
Check if two numbers differ with less than a tolerance.
Definition: TGeoShape.cxx:328
static Double_t Tolerance()
Definition: TGeoShape.h:93
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:591
static const double x2[5]
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString...
Definition: TString.cxx:2335
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:33
const int ny
Definition: kalman.C:17
static Double_t fgEpsMch
Definition: TGeoShape.h:31
static Bool_t IsCrossingSemiplane(const Double_t *point, const Double_t *dir, Double_t cphi, Double_t sphi, Double_t &snext, Double_t &rxy)
Compute distance from POINT to semiplane defined by PHI angle along DIR.
Definition: TGeoShape.cxx:299
static const double x4[22]
TGeoMedium * GetMedium() const
Definition: TGeoVolume.h:189
static TGeoMatrix * fgTransform
Definition: TGeoShape.h:30
void LocalToMaster(const Double_t *local, Double_t *master) const
Definition: TGeoManager.h:508
TGeoMaterial * GetMaterial() const
Definition: TGeoVolume.h:188
void LocalToMasterBomb(const Double_t *local, Double_t *master) const
Definition: TGeoManager.h:510
Double_t ATan2(Double_t, Double_t)
Definition: TMath.h:454
static void SetTransform(TGeoMatrix *matrix)
Set current transformation matrix that applies to shape.
Definition: TGeoShape.cxx:546
virtual Int_t ShapeDistancetoPrimitive(const TGeoShape *shape, Int_t numpoints, Int_t px, Int_t py) const =0
void SetSectionsValid(UInt_t mask)
Definition: TBuffer3D.h:67
const char * GetPointerName() const
Provide a pointer name containing uid.
Definition: TGeoShape.cxx:701
static Double_t SafetySeg(Double_t r, Double_t z, Double_t r1, Double_t z1, Double_t r2, Double_t z2, Bool_t outer)
Compute distance from point of coordinates (r,z) to segment (r1,z1):(r2,z2)
Definition: TGeoShape.cxx:496
Bool_t fLocalFrame
Definition: TBuffer3D.h:92
point * points
Definition: X3DBuffer.c:20
static Bool_t IsSegCrossing(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Double_t x3, Double_t y3, Double_t x4, Double_t y4)
Check if segments (A,B) and (C,D) are crossing, where: A(x1,y1), B(x2,y2), C(x3,y3), D(x4,y4)
Definition: TGeoShape.cxx:338
virtual const char * GetName() const
Get the shape name.
Definition: TGeoShape.cxx:250
TGeoShape()
Default constructor.
Definition: TGeoShape.cxx:167
Char_t GetTransparency() const
Definition: TGeoVolume.h:203
Int_t ShapeDistancetoPrimitive(Int_t numpoints, Int_t px, Int_t py) const
Returns distance to shape primitive mesh.
Definition: TGeoShape.cxx:261
Bool_t SectionsValid(UInt_t mask) const
Definition: TBuffer3D.h:69
Base abstract class for all shapes.
Definition: TGeoShape.h:27
TRandom2 r(17)
virtual Int_t GetValue(const char *name, Int_t dflt)
Returns the integer value for a resource.
Definition: TEnv.cxx:496
unsigned int r1[N_CITIES]
Definition: simanTSP.cxx:321
virtual Bool_t IsComposite() const
Definition: TGeoShape.h:132
virtual const TBuffer3D & GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
Stub implementation to avoid forcing implementation at this stage.
Definition: TGeoShape.cxx:691
void ResetShapeBit(UInt_t f)
Definition: TGeoShape.h:164
TObjArray * GetListOfShapes() const
Definition: TGeoManager.h:465
static Double_t DistToPhiMin(const Double_t *point, const Double_t *dir, Double_t s1, Double_t c1, Double_t s2, Double_t c2, Double_t sm, Double_t cm, Bool_t in=kTRUE)
compute distance from point (inside phi) to both phi planes. Return minimum.
Definition: TGeoShape.cxx:407
unsigned int UInt_t
Definition: RtypesCore.h:42
void TransformPoints(Double_t *points, UInt_t NbPoints) const
Tranform a set of points (LocalToMaster)
Definition: TGeoShape.cxx:554
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:925
TGeoHMatrix * GetCurrentMatrix() const
Definition: TGeoManager.h:481
Double_t E()
Definition: TMath.h:54
static Bool_t IsInPhiRange(const Double_t *point, Double_t phi1, Double_t phi2)
Static method to check if a point is in the phi range (phi1, phi2) [degrees].
Definition: TGeoShape.cxx:285
Generic 3D primitive description class.
Definition: TBuffer3D.h:19
TObject * fID
Definition: TBuffer3D.h:89
void GetHomogenousMatrix(Double_t *hmat) const
The homogenous matrix associated with the transformation is used for piling up&#39;s and visualization...
Definition: TGeoMatrix.cxx:364
TString fName
Definition: TNamed.h:36
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:389
Double_t Cos(Double_t)
Definition: TMath.h:424
virtual Color_t GetLineColor() const
Return the line color.
Definition: TAttLine.h:39
Bool_t fReflection
Definition: TBuffer3D.h:93
static Bool_t IsCloseToPhi(Double_t epsil, const Double_t *point, Double_t c1, Double_t s1, Double_t c2, Double_t s2)
True if point is closer than epsil to one of the phi planes defined by c1,s1 or c2,s2.
Definition: TGeoShape.cxx:271
return c2
Definition: legend2.C:14
static const double x1[5]
#define ClassImp(name)
Definition: Rtypes.h:279
double f(double x)
R__EXTERN TGeoManager * gGeoManager
Definition: TGeoManager.h:554
Int_t fShapeId
Definition: TGeoShape.h:74
double Double_t
Definition: RtypesCore.h:55
TGeoHMatrix * GetGLMatrix() const
Definition: TGeoManager.h:482
Int_t Type() const
Definition: TBuffer3D.h:87
R__EXTERN TEnv * gEnv
Definition: TEnv.h:174
Int_t fColor
Definition: TBuffer3D.h:90
Bool_t IsMatrixTransform() const
Definition: TGeoManager.h:386
static Double_t Big()
Definition: TGeoShape.h:90
Int_t AddShape(const TGeoShape *shape)
Add a shape to the list. Returns index of the shape in list.
static TGeoMatrix * GetTransform()
Returns current transformation matrix that applies to shape.
Definition: TGeoShape.cxx:538
void SetShapeBit(UInt_t f, Bool_t set)
Equivalent of TObject::SetBit.
Definition: TGeoShape.cxx:526
Mother of all ROOT objects.
Definition: TObject.h:37
you should not use this method at all Int_t Int_t z
Definition: TRolke.cxx:630
Int_t GetBombMode() const
Definition: TGeoManager.h:188
Abstract class for geometry painters.
virtual void DrawShape(TGeoShape *shape, Option_t *option="")=0
virtual void ExecuteShapeEvent(TGeoShape *shape, Int_t event, Int_t px, Int_t py)=0
Bool_t IsCleaning() const
Definition: TGeoManager.h:454
Double_t Sin(Double_t)
Definition: TMath.h:421
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t b
Definition: TRolke.cxx:630
#define NULL
Definition: Rtypes.h:82
static void NormalPhi(const Double_t *point, const Double_t *dir, Double_t *norm, Double_t c1, Double_t s1, Double_t c2, Double_t s2)
Static method to compute normal to phi planes.
Definition: TGeoShape.cxx:439
Short_t fTransparency
Definition: TBuffer3D.h:91
virtual void LocalToMasterBomb(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:430
virtual void Paint(Option_t *option="")
Paint this shape.
Definition: TGeoShape.cxx:736
Double_t Sqrt(Double_t x)
Definition: TMath.h:464
TGeoVolume * GetTopVolume() const
Definition: TGeoManager.h:497
TGeoShape * GetShape() const
Definition: TGeoVolume.h:204
virtual Int_t GetSize() const
Definition: TCollection.h:95
virtual void Draw(Option_t *option="")
Draw this shape.
Definition: TGeoShape.cxx:723
const Bool_t kTRUE
Definition: Rtypes.h:91
virtual void PaintShape(TGeoShape *shape, Option_t *option="")=0
double norm(double *x, double *p)
Definition: unuranDistr.cxx:40
Int_t GetBasicColor() const
Get the basic color (0-7).
Definition: TGeoShape.cxx:675
TVirtualGeoPainter * GetGeomPainter()
Make a default painter if none present. Returns pointer to it.
unsigned int r2[N_CITIES]
Definition: simanTSP.cxx:322
char name[80]
Definition: TGX11.cxx:109
UInt_t fShapeBits
Definition: TGeoShape.h:75
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:911
static const double x3[11]
const char * Data() const
Definition: TString.h:349