#include <G4TriangularFacet.hh>

Inheritance diagram for G4TriangularFacet:

Public Member Functions
G4int	AllocatedMemory ()

void	ApplyTranslation (const G4ThreeVector v)

G4ThreeVector	Distance (const G4ThreeVector &p)

G4double	Distance (const G4ThreeVector &p, G4double minDist)

G4double	Distance (const G4ThreeVector &p, G4double minDist, const G4bool outgoing)

G4double	Extent (const G4ThreeVector axis)

	G4TriangularFacet ()

	G4TriangularFacet (const G4ThreeVector &vt0, const G4ThreeVector &vt1, const G4ThreeVector &vt2, G4FacetVertexType)

	G4TriangularFacet (const G4TriangularFacet &right)

	G4TriangularFacet (G4TriangularFacet &&right)

G4double	GetArea () const

G4ThreeVector	GetCircumcentre () const

G4VFacet *	GetClone ()

G4GeometryType	GetEntityType () const

G4TriangularFacet *	GetFlippedFacet ()

G4int	GetNumberOfVertices () const

G4ThreeVector	GetPointOnFace () const

G4double	GetRadius () const

G4ThreeVector	GetSurfaceNormal () const

G4ThreeVector	GetVertex (G4int i) const

G4int	GetVertexIndex (G4int i) const

G4bool	Intersect (const G4ThreeVector &p, const G4ThreeVector &v, const G4bool outgoing, G4double &distance, G4double &distFromSurface, G4ThreeVector &normal)

G4bool	IsDefined () const

G4bool	IsInside (const G4ThreeVector &p) const

G4TriangularFacet &	operator= (const G4TriangularFacet &right)

G4TriangularFacet &	operator= (G4TriangularFacet &&right)

G4bool	operator== (const G4VFacet &right) const

void	SetSurfaceNormal (G4ThreeVector normal)

void	SetVertex (G4int i, const G4ThreeVector &val)

void	SetVertexIndex (G4int i, G4int j)

void	SetVertices (std::vector< G4ThreeVector > *v)

std::ostream &	StreamInfo (std::ostream &os) const

	~G4TriangularFacet ()

Protected Attributes
G4double	kCarTolerance

Static Protected Attributes
static const G4double	dirTolerance = 1.0E-14

Private Member Functions
void	CopyFrom (const G4TriangularFacet &rhs)

void	MoveFrom (G4TriangularFacet &rhs)

Private Attributes
G4double	fA

G4double	fArea = 0.0

G4double	fB

G4double	fC

G4ThreeVector	fCircumcentre

G4double	fDet

G4ThreeVector	fE1

G4ThreeVector	fE2

std::array< G4int, 3 >	fIndices

G4bool	fIsDefined = false

G4double	fRadius = 0.0

G4double	fSqrDist = 0.0

G4ThreeVector	fSurfaceNormal

std::vector< G4ThreeVector > *	fVertices = nullptr

Detailed Description

Definition at line 60 of file G4TriangularFacet.hh.

Constructor & Destructor Documentation

◆ G4TriangularFacet() [1/4]

G4TriangularFacet::G4TriangularFacet ( )

Definition at line 144 of file G4TriangularFacet.cc.

  : fSqrDist(0.)
{
  fVertices = new vector<G4ThreeVector>(3);
  G4ThreeVector zero(0,0,0);
  SetVertex(0, zero);
  SetVertex(1, zero);
  SetVertex(2, zero);
  for (G4int i = 0; i < 3; ++i) fIndices[i] = -1;
  fIsDefined = false;
  fSurfaceNormal.set(0,0,0);
  fA = fB = fC = 0;
  fE1 = zero;
  fE2 = zero;
  fDet = 0.0;
  fArea = fRadius = 0.0;
}

References fA, fArea, fB, fC, fDet, fE1, fE2, fIndices, fIsDefined, fRadius, fSurfaceNormal, fVertices, CLHEP::Hep3Vector::set(), SetVertex(), and anonymous_namespace{G4CascadeDeexciteBase.cc}::zero.

Referenced by GetClone(), and GetFlippedFacet().

◆ ~G4TriangularFacet()

G4TriangularFacet::~G4TriangularFacet ( )

Definition at line 164 of file G4TriangularFacet.cc.

{
  SetVertices(nullptr);
}

References SetVertices().

◆ G4TriangularFacet() [2/4]

G4TriangularFacet::G4TriangularFacet	(	const G4ThreeVector &	vt0,
		const G4ThreeVector &	vt1,
		const G4ThreeVector &	vt2,
		G4FacetVertexType	vertexType
	)

