G4GeomTools Class Reference

G4GeomTools is a collecting utilities which can be helpful for a wide range of geometry-related tasks. More...

#include <G4GeomTools.hh>

Static Public Member Functions
static G4double	TriangleArea (G4double Ax, G4double Ay, G4double Bx, G4double By, G4double Cx, G4double Cy)
static G4double	TriangleArea (const G4TwoVector &A, const G4TwoVector &B, const G4TwoVector &C)
static G4double	QuadArea (const G4TwoVector &A, const G4TwoVector &B, const G4TwoVector &C, const G4TwoVector &D)
static G4double	PolygonArea (const G4TwoVectorList &polygon)
static G4bool	PointInTriangle (G4double Px, G4double Py, G4double Ax, G4double Ay, G4double Bx, G4double By, G4double Cx, G4double Cy)
static G4bool	PointInTriangle (const G4TwoVector &P, const G4TwoVector &A, const G4TwoVector &B, const G4TwoVector &C)
static G4bool	PointInPolygon (const G4TwoVector &P, const G4TwoVectorList &Polygon)
static G4bool	IsConvex (const G4TwoVectorList &polygon)
static G4bool	TriangulatePolygon (const G4TwoVectorList &polygon, std::vector< G4int > &result)
static G4bool	TriangulatePolygon (const G4TwoVectorList &polygon, G4TwoVectorList &result)
static void	RemoveRedundantVertices (G4TwoVectorList &polygon, std::vector< G4int > &iout, G4double tolerance=0.0)
static G4bool	DiskExtent (G4double rmin, G4double rmax, G4double startPhi, G4double delPhi, G4TwoVector &pmin, G4TwoVector &pmax)
static void	DiskExtent (G4double rmin, G4double rmax, G4double sinPhiStart, G4double cosPhiStart, G4double sinPhiEnd, G4double cosPhiEnd, G4TwoVector &pmin, G4TwoVector &pmax)
static G4double	EllipsePerimeter (G4double a, G4double b)
static G4double	EllipticConeLateralArea (G4double a, G4double b, G4double h)
static G4ThreeVector	TriangleAreaNormal (const G4ThreeVector &A, const G4ThreeVector &B, const G4ThreeVector &C)
static G4ThreeVector	QuadAreaNormal (const G4ThreeVector &A, const G4ThreeVector &B, const G4ThreeVector &C, const G4ThreeVector &D)
static G4ThreeVector	PolygonAreaNormal (const G4ThreeVectorList &polygon)
static G4double	DistancePointSegment (const G4ThreeVector &P, const G4ThreeVector &A, const G4ThreeVector &B)
static G4ThreeVector	ClosestPointOnSegment (const G4ThreeVector &P, const G4ThreeVector &A, const G4ThreeVector &B)
static G4ThreeVector	ClosestPointOnTriangle (const G4ThreeVector &P, const G4ThreeVector &A, const G4ThreeVector &B, const G4ThreeVector &C)
static G4bool	SphereExtent (G4double rmin, G4double rmax, G4double startTheta, G4double delTheta, G4double startPhi, G4double delPhi, G4ThreeVector &pmin, G4ThreeVector &pmax)
static G4double	HypeStereo (G4double r0, G4double r, G4double h)
static void	TwistedTubeBoundingTrap (G4double twistAng, G4double endInnerRad, G4double endOuterRad, G4double dPhi, G4TwoVectorList &vertices)
static G4double	HyperboloidSurfaceArea (G4double dphi, G4double r0, G4double tanstereo, G4double zmin, G4double zmax)

Detailed Description

G4GeomTools is a collecting utilities which can be helpful for a wide range of geometry-related tasks.

Definition at line 50 of file G4GeomTools.hh.

Member Function Documentation

ClosestPointOnSegment()

G4ThreeVector G4GeomTools::ClosestPointOnSegment ( const G4ThreeVector & P, const G4ThreeVector & A, const G4ThreeVector & B )

static

Finds a point on a 3D line segment AB closest to point 'P'.

Definition at line 667 of file G4GeomTools.cc.

{
  G4ThreeVector AP = P - A;
  G4ThreeVector AB = B - A;
 
  G4double u = AP.dot(AB);
  if (u <= 0) { return A; }      // closest point is A
 
  G4double len2 = AB.mag2();
  if (u >= len2) { return B; }   // closest point is B
 
  G4double t = u/len2;
  return A + t*AB;           // closest point on segment
}

ClosestPointOnTriangle()

G4ThreeVector G4GeomTools::ClosestPointOnTriangle ( const G4ThreeVector & P, const G4ThreeVector & A, const G4ThreeVector & B, const G4ThreeVector & C )

static

Finds a point on a 3D triangle ABC closest to point 'P'.

Definition at line 696 of file G4GeomTools.cc.

