#include <G4BoundingEnvelope.hh>

Public Member Functions
G4bool	BoundingBoxVsVoxelLimits (const EAxis pAxis, const G4VoxelLimits &pVoxelLimits, const G4Transform3D &pTransform3D, G4double &pMin, G4double &pMax) const

G4bool	CalculateExtent (const EAxis pAxis, const G4VoxelLimits &pVoxelLimits, const G4Transform3D &pTransform3D, G4double &pMin, G4double &pMax) const

	G4BoundingEnvelope (const G4ThreeVector &pMin, const G4ThreeVector &pMax)

	G4BoundingEnvelope (const G4ThreeVector &pMin, const G4ThreeVector &pMax, const std::vector< const G4ThreeVectorList * > &polygons)

	G4BoundingEnvelope (const std::vector< const G4ThreeVectorList * > &polygons)

	~G4BoundingEnvelope ()=default

Private Member Functions
void	CheckBoundingBox ()

void	CheckBoundingPolygons ()

G4bool	ClipEdgesByVoxel (const std::vector< G4Segment3D > &pEdges, const G4VoxelLimits &pLimits, G4Segment3D &pExtent) const

void	ClipVoxelByPlanes (G4int pBits, const G4VoxelLimits &pLimits, const std::vector< G4Plane3D > &pPlanes, const G4Segment3D &pAABB, G4Segment3D &pExtent) const

void	CreateListOfEdges (const G4Polygon3D &baseA, const G4Polygon3D &baseB, std::vector< G4Segment3D > &pEdges) const

void	CreateListOfPlanes (const G4Polygon3D &baseA, const G4Polygon3D &baseB, std::vector< G4Plane3D > &pPlanes) const

G4double	FindScaleFactor (const G4Transform3D &pTransform3D) const

void	GetPrismAABB (const G4Polygon3D &pBaseA, const G4Polygon3D &pBaseB, G4Segment3D &pAABB) const

void	TransformVertices (const G4Transform3D &pTransform3D, std::vector< G4Polygon3D * > &pBases) const

Private Attributes
G4ThreeVector	fMax

G4ThreeVector	fMin

const std::vector< const G4ThreeVectorList * > *	fPolygons = nullptr

Detailed Description

Definition at line 58 of file G4BoundingEnvelope.hh.

Constructor & Destructor Documentation

◆ G4BoundingEnvelope() [1/3]

G4BoundingEnvelope::G4BoundingEnvelope	(	const G4ThreeVector &	pMin,
		const G4ThreeVector &	pMax
	)

Definition at line 44 of file G4BoundingEnvelope.cc.

  : fMin(pMin), fMax(pMax)
{
  // Check correctness of bounding box
  //
  CheckBoundingBox();
}

References CheckBoundingBox().

◆ G4BoundingEnvelope() [2/3]

G4BoundingEnvelope::G4BoundingEnvelope ( const std::vector< const G4ThreeVectorList * > & polygons )

Definition at line 57 of file G4BoundingEnvelope.cc.

  : fPolygons(&polygons)
{
  // Check correctness of polygons
  //
  CheckBoundingPolygons();
 
  // Set bounding box
  //
  G4double xmin =  kInfinity, ymin =  kInfinity, zmin =  kInfinity;
  G4double xmax = -kInfinity, ymax = -kInfinity, zmax = -kInfinity;
  for (auto ibase = fPolygons->cbegin(); ibase != fPolygons->cend(); ++ibase)
  { 
    for (auto ipoint = (*ibase)->cbegin(); ipoint != (*ibase)->cend(); ++ipoint)
    {
      G4double x = ipoint->x(); 
      if (x < xmin) xmin = x;
      if (x > xmax) xmax = x;
      G4double y = ipoint->y(); 
      if (y < ymin) ymin = y;
      if (y > ymax) ymax = y;
      G4double z = ipoint->z(); 
      if (z < zmin) zmin = z;
      if (z > zmax) zmax = z;
    }
  }
  fMin.set(xmin,ymin,zmin);
  fMax.set(xmax,ymax,zmax);
 
  // Check correctness of bounding box
  //
  CheckBoundingBox();
}

References CheckBoundingBox(), CheckBoundingPolygons(), fMax, fMin, fPolygons, kInfinity, and CLHEP::Hep3Vector::set().

◆ G4BoundingEnvelope() [3/3]

G4BoundingEnvelope::G4BoundingEnvelope	(	const G4ThreeVector &	pMin,
		const G4ThreeVector &	pMax,
		const std::vector< const G4ThreeVectorList * > &	polygons
	)

Definition at line 96 of file G4BoundingEnvelope.cc.

  : fMin(pMin), fMax(pMax), fPolygons(&polygons)
{
  // Check correctness of bounding box and polygons
  // 
  CheckBoundingBox();
  CheckBoundingPolygons();
}

References CheckBoundingBox(), and CheckBoundingPolygons().

◆ ~G4BoundingEnvelope()

G4BoundingEnvelope::~G4BoundingEnvelope ( )

default

Member Function Documentation

◆ BoundingBoxVsVoxelLimits()

G4bool G4BoundingEnvelope::BoundingBoxVsVoxelLimits	(	const EAxis	pAxis,
		const G4VoxelLimits &	pVoxelLimits,
		const G4Transform3D &	pTransform3D,
		G4double &	pMin,
		G4double &	pMax
	)		const

Definition at line 180 of file G4BoundingEnvelope.cc.