Definition at line 56 of file G4TriangularFacet.cc.

  : G4VFacet()
{
  fVertices = new vector<G4ThreeVector>(3);
 
  SetVertex(0, vt0);
  if (vertexType == ABSOLUTE)
  {
    SetVertex(1, vt1);
    SetVertex(2, vt2);
    fE1 = vt1 - vt0;
    fE2 = vt2 - vt0;
  }
  else
  {
    SetVertex(1, vt0 + vt1);
    SetVertex(2, vt0 + vt2);
    fE1 = vt1;
    fE2 = vt2;
  }
 
  G4ThreeVector E1xE2 = fE1.cross(fE2);
  fArea = 0.5 * E1xE2.mag();
  for (G4int i = 0; i < 3; ++i) fIndices[i] = -1;
 
  fIsDefined = true;
  G4double delta = kCarTolerance; // Set tolerance for checking
 
  // Check length of edges
  //
  G4double leng1 = fE1.mag();
  G4double leng2 = (fE2-fE1).mag();
  G4double leng3 = fE2.mag();
  if (leng1 <= delta || leng2 <= delta || leng3 <= delta) 
  {
    fIsDefined = false;
  }
 
  // Check min height of triangle
  //
  if (fIsDefined)
  {
    if (2.*fArea/std::max(std::max(leng1,leng2),leng3) <= delta)
    {
      fIsDefined = false;
    } 
  }
 
  // Define facet
  //
  if (!fIsDefined)
  {
    ostringstream message;
    message << "Facet is too small or too narrow." << G4endl
            << "Triangle area = " << fArea << G4endl
            << "P0 = " << GetVertex(0) << G4endl
            << "P1 = " << GetVertex(1) << G4endl
            << "P2 = " << GetVertex(2) << G4endl
            << "Side1 length (P0->P1) = " << leng1 << G4endl
            << "Side2 length (P1->P2) = " << leng2 << G4endl
            << "Side3 length (P2->P0) = " << leng3;
    G4Exception("G4TriangularFacet::G4TriangularFacet()",
        "GeomSolids1001", JustWarning, message);
    fSurfaceNormal.set(0,0,0);
    fA = fB = fC = 0.0;
    fDet = 0.0;
    fCircumcentre = vt0 + 0.5*fE1 + 0.5*fE2;
    fArea = fRadius = 0.0;
  }
  else
  { 
    fSurfaceNormal = E1xE2.unit();
    fA   = fE1.mag2();
    fB   = fE1.dot(fE2);
    fC   = fE2.mag2();
    fDet = std::fabs(fA*fC - fB*fB);
 
    fCircumcentre = 
      vt0 + (E1xE2.cross(fE1)*fC + fE2.cross(E1xE2)*fA) / (2.*E1xE2.mag2());
    fRadius = (fCircumcentre - vt0).mag();
  }
}

References ABSOLUTE, CLHEP::Hep3Vector::cross(), CLHEP::Hep3Vector::dot(), fA, fArea, fB, fC, fCircumcentre, fDet, fE1, fE2, fIndices, fIsDefined, fRadius, fSurfaceNormal, fVertices, G4endl, G4Exception(), GetVertex(), JustWarning, G4VFacet::kCarTolerance, CLHEP::Hep3Vector::mag(), CLHEP::Hep3Vector::mag2(), G4INCL::Math::max(), CLHEP::Hep3Vector::set(), SetVertex(), and CLHEP::Hep3Vector::unit().

◆ G4TriangularFacet() [3/4]

G4TriangularFacet::G4TriangularFacet ( const G4TriangularFacet & right )

Definition at line 203 of file G4TriangularFacet.cc.

  : G4VFacet(rhs)
{
  CopyFrom(rhs);
}

References CopyFrom().

◆ G4TriangularFacet() [4/4]

G4TriangularFacet::G4TriangularFacet ( G4TriangularFacet && right )

Definition at line 211 of file G4TriangularFacet.cc.

  : G4VFacet(rhs)
{
  MoveFrom(rhs);
}

References MoveFrom().

Member Function Documentation

◆ AllocatedMemory()

G4int G4TriangularFacet::AllocatedMemory ( )

inlinevirtual

Implements G4VFacet.

Definition at line 162 of file G4TriangularFacet.hh.

{
  G4int size = sizeof(*this);
  size += GetNumberOfVertices() * sizeof(G4ThreeVector);
  return size;
}

References GetNumberOfVertices().

Referenced by G4QuadrangularFacet::AllocatedMemory().

◆ ApplyTranslation()

void G4VFacet::ApplyTranslation ( const G4ThreeVector v )

inherited

Definition at line 85 of file G4VFacet.cc.

{
  G4int n = GetNumberOfVertices();
  for (G4int i = 0; i < n; ++i)
  {
    SetVertex(i, GetVertex(i) + v);
  }
}

References G4VFacet::GetNumberOfVertices(), G4VFacet::GetVertex(), CLHEP::detail::n, and G4VFacet::SetVertex().

◆ CopyFrom()

void G4TriangularFacet::CopyFrom ( const G4TriangularFacet & rhs )

private

Definition at line 171 of file G4TriangularFacet.cc.

{
  char *p = (char *) &rhs;
  copy(p, p + sizeof(*this), (char *)this);
 
  if (fIndices[0] < 0 && fVertices == nullptr)
  {
    fVertices = new vector<G4ThreeVector>(3);
    for (G4int i = 0; i < 3; ++i) (*fVertices)[i] = (*rhs.fVertices)[i];
  }
}

References field_utils::copy(), fIndices, and fVertices.

Referenced by G4TriangularFacet(), and operator=().

◆ Distance() [1/3]

G4ThreeVector G4TriangularFacet::Distance ( const G4ThreeVector & p )

Definition at line 290 of file G4TriangularFacet.cc.

