G4PolyPhiFace Class Reference

#include <G4PolyPhiFace.hh>

Inheritance diagram for G4PolyPhiFace:

Collaboration diagram for G4PolyPhiFace:

Public Member Functions
	G4PolyPhiFace (const G4ReduciblePolygon *rz, G4double phi, G4double deltaPhi, G4double phiOther)
virtual	~G4PolyPhiFace ()
	G4PolyPhiFace (const G4PolyPhiFace &source)
G4PolyPhiFace &	operator= (const G4PolyPhiFace &source)
G4bool	Intersect (const G4ThreeVector &p, const G4ThreeVector &v, G4bool outgoing, G4double surfTolerance, G4double &distance, G4double &distFromSurface, G4ThreeVector &normal, G4bool &allBehind)
G4double	Distance (const G4ThreeVector &p, G4bool outgoing)
EInside	Inside (const G4ThreeVector &p, G4double tolerance, G4double *bestDistance)
G4ThreeVector	Normal (const G4ThreeVector &p, G4double *bestDistance)
G4double	Extent (const G4ThreeVector axis)
void	CalculateExtent (const EAxis axis, const G4VoxelLimits &voxelLimit, const G4AffineTransform &tranform, G4SolidExtentList &extentList)
G4VCSGface *	Clone ()
G4double	SurfaceArea ()
G4double	SurfaceTriangle (G4ThreeVector p1, G4ThreeVector p2, G4ThreeVector p3, G4ThreeVector *p4)
G4ThreeVector	GetPointOnFace ()
	G4PolyPhiFace (__void__ &)
void	Diagnose (G4VSolid *solid)
Public Member Functions inherited from G4VCSGface
	G4VCSGface ()
virtual	~G4VCSGface ()

Protected Member Functions
G4bool	InsideEdgesExact (G4double r, G4double z, G4double normSign, const G4ThreeVector &p, const G4ThreeVector &v)
G4bool	InsideEdges (G4double r, G4double z)
G4bool	InsideEdges (G4double r, G4double z, G4double *distRZ2, G4PolyPhiFaceVertex **base3Dnorm=0, G4ThreeVector **head3Dnorm=0)
G4double	ExactZOrder (G4double z, G4double qx, G4double qy, G4double qz, const G4ThreeVector &v, G4double normSign, const G4PolyPhiFaceVertex *vert) const
void	CopyStuff (const G4PolyPhiFace &source)
G4double	Area2 (G4TwoVector a, G4TwoVector b, G4TwoVector c)
G4bool	Left (G4TwoVector a, G4TwoVector b, G4TwoVector c)
G4bool	LeftOn (G4TwoVector a, G4TwoVector b, G4TwoVector c)
G4bool	Collinear (G4TwoVector a, G4TwoVector b, G4TwoVector c)
G4bool	IntersectProp (G4TwoVector a, G4TwoVector b, G4TwoVector c, G4TwoVector d)
G4bool	Between (G4TwoVector a, G4TwoVector b, G4TwoVector c)
G4bool	Intersect (G4TwoVector a, G4TwoVector b, G4TwoVector c, G4TwoVector d)
G4bool	Diagonalie (G4PolyPhiFaceVertex *a, G4PolyPhiFaceVertex *b)
G4bool	InCone (G4PolyPhiFaceVertex *a, G4PolyPhiFaceVertex *b)
G4bool	Diagonal (G4PolyPhiFaceVertex *a, G4PolyPhiFaceVertex *b)
void	EarInit ()
void	Triangulate ()

Protected Attributes
G4int	numEdges
G4PolyPhiFaceEdge *	edges
G4PolyPhiFaceVertex *	corners
G4ThreeVector	normal
G4ThreeVector	radial
G4ThreeVector	surface
G4ThreeVector	surface_point
G4double	rMin
G4double	rMax
G4double	zMin
G4double	zMax
G4bool	allBehind
G4double	kCarTolerance
G4double	fSurfaceArea
G4PolyPhiFaceVertex *	triangles

Detailed Description

Definition at line 85 of file G4PolyPhiFace.hh.

Constructor & Destructor Documentation

G4PolyPhiFace::G4PolyPhiFace	(	const G4ReduciblePolygon *	rz,
		G4double	phi,
		G4double	deltaPhi,
		G4double	phiOther
	)