{
  pMin =  kInfinity;
  pMax = -kInfinity;
  G4double xminlim = pVoxelLimits.GetMinXExtent();
  G4double xmaxlim = pVoxelLimits.GetMaxXExtent();
  G4double yminlim = pVoxelLimits.GetMinYExtent();
  G4double ymaxlim = pVoxelLimits.GetMaxYExtent();
  G4double zminlim = pVoxelLimits.GetMinZExtent();
  G4double zmaxlim = pVoxelLimits.GetMaxZExtent();
 
  // Special case of pure translation
  //
  if (pTransform3D.xx()==1 && pTransform3D.yy()==1 && pTransform3D.zz()==1)
  {
    G4double xmin = fMin.x() + pTransform3D.dx();
    G4double xmax = fMax.x() + pTransform3D.dx();
    G4double ymin = fMin.y() + pTransform3D.dy();
    G4double ymax = fMax.y() + pTransform3D.dy();
    G4double zmin = fMin.z() + pTransform3D.dz();
    G4double zmax = fMax.z() + pTransform3D.dz();
 
    if (xmin-kCarTolerance > xmaxlim) return true;
    if (xmax+kCarTolerance < xminlim) return true;
    if (ymin-kCarTolerance > ymaxlim) return true;
    if (ymax+kCarTolerance < yminlim) return true;
    if (zmin-kCarTolerance > zmaxlim) return true;
    if (zmax+kCarTolerance < zminlim) return true;
 
    if (xmin >= xminlim && xmax <= xmaxlim &&
        ymin >= yminlim && ymax <= ymaxlim &&
        zmin >= zminlim && zmax <= zmaxlim)
    {
      if (pAxis == kXAxis)
      {
        pMin = (xmin-kCarTolerance < xminlim) ? xminlim : xmin;
        pMax = (xmax+kCarTolerance > xmaxlim) ? xmaxlim : xmax;
      } 
      else if (pAxis == kYAxis)
      {
        pMin = (ymin-kCarTolerance < yminlim) ? yminlim : ymin;
        pMax = (ymax+kCarTolerance > ymaxlim) ? ymaxlim : ymax;
      }
      else if (pAxis == kZAxis)
      {
        pMin = (zmin-kCarTolerance < zminlim) ? zminlim : zmin;
        pMax = (zmax+kCarTolerance > zmaxlim) ? zmaxlim : zmax;
      }
      pMin -= kCarTolerance;
      pMax += kCarTolerance;
      return true;
    }
  }
 
  // Find max scale factor of the transformation, set delta 
  // equal to kCarTolerance multiplied by the scale factor
  //
  G4double scale = FindScaleFactor(pTransform3D);
  G4double delta = kCarTolerance*scale;  
 
  // Set the sphere surrounding the bounding box
  //
  G4Point3D center = pTransform3D*G4Point3D(0.5*(fMin+fMax));
  G4double  radius = 0.5*scale*(fMax-fMin).mag() + delta;
 
  // Check if the sphere surrounding the bounding box is outside
  // the voxel limits
  //
  if (center.x()-radius > xmaxlim) return true;
  if (center.y()-radius > ymaxlim) return true;
  if (center.z()-radius > zmaxlim) return true;
  if (center.x()+radius < xminlim) return true;
  if (center.y()+radius < yminlim) return true;
  if (center.z()+radius < zminlim) return true;
  return false;
}

References HepGeom::Transform3D::dx(), HepGeom::Transform3D::dy(), HepGeom::Transform3D::dz(), FindScaleFactor(), fMax, fMin, G4VoxelLimits::GetMaxXExtent(), G4VoxelLimits::GetMaxYExtent(), G4VoxelLimits::GetMaxZExtent(), G4VoxelLimits::GetMinXExtent(), G4VoxelLimits::GetMinYExtent(), G4VoxelLimits::GetMinZExtent(), kCarTolerance, kInfinity, kXAxis, kYAxis, kZAxis, pMax, pMin, HepGeom::BasicVector3D< T >::x(), CLHEP::Hep3Vector::x(), HepGeom::Transform3D::xx(), HepGeom::BasicVector3D< T >::y(), CLHEP::Hep3Vector::y(), HepGeom::Transform3D::yy(), HepGeom::BasicVector3D< T >::z(), CLHEP::Hep3Vector::z(), and HepGeom::Transform3D::zz().

◆ CalculateExtent()

G4bool G4BoundingEnvelope::CalculateExtent	(	const EAxis	pAxis,
		const G4VoxelLimits &	pVoxelLimits,
		const G4Transform3D &	pTransform3D,
		G4double &	pMin,
		G4double &	pMax
	)		const

Definition at line 266 of file G4BoundingEnvelope.cc.

