TGeoShape.cxx

Go to the documentation of this file.

// @(#)root/geom:$Id$
// Author: Andrei Gheata   31/01/02
 
/*************************************************************************
 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers.               *
 * All rights reserved.                                                  *
 *                                                                       *
 * For the licensing terms see $ROOTSYS/LICENSE.                         *
 * For the list of contributors see $ROOTSYS/README/CREDITS.             *
 *************************************************************************/
 
/** \class TGeoShape
\ingroup Shapes_classes
Base abstract class for all shapes.
 
  Shapes are geometrical objects that provide the basic modelling
functionality. They provide the definition of the LOCAL frame of coordinates,
with respect to which they are defined. Any implementation of a shape deriving
from the base TGeoShape class has to provide methods for :
 
  - finding out if a point defined in their local frame is or not contained
    inside;
  - computing the distance from a local point to getting outside/entering the
    shape, given a known direction;
  - computing the maximum distance in any direction from a local point that
    does NOT result in a boundary crossing of the shape (safe distance);
  - computing the cosines of the normal vector to the crossed shape surface,
    given a starting local point and an ongoing direction.
    All the features above are globally managed by the modeller in order to
    provide navigation functionality. In addition to those, shapes have also to
    implement additional specific abstract methods :
  - computation of the minimal box bounding the shape, given that this box have
    to be aligned with the local coordinates;
  - algorithms for dividing the shape along a given axis and producing resulting
    divisions volumes.
 
  The modeler currently provides a set of 16 basic shapes, which we will call
primitives. It also provides a special class allowing the creation of shapes
made as a result of boolean operations between primitives. These are called
composite shapes and the composition operation can be recursive (composition
of composites). This allows the creation of a quite large number of different
shape topologies and combinations.
 
  Named shapes register themselves to the manager class at creation time. The
manager is responsible for their final deletion. Shapes can be created using their
default constructor if their retrieval by name is not needed, but in this case
they are owned by the user. A shape may be referenced by several volumes,
therefore its deletion is not possible once volumes were defined based on it.
 
### Creating shapes
 
  Shape objects embed only the minimum set of parameters that are fully
describing a valid physical shape. For instance, a tube is represented by
its half length, the minimum radius and the maximum radius. Shapes are used
together with media in order to create volumes, which in their turn
are the main components of the geometrical tree. A specific shape can be created
stand-alone :
 
~~~ {.cpp}
  TGeoBBox *box = new TGeoBBox("s_box", halfX, halfY, halfZ); // named
  TGeoTube *tub = new TGeoTube(rmin, rmax, halfZ);            // no name
  ...  (see each specific shape constructors)
~~~
 
  Sometimes it is much easier to create a volume having a given shape in one
step, since shapes are not directly linked in the geometrical tree but volumes
are :
 
~~~ {.cpp}
  TGeoVolume *vol_box = gGeoManager->MakeBox("BOX_VOL", "mat1", halfX, halfY, halfZ);
  TGeoVolume *vol_tub = gGeoManager->MakeTube("TUB_VOL", "mat2", rmin, rmax, halfZ);
  ...  (see MakeXXX() utilities in TGeoManager class)
~~~
 
### Shape queries
 
Note that global queries related to a geometry are handled by the manager class.
However, shape-related queries might be sometimes useful.
 
#### `Bool_t TGeoShape::Contains(const Double_t *point[3])`
 
this method returns true if POINT is actually inside the shape. The point
has to be defined in the local shape reference. For instance, for a box having
DX, DY and DZ half-lengths a point will be considered inside if :
 
~~~ {.cpp}
  | -DX <= point[0] <= DX
  | -DY <= point[1] <= DY
  | -DZ <= point[2] <= DZ
~~~
 
#### `Double_t TGeoShape::DistFromInside(Double_t *point[3], Double_t *dir[3], Int_t iact, Double_t step, Double_t
*safe)`
 
computes the distance to exiting a shape from a given point INSIDE, along
a given direction. The direction is given by its director cosines with respect
to the local shape coordinate system. This method provides additional
information according the value of IACT input parameter :
 
  - IACT = 0     => compute only safe distance and fill it at the location
                    given by SAFE
  - IACT = 1     => a proposed STEP is supplied. The safe distance is computed
                    first. If this is bigger than STEP than the proposed step
                    is approved and returned by the method since it does not
                    cross the shape boundaries. Otherwise, the distance to
                    exiting the shape is computed and returned.
  - IACT = 2     => compute both safe distance and distance to exiting, ignoring
                    the proposed step.
  - IACT > 2     => compute only the distance to exiting, ignoring anything else.
 
#### `Double_t TGeoShape::DistFromOutside(Double_t *point[3], Double_t *dir[3], Int_t iact, Double_t step, Double_t
*safe)`
 
computes the distance to entering a shape from a given point OUTSIDE. Acts
in the same way as B).
 
#### `Double_t Safety(const Double_t *point[3], Bool_t inside)`
 
compute maximum shift of a point in any direction that does not change its
INSIDE/OUTSIDE state (does not cross shape boundaries). The state of the point
have to be properly supplied.
 
#### `Double_t *Normal(Double_t *point[3], Double_t *dir[3], Bool_t inside)`
 
returns director cosines of normal to the crossed shape surface from a
given point towards a direction. One has to specify if the point is inside
or outside shape. According to this, the normal will be outwards or inwards
shape respectively. Normal components are statically stored by shape class,
so it has to be copied after retrieval in a different array.
 
### Dividing shapes
 
Shapes can generally be divided along a given axis. Supported axis are
X, Y, Z, Rxy, Phi, Rxyz. A given shape cannot be divided however on any axis.
The general rule is that that divisions are possible on whatever axis that
produces still known shapes as slices. The division of shapes should not be
performed by TGeoShape::Divide() calls, but rather by TGeoVolume::Divide().
The algorithm for dividing a specific shape is known by the shape object, but
is always invoked in a generic way from the volume level. Details on how to
do that can be found in TGeoVolume class. One can see how all division options
are interpreted and which is their result inside specific shape classes.
 
\image html geom_t_shape.png width=600px
*/
 