{
  G4ThreeVector diff  = A - P;
  G4ThreeVector edge0 = B - A;
  G4ThreeVector edge1 = C - A;
 
  G4double a = edge0.mag2();
  G4double b = edge0.dot(edge1);
  G4double c = edge1.mag2();
  G4double d = diff.dot(edge0);
  G4double e = diff.dot(edge1);
 
  G4double det = a*c - b*b;
  G4double t0  = b*e - c*d;
  G4double t1  = b*d - a*e;
 
  /*
             ^ t1
         \ 2 |
          \  |
           \ |     regions
            \|
             C
             |\
         3   | \   1
             |  \
             | 0 \
             |    \
        ---- A --- B ----> t0
             |      \
         4   |   5   \   6
             |        \
  */
 
  G4int region = -1;
  if (t0+t1 <= det)
  {
    region = (t0 < 0) ? ((t1 < 0) ? 4 : 3) : ((t1 < 0) ? 5 : 0);
  }
  else
  {
    region = (t0 < 0) ? 2 : ((t1 < 0) ? 6 : 1);
  }
 
  switch (region)
  {
  case 0:  // interior of triangle
    {
      G4double invDet = 1./det;
      return A + (t0*invDet)*edge0 + (t1*invDet)*edge1;
    }
  case 1:  // edge BC
    {
      G4double numer = c + e - b - d;
      if (numer <= 0) { return C; }
      G4double denom = a - 2*b + c;
      return (numer >= denom) ? B : C + (numer/denom)*(edge0-edge1);
    }
  case 2:  // edge AC or BC
    {
      G4double tmp0 = b + d;
      G4double tmp1 = c + e;
      if (tmp1 > tmp0)
      {
        G4double numer = tmp1 - tmp0;
        G4double denom = a - 2*b + c;
        return (numer >= denom) ? B : C + (numer/denom)*(edge0-edge1);
      }
      // same:  (e >= 0) ? A : ((-e >= c) ? C : A + (-e/c)*edge1)
      return (tmp1 <= 0) ? C : (( e >= 0) ? A : A + (-e/c)*edge1);
    }
  case 3:  // edge AC
    return (e >= 0) ? A : ((-e >= c) ? C : A + (-e/c)*edge1);
 
  case 4:  // edge AB or AC
    if (d < 0) { return (-d >= a) ? B : A + (-d/a)*edge0; }
    return (e >= 0) ? A : ((-e >= c) ? C : A + (-e/c)*edge1);
 
  case 5:  // edge AB
    return (d >= 0) ? A : ((-d >= a) ? B : A + (-d/a)*edge0);
 
  case 6:  // edge AB or BC
    {
      G4double tmp0 = b + e;
      G4double tmp1 = a + d;
      if (tmp1 > tmp0)
      {
        G4double numer = tmp1 - tmp0;
        G4double denom = a - 2*b + c;
        return (numer >= denom) ? C : B + (numer/denom)*(edge1-edge0);
      }
      // same:  (d >= 0) ? A : ((-d >= a) ? B : A + (-d/a)*edge0)
      return (tmp1 <= 0) ? B : (( d >= 0) ? A : A + (-d/a)*edge0);
    }
  default: // impossible case
    return {kInfinity,kInfinity,kInfinity};
  }
}

DiskExtent() [1/2]

void G4GeomTools::DiskExtent ( G4double rmin, G4double rmax, G4double sinPhiStart, G4double cosPhiStart, G4double sinPhiEnd, G4double cosPhiEnd, G4TwoVector & pmin, G4TwoVector & pmax )

static

Calculates the bounding rectangle of a disk sector. Faster version without check of parameters.

Definition at line 432 of file G4GeomTools.cc.