{
  pMin =  kInfinity;
  pMax = -kInfinity;
  G4double xminlim = pVoxelLimits.GetMinXExtent();
  G4double xmaxlim = pVoxelLimits.GetMaxXExtent();
  G4double yminlim = pVoxelLimits.GetMinYExtent();
  G4double ymaxlim = pVoxelLimits.GetMaxYExtent();
  G4double zminlim = pVoxelLimits.GetMinZExtent();
  G4double zmaxlim = pVoxelLimits.GetMaxZExtent();
 
  // Special case of pure translation
  //
  if (pTransform3D.xx()==1 && pTransform3D.yy()==1 && pTransform3D.zz()==1)
  {
    G4double xmin = fMin.x() + pTransform3D.dx();
    G4double xmax = fMax.x() + pTransform3D.dx();
    G4double ymin = fMin.y() + pTransform3D.dy();
    G4double ymax = fMax.y() + pTransform3D.dy();
    G4double zmin = fMin.z() + pTransform3D.dz();
    G4double zmax = fMax.z() + pTransform3D.dz();
 
    if (xmin-kCarTolerance > xmaxlim) return false;
    if (xmax+kCarTolerance < xminlim) return false;
    if (ymin-kCarTolerance > ymaxlim) return false;
    if (ymax+kCarTolerance < yminlim) return false;
    if (zmin-kCarTolerance > zmaxlim) return false;
    if (zmax+kCarTolerance < zminlim) return false;
 
    if (fPolygons == nullptr)
    {
      if (pAxis == kXAxis)
      {
        pMin = (xmin-kCarTolerance < xminlim) ? xminlim : xmin;
        pMax = (xmax+kCarTolerance > xmaxlim) ? xmaxlim : xmax;
      } 
      else if (pAxis == kYAxis)
      {
        pMin = (ymin-kCarTolerance < yminlim) ? yminlim : ymin;
        pMax = (ymax+kCarTolerance > ymaxlim) ? ymaxlim : ymax;
      }
      else if (pAxis == kZAxis)
      {
        pMin = (zmin-kCarTolerance < zminlim) ? zminlim : zmin;
        pMax = (zmax+kCarTolerance > zmaxlim) ? zmaxlim : zmax;
      }
      pMin -= kCarTolerance;
      pMax += kCarTolerance;
      return true;
    }
  }
 
  // Find max scale factor of the transformation, set delta 
  // equal to kCarTolerance multiplied by the scale factor
  //
  G4double scale = FindScaleFactor(pTransform3D);
  G4double delta = kCarTolerance*scale;  
 
  // Set the sphere surrounding the bounding box
  //
  G4Point3D center = pTransform3D*G4Point3D(0.5*(fMin+fMax));
  G4double  radius = 0.5*scale*(fMax-fMin).mag() + delta;
 
  // Check if the sphere surrounding the bounding box is within
  // the voxel limits, if so then transform only one coordinate
  //
  if (center.x()-radius >= xminlim && center.x()+radius <= xmaxlim &&
      center.y()-radius >= yminlim && center.y()+radius <= ymaxlim &&
      center.z()-radius >= zminlim && center.z()+radius <= zmaxlim )
  {
    G4double cx, cy, cz, cd;
    if (pAxis == kXAxis)
    { 
      cx = pTransform3D.xx(); 
      cy = pTransform3D.xy();
      cz = pTransform3D.xz();
      cd = pTransform3D.dx();
    }
    else if (pAxis == kYAxis)
    { 
      cx = pTransform3D.yx(); 
      cy = pTransform3D.yy();
      cz = pTransform3D.yz();
      cd = pTransform3D.dy();
    }
    else if (pAxis == kZAxis)
    { 
      cx = pTransform3D.zx(); 
      cy = pTransform3D.zy();
      cz = pTransform3D.zz();
      cd = pTransform3D.dz();
    }
    else
    { 
      cx = cy = cz = cd = kInfinity;
    }
    G4double emin = kInfinity, emax = -kInfinity;
    if (fPolygons == nullptr)
    {
      G4double coor;
      coor = cx*fMin.x() + cy*fMin.y() + cz*fMin.z() + cd;
      if (coor < emin) emin = coor;
      if (coor > emax) emax = coor;
      coor = cx*fMax.x() + cy*fMin.y() + cz*fMin.z() + cd;
      if (coor < emin) emin = coor;
      if (coor > emax) emax = coor;
      coor = cx*fMax.x() + cy*fMax.y() + cz*fMin.z() + cd;
      if (coor < emin) emin = coor;
      if (coor > emax) emax = coor;
      coor = cx*fMin.x() + cy*fMax.y() + cz*fMin.z() + cd;
      if (coor < emin) emin = coor;
      if (coor > emax) emax = coor;
      coor = cx*fMin.x() + cy*fMin.y() + cz*fMax.z() + cd;
      if (coor < emin) emin = coor;
      if (coor > emax) emax = coor;
      coor = cx*fMax.x() + cy*fMin.y() + cz*fMax.z() + cd;
      if (coor < emin) emin = coor;
      if (coor > emax) emax = coor;
      coor = cx*fMax.x() + cy*fMax.y() + cz*fMax.z() + cd;
      if (coor < emin) emin = coor;
      if (coor > emax) emax = coor;
      coor = cx*fMin.x() + cy*fMax.y() + cz*fMax.z() + cd;
      if (coor < emin) emin = coor;
      if (coor > emax) emax = coor;
    }
    else
    {
      for (auto ibase=fPolygons->cbegin(); ibase!=fPolygons->cend(); ++ibase)
      { 
        for (auto ipoint=(*ibase)->cbegin(); ipoint!=(*ibase)->cend(); ++ipoint)
        {
          G4double coor = ipoint->x()*cx + ipoint->y()*cy + ipoint->z()*cz + cd;
          if (coor < emin) emin = coor;
          if (coor > emax) emax = coor;
        }
      }
    }
    pMin = emin - delta;
    pMax = emax + delta;
    return true;
  } 
 
  // Check if the sphere surrounding the bounding box is outside
  // the voxel limits
  //
  if (center.x()-radius > xmaxlim) return false;
  if (center.y()-radius > ymaxlim) return false;
  if (center.z()-radius > zmaxlim) return false;
  if (center.x()+radius < xminlim) return false;
  if (center.y()+radius < yminlim) return false;
  if (center.z()+radius < zminlim) return false;
 
  // Allocate memory for transformed polygons
  //
  G4int nbases = (fPolygons == nullptr) ? 2 : fPolygons->size();
  std::vector<G4Polygon3D*> bases(nbases);
  if (fPolygons == nullptr)
  {
    bases[0] = new G4Polygon3D(4); 
    bases[1] = new G4Polygon3D(4);
  }
  else
  {
    for (G4int i=0; i<nbases; ++i)
    { 
      bases[i] = new G4Polygon3D((*fPolygons)[i]->size());
    }
  }
 
  //  Transform vertices 
  //
  TransformVertices(pTransform3D, bases);
 
  // Create adjusted G4VoxelLimits box. New limits are extended by 
  // delta, kCarTolerance multiplied by max scale factor of
  // the transformation
  //
  EAxis axis[] = { kXAxis,kYAxis,kZAxis };
  G4VoxelLimits limits; // default is unlimited 
  for (auto i=0; i<3; ++i)
  {
    if (pVoxelLimits.IsLimited(axis[i]))
    {
      G4double emin = pVoxelLimits.GetMinExtent(axis[i]) - delta;
      G4double emax = pVoxelLimits.GetMaxExtent(axis[i]) + delta;
      limits.AddLimit(axis[i], emin, emax);
    }
  }
 
  // Main loop along the set of prisms
  //
  G4Segment3D extent;
  extent.first  = G4Point3D( kInfinity, kInfinity, kInfinity);
  extent.second = G4Point3D(-kInfinity,-kInfinity,-kInfinity);
  for (G4int k=0; k<nbases-1; ++k)
  {
    // Find bounding box of current prism
    G4Polygon3D* baseA = bases[k];
    G4Polygon3D* baseB = bases[k+1];
    G4Segment3D  prismAABB;
    GetPrismAABB(*baseA, *baseB, prismAABB);
 
    // Check if prismAABB is completely within the voxel limits
    if (prismAABB.first.x() >= limits.GetMinXExtent() &&
        prismAABB.first.y() >= limits.GetMinYExtent() &&
        prismAABB.first.z() >= limits.GetMinZExtent() &&
        prismAABB.second.x()<= limits.GetMaxXExtent() &&
        prismAABB.second.y()<= limits.GetMaxYExtent() &&
        prismAABB.second.z()<= limits.GetMaxZExtent()) 
    {
      if (extent.first.x()  > prismAABB.first.x())
        extent.first.setX( prismAABB.first.x() ); 
      if (extent.first.y()  > prismAABB.first.y())
        extent.first.setY( prismAABB.first.y() ); 
      if (extent.first.z()  > prismAABB.first.z())
        extent.first.setZ( prismAABB.first.z() ); 
      if (extent.second.x() < prismAABB.second.x())
        extent.second.setX(prismAABB.second.x()); 
      if (extent.second.y() < prismAABB.second.y())
        extent.second.setY(prismAABB.second.y()); 
      if (extent.second.z() < prismAABB.second.z())
        extent.second.setZ(prismAABB.second.z());
      continue; 
    }
 
    // Check if prismAABB is outside the voxel limits
    if (prismAABB.first.x()  > limits.GetMaxXExtent()) continue;
    if (prismAABB.first.y()  > limits.GetMaxYExtent()) continue;
    if (prismAABB.first.z()  > limits.GetMaxZExtent()) continue;
    if (prismAABB.second.x() < limits.GetMinXExtent()) continue;
    if (prismAABB.second.y() < limits.GetMinYExtent()) continue;
    if (prismAABB.second.z() < limits.GetMinZExtent()) continue;
 
    // Clip edges of the prism by adjusted G4VoxelLimits box
    std::vector<G4Segment3D> vecEdges; 
    CreateListOfEdges(*baseA, *baseB, vecEdges);
    if (ClipEdgesByVoxel(vecEdges, limits, extent)) continue;
 
    // Some edges of the prism are completely outside of the voxel
    // limits, clip selected edges (see bits) of adjusted G4VoxelLimits
    // by the prism 
    G4int bits = 0x000;
    if (limits.GetMinXExtent() < prismAABB.first.x())
      bits |= 0x988; // 1001 1000 1000
    if (limits.GetMaxXExtent() > prismAABB.second.x())
      bits |= 0x622; // 0110 0010 0010
 
    if (limits.GetMinYExtent() < prismAABB.first.y())
      bits |= 0x311; // 0011 0001 0001
    if (limits.GetMaxYExtent() > prismAABB.second.y())
      bits |= 0xC44; // 1100 0100 0100
 
    if (limits.GetMinZExtent() < prismAABB.first.z())
      bits |= 0x00F; // 0000 0000 1111
    if (limits.GetMaxZExtent() > prismAABB.second.z())
      bits |= 0x0F0; // 0000 1111 0000
    if (bits == 0xFFF) continue;
 
    std::vector<G4Plane3D> vecPlanes; 
    CreateListOfPlanes(*baseA, *baseB, vecPlanes);
    ClipVoxelByPlanes(bits, limits, vecPlanes, prismAABB, extent);
  } // End of the main loop
 
  // Free memory
  //
  for (G4int i=0; i<nbases; ++i) { delete bases[i]; bases[i] = nullptr; }
 
  // Final adjustment of the extent
  // 
  G4double emin = 0, emax = 0;
  if (pAxis == kXAxis) { emin = extent.first.x(); emax = extent.second.x(); }
  if (pAxis == kYAxis) { emin = extent.first.y(); emax = extent.second.y(); }
  if (pAxis == kZAxis) { emin = extent.first.z(); emax = extent.second.z(); }
 
  if (emin > emax) return false;
  emin -= delta;
  emax += delta;
  G4double minlim = pVoxelLimits.GetMinExtent(pAxis);
  G4double maxlim = pVoxelLimits.GetMaxExtent(pAxis);
  pMin = (emin < minlim) ? minlim-kCarTolerance : emin;
  pMax = (emax > maxlim) ? maxlim+kCarTolerance : emax;
  return true;
}