Definition at line 61 of file G4PolyPhiFace.cc.

  : fSurfaceArea(0.), triangles(0)  
{
  kCarTolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();

  numEdges = rz->NumVertices();
  
  rMin = rz->Amin();
  rMax = rz->Amax();
  zMin = rz->Bmin();
  zMax = rz->Bmax();

  //
  // Is this the "starting" phi edge of the two?
  //
  G4bool start = (phiOther > phi);
  
  //
  // Build radial vector
  //
  radial = G4ThreeVector( std::cos(phi), std::sin(phi), 0.0 );
  
  //
  // Build normal
  //
  G4double zSign = start ? 1 : -1;
  normal = G4ThreeVector( zSign*radial.y(), -zSign*radial.x(), 0 );
  
  //
  // Is allBehind?
  //
  allBehind = (zSign*(std::cos(phiOther)*radial.y() - std::sin(phiOther)*radial.x()) < 0);
  
  //
  // Adjacent edges
  //
  G4double midPhi = phi + (start ? +0.5 : -0.5)*deltaPhi;
  G4double cosMid = std::cos(midPhi), 
                 sinMid = std::sin(midPhi);

  //
  // Allocate corners
  //
  corners = new G4PolyPhiFaceVertex[numEdges];
  //
  // Fill them
  //
  G4ReduciblePolygonIterator iterRZ(rz);
  
  G4PolyPhiFaceVertex *corn = corners;
  G4PolyPhiFaceVertex *helper=corners;

  iterRZ.Begin();
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    corn->r = iterRZ.GetA();
    corn->z = iterRZ.GetB();
    corn->x = corn->r*radial.x();
    corn->y = corn->r*radial.y();

    // Add pointer on prev corner
    //
    if( corn == corners )
      { corn->prev = corners+numEdges-1;}
    else
      { corn->prev = helper; }

    // Add pointer on next corner
    //
    if( corn < corners+numEdges-1 )
      { corn->next = corn+1;}
    else
      { corn->next = corners; }

    helper = corn;
  } while( ++corn, iterRZ.Next() );

  //
  // Allocate edges
  //
  edges = new G4PolyPhiFaceEdge[numEdges];

  //
  // Fill them
  //
  G4double rFact = std::cos(0.5*deltaPhi);
  G4double rFactNormalize = 1.0/std::sqrt(1.0+rFact*rFact);

  G4PolyPhiFaceVertex *prev = corners+numEdges-1,
                      *here = corners;
  G4PolyPhiFaceEdge   *edge = edges;
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    G4ThreeVector sideNorm;
    
    edge->v0 = prev;
    edge->v1 = here;

    G4double dr = here->r - prev->r,
             dz = here->z - prev->z;
                         
    edge->length = std::sqrt( dr*dr + dz*dz );

    edge->tr = dr/edge->length;
    edge->tz = dz/edge->length;
    
    if ((here->r < DBL_MIN) && (prev->r < DBL_MIN))
    {
      //
      // Sigh! Always exceptions!
      // This edge runs at r==0, so its adjoing surface is not a
      // PolyconeSide or PolyhedraSide, but the opposite PolyPhiFace.
      //
      G4double zSignOther = start ? -1 : 1;
      sideNorm = G4ThreeVector(  zSignOther*std::sin(phiOther), 
                                -zSignOther*std::cos(phiOther), 0 );
    }
    else
    {
      sideNorm = G4ThreeVector( edge->tz*cosMid,
                                edge->tz*sinMid,
                               -edge->tr*rFact );
      sideNorm *= rFactNormalize;
    }
    sideNorm += normal;
    
    edge->norm3D = sideNorm.unit();
  } while( edge++, prev=here, ++here < corners+numEdges );

  //
  // Go back and fill in corner "normals"
  //
  G4PolyPhiFaceEdge *prevEdge = edges+numEdges-1;
  edge = edges;
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    //
    // Calculate vertex 2D normals (on the phi surface)
    //
    G4double rPart = prevEdge->tr + edge->tr;
    G4double zPart = prevEdge->tz + edge->tz;
    G4double norm = std::sqrt( rPart*rPart + zPart*zPart );
    G4double rNorm = +zPart/norm;
    G4double zNorm = -rPart/norm;
    
    edge->v0->rNorm = rNorm;
    edge->v0->zNorm = zNorm;
    
    //
    // Calculate the 3D normals.
    //
    // Find the vector perpendicular to the z axis
    // that defines the plane that contains the vertex normal
    //
    G4ThreeVector xyVector;
    
    if (edge->v0->r < DBL_MIN)
    {
      //
      // This is a vertex at r==0, which is a special
      // case. The normal we will construct lays in the
      // plane at the center of the phi opening.
      //
      // We also know that rNorm < 0
      //
      G4double zSignOther = start ? -1 : 1;
      G4ThreeVector normalOther(  zSignOther*std::sin(phiOther), 
                                 -zSignOther*std::cos(phiOther), 0 );
                
      xyVector = - normal - normalOther;
    }
    else
    {
      //
      // This is a vertex at r > 0. The plane
      // is the average of the normal and the
      // normal of the adjacent phi face
      //
      xyVector = G4ThreeVector( cosMid, sinMid, 0 );
      if (rNorm < 0)
        xyVector -= normal;
      else
        xyVector += normal;
    }
    
    //
    // Combine it with the r/z direction from the face
    //
    edge->v0->norm3D = rNorm*xyVector.unit() + G4ThreeVector( 0, 0, zNorm );
  } while(  prevEdge=edge, ++edge < edges+numEdges );
  
  //
  // Build point on surface
  //
  G4double rAve = 0.5*(rMax-rMin),
           zAve = 0.5*(zMax-zMin);
  surface = G4ThreeVector( rAve*radial.x(), rAve*radial.y(), zAve );
}