#include "TObjArray.h"
#include "TEnv.h"
#include "TError.h"
 
#include "TGeoMatrix.h"
#include "TGeoManager.h"
#include "TGeoVolume.h"
#include "TGeoShape.h"
#include "TVirtualGeoPainter.h"
#include "TBuffer3D.h"
#include "TBuffer3DTypes.h"
#include "TMath.h"
 
ClassImp(TGeoShape);
 
TGeoMatrix *TGeoShape::fgTransform = nullptr;
Double_t TGeoShape::fgEpsMch = 2.220446049250313e-16;
 
////////////////////////////////////////////////////////////////////////////////
/// Default constructor
 
TGeoShape::TGeoShape()
{
   fShapeBits = 0;
   fShapeId = 0;
   if (!gGeoManager) {
      gGeoManager = new TGeoManager("Geometry", "default geometry");
      // gROOT->AddGeoManager(gGeoManager);
   }
   //   fShapeId = gGeoManager->GetListOfShapes()->GetSize();
   //   gGeoManager->AddShape(this);
}
 
////////////////////////////////////////////////////////////////////////////////
/// Default constructor
 
TGeoShape::TGeoShape(const char *name) : TNamed(name, "")
{
   fShapeBits = 0;
   fShapeId = 0;
   if (!gGeoManager) {
      gGeoManager = new TGeoManager("Geometry", "default geometry");
      // gROOT->AddGeoManager(gGeoManager);
   }
   fShapeId = gGeoManager->GetListOfShapes()->GetSize();
   gGeoManager->AddShape(this);
}
 
////////////////////////////////////////////////////////////////////////////////
/// Destructor
 
TGeoShape::~TGeoShape()
{
   if (gGeoManager && !gGeoManager->IsCleaning())
      gGeoManager->GetListOfShapes()->Remove(this);
}
 
////////////////////////////////////////////////////////////////////////////////
/// Test for shape navigation methods. Summary for test numbers:
///
///  - 1: DistFromInside/Outside. Sample points inside the shape. Generate
///    directions randomly in cos(theta). Compute DistFromInside and move the
///    point with bigger distance. Compute DistFromOutside back from new point.
///    Plot d-(d1+d2)
 
void TGeoShape::CheckShape(Int_t testNo, Int_t nsamples, Option_t *option)
{
   if (!gGeoManager) {
      Error("CheckShape", "No geometry manager");
      return;
   }
   TGeoShape *shape = (TGeoShape *)this;
   gGeoManager->CheckShape(shape, testNo, nsamples, option);
}
 