{
  static const G4double kCarTolerance =
         G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
 
  // check if 360 degrees
  //
  pmin.set(-rmax,-rmax);
  pmax.set( rmax, rmax);
 
  if (std::abs(sinEnd-sinStart) < kCarTolerance &&
      std::abs(cosEnd-cosStart) < kCarTolerance) { return; }
 
  // get start and end quadrants
  //
  //      1 | 0
  //     ---+---
  //      3 | 2
  //
  G4int icase = (cosEnd < 0) ? 1 : 0;
  if (sinEnd   < 0) { icase += 2; }
  if (cosStart < 0) { icase += 4; }
  if (sinStart < 0) { icase += 8; }
 
  switch (icase)
  {
  // start quadrant 0
  case  0:                                 // start->end : 0->0
    if (sinEnd < sinStart) { break; }
    pmin.set(rmin*cosEnd,rmin*sinStart);
    pmax.set(rmax*cosStart,rmax*sinEnd  );
    break;
  case  1:                                 // start->end : 0->1
    pmin.set(rmax*cosEnd,std::min(rmin*sinStart,rmin*sinEnd));
    pmax.set(rmax*cosStart,rmax  );
    break;
  case  2:                                 // start->end : 0->2
    pmin.set(-rmax,-rmax);
    pmax.set(std::max(rmax*cosStart,rmax*cosEnd),rmax);
    break;
  case  3:                                 // start->end : 0->3
    pmin.set(-rmax,rmax*sinEnd);
    pmax.set(rmax*cosStart,rmax);
    break;
  // start quadrant 1
  case  4:                                 // start->end : 1->0
    pmin.set(-rmax,-rmax);
    pmax.set(rmax,std::max(rmax*sinStart,rmax*sinEnd));
    break;
  case  5:                                 // start->end : 1->1
    if (sinEnd > sinStart) { break; }
    pmin.set(rmax*cosEnd,rmin*sinEnd  );
    pmax.set(rmin*cosStart,rmax*sinStart);
    break;
  case  6:                                 // start->end : 1->2
    pmin.set(-rmax,-rmax);
    pmax.set(rmax*cosEnd,rmax*sinStart);
    break;
  case  7:                                 // start->end : 1->3
    pmin.set(-rmax,rmax*sinEnd);
    pmax.set(std::max(rmin*cosStart,rmin*cosEnd),rmax*sinStart);
    break;
  // start quadrant 2
  case  8:                                 // start->end : 2->0
    pmin.set(std::min(rmin*cosStart,rmin*cosEnd),rmax*sinStart);
    pmax.set(rmax,rmax*sinEnd);
    break;
  case  9:                                 // start->end : 2->1
    pmin.set(rmax*cosEnd,rmax*sinStart);
    pmax.set(rmax,rmax);
    break;
  case 10:                                 // start->end : 2->2
    if (sinEnd < sinStart) { break; }
    pmin.set(rmin*cosStart,rmax*sinStart);
    pmax.set(rmax*cosEnd,rmin*sinEnd  );
    break;
  case 11:                                 // start->end : 2->3
    pmin.set(-rmax,std::min(rmax*sinStart,rmax*sinEnd));
    pmax.set(rmax,rmax);
    break;
  // start quadrant 3
  case 12:                                 // start->end : 3->0
    pmin.set(rmax*cosStart,-rmax);
    pmax.set(rmax,rmax*sinEnd);
    break;
  case 13:                                 // start->end : 3->1
    pmin.set(std::min(rmax*cosStart,rmax*cosEnd),-rmax);
    pmax.set(rmax,rmax);
    break;
  case 14:                                 // start->end : 3->2
    pmin.set(rmax*cosStart,-rmax);
    pmax.set(rmax*cosEnd,std::max(rmin*sinStart,rmin*sinEnd));
    break;
  case 15:                                 // start->end : 3->3
    if (sinEnd > sinStart) { break; }
    pmin.set(rmax*cosStart,rmax*sinEnd);
    pmax.set(rmin*cosEnd,rmin*sinStart);
    break;
  }
  return;
}

DiskExtent() [2/2]

G4bool G4GeomTools::DiskExtent ( G4double rmin, G4double rmax, G4double startPhi, G4double delPhi, G4TwoVector & pmin, G4TwoVector & pmax )

static

Calculates the bounding rectangle of a disk sector. It returns false if the input parameters do not meet the following criteria: rmin >= 0 rmax > rmin + kCarTolerance delPhi > 0 + kCarTolerance.

Definition at line 399 of file G4GeomTools.cc.

{
  static const G4double kCarTolerance =
         G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
 
  // check parameters
  //
  pmin.set(0,0);
  pmax.set(0,0);
  if (rmin   <  0) { return false; }
  if (rmax   <= rmin + kCarTolerance) { return false; }
  if (delPhi <= 0    + kCarTolerance) { return false; }
 
  // calculate extent
  //
  pmin.set(-rmax,-rmax);
  pmax.set( rmax, rmax);
  if (delPhi >= CLHEP::twopi) { return true; }
 
  DiskExtent(rmin,rmax,
             std::sin(startPhi),std::cos(startPhi),
             std::sin(startPhi+delPhi),std::cos(startPhi+delPhi),
             pmin,pmax);
  return true;
}

Referenced by G4Cons::BoundingLimits(), G4CutTubs::BoundingLimits(), G4GenericPolycone::BoundingLimits(), G4Polycone::BoundingLimits(), G4Sphere::BoundingLimits(), G4Torus::BoundingLimits(), G4Tubs::BoundingLimits(), G4TwistedTubs::BoundingLimits(), G4Torus::CalculateExtent(), DiskExtent(), and SphereExtent().

DistancePointSegment()

G4double G4GeomTools::DistancePointSegment ( const G4ThreeVector & P, const G4ThreeVector & A, const G4ThreeVector & B )

static

Calculates the distance between a point 'P' and line segment AB in 3D.