Here is the call graph for this function:

G4PolyPhiFace::~G4PolyPhiFace ( )

virtual

Definition at line 301 of file G4PolyPhiFace.cc.

{
  delete [] edges;
  delete [] corners;
}

G4PolyPhiFace::G4PolyPhiFace

(

const G4PolyPhiFace &

source

)

Definition at line 311 of file G4PolyPhiFace.cc.

  : G4VCSGface()
{
  CopyStuff( source );
}

Here is the call graph for this function:

G4PolyPhiFace::G4PolyPhiFace

(

__void__ &

)

Definition at line 291 of file G4PolyPhiFace.cc.

  : numEdges(0), edges(0), corners(0), rMin(0.), rMax(0.), zMin(0.), zMax(0.),
    allBehind(false), kCarTolerance(0.), fSurfaceArea(0.), triangles(0)
{
}

Member Function Documentation

G4double G4PolyPhiFace::Area2 ( G4TwoVector  a, G4TwoVector  b, G4TwoVector  c  )

protected

Definition at line 960 of file G4PolyPhiFace.cc.

{
  return ((b.x()-a.x())*(c.y()-a.y())-
          (c.x()-a.x())*(b.y()-a.y()));
}

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4PolyPhiFace::Between ( G4TwoVector  a, G4TwoVector  b, G4TwoVector  c  )

protected

Definition at line 1018 of file G4PolyPhiFace.cc.

{
  if( !Collinear(a,b,c) ) { return false; }

  if(a.x()!=b.x())
  {
    return ((a.x()<=c.x())&&(c.x()<=b.x()))||
           ((a.x()>=c.x())&&(c.x()>=b.x()));
  }
  else
  {
    return ((a.y()<=c.y())&&(c.y()<=b.y()))||
           ((a.y()>=c.y())&&(c.y()>=b.y()));
  }
}

Here is the call graph for this function:

Here is the caller graph for this function:

void G4PolyPhiFace::CalculateExtent ( const EAxis  axis, const G4VoxelLimits &  voxelLimit, const G4AffineTransform &  tranform, G4SolidExtentList &  extentList  )

virtual

Implements G4VCSGface.

Definition at line 628 of file G4PolyPhiFace.cc.

{
  //
  // Construct a (sometimes big) clippable polygon, 
  //
  // Perform the necessary transformations while doing so
  //
  G4ClippablePolygon polygon;
  
  G4PolyPhiFaceVertex *corner = corners;
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    G4ThreeVector point( 0, 0, corner->z );
    point += radial*corner->r;
    
    polygon.AddVertexInOrder( transform.TransformPoint( point ) );
  } while( ++corner < corners + numEdges );
  
  //
  // Clip away
  //
  if (polygon.PartialClip( voxelLimit, axis ))
  {
    //
    // Add it to the list
    //
    polygon.SetNormal( transform.TransformAxis(normal) );
    extentList.AddSurface( polygon );
  }
}

Here is the call graph for this function:

G4VCSGface* G4PolyPhiFace::Clone ( )

inlinevirtual

Implements G4VCSGface.

G4bool G4PolyPhiFace::Collinear ( G4TwoVector  a, G4TwoVector  b, G4TwoVector  c  )

protected

Definition at line 991 of file G4PolyPhiFace.cc.

{
  return Area2(a,b,c)==0;
}

Here is the call graph for this function:

Here is the caller graph for this function:

void G4PolyPhiFace::CopyStuff ( const G4PolyPhiFace &  source )

protected

Definition at line 337 of file G4PolyPhiFace.cc.

{
  //
  // The simple stuff
  //
  numEdges  = source.numEdges;
  normal    = source.normal;
  radial    = source.radial;
  surface   = source.surface;
  rMin    = source.rMin;
  rMax    = source.rMax;
  zMin    = source.zMin;
  zMax    = source.zMax;
  allBehind  = source.allBehind;
  triangles  = 0;

  kCarTolerance = source.kCarTolerance;
  fSurfaceArea = source.fSurfaceArea;

  //
  // Corner dynamic array
  //
  corners = new G4PolyPhiFaceVertex[numEdges];
  G4PolyPhiFaceVertex *corn = corners,
                      *sourceCorn = source.corners;
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    *corn = *sourceCorn;
  } while( ++sourceCorn, ++corn < corners+numEdges );
  
  //
  // Edge dynamic array
  //
  edges = new G4PolyPhiFaceEdge[numEdges];

  G4PolyPhiFaceVertex *prev = corners+numEdges-1,
                      *here = corners;
  G4PolyPhiFaceEdge   *edge = edges,
                      *sourceEdge = source.edges;
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    *edge = *sourceEdge;
    edge->v0 = prev;
    edge->v1 = here;
  } while( ++sourceEdge, ++edge, prev=here, ++here < corners+numEdges );
}

Here is the caller graph for this function:

void G4PolyPhiFace::Diagnose

(

G4VSolid *

solid

)