////////////////////////////////////////////////////////////////////////////////
/// Compute machine round-off double precision error as the smallest number that
/// if added to 1.0 is different than 1.0.
 
Double_t TGeoShape::ComputeEpsMch()
{
   Double_t temp1 = 1.0;
   Double_t temp2 = 1.0 + temp1;
   Double_t mchEps = 0.;
   while (temp2 > 1.0) {
      mchEps = temp1;
      temp1 /= 2;
      temp2 = 1.0 + temp1;
   }
   fgEpsMch = mchEps;
   return fgEpsMch;
}
 
////////////////////////////////////////////////////////////////////////////////
/// static function returning the machine round-off error
 
Double_t TGeoShape::EpsMch()
{
   return fgEpsMch;
}
 
////////////////////////////////////////////////////////////////////////////////
/// Get the shape name.
 
const char *TGeoShape::GetName() const
{
   if (!fName[0]) {
      return ((TObject *)this)->ClassName();
   }
   return TNamed::GetName();
}
 
////////////////////////////////////////////////////////////////////////////////
/// Returns distance to shape primitive mesh.
 
Int_t TGeoShape::ShapeDistancetoPrimitive(Int_t numpoints, Int_t px, Int_t py) const
{
   TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter();
   if (!painter)
      return 9999;
   return painter->ShapeDistancetoPrimitive(this, numpoints, px, py);
}
 
////////////////////////////////////////////////////////////////////////////////
/// Implementation of the inside function using just Contains and GetNormal
 
TGeoShape::EInside TGeoShape::Inside(const Double_t *point) const
{
   return tgeo_impl::Inside(point, this);
}
 
////////////////////////////////////////////////////////////////////////////////
/// True if point is closer than epsil to one of the phi planes defined by c1,s1 or c2,s2
 
Bool_t
TGeoShape::IsCloseToPhi(Double_t epsil, const Double_t *point, Double_t c1, Double_t s1, Double_t c2, Double_t s2)
{
   Double_t saf1 = TGeoShape::Big();
   Double_t saf2 = TGeoShape::Big();
   if (point[0] * c1 + point[1] * s1 >= 0)
      saf1 = TMath::Abs(-point[0] * s1 + point[1] * c1);
   if (point[0] * c2 + point[1] * s2 >= 0)
      saf2 = TMath::Abs(point[0] * s2 - point[1] * c2);
   Double_t saf = TMath::Min(saf1, saf2);
   if (saf < epsil)
      return kTRUE;
   return kFALSE;
}
 
////////////////////////////////////////////////////////////////////////////////
/// Static method to check if a point is in the phi range (phi1, phi2) [degrees]
 
Bool_t TGeoShape::IsInPhiRange(const Double_t *point, Double_t phi1, Double_t phi2)
{
   Double_t phi = TMath::ATan2(point[1], point[0]) * TMath::RadToDeg();
   while (phi < phi1)
      phi += 360.;
   Double_t ddp = phi - phi1;
   if (ddp > phi2 - phi1)
      return kFALSE;
   return kTRUE;
}
 
////////////////////////////////////////////////////////////////////////////////
/// Compute distance from POINT to semiplane defined by PHI angle along DIR. Computes
/// also radius at crossing point. This might be negative in case the crossing is
/// on the other side of the semiplane.
 
Bool_t TGeoShape::IsCrossingSemiplane(const Double_t *point, const Double_t *dir, Double_t cphi, Double_t sphi,
                                      Double_t &snext, Double_t &rxy)
{
   snext = rxy = TGeoShape::Big();
   Double_t nx = -sphi;
   Double_t ny = cphi;
   Double_t rxy0 = point[0] * cphi + point[1] * sphi;
   Double_t rdotn = point[0] * nx + point[1] * ny;
   if (TMath::Abs(rdotn) < TGeoShape::Tolerance()) {
      snext = 0.0;
      rxy = rxy0;
      return kTRUE;
   }
   if (rdotn < 0) {
      rdotn = -rdotn;
   } else {
      nx = -nx;
      ny = -ny;
   }
   Double_t ddotn = dir[0] * nx + dir[1] * ny;
   if (ddotn <= 0)
      return kFALSE;
   snext = rdotn / ddotn;
   rxy = rxy0 + snext * (dir[0] * cphi + dir[1] * sphi);
   if (rxy < 0)
      return kFALSE;
   return kTRUE;
}
 