Definition at line 646 of file G4GeomTools.cc.

{
  G4ThreeVector AP = P - A;
  G4ThreeVector AB = B - A;
 
  G4double u = AP.dot(AB);
  if (u <= 0) { return AP.mag(); }       // closest point is A
 
  G4double len2 = AB.mag2();
  if (u >= len2) { return (B-P).mag(); } // closest point is B
 
  return ((u/len2)*AB - AP).mag();   // distance to line
}

EllipsePerimeter()

G4double G4GeomTools::EllipsePerimeter ( G4double a, G4double b )

static

Computes the circumference (perimeter) of an ellipse.

Definition at line 541 of file G4GeomTools.cc.

{
  G4double x = std::abs(pA);
  G4double y = std::abs(pB);
  G4double a = std::max(x,y);
  G4double b = std::min(x,y);
  G4double e = std::sqrt((1. - b/a)*(1. + b/a));
  return 4. * a * comp_ellint_2(e);
}

EllipticConeLateralArea()

G4double G4GeomTools::EllipticConeLateralArea ( G4double a, G4double b, G4double h )

static

Computes the lateral surface area of an elliptic cone.

Definition at line 555 of file G4GeomTools.cc.

{
  G4double x = std::abs(pA);
  G4double y = std::abs(pB);
  G4double h = std::abs(pH);
  G4double a = std::max(x,y);
  G4double b = std::min(x,y);
  G4double e = std::sqrt((1. - b/a)*(1. + b/a)) / std::hypot(1.,b/h);
  return 2. * a * std::hypot(b,h) * comp_ellint_2(e);
}

Referenced by G4EllipticalCone::GetSurfaceArea().

HyperboloidSurfaceArea()

G4double G4GeomTools::HyperboloidSurfaceArea ( G4double dphi, G4double r0, G4double tanstereo, G4double zmin, G4double zmax )

static

Calculate surface area of the hyperboloid between 'zmin' and 'zmax'.

Parameters

[in]	dphi	Delta phi.
[in]	r0	The radius at z = 0.
[in]	tanstereo	The tangent of the stereo angle.
[in]	zmin	Minimum Z.
[in]	zmax	Maximum Z.

Definition at line 913 of file G4GeomTools.cc.

{
  static const G4double kCarTolerance =
         G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
 
  G4double a = std::abs(r0); // radius at z = 0
  G4double t = std::abs(tanstereo); // tan(stereo)
  G4double phi = std::abs(dphi); // delta phi
 
  // Check spesial cases: cylindrical and conical surfaces
  if (t < kCarTolerance) { return a*std::abs(zmax - zmin)*phi; } // cylinder
  G4double rmin = std::hypot(t*zmin, a); // radius at zmin
  G4double rmax = std::hypot(t*zmax, a); // radius at zmax
  if (a <  kCarTolerance) // cone
  {
    G4double smin = rmin*std::hypot(rmin, zmin);
    G4double smax = rmax*std::hypot(rmax, zmax);
    return (zmin*zmax < 0.) ? (smin + smax)*phi/2. : std::abs(smax - smin)*phi/2.;
  }
  // Find surface area
  G4double tt = t*t;
  G4double aa = a*a;
  G4double cc = aa/tt;
  G4double k = std::sqrt(aa + cc)/cc;
 
  G4double hmin = std::abs(zmin);
  G4double smin = a*(hmin*std::hypot(1., k*hmin) + std::asinh(k*hmin)/k);
  if (zmax == -zmin) { return smin*phi; }
 
  G4double hmax = std::abs(zmax);
  G4double smax = a*(hmax*std::hypot(1., k*hmax) + std::asinh(k*hmax)/k);
  return (zmin*zmax < 0.) ? (smin + smax)*phi/2. :std::abs(smax - smin)*phi/2.;
}

Referenced by G4Hype::SetInnerRadius(), G4Hype::SetInnerStereo(), G4Hype::SetOuterRadius(), G4Hype::SetOuterStereo(), and G4Hype::SetZHalfLength().

HypeStereo()

G4double G4GeomTools::HypeStereo ( G4double r0, G4double r, G4double h )

static

Calculates the hyperbolic surface stereo. Stereo is a half angle at the intersection point of the two lines in the tangent plane cross-section.

Definition at line 861 of file G4GeomTools.cc.

{
  static const G4double kCarTolerance =
         G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
  if (std::abs(r - r0) < kCarTolerance) { return 0.; }
  return std::atan(std::sqrt((r - r0)*(r + r0))/std::abs(h));
}

IsConvex()

G4bool G4GeomTools::IsConvex ( const G4TwoVectorList & polygon )

static

Decides if a 2D 'polygon' is convex, i.e. if all internal angles are less than pi.

Definition at line 166 of file G4GeomTools.cc.