{
  G4ThreeVector D  = GetVertex(0) - p;
  G4double d = fE1.dot(D);
  G4double e = fE2.dot(D);
  G4double f = D.mag2();
  G4double q = fB*e - fC*d;
  G4double t = fB*d - fA*e;
  fSqrDist = 0.;
 
  if (q+t <= fDet)
  {
    if (q < 0.0)
    {
      if (t < 0.0)
      {
        //
        // We are in region 4.
        //
        if (d < 0.0)
        {
          t = 0.0;
          if (-d >= fA) {q = 1.0; fSqrDist = fA + 2.0*d + f;}
          else         {q = -d/fA; fSqrDist = d*q + f;}
        }
        else
        {
          q = 0.0;
          if       (e >= 0.0) {t = 0.0; fSqrDist = f;}
          else if (-e >= fC)   {t = 1.0; fSqrDist = fC + 2.0*e + f;}
          else                {t = -e/fC; fSqrDist = e*t + f;}
        }
      }
      else
      {
        //
        // We are in region 3.
        //
        q = 0.0;
        if      (e >= 0.0) {t = 0.0; fSqrDist = f;}
        else if (-e >= fC)  {t = 1.0; fSqrDist = fC + 2.0*e + f;}
        else               {t = -e/fC; fSqrDist = e*t + f;}
      }
    }
    else if (t < 0.0)
    {
      //
      // We are in region 5.
      //
      t = 0.0;
      if      (d >= 0.0) {q = 0.0; fSqrDist = f;}
      else if (-d >= fA)  {q = 1.0; fSqrDist = fA + 2.0*d + f;}
      else               {q = -d/fA; fSqrDist = d*q + f;}
    }
    else
    {
      //
      // We are in region 0.
      //
      G4double dist = fSurfaceNormal.dot(D);
      fSqrDist = dist*dist;
      return fSurfaceNormal*dist;
    }
  }
  else
  {
    if (q < 0.0)
    {
      //
      // We are in region 2.
      //
      G4double tmp0 = fB + d;
      G4double tmp1 = fC + e;
      if (tmp1 > tmp0)
      {
        G4double numer = tmp1 - tmp0;
        G4double denom = fA - 2.0*fB + fC;
        if (numer >= denom) {q = 1.0; t = 0.0; fSqrDist = fA + 2.0*d + f;}
        else
        {
          q       = numer/denom;
          t       = 1.0 - q;
          fSqrDist = q*(fA*q + fB*t +2.0*d) + t*(fB*q + fC*t + 2.0*e) + f;
        }
      }
      else
      {
        q = 0.0;
        if      (tmp1 <= 0.0) {t = 1.0; fSqrDist = fC + 2.0*e + f;}
        else if (e >= 0.0)    {t = 0.0; fSqrDist = f;}
        else                  {t = -e/fC; fSqrDist = e*t + f;}
      }
    }
    else if (t < 0.0)
    {
      //
      // We are in region 6.
      //
      G4double tmp0 = fB + e;
      G4double tmp1 = fA + d;
      if (tmp1 > tmp0)
      {
        G4double numer = tmp1 - tmp0;
        G4double denom = fA - 2.0*fB + fC;
        if (numer >= denom) {t = 1.0; q = 0.0; fSqrDist = fC + 2.0*e + f;}
        else
        {
          t       = numer/denom;
          q       = 1.0 - t;
          fSqrDist = q*(fA*q + fB*t +2.0*d) + t*(fB*q + fC*t + 2.0*e) + f;
        }
      }
      else
      {
        t = 0.0;
        if      (tmp1 <= 0.0) {q = 1.0; fSqrDist = fA + 2.0*d + f;}
        else if (d >= 0.0)    {q = 0.0; fSqrDist = f;}
        else                  {q = -d/fA; fSqrDist = d*q + f;}
      }
    }
    else
      //
      // We are in region 1.
      //
    {
      G4double numer = fC + e - fB - d;
      if (numer <= 0.0)
      {
        q       = 0.0;
        t       = 1.0;
        fSqrDist = fC + 2.0*e + f;
      }
      else
      {
        G4double denom = fA - 2.0*fB + fC;
        if (numer >= denom) {q = 1.0; t = 0.0; fSqrDist = fA + 2.0*d + f;}
        else
        {
          q       = numer/denom;
          t       = 1.0 - q;
          fSqrDist = q*(fA*q + fB*t + 2.0*d) + t*(fB*q + fC*t + 2.0*e) + f;
        }
      }
    }
  } 
  //
  //
  // Do fA check for rounding errors in the distance-squared.  It appears that
  // the conventional methods for calculating fSqrDist breaks down when very
  // near to or at the surface (as required by transport).
  // We'll therefore also use the magnitude-squared of the vector displacement.
  // (Note that I've also tried to get around this problem by using the
  // existing equations for
  //
  //    fSqrDist = function(fA,fB,fC,d,q,t)
  //
  // and use fA more accurate addition process which minimises errors and
  // breakdown of cummutitivity [where (A+B)+C != A+(B+C)] but this still
  // doesn't work.
  // Calculation from u = D + q*fE1 + t*fE2 is less efficient, but appears
  // more robust.
  //
  if (fSqrDist < 0.0) fSqrDist = 0.;
  G4ThreeVector u = D + q*fE1 + t*fE2;
  G4double u2 = u.mag2();
  //
  // The following (part of the roundoff error check) is from Oliver Merle'q
  // updates.
  //
  if (fSqrDist > u2) fSqrDist = u2;
 
  return u;
}

References D(), CLHEP::Hep3Vector::dot(), fA, fB, fC, fDet, fE1, fE2, fSqrDist, fSurfaceNormal, GetVertex(), and CLHEP::Hep3Vector::mag2().