◆ CheckBoundingBox()

void G4BoundingEnvelope::CheckBoundingBox ( )

private

Definition at line 112 of file G4BoundingEnvelope.cc.

{
  if (fMin.x() >= fMax.x() || fMin.y() >= fMax.y() || fMin.z() >= fMax.z())
  {
    std::ostringstream message;
    message << "Badly defined bounding box (min >= max)!"
            << "\npMin = " << fMin 
            << "\npMax = " << fMax;
    G4Exception("G4BoundingEnvelope::CheckBoundingBox()",
                "GeomMgt0001", JustWarning, message);
  }
}

References fMax, fMin, G4Exception(), JustWarning, CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

Referenced by G4BoundingEnvelope().

◆ CheckBoundingPolygons()

void G4BoundingEnvelope::CheckBoundingPolygons ( )

private

Definition at line 130 of file G4BoundingEnvelope.cc.

{
  G4int nbases = fPolygons->size();
  if (nbases < 2)
  {
    std::ostringstream message;
    message << "Wrong number of polygons in the sequence: " << nbases
            << "\nShould be at least two!";
    G4Exception("G4BoundingEnvelope::CheckBoundingPolygons()", 
                "GeomMgt0001", FatalException, message);
    return;
  }
 
  G4int nsize  = std::max((*fPolygons)[0]->size(),(*fPolygons)[1]->size());
  if (nsize < 3)
  {
    std::ostringstream message;
    message << "Badly constructed polygons!"
            << "\nNumber of polygons: " << nbases
            << "\nPolygon #0 size: " << (*fPolygons)[0]->size()
            << "\nPolygon #1 size: " << (*fPolygons)[1]->size()
            << "\n...";
    G4Exception("G4BoundingEnvelope::CheckBoundingPolygons()", 
                "GeomMgt0001", FatalException, message);
    return;
  }
 
  for (G4int k=0; k<nbases; ++k)
  {
    G4int np = (*fPolygons)[k]->size();
    if (np == nsize)            continue;
    if (np == 1 && k==0)        continue;
    if (np == 1 && k==nbases-1) continue;
    std::ostringstream message;
    message << "Badly constructed polygons!"
            << "\nNumber of polygons: " << nbases
            << "\nPolygon #" << k << " size: " << np
            << "\nexpected size: " << nsize;
    G4Exception("G4BoundingEnvelope::SetBoundingPolygons()", 
                "GeomMgt0001", FatalException, message);
    return;
  }  
}