////////////////////////////////////////////////////////////////////////////////
/// Check if two numbers differ with less than a tolerance.
 
Bool_t TGeoShape::IsSameWithinTolerance(Double_t a, Double_t b)
{
   if (TMath::Abs(a - b) < 1.E-10)
      return kTRUE;
   return kFALSE;
}
 
////////////////////////////////////////////////////////////////////////////////
/// Check if segments (A,B) and (C,D) are crossing,
/// where: A(x1,y1), B(x2,y2), C(x3,y3), D(x4,y4)
 
Bool_t TGeoShape::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)
{
   Double_t eps = TGeoShape::Tolerance();
   Bool_t stand1 = kFALSE;
   Double_t dx1 = x2 - x1;
   Bool_t stand2 = kFALSE;
   Double_t dx2 = x4 - x3;
   Double_t xm = 0.;
   Double_t ym = 0.;
   Double_t a1 = 0.;
   Double_t b1 = 0.;
   Double_t a2 = 0.;
   Double_t b2 = 0.;
   if (TMath::Abs(dx1) < eps)
      stand1 = kTRUE;
   if (TMath::Abs(dx2) < eps)
      stand2 = kTRUE;
   if (!stand1) {
      a1 = (x2 * y1 - x1 * y2) / dx1;
      b1 = (y2 - y1) / dx1;
   }
   if (!stand2) {
      a2 = (x4 * y3 - x3 * y4) / dx2;
      b2 = (y4 - y3) / dx2;
   }
   if (stand1 && stand2) {
      // Segments parallel and vertical
      if (TMath::Abs(x1 - x3) < eps) {
         // Check if segments are overlapping
         if ((y3 - y1) * (y3 - y2) < -eps || (y4 - y1) * (y4 - y2) < -eps || (y1 - y3) * (y1 - y4) < -eps ||
             (y2 - y3) * (y2 - y4) < -eps)
            return kTRUE;
         return kFALSE;
      }
      // Different x values
      return kFALSE;
   }
 
   if (stand1) {
      // First segment vertical
      xm = x1;
      ym = a2 + b2 * xm;
   } else {
      if (stand2) {
         // Second segment vertical
         xm = x3;
         ym = a1 + b1 * xm;
      } else {
         // Normal crossing
         if (TMath::Abs(b1 - b2) < eps) {
            // Parallel segments, are they aligned
            if (TMath::Abs(y3 - (a1 + b1 * x3)) > eps)
               return kFALSE;
            // Aligned segments, are they overlapping
            if ((x3 - x1) * (x3 - x2) < -eps || (x4 - x1) * (x4 - x2) < -eps || (x1 - x3) * (x1 - x4) < -eps ||
                (x2 - x3) * (x2 - x4) < -eps)
               return kTRUE;
            return kFALSE;
         }
         xm = (a1 - a2) / (b2 - b1);
         ym = (a1 * b2 - a2 * b1) / (b2 - b1);
      }
   }
   // Check if crossing point is both between A,B and C,D
   Double_t check = (xm - x1) * (xm - x2) + (ym - y1) * (ym - y2);
   if (check > -eps)
      return kFALSE;
   check = (xm - x3) * (xm - x4) + (ym - y3) * (ym - y4);
   if (check > -eps)
      return kFALSE;
   return kTRUE;
}
 
////////////////////////////////////////////////////////////////////////////////
/// compute distance from point (inside phi) to both phi planes. Return minimum.
 
Double_t TGeoShape::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)
{
   Double_t sfi1 = TGeoShape::Big();
   Double_t sfi2 = TGeoShape::Big();
   Double_t s = 0;
   Double_t un = dir[0] * s1 - dir[1] * c1;
   if (!in)
      un = -un;
   if (un > 0) {
      s = -point[0] * s1 + point[1] * c1;
      if (!in)
         s = -s;
      if (s >= 0) {
         s /= un;
         if (((point[0] + s * dir[0]) * sm - (point[1] + s * dir[1]) * cm) >= 0)
            sfi1 = s;
      }
   }
   un = -dir[0] * s2 + dir[1] * c2;
   if (!in)
      un = -un;
   if (un > 0) {
      s = point[0] * s2 - point[1] * c2;
      if (!in)
         s = -s;
      if (s >= 0) {
         s /= un;
         if ((-(point[0] + s * dir[0]) * sm + (point[1] + s * dir[1]) * cm) >= 0)
            sfi2 = s;
      }
   }
   return TMath::Min(sfi1, sfi2);
}
 