Definition at line 272 of file G4PolyPhiFace.cc.

{
  G4PolyPhiFaceVertex   *corner = corners;
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    G4ThreeVector test(corner->x, corner->y, corner->z);
    test -= 1E-6*corner->norm3D;
    
    if (owner->Inside(test) != kInside) 
      G4Exception( "G4PolyPhiFace::Diagnose()", "GeomSolids0002",
                   FatalException, "Bad vertex normal found." );
  } while( ++corner < corners+numEdges );
}

Here is the call graph for this function:

G4bool G4PolyPhiFace::Diagonal ( G4PolyPhiFaceVertex *  a, G4PolyPhiFaceVertex *  b  )

protected

Definition at line 1118 of file G4PolyPhiFace.cc.

{ 
  return InCone(a,b) && InCone(b,a) && Diagonalie(a,b);
}

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4PolyPhiFace::Diagonalie ( G4PolyPhiFaceVertex *  a, G4PolyPhiFaceVertex *  b  )

protected

Definition at line 1057 of file G4PolyPhiFace.cc.

{
  G4PolyPhiFaceVertex   *corner = triangles;
  G4PolyPhiFaceVertex   *corner_next=triangles;

  // For each Edge (corner,corner_next) 
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    corner_next=corner->next;

    // Skip edges incident to a of b
    //
    if( (corner!=a)&&(corner_next!=a)
      &&(corner!=b)&&(corner_next!=b) )
    {
       G4TwoVector rz1,rz2,rz3,rz4;
       rz1 = G4TwoVector(a->r,a->z);
       rz2 = G4TwoVector(b->r,b->z);
       rz3 = G4TwoVector(corner->r,corner->z);
       rz4 = G4TwoVector(corner_next->r,corner_next->z);
       if( Intersect(rz1,rz2,rz3,rz4) ) { return false; }
    }
    corner=corner->next;   
   
  } while( corner != triangles );

  return true;
}

Here is the call graph for this function:

Here is the caller graph for this function:

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

virtual

Implements G4VCSGface.

Definition at line 445 of file G4PolyPhiFace.cc.

{
  G4double normSign = outgoing ? +1 : -1;
  //
  // Correct normal? 
  //
  G4ThreeVector ps = p - surface;
  G4double distPhi = -normSign*normal.dot(ps);
  
  if (distPhi < -0.5*kCarTolerance) 
    return kInfinity;
  else if (distPhi < 0)
    distPhi = 0.0;
  
  //
  // Calculate projected point in r,z
  //
  G4double r = radial.dot(p);
  
  //
  // Are we inside the face?
  //
  G4double distRZ2;
  
  if (InsideEdges( r, p.z(), &distRZ2, 0 ))
  {
    //
    // Yup, answer is just distPhi
    //
    return distPhi;
  }
  else
  {
    //
    // Nope. Penalize by distance out
    //
    return std::sqrt( distPhi*distPhi + distRZ2 );
  }
}  

Here is the call graph for this function:

void G4PolyPhiFace::EarInit ( )

protected

Definition at line 1127 of file G4PolyPhiFace.cc.

{
  G4PolyPhiFaceVertex *corner = triangles;
  G4PolyPhiFaceVertex *c_prev,*c_next;
  
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
     // We need to determine three consecutive vertices
     //
     c_next=corner->next;
     c_prev=corner->prev; 

     // Calculation of ears
     //
     corner->ear=Diagonal(c_prev,c_next);   
     corner=corner->next;

  } while( corner!=triangles );
}

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4PolyPhiFace::ExactZOrder ( G4double  z, G4double  qx, G4double  qy, G4double  qz, const G4ThreeVector &  v, G4double  normSign, const G4PolyPhiFaceVertex *  vert  ) const

inlineprotected

Here is the caller graph for this function:

G4double G4PolyPhiFace::Extent ( const G4ThreeVector  axis )

virtual

Implements G4VCSGface.

Definition at line 606 of file G4PolyPhiFace.cc.

{
  G4double max = -kInfinity;
  
  G4PolyPhiFaceVertex *corner = corners;
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    G4double here = axis.x()*corner->r*radial.x()
            + axis.y()*corner->r*radial.y()
            + axis.z()*corner->z;
    if (here > max) max = here;
  } while( ++corner < corners + numEdges );
  
  return max;
}  

Here is the call graph for this function:

G4ThreeVector G4PolyPhiFace::GetPointOnFace ( )

virtual

Implements G4VCSGface.

Definition at line 947 of file G4PolyPhiFace.cc.

{
  Triangulate();
  return surface_point;
}

Here is the call graph for this function:

G4bool G4PolyPhiFace::InCone ( G4PolyPhiFaceVertex *  a, G4PolyPhiFaceVertex *  b  )

protected

Definition at line 1091 of file G4PolyPhiFace.cc.