Referenced by G4QuadrangularFacet::Distance(), Distance(), and Intersect().

◆ Distance() [2/3]

G4double G4TriangularFacet::Distance	(	const G4ThreeVector &	p,
		G4double	minDist
	)

virtual

Implements G4VFacet.

Definition at line 474 of file G4TriangularFacet.cc.

{
  //
  // Start with quicky test to determine if the surface of the sphere enclosing
  // the triangle is any closer to p than minDist.  If not, then don't bother
  // about more accurate test.
  //
  G4double dist = kInfinity;
  if ((p-fCircumcentre).mag()-fRadius < minDist)
  {
    //
    // It's possible that the triangle is closer than minDist,
    // so do more accurate assessment.
    //
    dist = Distance(p).mag();
  }
  return dist;
}

References Distance(), fCircumcentre, fRadius, kInfinity, and CLHEP::Hep3Vector::mag().

◆ Distance() [3/3]

G4double G4TriangularFacet::Distance	(	const G4ThreeVector &	p,
		G4double	minDist,
		const G4bool	outgoing
	)

virtual

Implements G4VFacet.

Definition at line 509 of file G4TriangularFacet.cc.

{
  //
  // Start with quicky test to determine if the surface of the sphere enclosing
  // the triangle is any closer to p than minDist.  If not, then don't bother
  // about more accurate test.
  //
  G4double dist = kInfinity;
  if ((p-fCircumcentre).mag()-fRadius < minDist)
  {
    //
    // It's possible that the triangle is closer than minDist,
    // so do more accurate assessment.
    //
    G4ThreeVector v  = Distance(p);
    G4double dist1 = sqrt(fSqrDist);
    G4double dir = v.dot(fSurfaceNormal);
    G4bool wrongSide = (dir > 0.0 && !outgoing) || (dir < 0.0 && outgoing);
    if (dist1 <= kCarTolerance)
    {
      //
      // Point p is very close to triangle.  Check if it's on the wrong side,
      // in which case return distance of 0.0 otherwise .
      //
      if (wrongSide) dist = 0.0;
      else dist = dist1;
    }
    else if (!wrongSide) dist = dist1;
  }
  return dist;
}

References Distance(), CLHEP::Hep3Vector::dot(), fCircumcentre, fRadius, fSqrDist, fSurfaceNormal, G4VFacet::kCarTolerance, and kInfinity.

◆ Extent()

G4double G4TriangularFacet::Extent ( const G4ThreeVector axis )

virtual

Implements G4VFacet.

Definition at line 550 of file G4TriangularFacet.cc.

{
  G4double ss = GetVertex(0).dot(axis);
  G4double sp = GetVertex(1).dot(axis);
  if (sp > ss) ss = sp;
  sp = GetVertex(2).dot(axis);
  if (sp > ss) ss = sp;
  return ss;
}

References CLHEP::Hep3Vector::dot(), GetVertex(), and G4InuclParticleNames::sp.

◆ GetArea()

G4double G4TriangularFacet::GetArea ( ) const

virtual

Implements G4VFacet.

Definition at line 793 of file G4TriangularFacet.cc.

{
  return fArea;
}

References fArea.

Referenced by G4QuadrangularFacet::GetArea(), and G4QuadrangularFacet::GetPointOnFace().

◆ GetCircumcentre()

G4ThreeVector G4TriangularFacet::GetCircumcentre ( ) const

inlinevirtual

Implements G4VFacet.

Definition at line 152 of file G4TriangularFacet.hh.

{
  return fCircumcentre;
}

References fCircumcentre.

◆ GetClone()

G4VFacet * G4TriangularFacet::GetClone ( )

virtual

Implements G4VFacet.

Definition at line 255 of file G4TriangularFacet.cc.

{
  G4TriangularFacet* fc =
    new G4TriangularFacet (GetVertex(0), GetVertex(1), GetVertex(2), ABSOLUTE);
  return fc;
}

References ABSOLUTE, G4TriangularFacet(), and GetVertex().

◆ GetEntityType()

G4GeometryType G4TriangularFacet::GetEntityType ( ) const

virtual

Implements G4VFacet.

Definition at line 800 of file G4TriangularFacet.cc.

{
  return "G4TriangularFacet";
}

◆ GetFlippedFacet()

G4TriangularFacet * G4TriangularFacet::GetFlippedFacet ( )

Definition at line 269 of file G4TriangularFacet.cc.

{
  G4TriangularFacet* flipped =
    new G4TriangularFacet (GetVertex(0), GetVertex(1), GetVertex(2), ABSOLUTE);
  return flipped;
}

References ABSOLUTE, G4TriangularFacet(), and GetVertex().

◆ GetNumberOfVertices()

G4int G4TriangularFacet::GetNumberOfVertices ( ) const

inlinevirtual

Implements G4VFacet.

Definition at line 136 of file G4TriangularFacet.hh.

{
  return 3;
}

Referenced by AllocatedMemory().

◆ GetPointOnFace()

G4ThreeVector G4TriangularFacet::GetPointOnFace ( ) const

virtual

Implements G4VFacet.

Definition at line 779 of file G4TriangularFacet.cc.

{
  G4double u = G4UniformRand();
  G4double v = G4UniformRand();
  if (u+v > 1.) { u = 1. - u; v = 1. - v; }
  return GetVertex(0) + u*fE1 + v*fE2;
}