{
  static const G4double kCarTolerance =
         G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
 
  G4bool gotNegative = false;
  G4bool gotPositive = false;
  auto  n = (G4int)polygon.size();
  if (n <= 0) { return false; }
  for (G4int icur=0; icur<n; ++icur)
  {
    G4int iprev = (icur ==   0) ? n-1 : icur-1;
    G4int inext = (icur == n-1) ?   0 : icur+1;
    G4TwoVector e1 = polygon[icur]  - polygon[iprev];
    G4TwoVector e2 = polygon[inext] - polygon[icur];
    G4double cross = e1.x()*e2.y() - e1.y()*e2.x();
    if (std::abs(cross) < kCarTolerance) { return false; }
    if (cross <  0) { gotNegative = true; }
    if (cross >  0) { gotPositive = true; }
    if (gotNegative && gotPositive) { return false; }
  }
  return true;
}

Referenced by G4ExtrudedSolid::G4ExtrudedSolid(), and G4ExtrudedSolid::G4ExtrudedSolid().

PointInPolygon()

G4bool G4GeomTools::PointInPolygon ( const G4TwoVector & P, const G4TwoVectorList & Polygon )

static

Decides if a point P is inside the 'Polygon'.

Definition at line 146 of file G4GeomTools.cc.

{
  auto  Nv = (G4int)v.size();
  G4bool in = false;
  for (G4int i = 0, k = Nv - 1; i < Nv; k = i++)
  {
    if ((v[i].y() > p.y()) != (v[k].y() > p.y()))
    {
      G4double ctg = (v[k].x()-v[i].x())/(v[k].y()-v[i].y());
      in ^= static_cast<int>(p.x() < (p.y()-v[i].y())*ctg + v[i].x());
    }
  }
  return in;
}

PointInTriangle() [1/2]

G4bool G4GeomTools::PointInTriangle ( const G4TwoVector & P, const G4TwoVector & A, const G4TwoVector & B, const G4TwoVector & C )

static

Decides if a point P is inside the triangle ABC.

Definition at line 118 of file G4GeomTools.cc.

{
  G4double Ax = A.x(), Ay = A.y();
  G4double Bx = B.x(), By = B.y();
  G4double Cx = C.x(), Cy = C.y();
  G4double Px = P.x(), Py = P.y();
  if ((Bx-Ax)*(Cy-Ay) - (By-Ay)*(Cx-Ax) > 0.)
  {
    if ((Ax-Cx)*(Py-Cy) - (Ay-Cy)*(Px-Cx) < 0.) { return false; }
    if ((Bx-Ax)*(Py-Ay) - (By-Ay)*(Px-Ax) < 0.) { return false; }
    if ((Cx-Bx)*(Py-By) - (Cy-By)*(Px-Bx) < 0.) { return false; }
  }
  else
  {
    if ((Ax-Cx)*(Py-Cy) - (Ay-Cy)*(Px-Cx) > 0.) { return false; }
    if ((Bx-Ax)*(Py-Ay) - (By-Ay)*(Px-Ax) > 0.) { return false; }
    if ((Cx-Bx)*(Py-By) - (Cy-By)*(Px-Bx) > 0.) { return false; }
  }
  return true;
}

PointInTriangle() [2/2]

G4bool G4GeomTools::PointInTriangle ( G4double Px, G4double Py, G4double Ax, G4double Ay, G4double Bx, G4double By, G4double Cx, G4double Cy )

static

Decides if a point (Px,Py) is inside the triangle (Ax,Ay)(Bx,By)(Cx,Cy).

Definition at line 93 of file G4GeomTools.cc.

{
  if ((Bx-Ax)*(Cy-Ay) - (By-Ay)*(Cx-Ax) > 0.)
  {
    if ((Ax-Cx)*(Py-Cy) - (Ay-Cy)*(Px-Cx) < 0.) { return false; }
    if ((Bx-Ax)*(Py-Ay) - (By-Ay)*(Px-Ax) < 0.) { return false; }
    if ((Cx-Bx)*(Py-By) - (Cy-By)*(Px-Bx) < 0.) { return false; }
  }
  else
  {
    if ((Ax-Cx)*(Py-Cy) - (Ay-Cy)*(Px-Cx) > 0.) { return false; }
    if ((Bx-Ax)*(Py-Ay) - (By-Ay)*(Px-Ax) > 0.) { return false; }
    if ((Cx-Bx)*(Py-By) - (Cy-By)*(Px-Bx) > 0.) { return false; }
  }
  return true;
}

PolygonArea()

G4double G4GeomTools::PolygonArea ( const G4TwoVectorList & polygon )

static

Calculates the area of a 2D polygon, returned value is positive if the vertices of the polygon are defined in anticlockwise order, otherwise it is negative.

Definition at line 76 of file G4GeomTools.cc.