{
  // a0,a and a1 are consecutive vertices
  //
  G4PolyPhiFaceVertex *a0,*a1;
  a1=a->next;
  a0=a->prev;

  G4TwoVector arz,arz0,arz1,brz;
  arz=G4TwoVector(a->r,a->z);arz0=G4TwoVector(a0->r,a0->z);
  arz1=G4TwoVector(a1->r,a1->z);brz=G4TwoVector(b->r,b->z);
    
  
  if(LeftOn(arz,arz1,arz0))  // If a is convex vertex
  {
    return Left(arz,brz,arz0)&&Left(brz,arz,arz1);
  }
  else                       // Else a is reflex
  {
    return !( LeftOn(arz,brz,arz1)&&LeftOn(brz,arz,arz0));
  }
}

Here is the call graph for this function:

Here is the caller graph for this function:

EInside G4PolyPhiFace::Inside ( const G4ThreeVector &  p, G4double  tolerance, G4double *  bestDistance  )

virtual

Implements G4VCSGface.

Definition at line 489 of file G4PolyPhiFace.cc.

{
  //
  // Get distance along phi, which if negative means the point
  // is nominally inside the shape.
  //
  G4ThreeVector ps = p - surface;
  G4double distPhi = normal.dot(ps);
  
  //
  // Calculate projected point in r,z
  //
  G4double r = radial.dot(p);
  
  //
  // Are we inside the face?
  //
  G4double distRZ2;
  G4PolyPhiFaceVertex *base3Dnorm=0;
  G4ThreeVector       *head3Dnorm=0;
  
  if (InsideEdges( r, p.z(), &distRZ2, &base3Dnorm, &head3Dnorm ))
  {
    //
    // Looks like we're inside. Distance is distance in phi.
    //
    *bestDistance = std::fabs(distPhi);
    
    //
    // Use distPhi to decide fate
    //
    if (distPhi < -tolerance) return kInside;
    if (distPhi <  tolerance) return kSurface;
    return kOutside;
  }
  else
  {
    //
    // We're outside the extent of the face,
    // so the distance is penalized by distance from edges in RZ
    //
    *bestDistance = std::sqrt( distPhi*distPhi + distRZ2 );
    
    //
    // Use edge normal to decide fate
    //
    G4ThreeVector cc( base3Dnorm->r*radial.x(),
                      base3Dnorm->r*radial.y(),
                      base3Dnorm->z );
    cc = p - cc;
    G4double normDist = head3Dnorm->dot(cc);
    if ( distRZ2 > tolerance*tolerance )
    {
      //
      // We're far enough away that kSurface is not possible
      //
      return normDist < 0 ? kInside : kOutside;
    }
    
    if (normDist < -tolerance) return kInside;
    if (normDist <  tolerance) return kSurface;
    return kOutside;
  }
}  

Here is the call graph for this function:

G4bool G4PolyPhiFace::InsideEdges ( G4double  r, G4double  z  )

protected

Definition at line 816 of file G4PolyPhiFace.cc.

{
  //
  // Quick check of extent
  //
  if ( r < rMin || r > rMax ) return false;
  if ( z < zMin || z > zMax ) return false;
  
  //
  // More thorough check
  //
  G4double notUsed;
  
  return InsideEdges( r, z, &notUsed, 0 );
}

Here is the caller graph for this function:

G4bool G4PolyPhiFace::InsideEdges ( G4double  r, G4double  z, G4double *  distRZ2, G4PolyPhiFaceVertex **  base3Dnorm = 0, G4ThreeVector **  head3Dnorm = 0  )

protected

Definition at line 838 of file G4PolyPhiFace.cc.

{
  G4double bestDistance2 = kInfinity;
  G4bool   answer = 0;
  
  G4PolyPhiFaceEdge *edge = edges;
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    G4PolyPhiFaceVertex *testMe;
    //
    // Get distance perpendicular to the edge
    //
    G4double dr = (r-edge->v0->r), dz = (z-edge->v0->z);

    G4double distOut = dr*edge->tz - dz*edge->tr;
    G4double distance2 = distOut*distOut;
    if (distance2 > bestDistance2) continue;        // No hope!

    //
    // Check to see if normal intersects edge within the edge's boundary
    //
    G4double q = dr*edge->tr + dz*edge->tz;

    //
    // If it doesn't, penalize distance2 appropriately
    //
    if (q < 0)
    {
      distance2 += q*q;
      testMe = edge->v0;
    }
    else if (q > edge->length)
    {
      G4double s2 = q-edge->length;
      distance2 += s2*s2;
      testMe = edge->v1;
    }
    else
    {
      testMe = 0;
    }
    
    //
    // Closest edge so far?
    //
    if (distance2 < bestDistance2)
    {
      bestDistance2 = distance2;
      if (testMe)
      {
        G4double distNorm = dr*testMe->rNorm + dz*testMe->zNorm;
        answer = (distNorm <= 0);
        if (base3Dnorm)
        {
          *base3Dnorm = testMe;
          *head3Dnorm = &testMe->norm3D;
        }
      }
      else
      {
        answer = (distOut <= 0);                        
        if (base3Dnorm)
        {
          *base3Dnorm = edge->v0;
          *head3Dnorm = &edge->norm3D;
        }
      }
    }
  } while( ++edge < edges + numEdges );
  
  *bestDist2 = bestDistance2;
  return answer;
}