References fE1, fE2, G4UniformRand, and GetVertex().

Referenced by G4QuadrangularFacet::GetPointOnFace().

◆ GetRadius()

G4double G4TriangularFacet::GetRadius ( ) const

inlinevirtual

Implements G4VFacet.

Definition at line 157 of file G4TriangularFacet.hh.

{
  return fRadius;
}

References fRadius.

◆ GetSurfaceNormal()

G4ThreeVector G4TriangularFacet::GetSurfaceNormal ( ) const

virtual

Implements G4VFacet.

Definition at line 807 of file G4TriangularFacet.cc.

{
  return fSurfaceNormal;
}

References fSurfaceNormal.

Referenced by G4QuadrangularFacet::GetSurfaceNormal().

◆ GetVertex()

G4ThreeVector G4TriangularFacet::GetVertex ( G4int i ) const

inlinevirtual

Implements G4VFacet.

Definition at line 141 of file G4TriangularFacet.hh.

{      
  G4int indice = fIndices[i];
  return indice < 0 ? (*fVertices)[i] : (*fVertices)[indice];
}

References fIndices.

Referenced by Distance(), Extent(), G4TriangularFacet(), GetClone(), GetFlippedFacet(), GetPointOnFace(), G4QuadrangularFacet::GetVertex(), and Intersect().

◆ GetVertexIndex()

G4int G4TriangularFacet::GetVertexIndex ( G4int i ) const

inlinevirtual

Implements G4VFacet.

Definition at line 169 of file G4TriangularFacet.hh.

{
  if (i < 3) return fIndices[i];
  else       return 999999999;
}

References fIndices.

Referenced by G4QuadrangularFacet::GetVertexIndex().

◆ Intersect()

G4bool G4TriangularFacet::Intersect	(	const G4ThreeVector &	p,
		const G4ThreeVector &	v,
		const G4bool	outgoing,
		G4double &	distance,
		G4double &	distFromSurface,
		G4ThreeVector &	normal
	)

virtual

Implements G4VFacet.

Definition at line 587 of file G4TriangularFacet.cc.

{
  //
  // Check whether the direction of the facet is consistent with the vector v
  // and the need to be outgoing or ingoing.  If inconsistent, disregard and
  // return false.
  //
  G4double w = v.dot(fSurfaceNormal);
  if ((outgoing && w < -dirTolerance) || (!outgoing && w > dirTolerance))
  {
    distance = kInfinity;
    distFromSurface = kInfinity;
    normal.set(0,0,0);
    return false;
  } 
  //
  // Calculate the orthogonal distance from p to the surface containing the
  // triangle.  Then determine if we're on the right or wrong side of the
  // surface (at fA distance greater than kCarTolerance to be consistent with
  // "outgoing".
  //
  const G4ThreeVector& p0 = GetVertex(0);
  G4ThreeVector D  = p0 - p;
  distFromSurface  = D.dot(fSurfaceNormal);
  G4bool wrongSide = (outgoing && distFromSurface < -0.5*kCarTolerance) ||
    (!outgoing && distFromSurface >  0.5*kCarTolerance);
 
  if (wrongSide)
  {
    distance = kInfinity;
    distFromSurface = kInfinity;
    normal.set(0,0,0);
    return false;
  }
 
  wrongSide = (outgoing && distFromSurface < 0.0)
           || (!outgoing && distFromSurface > 0.0);
  if (wrongSide)
  {
    //
    // We're slightly on the wrong side of the surface.  Check if we're close
    // enough using fA precise distance calculation.
    //
    G4ThreeVector u = Distance(p);
    if (fSqrDist <= kCarTolerance*kCarTolerance)
    {
      //
      // We're very close.  Therefore return fA small negative number
      // to pretend we intersect.
      //
      // distance = -0.5*kCarTolerance
      distance = 0.0;
      normal = fSurfaceNormal;
      return true;
    }
    else
    {
      //
      // We're close to the surface containing the triangle, but sufficiently
      // far from the triangle, and on the wrong side compared to the directions
      // of the surface normal and v.  There is no intersection.
      //
      distance = kInfinity;
      distFromSurface = kInfinity;
      normal.set(0,0,0);
      return false;
    }
  }
  if (w < dirTolerance && w > -dirTolerance)
  {
    //
    // The ray is within the plane of the triangle. Project the problem into 2D
    // in the plane of the triangle. First try to create orthogonal unit vectors
    // mu and nu, where mu is fE1/|fE1|.  This is kinda like
    // the original algorithm due to Rickard Holmberg, but with better
    // mathematical justification than the original method ... however,
    // beware Rickard's was less time-consuming.
    //
    // Note that vprime is not fA unit vector.  We need to keep it unnormalised
    // since the values of distance along vprime (s0 and s1) for intersection
    // with the triangle will be used to determine if we cut the plane at the
    // same time.
    //
    G4ThreeVector mu = fE1.unit();
    G4ThreeVector nu = fSurfaceNormal.cross(mu);
    G4TwoVector pprime(p.dot(mu), p.dot(nu));
    G4TwoVector vprime(v.dot(mu), v.dot(nu));
    G4TwoVector P0prime(p0.dot(mu), p0.dot(nu));
    G4TwoVector E0prime(fE1.mag(), 0.0);
    G4TwoVector E1prime(fE2.dot(mu), fE2.dot(nu));
    G4TwoVector loc[2];
    if (G4TessellatedGeometryAlgorithms::IntersectLineAndTriangle2D(pprime,
                                    vprime, P0prime, E0prime, E1prime, loc))
    {
      //
      // There is an intersection between the line and triangle in 2D.
      // Now check which part of the line intersects with the plane
      // containing the triangle in 3D.
      //
      G4double vprimemag = vprime.mag();
      G4double s0        = (loc[0] - pprime).mag()/vprimemag;
      G4double s1        = (loc[1] - pprime).mag()/vprimemag;
      G4double normDist0 = fSurfaceNormal.dot(s0*v) - distFromSurface;
      G4double normDist1 = fSurfaceNormal.dot(s1*v) - distFromSurface;
 
      if ((normDist0 < 0.0 && normDist1 < 0.0)
       || (normDist0 > 0.0 && normDist1 > 0.0)
       || (normDist0 == 0.0 && normDist1 == 0.0) ) 
      {
        distance        = kInfinity;
        distFromSurface = kInfinity;
        normal.set(0,0,0);
        return false;
      }
      else
      {
        G4double dnormDist = normDist1 - normDist0;
        if (fabs(dnormDist) < DBL_EPSILON)
        {
          distance = s0;
          normal   = fSurfaceNormal;
          if (!outgoing) distFromSurface = -distFromSurface;
          return true;
        }
        else
        {
          distance = s0 - normDist0*(s1-s0)/dnormDist;
          normal   = fSurfaceNormal;
          if (!outgoing) distFromSurface = -distFromSurface;
          return true;
        }
      }
    }
    else
    {
      distance = kInfinity;
      distFromSurface = kInfinity;
      normal.set(0,0,0);
      return false;
    }
  }
  //
  //
  // Use conventional algorithm to determine the whether there is an
  // intersection.  This involves determining the point of intersection of the
  // line with the plane containing the triangle, and then calculating if the
  // point is within the triangle.
  //
  distance = distFromSurface / w;
  G4ThreeVector pp = p + v*distance;
  G4ThreeVector DD = p0 - pp;
  G4double d = fE1.dot(DD);
  G4double e = fE2.dot(DD);
  G4double ss = fB*e - fC*d;
  G4double t = fB*d - fA*e;
 
  G4double sTolerance = (fabs(fB)+ fabs(fC) + fabs(d) + fabs(e))*kCarTolerance;
  G4double tTolerance = (fabs(fA)+ fabs(fB) + fabs(d) + fabs(e))*kCarTolerance;
  G4double detTolerance = (fabs(fA)+ fabs(fC) + 2*fabs(fB) )*kCarTolerance;
 
  //if (ss < 0.0 || t < 0.0 || ss+t > fDet)
  if (ss < -sTolerance || t < -tTolerance || ( ss+t - fDet ) > detTolerance)
  {
    //
    // The intersection is outside of the triangle.
    //
    distance = distFromSurface = kInfinity;
    normal.set(0,0,0);
    return false;
  }
  else
  {
    //
    // There is an intersection.  Now we only need to set the surface normal.
    //
    normal = fSurfaceNormal;
    if (!outgoing) distFromSurface = -distFromSurface;
    return true;
  }
}