{
  auto  n = (G4int)p.size();
  if (n < 3) { return 0.0; } // degenerate polygon
 
  G4double area = p[n-1].x()*p[0].y() - p[0].x()*p[n-1].y();
  for(G4int i=1; i<n; ++i)
  {
    area += p[i-1].x()*p[i].y() - p[i].x()*p[i-1].y();
  }
  return area*0.5;
}

Referenced by G4GenericPolycone::CalculateExtent(), G4Polycone::CalculateExtent(), G4Polyhedra::CalculateExtent(), G4ExtrudedSolid::G4ExtrudedSolid(), G4ExtrudedSolid::G4ExtrudedSolid(), and TriangulatePolygon().

PolygonAreaNormal()

G4ThreeVector G4GeomTools::PolygonAreaNormal ( const G4ThreeVectorList & polygon )

static

Finds the normal to the plane of a 3D polygon; the length of the normal is equal to the area of the polygon.

Definition at line 630 of file G4GeomTools.cc.

{
  auto  n = (G4int)p.size();
  if (n < 3) { return {0,0,0}; } // degerate polygon
  G4ThreeVector normal = p[n-1].cross(p[0]);
  for(G4int i=1; i<n; ++i)
  {
    normal += p[i-1].cross(p[i]);
  }
  return normal*0.5;
}

QuadArea()

G4double G4GeomTools::QuadArea ( const G4TwoVector & A, const G4TwoVector & B, const G4TwoVector & C, const G4TwoVector & D )

static

Calculates the area of a 2D quadrilateral, returned value is positive if the vertices of the quadrilateral are given in anticlockwise order, otherwise it is negative.

Definition at line 64 of file G4GeomTools.cc.

{
  return ((C.x()-A.x())*(D.y()-B.y()) - (C.y()-A.y())*(D.x()-B.x()))*0.5;
}

QuadAreaNormal()

G4ThreeVector G4GeomTools::QuadAreaNormal ( const G4ThreeVector & A, const G4ThreeVector & B, const G4ThreeVector & C, const G4ThreeVector & D )

static

Finds the normal to the plane of a 3D quadrilateral ABCD; the length of the normal is equal to the area of the quadrilateral.

Definition at line 618 of file G4GeomTools.cc.

{
  return ((C-A).cross(D-B))*0.5;
}

Referenced by G4Polyhedra::GetSurfaceArea(), G4Trap::GetSurfaceArea(), and G4Trap::SetCachedValues().

RemoveRedundantVertices()

void G4GeomTools::RemoveRedundantVertices ( G4TwoVectorList & polygon, std::vector< G4int > & iout, G4double tolerance = 0.0 )

static

Removes collinear and coincident points from a 2D 'polygon'. Indices of removed points are available in 'iout'. Allows to specify a 'tolerance'.

Definition at line 310 of file G4GeomTools.cc.

{
  iout.resize(0);
  // set tolerance squared
  G4double delta = sqr(tolerance);
  // set special value to mark vertices for removal
  G4double removeIt = kInfinity;
 
  auto  nv = (G4int)polygon.size();
 
  // Main loop: check every three consecutive points, if the points
  // are collinear then mark middle point for removal
  //
  G4int icur = 0, iprev = 0, inext = 0, nout = 0;
  for (G4int i=0; i<nv; ++i)
  {
    icur = i;                    // index of current point
 
    for (G4int k=1; k<nv+1; ++k) // set index of previous point
    {
      iprev = icur - k;
      if (iprev < 0) { iprev += nv; }
      if (polygon[iprev].x() != removeIt) { break; }
    }
 
    for (G4int k=1; k<nv+1; ++k) // set index of next point
    {
      inext = icur + k;
      if (inext >= nv) { inext -= nv; }
      if (polygon[inext].x() != removeIt) { break; }
    }
 
    if (iprev == inext) { break; }   // degenerate polygon, stop
 
    // Calculate parameters of triangle (iprev->icur->inext),
    // if triangle is too small or too narrow then mark current
    // point for removal
    G4TwoVector e1 = polygon[iprev] - polygon[icur];
    G4TwoVector e2 = polygon[inext] - polygon[icur];
 
    // Check length of edges, then check height of the triangle
    G4double leng1 = e1.mag2();
    G4double leng2 = e2.mag2();
    G4double leng3 = (e2-e1).mag2();
    if (leng1 <= delta || leng2 <= delta || leng3 <= delta)
    {
      polygon[icur].setX(removeIt); ++nout;
    }
    else
    {
      G4double lmax = std::max(std::max(leng1,leng2),leng3);
      G4double area = std::abs(e1.x()*e2.y()-e1.y()*e2.x())*0.5;
      if (area/std::sqrt(lmax) <= std::abs(tolerance))
      {
        polygon[icur].setX(removeIt); ++nout;
      }
    }
  }
 
  // Remove marked points
  //
  icur = 0;
  if (nv - nout < 3)           // degenerate polygon, remove all points
  {
    for (G4int i=0; i<nv; ++i) { iout.push_back(i); }
    polygon.resize(0);
    nv = 0;
  }
  for (G4int i=0; i<nv; ++i) // move points, if required
  {
    if (polygon[i].x() != removeIt)
    {
      polygon[icur++] = polygon[i];
    }
    else
    {
      iout.push_back(i);
    }
  }
  if (icur < nv) { polygon.resize(icur); }
  return;
}