G4bool G4PolyPhiFace::InsideEdgesExact ( G4double  r, G4double  z, G4double  normSign, const G4ThreeVector &  p, const G4ThreeVector &  v  )

protected

Definition at line 683 of file G4PolyPhiFace.cc.

{
  //
  // Quick check of extent
  //
  if ( (r < rMin-kCarTolerance)
    || (r > rMax+kCarTolerance) ) return false;
       
  if ( (z < zMin-kCarTolerance)
    || (z > zMax+kCarTolerance) ) return false;
  
  //
  // Exact check: loop over all vertices
  //
  G4double qx = p.x() + v.x(),
           qy = p.y() + v.y(),
           qz = p.z() + v.z();

  G4int answer = 0;
  G4PolyPhiFaceVertex *corn = corners, 
                      *prev = corners+numEdges-1;

  G4double cornZ, prevZ;
  
  prevZ = ExactZOrder( z, qx, qy, qz, v, normSign, prev );
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    //
    // Get z order of this vertex, and compare to previous vertex
    //
    cornZ = ExactZOrder( z, qx, qy, qz, v, normSign, corn );
    
    if (cornZ < 0)
    {
      if (prevZ < 0) continue;
    }
    else if (cornZ > 0)
    {
      if (prevZ > 0) continue;
    }
    else
    {
      //
      // By chance, we overlap exactly (within precision) with 
      // the current vertex. Continue if the same happened previously
      // (e.g. the previous vertex had the same z value)
      //
      if (prevZ == 0) continue;
      
      //
      // Otherwise, to decide what to do, we need to know what is
      // coming up next. Specifically, we need to find the next vertex
      // with a non-zero z order.
      //
      // One might worry about infinite loops, but the above conditional
      // should prevent it
      //
      G4PolyPhiFaceVertex *next = corn;
      G4double nextZ;
      do    // Loop checking, 13.08.2015, G.Cosmo
      {
        next++;
        if (next == corners+numEdges) next = corners;

        nextZ = ExactZOrder( z, qx, qy, qz, v, normSign, next );
      } while( nextZ == 0 );
      
      //
      // If we won't be changing direction, go to the next vertex
      //
      if (nextZ*prevZ < 0) continue;
    }
  
      
    //
    // We overlap in z with the side of the face that stretches from
    // vertex "prev" to "corn". On which side (left or right) do
    // we lay with respect to this segment?
    //  
    G4ThreeVector qa( qx - prev->x, qy - prev->y, qz - prev->z ),
                  qb( qx - corn->x, qy - corn->y, qz - corn->z );

    G4double aboveOrBelow = normSign*qa.cross(qb).dot(v);
    
    if (aboveOrBelow > 0) 
      answer++;
    else if (aboveOrBelow < 0)
      answer--;
    else
    {
      //
      // A precisely zero answer here means we exactly
      // intersect (within roundoff) the edge of the face.
      // Return true in this case.
      //
      return true;
    }
  } while( prevZ = cornZ, prev=corn, ++corn < corners+numEdges );
  
  return answer!=0;
}

Here is the call graph for this function:

Here is the caller graph for this function:

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

virtual

Implements G4VCSGface.

Definition at line 388 of file G4PolyPhiFace.cc.

{
  G4double normSign = outgoing ? +1 : -1;
  
  //
  // These don't change
  //
  isAllBehind = allBehind;
  aNormal = normal;

  //
  // Correct normal? Here we have straight sides, and can safely ignore
  // intersections where the dot product with the normal is zero.
  //
  G4double dotProd = normSign*normal.dot(v);
  
  if (dotProd <= 0) return false;

  //
  // Calculate distance to surface. If the side is too far
  // behind the point, we must reject it.
  //
  G4ThreeVector ps = p - surface;
  distFromSurface = -normSign*ps.dot(normal);
    
  if (distFromSurface < -surfTolerance) return false;

  //
  // Calculate precise distance to intersection with the side
  // (along the trajectory, not normal to the surface)
  //
  distance = distFromSurface/dotProd;

  //
  // Calculate intersection point in r,z
  //
  G4ThreeVector ip = p + distance*v;
  
  G4double r = radial.dot(ip);
  
  //
  // And is it inside the r/z extent?
  //
  return InsideEdgesExact( r, ip.z(), normSign, p, v );
}

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4PolyPhiFace::Intersect ( G4TwoVector  a, G4TwoVector  b, G4TwoVector  c, G4TwoVector  d  )

protected

Definition at line 1038 of file G4PolyPhiFace.cc.

{
 if( IntersectProp(a,b,c,d) )
   { return true; }
 else if( Between(a,b,c)||
          Between(a,b,d)||
          Between(c,d,a)||
          Between(c,d,b) )
   { return true; }
 else
   { return false; }
}

Here is the call graph for this function:

G4bool G4PolyPhiFace::IntersectProp ( G4TwoVector  a, G4TwoVector  b, G4TwoVector  c, G4TwoVector  d  )