References FatalException, fPolygons, G4Exception(), and G4INCL::Math::max().

Referenced by G4BoundingEnvelope().

◆ ClipEdgesByVoxel()

G4bool G4BoundingEnvelope::ClipEdgesByVoxel	(	const std::vector< G4Segment3D > &	pEdges,
		const G4VoxelLimits &	pLimits,
		G4Segment3D &	pExtent
	)		const

private

Definition at line 831 of file G4BoundingEnvelope.cc.

{
  G4bool    done = true;
  G4Point3D emin = pExtent.first;
  G4Point3D emax = pExtent.second;
 
  G4int nedges = pEdges.size();
  for (G4int k=0; k<nedges; ++k)
  {
    G4Point3D p1 = pEdges[k].first;
    G4Point3D p2 = pEdges[k].second;
    if (std::abs(p1.x()-p2.x())+
        std::abs(p1.y()-p2.y())+
        std::abs(p1.z()-p2.z()) < kCarTolerance) continue;
    G4double  d1, d2;
    // Clip current edge by X min
    d1 = pBox.GetMinXExtent() - p1.x(); 
    d2 = pBox.GetMinXExtent() - p2.x(); 
    if (d1 > 0.0)
    {
      if (d2 > 0.0) { done = false; continue; } // go to next edge
      p1 = (p2*d1-p1*d2)/(d1-d2);                   // move p1
    }
    else
    {
      if (d2 > 0.0) { p2 = (p1*d2-p2*d1)/(d2-d1); } // move p2
    }
 
    // Clip current edge by X max
    d1 = p1.x() - pBox.GetMaxXExtent(); 
    d2 = p2.x() - pBox.GetMaxXExtent(); 
    if (d1 > 0.)
    {
      if (d2 > 0.) { done = false; continue; } // go to next edge
      p1 = (p2*d1-p1*d2)/(d1-d2);
    }
    else
    {
      if (d2 > 0.) { p2 = (p1*d2-p2*d1)/(d2-d1); }
    }
 
    // Clip current edge by Y min
    d1 = pBox.GetMinYExtent() - p1.y(); 
    d2 = pBox.GetMinYExtent() - p2.y(); 
    if (d1 > 0.)
    {
      if (d2 > 0.) { done = false; continue; } // go to next edge
      p1 = (p2*d1-p1*d2)/(d1-d2);
    }
    else
    {
      if (d2 > 0.) { p2 = (p1*d2-p2*d1)/(d2-d1); }
    }
 
    // Clip current edge by Y max
    d1 = p1.y() - pBox.GetMaxYExtent(); 
    d2 = p2.y() - pBox.GetMaxYExtent(); 
    if (d1 > 0.)
    {
      if (d2 > 0.) { done = false; continue; } // go to next edge
      p1 = (p2*d1-p1*d2)/(d1-d2);
    }
    else
    {
      if (d2 > 0.) { p2 = (p1*d2-p2*d1)/(d2-d1); }
    }
 
    // Clip current edge by Z min
    d1 = pBox.GetMinZExtent() - p1.z(); 
    d2 = pBox.GetMinZExtent() - p2.z(); 
    if (d1 > 0.)
    {
      if (d2 > 0.) { done = false; continue; } // go to next edge
      p1 = (p2*d1-p1*d2)/(d1-d2);
    }
    else
    {
      if (d2 > 0.) { p2 = (p1*d2-p2*d1)/(d2-d1); }
    }
 
    // Clip current edge by Z max
    d1 = p1.z() - pBox.GetMaxZExtent(); 
    d2 = p2.z() - pBox.GetMaxZExtent(); 
    if (d1 > 0.)
    {
      if (d2 > 0.) { done = false; continue; } // go to next edge
      p1 = (p2*d1-p1*d2)/(d1-d2);
    }
    else
    {
      if (d2 > 0.) { p2 = (p1*d2-p2*d1)/(d2-d1); }
    }
 
    // Adjust current extent
    emin.setX(std::min(std::min(p1.x(),p2.x()),emin.x())); 
    emin.setY(std::min(std::min(p1.y(),p2.y()),emin.y())); 
    emin.setZ(std::min(std::min(p1.z(),p2.z()),emin.z())); 
 
    emax.setX(std::max(std::max(p1.x(),p2.x()),emax.x())); 
    emax.setY(std::max(std::max(p1.y(),p2.y()),emax.y())); 
    emax.setZ(std::max(std::max(p1.z(),p2.z()),emax.z())); 
  }
 
  // Return true if all edges (at least partially) are inside
  // the voxel limits, otherwise return false
  pExtent.first  = emin;
  pExtent.second = emax;
 
  return done;
}

References d1, d2, emax, G4VoxelLimits::GetMaxXExtent(), G4VoxelLimits::GetMaxYExtent(), G4VoxelLimits::GetMaxZExtent(), G4VoxelLimits::GetMinXExtent(), G4VoxelLimits::GetMinYExtent(), G4VoxelLimits::GetMinZExtent(), kCarTolerance, G4INCL::Math::max(), G4INCL::Math::min(), HepGeom::BasicVector3D< T >::setX(), HepGeom::BasicVector3D< T >::setY(), HepGeom::BasicVector3D< T >::setZ(), HepGeom::BasicVector3D< T >::x(), HepGeom::BasicVector3D< T >::y(), and HepGeom::BasicVector3D< T >::z().

Referenced by CalculateExtent().

◆ ClipVoxelByPlanes()

void G4BoundingEnvelope::ClipVoxelByPlanes	(	G4int	pBits,
		const G4VoxelLimits &	pLimits,
		const std::vector< G4Plane3D > &	pPlanes,
		const G4Segment3D &	pAABB,
		G4Segment3D &	pExtent
	)		const