Referenced by G4Polycone::CalculateExtent(), G4Polyhedra::CalculateExtent(), G4ExtrudedSolid::G4ExtrudedSolid(), and G4ExtrudedSolid::G4ExtrudedSolid().

SphereExtent()

G4bool G4GeomTools::SphereExtent ( G4double rmin, G4double rmax, G4double startTheta, G4double delTheta, G4double startPhi, G4double delPhi, G4ThreeVector & pmin, G4ThreeVector & pmax )

static

Calculates the bounding box of a spherical sector,

Returns

false if input parameters do not meet the following criteria: rmin >= 0 rmax > rmin + kCarTolerance startTheta >= 0 && <= pi; delTheta > 0 + kCarTolerance delPhi > 0 + kCarTolerance.

Definition at line 803 of file G4GeomTools.cc.

{
  static const G4double kCarTolerance =
         G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
 
  // check parameters
  //
  pmin.set(0,0,0);
  pmax.set(0,0,0);
  if (rmin     <  0) { return false; }
  if (rmax     <= rmin + kCarTolerance) { return false; }
  if (delTheta <= 0    + kCarTolerance) { return false; }
  if (delPhi   <= 0    + kCarTolerance) { return false; }
 
  G4double stheta = startTheta;
  G4double dtheta = delTheta;
  if (stheta < 0 && stheta > CLHEP::pi) { return false; }
  if (stheta + dtheta > CLHEP::pi) { dtheta = CLHEP::pi - stheta; }
  if (dtheta <= 0 + kCarTolerance) { return false; }
 
  // calculate extent
  //
  pmin.set(-rmax,-rmax,-rmax);
  pmax.set( rmax, rmax, rmax);
  if (dtheta >= CLHEP::pi && delPhi >= CLHEP::twopi) { return true; }
 
  G4double etheta   = stheta + dtheta;
  G4double sinStart = std::sin(stheta);
  G4double cosStart = std::cos(stheta);
  G4double sinEnd   = std::sin(etheta);
  G4double cosEnd   = std::cos(etheta);
 
  G4double rhomin = rmin*std::min(sinStart,sinEnd);
  G4double rhomax = rmax;
  if (stheta > CLHEP::halfpi) { rhomax = rmax*sinStart; }
  if (etheta < CLHEP::halfpi) { rhomax = rmax*sinEnd; }
 
  G4TwoVector xymin,xymax;
  DiskExtent(rhomin,rhomax,
             std::sin(startPhi),std::cos(startPhi),
             std::sin(startPhi+delPhi),std::cos(startPhi+delPhi),
             xymin,xymax);
 
  G4double zmin = std::min(rmin*cosEnd,rmax*cosEnd);
  G4double zmax = std::max(rmin*cosStart,rmax*cosStart);
  pmin.set(xymin.x(),xymin.y(),zmin);
  pmax.set(xymax.x(),xymax.y(),zmax);
  return true;
}

TriangleArea() [1/2]

G4double G4GeomTools::TriangleArea ( const G4TwoVector & A, const G4TwoVector & B, const G4TwoVector & C )

static

Definition at line 52 of file G4GeomTools.cc.

{
  G4double Ax = A.x(), Ay = A.y();
  return ((B.x()-Ax)*(C.y()-Ay) - (B.y()-Ay)*(C.x()-Ax))*0.5;
}

TriangleArea() [2/2]

G4double G4GeomTools::TriangleArea ( G4double Ax, G4double Ay, G4double Bx, G4double By, G4double Cx, G4double Cy )

static

Functions to calculate the area of 2D triangle, returned value is positive if the vertices of the triangle are given in anticlockwise order, otherwise it is negative.

Definition at line 41 of file G4GeomTools.cc.

{
  return ((Bx-Ax)*(Cy-Ay) - (By-Ay)*(Cx-Ax))*0.5;
}

TriangleAreaNormal()

G4ThreeVector G4GeomTools::TriangleAreaNormal ( const G4ThreeVector & A, const G4ThreeVector & B, const G4ThreeVector & C )

static

Finds the normal to the plane of a 3D triangle ABC; the length of the normal is equal to the area of the triangle.

Definition at line 607 of file G4GeomTools.cc.

{
  return ((B-A).cross(C-A))*0.5;
}