protected

Definition at line 1002 of file G4PolyPhiFace.cc.

{
  if( Collinear(a,b,c) || Collinear(a,b,d)||
      Collinear(c,d,a) || Collinear(c,d,b) )  { return false; }

  G4bool Positive;
  Positive = !(Left(a,b,c))^!(Left(a,b,d));
  return Positive && (!Left(c,d,a)^!Left(c,d,b));
}

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4PolyPhiFace::Left ( G4TwoVector  a, G4TwoVector  b, G4TwoVector  c  )

protected

Definition at line 971 of file G4PolyPhiFace.cc.

{
  return Area2(a,b,c)>0;
}

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4PolyPhiFace::LeftOn ( G4TwoVector  a, G4TwoVector  b, G4TwoVector  c  )

protected

Definition at line 981 of file G4PolyPhiFace.cc.

{
  return Area2(a,b,c)>=0;
}

Here is the call graph for this function:

Here is the caller graph for this function:

G4ThreeVector G4PolyPhiFace::Normal ( const G4ThreeVector &  p, G4double *  bestDistance  )

virtual

Implements G4VCSGface.

Definition at line 563 of file G4PolyPhiFace.cc.

{
  //
  // Get distance along phi, which if negative means the point
  // is nominally inside the shape.
  //
  G4double distPhi = normal.dot(p);

  //
  // Calculate projected point in r,z
  //
  G4double r = radial.dot(p);
  
  //
  // Are we inside the face?
  //
  G4double distRZ2;
  
  if (InsideEdges( r, p.z(), &distRZ2, 0 ))
  {
    //
    // Yup, answer is just distPhi
    //
    *bestDistance = std::fabs(distPhi);
  }
  else
  {
    //
    // Nope. Penalize by distance out
    //
    *bestDistance = std::sqrt( distPhi*distPhi + distRZ2 );
  }
  
  return normal;
}

Here is the call graph for this function:

G4PolyPhiFace & G4PolyPhiFace::operator=

(

const G4PolyPhiFace &

source

)

Definition at line 321 of file G4PolyPhiFace.cc.

{
  if (this == &source)  { return *this; }

  delete [] edges;
  delete [] corners;
  
  CopyStuff( source );
  
  return *this;
}

Here is the call graph for this function:

G4double G4PolyPhiFace::SurfaceArea ( )

virtual

Implements G4VCSGface.

Definition at line 938 of file G4PolyPhiFace.cc.

{
  if ( fSurfaceArea==0. ) { Triangulate(); }
  return fSurfaceArea;
}

Here is the call graph for this function:

G4double G4PolyPhiFace::SurfaceTriangle	(	G4ThreeVector	p1,
		G4ThreeVector	p2,
		G4ThreeVector	p3,
		G4ThreeVector *	p4
	)

Definition at line 919 of file G4PolyPhiFace.cc.

{
  G4ThreeVector v, w;
  
  v = p3 - p1;
  w = p1 - p2;
  G4double lambda1 = G4UniformRand();
  G4double lambda2 = lambda1*G4UniformRand();
 
  *p4=p2 + lambda1*w + lambda2*v;
  return 0.5*(v.cross(w)).mag();
}

Here is the caller graph for this function:

void G4PolyPhiFace::Triangulate ( )

protected

Definition at line 1151 of file G4PolyPhiFace.cc.