////////////////////////////////////////////////////////////////////////////////
/// Static method to compute normal to phi planes.
 
void TGeoShape::NormalPhi(const Double_t *point, const Double_t *dir, Double_t *norm, Double_t c1, Double_t s1,
                          Double_t c2, Double_t s2)
{
   Double_t saf1 = TGeoShape::Big();
   Double_t saf2 = TGeoShape::Big();
   if (point[0] * c1 + point[1] * s1 >= 0)
      saf1 = TMath::Abs(-point[0] * s1 + point[1] * c1);
   if (point[0] * c2 + point[1] * s2 >= 0)
      saf2 = TMath::Abs(point[0] * s2 - point[1] * c2);
   Double_t c, s;
   if (saf1 < saf2) {
      c = c1;
      s = s1;
   } else {
      c = c2;
      s = s2;
   }
   norm[2] = 0;
   norm[0] = -s;
   norm[1] = c;
   if (dir[0] * norm[0] + dir[1] * norm[1] < 0) {
      norm[0] = s;
      norm[1] = -c;
   }
}
 
////////////////////////////////////////////////////////////////////////////////
/// Static method to compute safety w.r.t a phi corner defined by cosines/sines
/// of the angles phi1, phi2.
 
Double_t TGeoShape::SafetyPhi(const Double_t *point, Bool_t in, Double_t phi1, Double_t phi2)
{
   Bool_t inphi = TGeoShape::IsInPhiRange(point, phi1, phi2);
   if (inphi && !in)
      return -TGeoShape::Big();
   phi1 *= TMath::DegToRad();
   phi2 *= TMath::DegToRad();
   Double_t c1 = TMath::Cos(phi1);
   Double_t s1 = TMath::Sin(phi1);
   Double_t c2 = TMath::Cos(phi2);
   Double_t s2 = TMath::Sin(phi2);
   Double_t rsq = point[0] * point[0] + point[1] * point[1];
   Double_t rproj = point[0] * c1 + point[1] * s1;
   Double_t safsq = rsq - rproj * rproj;
   if (safsq < 0)
      return 0.;
   Double_t saf1 = (rproj < 0) ? TGeoShape::Big() : TMath::Sqrt(safsq);
   rproj = point[0] * c2 + point[1] * s2;
   safsq = rsq - rproj * rproj;
   if (safsq < 0)
      return 0.;
   Double_t saf2 = (rproj < 0) ? TGeoShape::Big() : TMath::Sqrt(safsq);
   Double_t safe = TMath::Min(saf1, saf2); // >0
   if (safe > 1E10) {
      if (in)
         return TGeoShape::Big();
      return -TGeoShape::Big();
   }
   return safe;
}
 
////////////////////////////////////////////////////////////////////////////////
/// Compute distance from point of coordinates (r,z) to segment (r1,z1):(r2,z2)
 
Double_t TGeoShape::SafetySeg(Double_t r, Double_t z, Double_t r1, Double_t z1, Double_t r2, Double_t z2, Bool_t outer)
{
   Double_t crossp = (z2 - z1) * (r - r1) - (z - z1) * (r2 - r1);
   crossp *= (outer) ? 1. : -1.;
   // Positive crossp means point on the requested side of the (1,2) segment
   if (crossp < -TGeoShape::Tolerance()) {
      //      if (((z-z1)*(z2-z)) > -1.E-10) return 0;
      if (outer)
         return TGeoShape::Big();
      else
         return 0.;
   }
   // Compute (1,P) dot (1,2)
   Double_t c1 = (z - z1) * (z2 - z1) + (r - r1) * (r2 - r1);
   // Negative c1 means point (1) is closest
   if (c1 < 1.E-10)
      return TMath::Sqrt((r - r1) * (r - r1) + (z - z1) * (z - z1));
   // Compute (2,P) dot (1,2)
   Double_t c2 = (z - z2) * (z2 - z1) + (r - r2) * (r2 - r1);
   // Positive c2 means point (2) is closest
   if (c2 > -1.E-10)
      return TMath::Sqrt((r - r2) * (r - r2) + (z - z2) * (z - z2));
   // The closest point is between (1) and (2)
   c2 = (z2 - z1) * (z2 - z1) + (r2 - r1) * (r2 - r1);
   // projected length factor with respect to (1,2) length
   Double_t alpha = c1 / c2;
   Double_t rp = r1 + alpha * (r2 - r1);
   Double_t zp = z1 + alpha * (z2 - z1);
   return TMath::Sqrt((r - rp) * (r - rp) + (z - zp) * (z - zp));
}
 