private

Definition at line 949 of file G4BoundingEnvelope.cc.

{
  G4Point3D emin = pExtent.first;
  G4Point3D emax = pExtent.second;
 
  // Create edges of the voxel limits box reducing them where
  // appropriate to avoid calculations with big numbers (kInfinity)
  //
  G4double xmin = std::max(pBox.GetMinXExtent(),pAABB.first.x() -1.);
  G4double xmax = std::min(pBox.GetMaxXExtent(),pAABB.second.x()+1.);
 
  G4double ymin = std::max(pBox.GetMinYExtent(),pAABB.first.y() -1.);
  G4double ymax = std::min(pBox.GetMaxYExtent(),pAABB.second.y()+1.);
 
  G4double zmin = std::max(pBox.GetMinZExtent(),pAABB.first.z() -1.);
  G4double zmax = std::min(pBox.GetMaxZExtent(),pAABB.second.z()+1.);
 
  std::vector<G4Segment3D> edges(12);
  G4int i = 0, bits = pBits;
  if (!(bits & 0x001))
  { 
    edges[i  ].first.set( xmin,ymin,zmin); 
    edges[i++].second.set(xmax,ymin,zmin);
  }
  if (!(bits & 0x002))
  {
    edges[i  ].first.set( xmax,ymin,zmin);
    edges[i++].second.set(xmax,ymax,zmin);
  }
  if (!(bits & 0x004))
  {
    edges[i  ].first.set( xmax,ymax,zmin);
    edges[i++].second.set(xmin,ymax,zmin);
  }
  if (!(bits & 0x008))
  {
    edges[i  ].first.set( xmin,ymax,zmin);
    edges[i++].second.set(xmin,ymin,zmin); 
  }
 
  if (!(bits & 0x010))
  { 
    edges[i  ].first.set( xmin,ymin,zmax); 
    edges[i++].second.set(xmax,ymin,zmax);
  }
  if (!(bits & 0x020))
  {
    edges[i  ].first.set( xmax,ymin,zmax);
    edges[i++].second.set(xmax,ymax,zmax);
  }
  if (!(bits & 0x040))
  {
    edges[i  ].first.set( xmax,ymax,zmax);
    edges[i++].second.set(xmin,ymax,zmax);
  }
  if (!(bits & 0x080))
  {
    edges[i  ].first.set( xmin,ymax,zmax);
    edges[i++].second.set(xmin,ymin,zmax); 
  }
 
  if (!(bits & 0x100))
  { 
    edges[i  ].first.set( xmin,ymin,zmin); 
    edges[i++].second.set(xmin,ymin,zmax); 
  }
  if (!(bits & 0x200))
  {
    edges[i  ].first.set( xmax,ymin,zmin);
    edges[i++].second.set(xmax,ymin,zmax);
  }
  if (!(bits & 0x400))
  {
    edges[i  ].first.set( xmax,ymax,zmin);
    edges[i++].second.set(xmax,ymax,zmax);
  }
  if (!(bits & 0x800))
  {
    edges[i  ].first.set( xmin,ymax,zmin);
    edges[i++].second.set(xmin,ymax,zmax);
  }
  edges.resize(i);
 
  // Clip the edges by the planes
  //
  for (auto iedge = edges.cbegin(); iedge != edges.cend(); ++iedge)
  {
    G4bool    exist = true;
    G4Point3D p1    = iedge->first;
    G4Point3D p2    = iedge->second;
    for (auto iplane = pPlanes.cbegin(); iplane != pPlanes.cend(); ++iplane)
    {
      // Clip current edge
      G4double d1 = iplane->distance(p1); 
      G4double d2 = iplane->distance(p2); 
      if (d1 > 0.0)
      {
        if (d2 > 0.0) { exist = false; break; } // go to next edge
        p1 = (p2*d1-p1*d2)/(d1-d2);                   // move p1
      }
      else
      {
        if (d2 > 0.0) { p2 = (p1*d2-p2*d1)/(d2-d1); } // move p2
      }
    }
    // Adjust the extent
    if (exist)
    {
      emin.setX(std::min(std::min(p1.x(),p2.x()),emin.x())); 
      emin.setY(std::min(std::min(p1.y(),p2.y()),emin.y())); 
      emin.setZ(std::min(std::min(p1.z(),p2.z()),emin.z())); 
 
      emax.setX(std::max(std::max(p1.x(),p2.x()),emax.x())); 
      emax.setY(std::max(std::max(p1.y(),p2.y()),emax.y())); 
      emax.setZ(std::max(std::max(p1.z(),p2.z()),emax.z())); 
    }
  }
 
  // Copy the extent back
  //
  pExtent.first  = emin;
  pExtent.second = emax;
}

References d1, d2, emax, G4VoxelLimits::GetMaxXExtent(), G4VoxelLimits::GetMaxYExtent(), G4VoxelLimits::GetMaxZExtent(), G4VoxelLimits::GetMinXExtent(), G4VoxelLimits::GetMinYExtent(), G4VoxelLimits::GetMinZExtent(), G4INCL::Math::max(), G4INCL::Math::min(), HepGeom::BasicVector3D< T >::setX(), HepGeom::BasicVector3D< T >::setY(), HepGeom::BasicVector3D< T >::setZ(), HepGeom::BasicVector3D< T >::x(), HepGeom::BasicVector3D< T >::y(), and HepGeom::BasicVector3D< T >::z().

Referenced by CalculateExtent().

◆ CreateListOfEdges()

void G4BoundingEnvelope::CreateListOfEdges	(	const G4Polygon3D &	baseA,
		const G4Polygon3D &	baseB,
		std::vector< G4Segment3D > &	pEdges
	)		const

private

Definition at line 685 of file G4BoundingEnvelope.cc.