References CLHEP::Hep3Vector::cross(), D(), DBL_EPSILON, G4VFacet::dirTolerance, Distance(), CLHEP::Hep3Vector::dot(), fA, fB, fC, fDet, fE1, fE2, fSqrDist, fSurfaceNormal, GetVertex(), G4TessellatedGeometryAlgorithms::IntersectLineAndTriangle2D(), G4VFacet::kCarTolerance, kInfinity, CLHEP::Hep3Vector::mag(), CLHEP::Hep2Vector::mag(), CLHEP::normal(), G4InuclParticleNames::pp, G4InuclParticleNames::s0, and CLHEP::Hep3Vector::unit().

Referenced by G4QuadrangularFacet::Intersect().

◆ IsDefined()

G4bool G4TriangularFacet::IsDefined ( ) const

inlinevirtual

Implements G4VFacet.

Definition at line 131 of file G4TriangularFacet.hh.

{
  return fIsDefined;
}

References fIsDefined.

Referenced by G4QuadrangularFacet::IsDefined().

◆ IsInside()

G4bool G4VFacet::IsInside ( const G4ThreeVector & p ) const

inherited

Definition at line 112 of file G4VFacet.cc.

{
  G4ThreeVector d =  p-GetVertex(0);
  G4double displacement = d.dot(GetSurfaceNormal());
  return displacement <= 0.0;
}

References CLHEP::Hep3Vector::dot(), G4VFacet::GetSurfaceNormal(), and G4VFacet::GetVertex().

Referenced by G4TessellatedSolid::SetExtremeFacets().

◆ MoveFrom()

void G4TriangularFacet::MoveFrom ( G4TriangularFacet & rhs )

private

Definition at line 185 of file G4TriangularFacet.cc.

{
  fSurfaceNormal = move(rhs.fSurfaceNormal);
  fArea = move(rhs.fArea);
  fCircumcentre = move(rhs.fCircumcentre);
  fRadius = move(rhs.fRadius);
  fIndices = move(rhs.fIndices);
  fA = move(rhs.fA); fB = move(rhs.fB); fC = move(rhs.fC);
  fDet = move(rhs.fDet);
  fSqrDist = move(rhs.fSqrDist);
  fE1 = move(rhs.fE1); fE2 = move(rhs.fE2);
  fIsDefined = move(rhs.fIsDefined);
  fVertices = move(rhs.fVertices);
  rhs.fVertices = nullptr;
}