////////////////////////////////////////////////////////////////////////////////
/// Equivalent of TObject::SetBit.
 
void TGeoShape::SetShapeBit(UInt_t f, Bool_t set)
{
   if (set) {
      SetShapeBit(f);
   } else {
      ResetShapeBit(f);
   }
}
 
////////////////////////////////////////////////////////////////////////////////
/// Returns current transformation matrix that applies to shape.
 
TGeoMatrix *TGeoShape::GetTransform()
{
   return fgTransform;
}
 
////////////////////////////////////////////////////////////////////////////////
/// Set current transformation matrix that applies to shape.
 
void TGeoShape::SetTransform(TGeoMatrix *matrix)
{
   fgTransform = matrix;
}
 
////////////////////////////////////////////////////////////////////////////////
/// Tranform a set of points (LocalToMaster)
 
void TGeoShape::TransformPoints(Double_t *points, UInt_t NbPnts) const
{
   UInt_t i, j;
   Double_t dmaster[3];
   if (fgTransform) {
      for (j = 0; j < NbPnts; j++) {
         i = 3 * j;
         fgTransform->LocalToMaster(&points[i], dmaster);
         points[i] = dmaster[0];
         points[i + 1] = dmaster[1];
         points[i + 2] = dmaster[2];
      }
      return;
   }
   if (!gGeoManager)
      return;
   Bool_t bomb = (gGeoManager->GetBombMode() == 0) ? kFALSE : kTRUE;
 
   for (j = 0; j < NbPnts; j++) {
      i = 3 * j;
      if (gGeoManager->IsMatrixTransform()) {
         TGeoHMatrix *glmat = gGeoManager->GetGLMatrix();
         if (bomb)
            glmat->LocalToMasterBomb(&points[i], dmaster);
         else
            glmat->LocalToMaster(&points[i], dmaster);
      } else {
         if (bomb)
            gGeoManager->LocalToMasterBomb(&points[i], dmaster);
         else
            gGeoManager->LocalToMaster(&points[i], dmaster);
      }
      points[i] = dmaster[0];
      points[i + 1] = dmaster[1];
      points[i + 2] = dmaster[2];
   }
}
 
////////////////////////////////////////////////////////////////////////////////
/// Fill the supplied buffer, with sections in desired frame
/// See TBuffer3D.h for explanation of sections, frame etc.
 