{
  G4int na = baseA.size();
  G4int nb = baseB.size();
  pEdges.clear();
  if (na == nb)
  {
    pEdges.resize(3*na);
    G4int k = na - 1;
    for (G4int i=0; i<na; ++i)
    {
      pEdges.push_back(G4Segment3D(baseA[i],baseB[i]));
      pEdges.push_back(G4Segment3D(baseA[i],baseA[k]));
      pEdges.push_back(G4Segment3D(baseB[i],baseB[k]));
      k = i;
    }
  }
  else if (nb == 1)
  {
    pEdges.resize(2*na);
    G4int k = na - 1;
    for (G4int i=0; i<na; ++i)
    {
      pEdges.push_back(G4Segment3D(baseA[i],baseA[k]));
      pEdges.push_back(G4Segment3D(baseA[i],baseB[0]));
      k = i;
    }
  }
  else if (na == 1)
  {
    pEdges.resize(2*nb);
    G4int k = nb - 1;
    for (G4int i=0; i<nb; ++i)
    {
      pEdges.push_back(G4Segment3D(baseB[i],baseB[k]));
      pEdges.push_back(G4Segment3D(baseB[i],baseA[0]));
      k = i;
    }
  }
}

Referenced by CalculateExtent().

◆ CreateListOfPlanes()

void G4BoundingEnvelope::CreateListOfPlanes	(	const G4Polygon3D &	baseA,
		const G4Polygon3D &	baseB,
		std::vector< G4Plane3D > &	pPlanes
	)		const

private

Definition at line 733 of file G4BoundingEnvelope.cc.

{
  // Find centers of the bases and internal point of the prism
  //
  G4int na = baseA.size();
  G4int nb = baseB.size();
  G4Point3D pa(0.,0.,0.), pb(0.,0.,0.), p0;
  G4Normal3D norm;
  for (G4int i=0; i<na; ++i) pa += baseA[i];
  for (G4int i=0; i<nb; ++i) pb += baseB[i];
  pa /= na; pb /= nb; p0 = (pa+pb)/2.;
 
  // Create list of planes
  //
  pPlanes.clear();
  if (na == nb)  // bases with equal number of vertices
  {
    G4int k = na - 1;
    for (G4int i=0; i<na; ++i)
    {
      norm = (baseB[k]-baseA[i]).cross(baseA[k]-baseB[i]);
      if (norm.mag2() > kCarTolerance)
      {
        pPlanes.push_back(G4Plane3D(norm,baseA[i]));
      }
      k = i;
    }
    norm = (baseA[2]-baseA[0]).cross(baseA[1]-pa);
    if (norm.mag2() > kCarTolerance)
    {
      pPlanes.push_back(G4Plane3D(norm,pa));
    }
    norm = (baseB[2]-baseB[0]).cross(baseB[1]-pb);
    if (norm.mag2() > kCarTolerance)
    {
      pPlanes.push_back(G4Plane3D(norm,pb));
    }
  }
  else if (nb == 1) // baseB has one vertex
  {
    G4int k = na - 1;
    for (G4int i=0; i<na; ++i)
    {
      norm = (baseA[i]-baseB[0]).cross(baseA[k]-baseB[0]);
      if (norm.mag2() > kCarTolerance)
      {
        pPlanes.push_back(G4Plane3D(norm,baseB[0]));
      }
      k = i;
    }
    norm = (baseA[2]-baseA[0]).cross(baseA[1]-pa);
    if (norm.mag2() > kCarTolerance)
    {
      pPlanes.push_back(G4Plane3D(norm,pa));
    }
  }
  else if (na == 1) // baseA has one vertex
  {
    G4int k = nb - 1;
    for (G4int i=0; i<nb; ++i)
    {
      norm = (baseB[i]-baseA[0]).cross(baseB[k]-baseA[0]);
      if (norm.mag2() > kCarTolerance)
      {
        pPlanes.push_back(G4Plane3D(norm,baseA[0]));
      }
      k = i;
    }
    norm = (baseB[2]-baseB[0]).cross(baseB[1]-pb);
    if (norm.mag2() > kCarTolerance)
    {
      pPlanes.push_back(G4Plane3D(norm,pb));
    }
  }
 
  // Ensure that normals of the planes point to outside
  //
  G4int nplanes = pPlanes.size();
  for (G4int i=0; i<nplanes; ++i)
  {
    pPlanes[i].normalize();
    if (pPlanes[i].distance(p0) > 0)
    {
      pPlanes[i] = G4Plane3D(-pPlanes[i].a(),-pPlanes[i].b(),
                             -pPlanes[i].c(),-pPlanes[i].d());
    } 
  }
}

References kCarTolerance, and HepGeom::BasicVector3D< T >::mag2().

Referenced by CalculateExtent().

◆ FindScaleFactor()

G4double G4BoundingEnvelope::FindScaleFactor ( const G4Transform3D & pTransform3D ) const

private

Definition at line 559 of file G4BoundingEnvelope.cc.

{
  if (pTransform3D.xx() == 1. &&
      pTransform3D.yy() == 1. &&
      pTransform3D.zz() == 1.) return 1.;
 
  G4double xx = pTransform3D.xx();
  G4double yx = pTransform3D.yx();
  G4double zx = pTransform3D.zx();
  G4double sxsx = xx*xx + yx*yx + zx*zx;
  G4double xy = pTransform3D.xy();
  G4double yy = pTransform3D.yy();
  G4double zy = pTransform3D.zy();
  G4double sysy = xy*xy + yy*yy + zy*zy;
  G4double xz = pTransform3D.xz();
  G4double yz = pTransform3D.yz();
  G4double zz = pTransform3D.zz();
  G4double szsz = xz*xz + yz*yz + zz*zz;
  G4double ss = std::max(std::max(sxsx,sysy),szsz);
  return (ss <= 1.) ? 1. : std::sqrt(ss);
}

References G4INCL::Math::max(), HepGeom::Transform3D::xx(), HepGeom::Transform3D::xy(), HepGeom::Transform3D::xz(), HepGeom::Transform3D::yx(), HepGeom::Transform3D::yy(), HepGeom::Transform3D::yz(), HepGeom::Transform3D::zx(), HepGeom::Transform3D::zy(), and HepGeom::Transform3D::zz().

Referenced by BoundingBoxVsVoxelLimits(), and CalculateExtent().

◆ GetPrismAABB()