References fA, fArea, fB, fC, fCircumcentre, fDet, fE1, fE2, fIndices, fIsDefined, fRadius, fSqrDist, fSurfaceNormal, and fVertices.

Referenced by G4TriangularFacet(), and operator=().

◆ operator=() [1/2]

G4TriangularFacet & G4TriangularFacet::operator= ( const G4TriangularFacet & right )

Definition at line 220 of file G4TriangularFacet.cc.

{
  SetVertices(nullptr);
 
  if (this != &rhs)
  {
    delete fVertices;
    CopyFrom(rhs);
  }
 
  return *this;
}

References CopyFrom(), fVertices, and SetVertices().

◆ operator=() [2/2]

G4TriangularFacet & G4TriangularFacet::operator= ( G4TriangularFacet && right )

Definition at line 236 of file G4TriangularFacet.cc.

{
  SetVertices(nullptr);
 
  if (this != &rhs)
  {
    delete fVertices;
    MoveFrom(rhs);
  }
 
  return *this;
}

References fVertices, MoveFrom(), and SetVertices().

◆ operator==()

G4bool G4VFacet::operator== ( const G4VFacet & right ) const

inherited

Definition at line 57 of file G4VFacet.cc.

{
  G4double tolerance = kCarTolerance*kCarTolerance/4.0;
 
  if (GetNumberOfVertices() != right.GetNumberOfVertices())
    return false;
  else if ((GetCircumcentre()-right.GetCircumcentre()).mag2() > tolerance)
    return false;
  else if (std::fabs((right.GetSurfaceNormal()).dot(GetSurfaceNormal())) < 0.9999999999)
    return false;
 
  G4bool coincident  = true;
  G4int i = 0;
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    coincident = false;
    G4int j   = 0; 
    do    // Loop checking, 13.08.2015, G.Cosmo
    {
      coincident = (GetVertex(i)-right.GetVertex(j)).mag2() < tolerance;
    } while (!coincident && ++j < GetNumberOfVertices());
  } while (coincident && ++i < GetNumberOfVertices());
 
  return coincident;
}

References G4VFacet::GetCircumcentre(), G4VFacet::GetNumberOfVertices(), G4VFacet::GetSurfaceNormal(), G4VFacet::GetVertex(), and G4VFacet::kCarTolerance.

◆ SetSurfaceNormal()

void G4TriangularFacet::SetSurfaceNormal ( G4ThreeVector normal )

Definition at line 814 of file G4TriangularFacet.cc.

{
  fSurfaceNormal = normal;
}

References fSurfaceNormal, and CLHEP::normal().

Referenced by G4QuadrangularFacet::G4QuadrangularFacet().

◆ SetVertex()

void G4TriangularFacet::SetVertex	(	G4int	i,
		const G4ThreeVector &	val
	)

inlinevirtual

Implements G4VFacet.

Definition at line 147 of file G4TriangularFacet.hh.

{
  (*fVertices)[i] = val;
}

Referenced by G4TriangularFacet(), and G4QuadrangularFacet::SetVertex().

◆ SetVertexIndex()

void G4TriangularFacet::SetVertexIndex	(	G4int	i,
		G4int	j
	)

inlinevirtual

Implements G4VFacet.

Definition at line 175 of file G4TriangularFacet.hh.

{
  fIndices[i] = j;
}

References fIndices.

Referenced by G4QuadrangularFacet::SetVertexIndex().

◆ SetVertices()

void G4TriangularFacet::SetVertices ( std::vector< G4ThreeVector > * v )

inlinevirtual

Implements G4VFacet.

Definition at line 180 of file G4TriangularFacet.hh.

{
  if (fIndices[0] < 0 && fVertices)
  {
    delete fVertices;
    fVertices = nullptr;
  }
  fVertices = v;
}

References fIndices, and fVertices.

Referenced by operator=(), G4QuadrangularFacet::SetVertices(), and ~G4TriangularFacet().

◆ StreamInfo()

std::ostream & G4VFacet::StreamInfo ( std::ostream & os ) const

inherited

Definition at line 96 of file G4VFacet.cc.

{
  os << G4endl;
  os << "*********************************************************************"
     << G4endl;
  os << "FACET TYPE       = " << GetEntityType() << G4endl;
  os << "ABSOLUTE VECTORS = " << G4endl;
  G4int n = GetNumberOfVertices();
  for (G4int i = 0; i < n; ++i)
    os << "P[" << i << "]      = " << GetVertex(i) << G4endl;
  os << "*********************************************************************"
     << G4endl;
 
  return os;
}

References G4endl, G4VFacet::GetEntityType(), G4VFacet::GetNumberOfVertices(), G4VFacet::GetVertex(), and CLHEP::detail::n.

Referenced by G4TessellatedSolid::AddFacet(), and G4TessellatedSolid::StreamInfo().

Field Documentation

◆ dirTolerance

const G4double G4VFacet::dirTolerance = 1.0E-14

staticprotectedinherited

Definition at line 92 of file G4VFacet.hh.

Referenced by G4QuadrangularFacet::Distance(), and Intersect().

◆ fA

G4double G4TriangularFacet::fA

private

Definition at line 120 of file G4TriangularFacet.hh.