void TGeoShape::FillBuffer3D(TBuffer3D &buffer, Int_t reqSections, Bool_t localFrame) const
{
   // Catch this common potential error here
   // We have to set kRawSize (unless already done) to allocate buffer space
   // before kRaw can be filled
   if (reqSections & TBuffer3D::kRaw) {
      if (!(reqSections & TBuffer3D::kRawSizes) && !buffer.SectionsValid(TBuffer3D::kRawSizes)) {
         R__ASSERT(kFALSE);
      }
   }
 
   if (reqSections & TBuffer3D::kCore) {
      // If writing core section all others will be invalid
      buffer.ClearSectionsValid();
 
      // Check/grab some objects we need
      if (!gGeoManager) {
         R__ASSERT(kFALSE);
         return;
      }
      const TGeoVolume *paintVolume = gGeoManager->GetPaintVolume();
      if (!paintVolume)
         paintVolume = gGeoManager->GetTopVolume();
      if (!paintVolume) {
         buffer.fID = const_cast<TGeoShape *>(this);
         buffer.fColor = 0;
         buffer.fTransparency = 0;
         //         R__ASSERT(kFALSE);
         //         return;
      } else {
         buffer.fID = const_cast<TGeoVolume *>(paintVolume);
         buffer.fColor = paintVolume->GetLineColor();
 
         buffer.fTransparency = paintVolume->GetTransparency();
         Double_t visdensity = gGeoManager->GetVisDensity();
         if (visdensity > 0 && paintVolume->GetMedium()) {
            if (paintVolume->GetMaterial()->GetDensity() < visdensity) {
               buffer.fTransparency = 90;
            }
         }
      }
 
      buffer.fLocalFrame = localFrame;
      Bool_t r1, r2 = kFALSE;
      r1 = gGeoManager->IsMatrixReflection();
      if (paintVolume && paintVolume->GetShape()) {
         if (paintVolume->GetShape()->IsReflected()) {
            // Temporary trick to deal with reflected shapes.
            // Still lighting gets wrong...
            if (buffer.Type() < TBuffer3DTypes::kTube)
               r2 = kTRUE;
         }
      }
      buffer.fReflection = ((r1 & (!r2)) | (r2 & !(r1)));
 
      // Set up local -> master translation matrix
      if (localFrame) {
         TGeoMatrix *localMasterMat = nullptr;
         if (TGeoShape::GetTransform()) {
            localMasterMat = TGeoShape::GetTransform();
         } else {
            localMasterMat = gGeoManager->GetCurrentMatrix();
 
            // For overlap drawing the correct matrix needs to obtained in
            // from GetGLMatrix() - this should not be applied in the case
            // of composite shapes
            if (gGeoManager->IsMatrixTransform() && !IsComposite()) {
               localMasterMat = gGeoManager->GetGLMatrix();
            }
         }
         if (!localMasterMat) {
            R__ASSERT(kFALSE);
            return;
         }
         localMasterMat->GetHomogenousMatrix(buffer.fLocalMaster);
      } else {
         buffer.SetLocalMasterIdentity();
      }
 
      buffer.SetSectionsValid(TBuffer3D::kCore);
   }
}
 
////////////////////////////////////////////////////////////////////////////////
/// Get the basic color (0-7).
 
Int_t TGeoShape::GetBasicColor() const
{
   Int_t basicColor = 0; // TODO: Check on sensible fallback
   if (gGeoManager) {
      const TGeoVolume *volume = gGeoManager->GetPaintVolume();
      if (volume) {
         basicColor = ((volume->GetLineColor() % 8) - 1) * 4;
         if (basicColor < 0)
            basicColor = 0;
      }
   }
   return basicColor;
}
 
////////////////////////////////////////////////////////////////////////////////
/// Stub implementation to avoid forcing implementation at this stage
 
const TBuffer3D &TGeoShape::GetBuffer3D(Int_t /*reqSections*/, Bool_t /*localFrame*/) const
{
   static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
   Warning("GetBuffer3D", "this must be implemented for shapes in a TGeoPainter hierarchy. This will be come a pure "
                          "virtual fn eventually.");
   return buffer;
}
 
////////////////////////////////////////////////////////////////////////////////
/// Provide a pointer name containing uid.
 
const char *TGeoShape::GetPointerName() const
{
   static TString name;
   Int_t uid = GetUniqueID();
   if (uid)
      name = TString::Format("p%s_%d", GetName(), uid);
   else
      name = TString::Format("p%s", GetName());
   return name.Data();
}
 
////////////////////////////////////////////////////////////////////////////////
/// Execute mouse actions on this shape.
 
void TGeoShape::ExecuteEvent(Int_t event, Int_t px, Int_t py)
{
   if (!gGeoManager)
      return;
   TVirtualGeoPainter *painter = gGeoManager->GetPainter();
   painter->ExecuteShapeEvent(this, event, px, py);
}
 
////////////////////////////////////////////////////////////////////////////////
/// Draw this shape.
 
void TGeoShape::Draw(Option_t *option)
{
   TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter();
   if (option && option[0]) {
      painter->DrawShape(this, option);
   } else {
      painter->DrawShape(this, gEnv->GetValue("Viewer3D.DefaultDrawOption", ""));
   }
}
 
////////////////////////////////////////////////////////////////////////////////
/// Paint this shape.
 
void TGeoShape::Paint(Option_t *option)
{
   TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter();
   if (option && option[0]) {
      painter->PaintShape(this, option);
   } else {
      painter->PaintShape(this, gEnv->GetValue("Viewer3D.DefaultDrawOption", ""));
   }
}