{ 
  // The copy of Polycone is made and this copy is reordered in order to 
  // have a list of triangles. This list is used for GetPointOnFace().

  G4PolyPhiFaceVertex *tri_help = new G4PolyPhiFaceVertex[numEdges];
  triangles = tri_help;
  G4PolyPhiFaceVertex *triang = triangles;

  std::vector<G4double> areas;
  std::vector<G4ThreeVector> points;
  G4double area=0.;
  G4PolyPhiFaceVertex *v0,*v1,*v2,*v3,*v4;
  v2=triangles;

  // Make copy for prev/next for triang=corners
  //
  G4PolyPhiFaceVertex *helper = corners;
  G4PolyPhiFaceVertex *helper2 = corners;
  do    // Loop checking, 13.08.2015, G.Cosmo
  {
    triang->r = helper->r;
    triang->z = helper->z;
    triang->x = helper->x;
    triang->y= helper->y;

    // add pointer on prev corner
    //
    if( helper==corners )
      { triang->prev=triangles+numEdges-1; }
    else
      { triang->prev=helper2; }

    // add pointer on next corner
    //
    if( helper<corners+numEdges-1 )
      { triang->next=triang+1; }
    else
      { triang->next=triangles; }
    helper2=triang;
    helper=helper->next;
    triang=triang->next;

  } while( helper!=corners );

  EarInit();
 
  G4int n=numEdges;
  G4int i=0;
  G4ThreeVector p1,p2,p3,p4;
  const G4int max_n_loops=numEdges*10000; // protection against infinite loop

  // Each step of outer loop removes one ear
  //
  while(n>3)      // Loop checking, 13.08.2015, G.Cosmo
  {               // Inner loop searches for one ear
    v2=triangles; 
    do    // Loop checking, 13.08.2015, G.Cosmo
    {
      if(v2->ear)  // Ear found. Fill variables
      {
        // (v1,v3) is diagonal
        //
        v3=v2->next; v4=v3->next;
        v1=v2->prev; v0=v1->prev;
        
        // Calculate areas and points

        p1=G4ThreeVector((v2)->x,(v2)->y,(v2)->z);
        p2=G4ThreeVector((v1)->x,(v1)->y,(v1)->z);
        p3=G4ThreeVector((v3)->x,(v3)->y,(v3)->z);

        G4double result1 = SurfaceTriangle(p1,p2,p3,&p4 );
        points.push_back(p4);
        areas.push_back(result1);
        area=area+result1;

        // Update earity of diagonal endpoints
        //
        v1->ear=Diagonal(v0,v3);
        v3->ear=Diagonal(v1,v4);

        // Cut off the ear v2 
        // Has to be done for a copy and not for real PolyPhiFace
        //
        v1->next=v3;
        v3->prev=v1;
        triangles=v3; // In case the head was v2
        n--;
 
        break; // out of inner loop
      }        // end if ear found

      v2=v2->next;
    
    } while( v2!=triangles );

    i++;
    if(i>=max_n_loops)
    {
      G4Exception( "G4PolyPhiFace::Triangulation()",
                   "GeomSolids0003", FatalException,
                   "Maximum number of steps is reached for triangulation!" );
    }
  }   // end outer while loop

  if(v2->next)
  {
     // add last triangle
     //
     v2=v2->next;
     p1=G4ThreeVector((v2)->x,(v2)->y,(v2)->z);
     p2=G4ThreeVector((v2->next)->x,(v2->next)->y,(v2->next)->z);
     p3=G4ThreeVector((v2->prev)->x,(v2->prev)->y,(v2->prev)->z);
     G4double result1 = SurfaceTriangle(p1,p2,p3,&p4 );
     points.push_back(p4);
     areas.push_back(result1);
     area=area+result1;
  }
 
  // Surface Area is stored
  //
  fSurfaceArea = area;
  
  // Second Step: choose randomly one surface
  //
  G4double chose = area*G4UniformRand();
   
  // Third Step: Get a point on choosen surface
  //
  G4double Achose1, Achose2;
  Achose1=0; Achose2=0.; 
  i=0;
  do     // Loop checking, 13.08.2015, G.Cosmo
  {
    Achose2+=areas[i];
    if(chose>=Achose1 && chose<Achose2)
    {
      G4ThreeVector point;
       point=points[i] ;
       surface_point=point;
       break;     
    }
    i++; Achose1=Achose2;
  } while( i<numEdges-2 );
 
  delete [] tri_help;
  tri_help = 0;
}

Here is the call graph for this function:

Here is the caller graph for this function:

Member Data Documentation

G4bool G4PolyPhiFace::allBehind

protected

Definition at line 231 of file G4PolyPhiFace.hh.

G4PolyPhiFaceVertex* G4PolyPhiFace::corners

protected

Definition at line 223 of file G4PolyPhiFace.hh.

G4PolyPhiFaceEdge* G4PolyPhiFace::edges

protected

Definition at line 222 of file G4PolyPhiFace.hh.

G4double G4PolyPhiFace::fSurfaceArea

protected

Definition at line 234 of file G4PolyPhiFace.hh.

G4double G4PolyPhiFace::kCarTolerance

protected

Definition at line 233 of file G4PolyPhiFace.hh.

G4ThreeVector G4PolyPhiFace::normal

protected

Definition at line 224 of file G4PolyPhiFace.hh.

G4int G4PolyPhiFace::numEdges

protected

Definition at line 221 of file G4PolyPhiFace.hh.

G4ThreeVector G4PolyPhiFace::radial

protected

Definition at line 225 of file G4PolyPhiFace.hh.

G4double G4PolyPhiFace::rMax

protected

Definition at line 229 of file G4PolyPhiFace.hh.

G4double G4PolyPhiFace::rMin

protected

Definition at line 229 of file G4PolyPhiFace.hh.

G4ThreeVector G4PolyPhiFace::surface

protected

Definition at line 226 of file G4PolyPhiFace.hh.

G4ThreeVector G4PolyPhiFace::surface_point

protected

Definition at line 227 of file G4PolyPhiFace.hh.

G4PolyPhiFaceVertex* G4PolyPhiFace::triangles

protected

Definition at line 235 of file G4PolyPhiFace.hh.

G4double G4PolyPhiFace::zMax

protected

Definition at line 229 of file G4PolyPhiFace.hh.

G4double G4PolyPhiFace::zMin

protected

Definition at line 229 of file G4PolyPhiFace.hh.

---

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

• geant4.10.03.p01/source/geometry/solids/specific/include/G4PolyPhiFace.hh
• geant4.10.03.p01/source/geometry/solids/specific/src/G4PolyPhiFace.cc