Referenced by G4Trap::GetPointOnSurface().

TriangulatePolygon() [1/2]

G4bool G4GeomTools::TriangulatePolygon ( const G4TwoVectorList & polygon, G4TwoVectorList & result )

static

Same using the function above and returning as 'result' a list of triangles.

Definition at line 194 of file G4GeomTools.cc.

{
  result.resize(0);
  std::vector<G4int> triangles;
  G4bool reply = TriangulatePolygon(polygon,triangles);
 
  auto  n = (G4int)triangles.size();
  for (G4int i=0; i<n; ++i) { result.push_back(polygon[triangles[i]]); }
  return reply;
}

TriangulatePolygon() [2/2]

G4bool G4GeomTools::TriangulatePolygon ( const G4TwoVectorList & polygon, std::vector< G4int > & result )

static

Simple implementation of "ear clipping" algorithm for triangulation of a simple contour/polygon, it places results in a std::vector as triplets of vertices. If triangulation is successful the function returns true, otherwise false.

Definition at line 210 of file G4GeomTools.cc.

{
  result.resize(0);
 
  // allocate and initialize list of Vertices in polygon
  //
  auto  n = (G4int)polygon.size();
  if (n < 3) { return false; }
 
  // we want a counter-clockwise polygon in V
  //
  G4double area = G4GeomTools::PolygonArea(polygon);
  auto  V = new G4int[n];
  if (area > 0.)
  {
    for (G4int i=0; i<n; ++i) { V[i] = i; }
  }
  else
  {
    for (G4int i=0; i<n; ++i) { V[i] = (n-1)-i; }
  }
 
  //  Triangulation: remove nv-2 Vertices, creating 1 triangle every time
  //
  G4int nv = n;
  G4int count = 2*nv; // error detection counter
  for(G4int b=nv-1; nv>2; )
  {
    // ERROR: if we loop, it is probably a non-simple polygon
    if ((count--) <= 0)
    {
      delete [] V;
      if (area < 0.) { std::reverse(result.begin(),result.end()); }
      return false;
    }
 
    // three consecutive vertices in current polygon, <a,b,c>
    G4int a = (b   < nv) ? b   : 0; // previous
          b = (a+1 < nv) ? a+1 : 0; // current
    G4int c = (b+1 < nv) ? b+1 : 0; // next
 
    if (CheckSnip(polygon, a,b,c, nv,V))
    {
      // output Triangle
      result.push_back(V[a]);
      result.push_back(V[b]);
      result.push_back(V[c]);
 
      // remove vertex b from remaining polygon
      nv--;
      for(G4int i=b; i<nv; ++i) { V[i] = V[i+1]; }
 
      count = 2*nv; // resest error detection counter
    }
  }
  delete [] V;
  if (area < 0.) { std::reverse(result.begin(),result.end()); }
  return true;
}

Referenced by G4ExtrudedSolid::CalculateExtent(), G4GenericPolycone::CalculateExtent(), G4Polycone::CalculateExtent(), G4Polyhedra::CalculateExtent(), and TriangulatePolygon().

TwistedTubeBoundingTrap()

void G4GeomTools::TwistedTubeBoundingTrap ( G4double twistAng, G4double endInnerRad, G4double endOuterRad, G4double dPhi, G4TwoVectorList & vertices )

static

Finds the XY-coordinates of the corners of a generic trap that bounds specified twisted tube.

Parameters

[in]	twistAng	The twist angle.
[in]	endInnerRad	The inner radius at z = halfZ.
[in]	endOuterRad	The outer radius at z = halfZ.
[in]	dPhi	Delta phi.
[in]	vertices	The corners of the generic trap.

Definition at line 875 of file G4GeomTools.cc.

{
  vertices.resize(8);
  G4double rmin = std::abs(endInnerRad);
  G4double rmax = std::abs(endOuterRad);
 
  // Set untwisted vertices
  G4double phi = dPhi/2.;
  G4double sinphi = std::sin(phi);
  G4double cosphi = std::cos(phi);
  G4double tanphi = std::tan(phi);
  vertices[0].set(rmin*cosphi, rmin*sinphi);
  vertices[1].set(rmax, rmax*tanphi);
  vertices[2].set(rmax,-rmax*tanphi);
  vertices[3].set(rmin*cosphi,-rmin*sinphi);
  vertices[4] = vertices[0];
  vertices[5] = vertices[1];
  vertices[6] = vertices[2];
  vertices[7] = vertices[3];
 
  // Twist vertices
  G4double ang = twistAng/2.;
  for(auto i = 0; i < 4; ++i)
  {
    vertices[i].rotate(-ang);    // vertices at -halfz
    vertices[i + 4].rotate(ang); // vertices at +halfz
  }
}

---

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

• G4GeomTools.hh
• G4GeomTools.cc