void G4BoundingEnvelope::GetPrismAABB	(	const G4Polygon3D &	pBaseA,
		const G4Polygon3D &	pBaseB,
		G4Segment3D &	pAABB
	)		const

private

Definition at line 637 of file G4BoundingEnvelope.cc.

{
  G4double xmin =  kInfinity, ymin =  kInfinity, zmin =  kInfinity;
  G4double xmax = -kInfinity, ymax = -kInfinity, zmax = -kInfinity;
 
  // First base
  //
  for (auto it1 = pBaseA.cbegin(); it1 != pBaseA.cend(); ++it1)
  { 
    G4double x = it1->x();
    if (x < xmin) xmin = x;
    if (x > xmax) xmax = x;
    G4double y = it1->y(); 
    if (y < ymin) ymin = y;
    if (y > ymax) ymax = y;
    G4double z = it1->z(); 
    if (z < zmin) zmin = z;
    if (z > zmax) zmax = z;
  }
 
  // Second base
  //
  for (auto it2 = pBaseB.cbegin(); it2 != pBaseB.cend(); ++it2)
  { 
    G4double x = it2->x(); 
    if (x < xmin) xmin = x;
    if (x > xmax) xmax = x;
    G4double y = it2->y(); 
    if (y < ymin) ymin = y;
    if (y > ymax) ymax = y;
    G4double z = it2->z(); 
    if (z < zmin) zmin = z;
    if (z > zmax) zmax = z;
  }
 
  // Set bounding box
  //
  pAABB.first  = G4Point3D(xmin,ymin,zmin);
  pAABB.second = G4Point3D(xmax,ymax,zmax);
}

References kInfinity.

Referenced by CalculateExtent().

◆ TransformVertices()

void G4BoundingEnvelope::TransformVertices	(	const G4Transform3D &	pTransform3D,
		std::vector< G4Polygon3D * > &	pBases
	)		const

private

Definition at line 586 of file G4BoundingEnvelope.cc.

{
  G4ThreeVectorList baseA(4), baseB(4);
  std::vector<const G4ThreeVectorList*> aabb(2);
  aabb[0] = &baseA; aabb[1] = &baseB;
  if (fPolygons == nullptr)
  {
    baseA[0].set(fMin.x(),fMin.y(),fMin.z());
    baseA[1].set(fMax.x(),fMin.y(),fMin.z());
    baseA[2].set(fMax.x(),fMax.y(),fMin.z());
    baseA[3].set(fMin.x(),fMax.y(),fMin.z());
    baseB[0].set(fMin.x(),fMin.y(),fMax.z());
    baseB[1].set(fMax.x(),fMin.y(),fMax.z());
    baseB[2].set(fMax.x(),fMax.y(),fMax.z());
    baseB[3].set(fMin.x(),fMax.y(),fMax.z());
  }
  std::vector<const G4ThreeVectorList*>::const_iterator ia, iaend;
  auto ib = pBases.begin();
  ia    = (fPolygons == nullptr) ? aabb.cbegin() : fPolygons->cbegin();
  iaend = (fPolygons == nullptr) ? aabb.cend()   : fPolygons->cend();
 
  if (pTransform3D.xx()==1 && pTransform3D.yy()==1 && pTransform3D.zz()==1)
  {
    G4ThreeVector offset = pTransform3D.getTranslation();
    for ( ; ia != iaend; ++ia, ++ib)
    { 
      auto ka = (*ia)->cbegin();
      auto kb = (*ib)->begin();
      for ( ; ka != (*ia)->cend(); ++ka, ++kb) { (*kb) = (*ka) + offset; }
    }
  }
  else
  {
    for ( ; ia != iaend; ++ia, ++ib)
    { 
      auto ka = (*ia)->cbegin();
      auto kb = (*ib)->begin();
      for ( ; ka != (*ia)->cend(); ++ka, ++kb)
      {
        (*kb) = pTransform3D*G4Point3D(*ka);
      }
    }
  }
}

References fMax, fMin, fPolygons, HepGeom::Transform3D::getTranslation(), CLHEP::Hep3Vector::x(), HepGeom::Transform3D::xx(), CLHEP::Hep3Vector::y(), HepGeom::Transform3D::yy(), CLHEP::Hep3Vector::z(), and HepGeom::Transform3D::zz().

Referenced by CalculateExtent().

Field Documentation

◆ fMax

G4ThreeVector G4BoundingEnvelope::fMax

private

Definition at line 139 of file G4BoundingEnvelope.hh.

Referenced by BoundingBoxVsVoxelLimits(), CalculateExtent(), CheckBoundingBox(), G4BoundingEnvelope(), and TransformVertices().

◆ fMin

G4ThreeVector G4BoundingEnvelope::fMin

private

Definition at line 139 of file G4BoundingEnvelope.hh.

Referenced by BoundingBoxVsVoxelLimits(), CalculateExtent(), CheckBoundingBox(), G4BoundingEnvelope(), and TransformVertices().

◆ fPolygons

const std::vector<const G4ThreeVectorList*>* G4BoundingEnvelope::fPolygons = nullptr

private

Definition at line 142 of file G4BoundingEnvelope.hh.

Referenced by CalculateExtent(), CheckBoundingPolygons(), G4BoundingEnvelope(), and TransformVertices().

The documentation for this class was generated from the following files:

source/geometry/management/include/G4BoundingEnvelope.hh
source/geometry/management/src/G4BoundingEnvelope.cc

Public Member Functions

Private Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

◆ G4BoundingEnvelope() [1/3]

◆ G4BoundingEnvelope() [2/3]

◆ G4BoundingEnvelope() [3/3]

◆ ~G4BoundingEnvelope()

Member Function Documentation

◆ BoundingBoxVsVoxelLimits()

◆ CalculateExtent()

◆ CheckBoundingBox()

◆ CheckBoundingPolygons()

◆ ClipEdgesByVoxel()

◆ ClipVoxelByPlanes()

◆ CreateListOfEdges()

◆ CreateListOfPlanes()

◆ FindScaleFactor()

◆ GetPrismAABB()

◆ TransformVertices()

Field Documentation

◆ fMax

◆ fMin

◆ fPolygons