Referenced by Distance(), G4TriangularFacet(), Intersect(), and MoveFrom().

◆ fArea

G4double G4TriangularFacet::fArea = 0.0

private

Definition at line 114 of file G4TriangularFacet.hh.

Referenced by G4TriangularFacet(), GetArea(), and MoveFrom().

◆ fB

G4double G4TriangularFacet::fB

private

Definition at line 120 of file G4TriangularFacet.hh.

Referenced by Distance(), G4TriangularFacet(), Intersect(), and MoveFrom().

◆ fC

G4double G4TriangularFacet::fC

private

Definition at line 120 of file G4TriangularFacet.hh.

Referenced by Distance(), G4TriangularFacet(), Intersect(), and MoveFrom().

◆ fCircumcentre

G4ThreeVector G4TriangularFacet::fCircumcentre

private

Definition at line 115 of file G4TriangularFacet.hh.

Referenced by Distance(), G4TriangularFacet(), GetCircumcentre(), and MoveFrom().

◆ fDet

G4double G4TriangularFacet::fDet

private

Definition at line 121 of file G4TriangularFacet.hh.

Referenced by Distance(), G4TriangularFacet(), Intersect(), and MoveFrom().

◆ fE1

G4ThreeVector G4TriangularFacet::fE1

private

Definition at line 123 of file G4TriangularFacet.hh.

Referenced by Distance(), G4TriangularFacet(), GetPointOnFace(), Intersect(), and MoveFrom().

◆ fE2

G4ThreeVector G4TriangularFacet::fE2

private

Definition at line 123 of file G4TriangularFacet.hh.

Referenced by Distance(), G4TriangularFacet(), GetPointOnFace(), Intersect(), and MoveFrom().

◆ fIndices

std::array<G4int, 3> G4TriangularFacet::fIndices

private

Definition at line 117 of file G4TriangularFacet.hh.

Referenced by CopyFrom(), G4TriangularFacet(), GetVertex(), GetVertexIndex(), MoveFrom(), SetVertexIndex(), and SetVertices().

◆ fIsDefined

G4bool G4TriangularFacet::fIsDefined = false

private

Definition at line 124 of file G4TriangularFacet.hh.

Referenced by G4TriangularFacet(), IsDefined(), and MoveFrom().

◆ fRadius

G4double G4TriangularFacet::fRadius = 0.0

private

Definition at line 116 of file G4TriangularFacet.hh.

Referenced by Distance(), G4TriangularFacet(), GetRadius(), and MoveFrom().

◆ fSqrDist

G4double G4TriangularFacet::fSqrDist = 0.0

private

Definition at line 122 of file G4TriangularFacet.hh.

Referenced by Distance(), Intersect(), and MoveFrom().

◆ fSurfaceNormal

G4ThreeVector G4TriangularFacet::fSurfaceNormal

private

Definition at line 113 of file G4TriangularFacet.hh.

Referenced by Distance(), G4TriangularFacet(), GetSurfaceNormal(), Intersect(), MoveFrom(), and SetSurfaceNormal().

◆ fVertices

std::vector<G4ThreeVector>* G4TriangularFacet::fVertices = nullptr

private

Definition at line 118 of file G4TriangularFacet.hh.

Referenced by CopyFrom(), G4TriangularFacet(), MoveFrom(), operator=(), and SetVertices().

◆ kCarTolerance

G4double G4VFacet::kCarTolerance

protectedinherited

Definition at line 93 of file G4VFacet.hh.

Referenced by Distance(), G4QuadrangularFacet::G4QuadrangularFacet(), G4TriangularFacet(), G4VFacet::G4VFacet(), Intersect(), and G4VFacet::operator==().

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

source/geometry/solids/specific/include/G4TriangularFacet.hh
source/geometry/solids/specific/src/G4TriangularFacet.cc

Public Member Functions

Protected Attributes

Static Protected Attributes

Private Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

◆ G4TriangularFacet() [1/4]

◆ ~G4TriangularFacet()

◆ G4TriangularFacet() [2/4]

◆ G4TriangularFacet() [3/4]

◆ G4TriangularFacet() [4/4]

Member Function Documentation

◆ AllocatedMemory()

◆ ApplyTranslation()

◆ CopyFrom()

◆ Distance() [1/3]

◆ Distance() [2/3]

◆ Distance() [3/3]

◆ Extent()

◆ GetArea()

◆ GetCircumcentre()

◆ GetClone()

◆ GetEntityType()

◆ GetFlippedFacet()

◆ GetNumberOfVertices()

◆ GetPointOnFace()

◆ GetRadius()

◆ GetSurfaceNormal()

◆ GetVertex()

◆ GetVertexIndex()

◆ Intersect()

◆ IsDefined()

◆ IsInside()

◆ MoveFrom()

◆ operator=() [1/2]

◆ operator=() [2/2]

◆ operator==()

◆ SetSurfaceNormal()

◆ SetVertex()

◆ SetVertexIndex()

◆ SetVertices()

◆ StreamInfo()

Field Documentation

◆ dirTolerance

◆ fA

◆ fArea

◆ fB

◆ fC

◆ fCircumcentre

◆ fDet

◆ fE1

◆ fE2

◆ fIndices

◆ fIsDefined

◆ fRadius

◆ fSqrDist

◆ fSurfaceNormal

◆ fVertices

◆ kCarTolerance