#include <G4TwistTubsHypeSide.hh>

Inheritance diagram for G4TwistTubsHypeSide:

Public Member Functions
	G4TwistTubsHypeSide (const G4String &name, const G4RotationMatrix &rot, const G4ThreeVector &tlate, const G4int handedness, const G4double kappa, const G4double tanstereo, const G4double r0, const EAxis axis0=kPhi, const EAxis axis1=kZAxis, G4double axis0min=-kInfinity, G4double axis1min=-kInfinity, G4double axis0max=kInfinity, G4double axis1max=kInfinity)

	G4TwistTubsHypeSide (const G4String &name, G4double EndInnerRadius[2], G4double EndOuterRadius[2], G4double DPhi, G4double EndPhi[2], G4double EndZ[2], G4double InnerRadius, G4double OuterRadius, G4double Kappa, G4double TanInnerStereo, G4double TanOuterStereo, G4int handedness)

virtual	~G4TwistTubsHypeSide ()

virtual G4int	DistanceToSurface (const G4ThreeVector &gp, const G4ThreeVector &gv, G4ThreeVector gxx[], G4double distance[], G4int areacode[], G4bool isvalid[], EValidate validate=kValidateWithTol)

virtual G4int	DistanceToSurface (const G4ThreeVector &gp, G4ThreeVector gxx[], G4double distance[], G4int areacode[])

virtual G4ThreeVector	GetNormal (const G4ThreeVector &xx, G4bool isGlobal=false)

virtual EInside	Inside (const G4ThreeVector &gp)

virtual G4double	GetRhoAtPZ (const G4ThreeVector &p, G4bool isglobal=false) const

virtual G4ThreeVector	SurfacePoint (G4double, G4double, G4bool isGlobal=false)

virtual G4double	GetBoundaryMin (G4double phi)

virtual G4double	GetBoundaryMax (G4double phi)

virtual G4double	GetSurfaceArea ()

virtual void	GetFacets (G4int m, G4int n, G4double xyz[][3], G4int faces[][4], G4int iside)

	G4TwistTubsHypeSide (__void__ &)

Public Member Functions inherited from G4VTwistSurface
	G4VTwistSurface (const G4String &name)

	G4VTwistSurface (const G4String &name, const G4RotationMatrix &rot, const G4ThreeVector &tlate, G4int handedness, const EAxis axis1, const EAxis axis2, G4double axis0min=-kInfinity, G4double axis1min=-kInfinity, G4double axis0max=kInfinity, G4double axis1max=kInfinity)

virtual	~G4VTwistSurface ()

virtual G4int	AmIOnLeftSide (const G4ThreeVector &me, const G4ThreeVector &vec, G4bool withTol=true)

virtual G4double	DistanceToBoundary (G4int areacode, G4ThreeVector &xx, const G4ThreeVector &p)

virtual G4double	DistanceToIn (const G4ThreeVector &gp, const G4ThreeVector &gv, G4ThreeVector &gxxbest)

virtual G4double	DistanceToOut (const G4ThreeVector &gp, const G4ThreeVector &gv, G4ThreeVector &gxxbest)

virtual G4double	DistanceTo (const G4ThreeVector &gp, G4ThreeVector &gxx)

virtual G4int	DistanceToSurface (const G4ThreeVector &gp, const G4ThreeVector &gv, G4ThreeVector gxx[], G4double distance[], G4int areacode[], G4bool isvalid[], EValidate validate=kValidateWithTol)=0

virtual G4int	DistanceToSurface (const G4ThreeVector &gp, G4ThreeVector gxx[], G4double distance[], G4int areacode[])=0

void	DebugPrint () const

virtual G4ThreeVector	GetNormal (const G4ThreeVector &xx, G4bool isGlobal)=0

virtual G4String	GetName () const

virtual void	GetBoundaryParameters (const G4int &areacode, G4ThreeVector &d, G4ThreeVector &x0, G4int &boundarytype) const

virtual G4ThreeVector	GetBoundaryAtPZ (G4int areacode, const G4ThreeVector &p) const

G4double	DistanceToPlaneWithV (const G4ThreeVector &p, const G4ThreeVector &v, const G4ThreeVector &x0, const G4ThreeVector &n0, G4ThreeVector &xx)

G4double	DistanceToPlane (const G4ThreeVector &p, const G4ThreeVector &x0, const G4ThreeVector &n0, G4ThreeVector &xx)

G4double	DistanceToPlane (const G4ThreeVector &p, const G4ThreeVector &x0, const G4ThreeVector &t1, const G4ThreeVector &t2, G4ThreeVector &xx, G4ThreeVector &n)

G4double	DistanceToLine (const G4ThreeVector &p, const G4ThreeVector &x0, const G4ThreeVector &d, G4ThreeVector &xx)

G4bool	IsAxis0 (G4int areacode) const

G4bool	IsAxis1 (G4int areacode) const

G4bool	IsOutside (G4int areacode) const

G4bool	IsInside (G4int areacode, G4bool testbitmode=false) const

G4bool	IsBoundary (G4int areacode, G4bool testbitmode=false) const

G4bool	IsCorner (G4int areacode, G4bool testbitmode=false) const

G4bool	IsValidNorm () const

G4bool	IsSameBoundary (G4VTwistSurface surface1, G4int areacode1, G4VTwistSurface surface2, G4int areacode2) const

G4int	GetAxisType (G4int areacode, G4int whichaxis) const

G4ThreeVector	ComputeGlobalPoint (const G4ThreeVector &lp) const

G4ThreeVector	ComputeLocalPoint (const G4ThreeVector &gp) const

G4ThreeVector	ComputeGlobalDirection (const G4ThreeVector &lp) const

G4ThreeVector	ComputeLocalDirection (const G4ThreeVector &gp) const

void	SetAxis (G4int i, const EAxis axis)

void	SetNeighbours (G4VTwistSurface ax0min, G4VTwistSurface ax1min, G4VTwistSurface ax0max, G4VTwistSurface ax1max)

virtual G4ThreeVector	SurfacePoint (G4double, G4double, G4bool isGlobal=false)=0

virtual G4double	GetBoundaryMin (G4double)=0

virtual G4double	GetBoundaryMax (G4double)=0

virtual G4double	GetSurfaceArea ()=0

virtual void	GetFacets (G4int m, G4int n, G4double xyz[][3], G4int faces[][4], G4int iside)=0

G4int	GetNode (G4int i, G4int j, G4int m, G4int n, G4int iside)

G4int	GetFace (G4int i, G4int j, G4int m, G4int n, G4int iside)

G4int	GetEdgeVisibility (G4int i, G4int j, G4int m, G4int n, G4int number, G4int orientation)

	G4VTwistSurface (__void__ &)

Additional Inherited Members
Public Types inherited from G4VTwistSurface
enum	EValidate { kDontValidate = 0 , kValidateWithTol = 1 , kValidateWithoutTol = 2 , kUninitialized = 3 }

Static Public Attributes inherited from G4VTwistSurface
static const G4int	sOutside = 0x00000000

static const G4int	sInside = 0x10000000

static const G4int	sBoundary = 0x20000000

static const G4int	sCorner = 0x40000000

static const G4int	sC0Min1Min = 0x40000101

static const G4int	sC0Max1Min = 0x40000201

static const G4int	sC0Max1Max = 0x40000202

static const G4int	sC0Min1Max = 0x40000102

static const G4int	sAxisMin = 0x00000101

static const G4int	sAxisMax = 0x00000202

static const G4int	sAxisX = 0x00000404

static const G4int	sAxisY = 0x00000808

static const G4int	sAxisZ = 0x00000C0C

static const G4int	sAxisRho = 0x00001010

static const G4int	sAxisPhi = 0x00001414

static const G4int	sAxis0 = 0x0000FF00

static const G4int	sAxis1 = 0x000000FF

static const G4int	sSizeMask = 0x00000303

static const G4int	sAxisMask = 0x0000FCFC

static const G4int	sAreaMask = 0XF0000000

Protected Member Functions inherited from G4VTwistSurface
G4VTwistSurface **	GetNeighbours ()

G4int	GetNeighbours (G4int areacode, G4VTwistSurface *surfaces[])

G4ThreeVector	GetCorner (G4int areacode) const

void	GetBoundaryAxis (G4int areacode, EAxis axis[]) const

void	GetBoundaryLimit (G4int areacode, G4double limit[]) const

virtual G4int	GetAreaCode (const G4ThreeVector &xx, G4bool withtol=true)=0

virtual void	SetBoundary (const G4int &axiscode, const G4ThreeVector &direction, const G4ThreeVector &x0, const G4int &boundarytype)

void	SetCorner (G4int areacode, G4double x, G4double y, G4double z)

Protected Attributes inherited from G4VTwistSurface
EAxis	fAxis [2]

G4double	fAxisMin [2]

G4double	fAxisMax [2]

CurrentStatus	fCurStatWithV

CurrentStatus	fCurStat

G4RotationMatrix	fRot

G4ThreeVector	fTrans

G4int	fHandedness

G4SurfCurNormal	fCurrentNormal

G4bool	fIsValidNorm

G4double	kCarTolerance

Detailed Description

Definition at line 43 of file G4TwistTubsHypeSide.hh.

Constructor & Destructor Documentation

◆ G4TwistTubsHypeSide() [1/3]

G4TwistTubsHypeSide::G4TwistTubsHypeSide	(	const G4String &	name,
		const G4RotationMatrix &	rot,
		const G4ThreeVector &	tlate,
		const G4int	handedness,
		const G4double	kappa,
		const G4double	tanstereo,
		const G4double	r0,
		const EAxis	axis0 = `kPhi`,
		const EAxis	axis1 = `kZAxis`,
		G4double	axis0min = `-kInfinity`,
		G4double	axis1min = `-kInfinity`,
		G4double	axis0max = `kInfinity`,
		G4double	axis1max = `kInfinity`
	)

Definition at line 40 of file G4TwistTubsHypeSide.cc.

  : G4VTwistSurface(name, rot, tlate, handedness, axis0, axis1,
                    axis0min, axis1min, axis0max, axis1max),
    fKappa(kappa), fTanStereo(tanstereo),
    fTan2Stereo(tanstereo*tanstereo), fR0(r0), fR02(r0*r0), fDPhi(twopi)
{
   if ( (axis0 == kZAxis) && (axis1 == kPhi) )
   {
      G4Exception("G4TwistTubsHypeSide::G4TwistTubsHypeSide()",
                  "GeomSolids0002", FatalErrorInArgument,
                  "Should swap axis0 and axis1!");
   }
   fInside.gp.set(kInfinity, kInfinity, kInfinity);
   fInside.inside = kOutside;
   fIsValidNorm = false;
   SetCorners();
   SetBoundaries();
}

◆ G4TwistTubsHypeSide() [2/3]

G4TwistTubsHypeSide::G4TwistTubsHypeSide	(	const G4String &	name,
		G4double	EndInnerRadius[2],
		G4double	EndOuterRadius[2],
		G4double	DPhi,
		G4double	EndPhi[2],
		G4double	EndZ[2],
		G4double	InnerRadius,
		G4double	OuterRadius,
		G4double	Kappa,
		G4double	TanInnerStereo,
		G4double	TanOuterStereo,
		G4int	handedness
	)

Definition at line 71 of file G4TwistTubsHypeSide.cc.

   : G4VTwistSurface(name)
{
 
   fHandedness = handedness;   // +z = +ve, -z = -ve
   fAxis[0]    = kPhi;
   fAxis[1]    = kZAxis;
   fAxisMin[0] = kInfinity;         // we cannot fix boundary min of Phi, 
   fAxisMax[0] = kInfinity;         // because it depends on z.
   fAxisMin[1] = EndZ[0];
   fAxisMax[1] = EndZ[1];
   fKappa      = Kappa;
   fDPhi       = DPhi ;
 
   if (handedness < 0)  // inner hyperbolic surface
   {
      fTanStereo  = TanInnerStereo;
      fR0         = InnerRadius;
   }
   else                 // outer hyperbolic surface
   {
      fTanStereo  = TanOuterStereo;
      fR0         = OuterRadius;
   }
   fTan2Stereo = fTanStereo * fTanStereo;
   fR02        = fR0 * fR0;
   
   fTrans.set(0, 0, 0);
   fIsValidNorm = false;
 
   fInside.gp.set(kInfinity, kInfinity, kInfinity);
   fInside.inside = kOutside;
   
   SetCorners(EndInnerRadius, EndOuterRadius, DPhi, EndPhi, EndZ) ; 
 
   SetBoundaries();
}

◆ ~G4TwistTubsHypeSide()

G4TwistTubsHypeSide::~G4TwistTubsHypeSide ( )

virtual

Definition at line 132 of file G4TwistTubsHypeSide.cc.

133{

134}

◆ G4TwistTubsHypeSide() [3/3]

G4TwistTubsHypeSide::G4TwistTubsHypeSide ( __void__ & a )

Definition at line 123 of file G4TwistTubsHypeSide.cc.

  : G4VTwistSurface(a), fKappa(0.), fTanStereo(0.), fTan2Stereo(0.),
    fR0(0.), fR02(0.), fDPhi(0.)
{
}

Member Function Documentation

◆ DistanceToSurface() [1/2]

G4int G4TwistTubsHypeSide::DistanceToSurface	(	const G4ThreeVector &	gp,
		const G4ThreeVector &	gv,
		G4ThreeVector	gxx[],
		G4double	distance[],
		G4int	areacode[],
		G4bool	isvalid[],
		EValidate	validate = `kValidateWithTol`
	)

virtual

Implements G4VTwistSurface.

Definition at line 251 of file G4TwistTubsHypeSide.cc.

{
   // Decide if and where a line intersects with a hyperbolic
   // surface (of infinite extent)
   //
   // Arguments:
   //     p       - (in) Point on trajectory
   //     v       - (in) Vector along trajectory
   //     r2      - (in) Square of radius at z = 0
   //     tan2phi - (in) std::tan(stereo)**2
   //     s       - (out) Up to two points of intersection, where the
   //                     intersection point is p + s*v, and if there are
   //                     two intersections, s[0] < s[1]. May be negative.
   // Returns:
   //     The number of intersections. If 0, the trajectory misses.
   //
   //
   // Equation of a line:
   //
   //       x = x0 + s*tx      y = y0 + s*ty      z = z0 + s*tz
   //
   // Equation of a hyperbolic surface:
   //
   //       x**2 + y**2 = r**2 + (z*tanPhi)**2
   //
   // Solution is quadratic:
   //
   //  a*s**2 + b*s + c = 0
   //
   // where:
   //
   //  a = tx**2 + ty**2 - (tz*tanPhi)**2
   //
   //  b = 2*( x0*tx + y0*ty - z0*tz*tanPhi**2 )
   //
   //  c = x0**2 + y0**2 - r**2 - (z0*tanPhi)**2
   //
      
   fCurStatWithV.ResetfDone(validate, &gp, &gv);
 
   if (fCurStatWithV.IsDone())
   {
      for (G4int i=0; i<fCurStatWithV.GetNXX(); ++i)
      {
         gxx[i] = fCurStatWithV.GetXX(i);
         distance[i] = fCurStatWithV.GetDistance(i);
         areacode[i] = fCurStatWithV.GetAreacode(i);
         isvalid[i]  = fCurStatWithV.IsValid(i);
      }
      return fCurStatWithV.GetNXX();
   }
   else   // initialize
   {
      for (auto i=0; i<2; ++i)
      {
         distance[i] = kInfinity;
         areacode[i] = sOutside;
         isvalid[i]  = false;
         gxx[i].set(kInfinity, kInfinity, kInfinity);
      }
   }
   
   G4ThreeVector p = ComputeLocalPoint(gp);
   G4ThreeVector v = ComputeLocalDirection(gv);
   G4ThreeVector xx[2]; 
   
   //
   // special case!  p is on origin.
   // 
 
   if (p.mag() == 0)
   {
      // p is origin. 
      // unique solution of 2-dimension question in r-z plane 
      // Equations:
      //    r^2 = fR02 + z^2*fTan2Stere0
      //    r = beta*z
      //        where 
      //        beta = vrho / vz
      // Solution (z value of intersection point):
      //    xxz = +- std::sqrt (fR02 / (beta^2 - fTan2Stereo))
      //
 
      G4double vz    = v.z();
      G4double absvz = std::fabs(vz);
      G4double vrho  = v.getRho();       
      G4double vslope = vrho/vz;
      G4double vslope2 = vslope * vslope;
      if (vrho == 0 || (vrho/absvz) <= (absvz*std::fabs(fTanStereo)/absvz))
      {
         // vz/vrho is bigger than slope of asymptonic line
         distance[0] = kInfinity;
         fCurStatWithV.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
                                        isvalid[0], 0, validate, &gp, &gv);
         return 0;
      }
       
      if (vz)
      { 
         G4double xxz  = std::sqrt(fR02 / (vslope2 - fTan2Stereo)) 
                        * (vz / std::fabs(vz)) ;
         G4double t = xxz / vz;
         xx[0].set(t*v.x(), t*v.y(), xxz);
      }
      else
      {
         // p.z = 0 && v.z =0
         xx[0].set(v.x()*fR0, v.y()*fR0, 0);  // v is a unit vector.
      }
      distance[0] = xx[0].mag();
      gxx[0]      = ComputeGlobalPoint(xx[0]);
 
      if (validate == kValidateWithTol)
      {
         areacode[0] = GetAreaCode(xx[0]);
         if (!IsOutside(areacode[0]))
         {
            if (distance[0] >= 0) isvalid[0] = true;
         }
      }
      else if (validate == kValidateWithoutTol)
      {
         areacode[0] = GetAreaCode(xx[0], false);
         if (IsInside(areacode[0]))
         {
            if (distance[0] >= 0) isvalid[0] = true;
         }
      }
      else  // kDontValidate                       
      {
         areacode[0] = sInside;
            if (distance[0] >= 0) isvalid[0] = true;
      }
                 
      fCurStatWithV.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
                                     isvalid[0], 1, validate, &gp, &gv);
      return 1;
   }
 
   //
   // special case end.
   // 
 
   G4double a = v.x()*v.x() + v.y()*v.y() - v.z()*v.z()*fTan2Stereo;
   G4double b = 2.0
              * ( p.x() * v.x() + p.y() * v.y() - p.z() * v.z() * fTan2Stereo );
   G4double c = p.x()*p.x() + p.y()*p.y() - fR02 - p.z()*p.z()*fTan2Stereo;
   G4double D = b*b - 4*a*c;          //discriminant
   G4int vout = 0;
   
   if (std::fabs(a) < DBL_MIN)
   {
      if (std::fabs(b) > DBL_MIN)            // single solution
      {
         distance[0] = -c/b;
         xx[0] = p + distance[0]*v;
         gxx[0] = ComputeGlobalPoint(xx[0]);
 
         if (validate == kValidateWithTol)
         {
            areacode[0] = GetAreaCode(xx[0]);
            if (!IsOutside(areacode[0]))
            {
               if (distance[0] >= 0) isvalid[0] = true;
            }
         }
         else if (validate == kValidateWithoutTol)
         {
            areacode[0] = GetAreaCode(xx[0], false);
            if (IsInside(areacode[0]))
            {
               if (distance[0] >= 0) isvalid[0] = true;
            }
         }
         else  // kDontValidate                       
         {
            areacode[0] = sInside;
               if (distance[0] >= 0) isvalid[0] = true;
         }
                 
         fCurStatWithV.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
                                        isvalid[0], 1, validate, &gp, &gv);
         vout = 1;
      }
      else
      {
        // if a=b=0 and c != 0, p is origin and v is parallel to asymptotic line
        // if a=b=c=0, p is on surface and v is paralell to stereo wire.
        // return distance = infinity
 
         fCurStatWithV.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
                                        isvalid[0], 0, validate, &gp, &gv);
         vout = 0;
      }
   }
   else if (D > DBL_MIN)          // double solutions
   {
      D = std::sqrt(D);
      G4double      factor = 0.5/a;
      G4double      tmpdist[2] = {kInfinity, kInfinity};
      G4ThreeVector tmpxx[2] ;
      G4int         tmpareacode[2] = {sOutside, sOutside};
      G4bool        tmpisvalid[2]  = {false, false};
 
      for (auto i=0; i<2; ++i)
      {
         tmpdist[i] = factor*(-b - D);
         D = -D;
         tmpxx[i] = p + tmpdist[i]*v;
        
         if (validate == kValidateWithTol)
         {
            tmpareacode[i] = GetAreaCode(tmpxx[i]);
            if (!IsOutside(tmpareacode[i]))
            {
               if (tmpdist[i] >= 0) tmpisvalid[i] = true;
               continue;
            }
         }
         else if (validate == kValidateWithoutTol)
         {
            tmpareacode[i] = GetAreaCode(tmpxx[i], false);
            if (IsInside(tmpareacode[i]))
            {
               if (tmpdist[i] >= 0) tmpisvalid[i] = true;
               continue;
            }
         }
         else  // kDontValidate
         {
            tmpareacode[i] = sInside;
               if (tmpdist[i] >= 0) tmpisvalid[i] = true;
            continue;
         }
      }      
 
      if (tmpdist[0] <= tmpdist[1])
      {
          distance[0] = tmpdist[0];
          distance[1] = tmpdist[1];
          xx[0]       = tmpxx[0];
          xx[1]       = tmpxx[1];
          gxx[0]      = ComputeGlobalPoint(tmpxx[0]);
          gxx[1]      = ComputeGlobalPoint(tmpxx[1]);
          areacode[0] = tmpareacode[0];
          areacode[1] = tmpareacode[1];
          isvalid[0]  = tmpisvalid[0];
          isvalid[1]  = tmpisvalid[1];
      }
      else
      {
          distance[0] = tmpdist[1];
          distance[1] = tmpdist[0];
          xx[0]       = tmpxx[1];
          xx[1]       = tmpxx[0];
          gxx[0]      = ComputeGlobalPoint(tmpxx[1]);
          gxx[1]      = ComputeGlobalPoint(tmpxx[0]);
          areacode[0] = tmpareacode[1];
          areacode[1] = tmpareacode[0];
          isvalid[0]  = tmpisvalid[1];
          isvalid[1]  = tmpisvalid[0];
      }
         
      fCurStatWithV.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
                                     isvalid[0], 2, validate, &gp, &gv);
      fCurStatWithV.SetCurrentStatus(1, gxx[1], distance[1], areacode[1],
                                     isvalid[1], 2, validate, &gp, &gv);
      vout = 2;
   }
   else
   {
      // if D<0, no solution
      // if D=0, just grazing the surfaces, return kInfinity
 
      fCurStatWithV.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
                                     isvalid[0], 0, validate, &gp, &gv);
      vout = 0;
   }
   return vout;
}

◆ DistanceToSurface() [2/2]

G4int G4TwistTubsHypeSide::DistanceToSurface	(	const G4ThreeVector &	gp,
		G4ThreeVector	gxx[],
		G4double	distance[],
		G4int	areacode[]
	)

virtual

Implements G4VTwistSurface.

Definition at line 541 of file G4TwistTubsHypeSide.cc.

{
    // Find the approximate distance of a point of a hyperbolic surface.
    // The distance must be an underestimate. 
    // It will also be nice (although not necessary) that the estimate is
    // always finite no matter how close the point is.
    //
    // We arranged G4Hype::ApproxDistOutside and G4Hype::ApproxDistInside
    // for this function. See these discriptions.
    
   const G4double halftol
     = 0.5 * G4GeometryTolerance::GetInstance()->GetRadialTolerance();
 
   fCurStat.ResetfDone(kDontValidate, &gp);
 
   if (fCurStat.IsDone())
   {
      for (G4int i=0; i<fCurStat.GetNXX(); ++i)
      {
         gxx[i] = fCurStat.GetXX(i);
         distance[i] = fCurStat.GetDistance(i);
         areacode[i] = fCurStat.GetAreacode(i);
      }
      return fCurStat.GetNXX();
   }
   else  // initialize
   {
      for (auto i=0; i<2; ++i)
      {
         distance[i] = kInfinity;
         areacode[i] = sOutside;
         gxx[i].set(kInfinity, kInfinity, kInfinity);
      }
   }
   
 
   G4ThreeVector p = ComputeLocalPoint(gp);
   G4ThreeVector xx;
 
   //
   // special case!
   // If p is on surface, return distance = 0 immediatery .
   //
   G4ThreeVector  lastgxx[2];
   for (auto i=0; i<2; ++i)
   {
      lastgxx[i] = fCurStatWithV.GetXX(i);
   }
 
   if ((gp - lastgxx[0]).mag() < halftol || (gp - lastgxx[1]).mag() < halftol)
   {
      // last winner, or last poststep point is on the surface.
      xx = p;             
      gxx[0] = gp;
      distance[0] = 0;      
 
      G4bool isvalid = true;
      fCurStat.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
                                isvalid, 1, kDontValidate, &gp);
 
      return 1;
   }
   //
   // special case end
   //
       
   G4double prho       = p.getRho();
   G4double pz         = std::fabs(p.z());           // use symmetry
   G4double r1         = std::sqrt(fR02 + pz * pz * fTan2Stereo);
   
   G4ThreeVector pabsz(p.x(), p.y(), pz);
    
   if (prho > r1 + halftol)   // p is outside of Hyperbolic surface
   {
      // First point xx1
      G4double t = r1 / prho;
      G4ThreeVector xx1(t * pabsz.x(), t * pabsz.y() , pz);
      
      // Second point xx2
      G4double z2 = (prho * fTanStereo + pz) / (1 + fTan2Stereo);
      G4double r2 = std::sqrt(fR02 + z2 * z2 * fTan2Stereo);
      t = r2 / prho;
      G4ThreeVector xx2(t * pabsz.x(), t * pabsz.y() , z2);
            
      G4double len = (xx2 - xx1).mag();
      if (len < DBL_MIN)
      {
         // xx2 = xx1?? I guess we
         // must have really bracketed the normal
         distance[0] = (pabsz - xx1).mag();
         xx = xx1;
      }
      else
      {
         distance[0] = DistanceToLine(pabsz, xx1, (xx2 - xx1) , xx);
      }
      
   }
   else if (prho < r1 - halftol)  // p is inside of Hyperbolic surface.
   {
      // First point xx1
      G4double t;
      G4ThreeVector xx1;
      if (prho < DBL_MIN)
      {
         xx1.set(r1, 0. , pz);
      }
      else
      {
         t = r1 / prho;
         xx1.set(t * pabsz.x(), t * pabsz.y() , pz);
      }
      
      // dr, dz is tangential vector of Hyparbolic surface at xx1
      // dr = r, dz = z*tan2stereo
      G4double dr  = pz * fTan2Stereo;
      G4double dz  = r1;
      G4double tanbeta   = dr / dz;
      G4double pztanbeta = pz * tanbeta;
      
      // Second point xx2 
      // xx2 is intersection between x-axis and tangential vector
      G4double r2 = r1 - pztanbeta;
      G4ThreeVector xx2;
      if (prho < DBL_MIN)
      {
         xx2.set(r2, 0. , 0.);
      }
      else
      {
         t  = r2 / prho;
         xx2.set(t * pabsz.x(), t * pabsz.y() , 0.);
      }
      
      G4ThreeVector d = xx2 - xx1;
      distance[0] = DistanceToLine(pabsz, xx1, d, xx);
          
   }
   else   // p is on Hyperbolic surface.
   {
      distance[0] = 0;
      xx.set(p.x(), p.y(), pz);
   }
 
   if (p.z() < 0)
   {
      G4ThreeVector tmpxx(xx.x(), xx.y(), -xx.z());
      xx = tmpxx;
   }
       
   gxx[0] = ComputeGlobalPoint(xx);
   areacode[0]    = sInside;
   G4bool isvalid = true;
   fCurStat.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
                             isvalid, 1, kDontValidate, &gp);
   return 1;
}

◆ GetBoundaryMax()

G4double G4TwistTubsHypeSide::GetBoundaryMax ( G4double phi )

inlinevirtual

Implements G4VTwistSurface.

Definition at line 182 of file G4TwistTubsHypeSide.hh.

{
  G4ThreeVector ptmp(0,0,z) ;  // temporary point with z Komponent only
  G4ThreeVector upperlimit;    // upper phi-boundary limit at z = ptmp.z()
  upperlimit = GetBoundaryAtPZ(sAxis0 & sAxisMax, ptmp);
  return   std::atan2( upperlimit.y(), upperlimit.x() ) ;
}

Referenced by GetFacets().

◆ GetBoundaryMin()

G4double G4TwistTubsHypeSide::GetBoundaryMin ( G4double phi )

inlinevirtual

Implements G4VTwistSurface.

Definition at line 173 of file G4TwistTubsHypeSide.hh.

{
  G4ThreeVector ptmp(0,0,z) ;  // temporary point with z Komponent only
  G4ThreeVector lowerlimit;    // lower phi-boundary limit at z = ptmp.z()
  lowerlimit = GetBoundaryAtPZ(sAxis0 & sAxisMin, ptmp);
  return  std::atan2( lowerlimit.y(), lowerlimit.x() ) ;  
}

Referenced by GetFacets().

◆ GetFacets()

void G4TwistTubsHypeSide::GetFacets	(	G4int	m,
		G4int	n,
		G4double	xyz[][3],
		G4int	faces[][4],
		G4int	iside
	)

virtual

Implements G4VTwistSurface.

Definition at line 996 of file G4TwistTubsHypeSide.cc.

{
  G4double z ;     // the two parameters for the surface equation
  G4double x,xmin,xmax ;
 
  G4ThreeVector p ;  // a point on the surface, given by (z,u)
 
  G4int nnode ;
  G4int nface ;
 
  // calculate the (n-1)*(k-1) vertices
 
  for ( G4int i = 0 ; i<n ; ++i )
  {
    z = fAxisMin[1] + i*(fAxisMax[1]-fAxisMin[1])/(n-1) ;
 
    for ( G4int j = 0 ; j<k ; ++j )
    {
      nnode = GetNode(i,j,k,n,iside) ;
 
      xmin = GetBoundaryMin(z) ; 
      xmax = GetBoundaryMax(z) ;
 
      if (fHandedness < 0)  // inner hyperbolic surface
      {
        x = xmin + j*(xmax-xmin)/(k-1) ;
      }
      else                  // outer hyperbolic surface
      {
        x = xmax - j*(xmax-xmin)/(k-1) ;
      }
 
      p = SurfacePoint(x,z,true) ;  // surface point in global coord.system
 
      xyz[nnode][0] = p.x() ;
      xyz[nnode][1] = p.y() ;
      xyz[nnode][2] = p.z() ;
 
      if ( i<n-1 && j<k-1 )    // clock wise filling
      {
        nface = GetFace(i,j,k,n,iside) ;
    faces[nface][0] = GetEdgeVisibility(i,j,k,n,0,1)
                        * ( GetNode(i  ,j  ,k,n,iside)+1) ;  
    faces[nface][1] = GetEdgeVisibility(i,j,k,n,1,1)
                        * ( GetNode(i+1,j  ,k,n,iside)+1) ;
    faces[nface][2] = GetEdgeVisibility(i,j,k,n,2,1)
                        * ( GetNode(i+1,j+1,k,n,iside)+1) ;
    faces[nface][3] = GetEdgeVisibility(i,j,k,n,3,1)
                        * ( GetNode(i  ,j+1,k,n,iside)+1) ;
      }
    }
  }
}

◆ GetNormal()

G4ThreeVector G4TwistTubsHypeSide::GetNormal	(	const G4ThreeVector &	xx,
		G4bool	isGlobal = `false`
	)

virtual

Implements G4VTwistSurface.

Definition at line 139 of file G4TwistTubsHypeSide.cc.

{
   // GetNormal returns a normal vector at a surface (or very close
   // to surface) point at tmpxx.
   // If isGlobal=true, it returns the normal in global coordinate.
   
   G4ThreeVector xx;
   if (isGlobal)
   {
      xx = ComputeLocalPoint(tmpxx);
      if ((xx - fCurrentNormal.p).mag() < 0.5 * kCarTolerance)
      {
         return ComputeGlobalDirection(fCurrentNormal.normal);
      }
   }
   else
   {
      xx = tmpxx;
      if (xx == fCurrentNormal.p)
      {
         return fCurrentNormal.normal;
      }
   }
   
   fCurrentNormal.p = xx;
 
   G4ThreeVector normal( xx.x(), xx.y(), -xx.z() * fTan2Stereo);
   normal *= fHandedness;
   normal = normal.unit();
 
   if (isGlobal)
   {
      fCurrentNormal.normal = ComputeLocalDirection(normal);
   }
   else
   {
      fCurrentNormal.normal = normal;
   }
   return fCurrentNormal.normal;
}

◆ GetRhoAtPZ()

G4double G4TwistTubsHypeSide::GetRhoAtPZ	(	const G4ThreeVector &	p,
		G4bool	isglobal = `false`
	)		const

inlinevirtual

Definition at line 149 of file G4TwistTubsHypeSide.hh.

{
  // Get Rho at p.z() on Hyperbolic Surface.
  G4ThreeVector tmpp;
  if (isglobal) { tmpp = fRot.inverse()*p - fTrans; }
  else          { tmpp = p; }
 
  return std::sqrt(fR02 + tmpp.z() * tmpp.z() * fTan2Stereo); 
}

Referenced by Inside().

◆ GetSurfaceArea()

G4double G4TwistTubsHypeSide::GetSurfaceArea ( )

inlinevirtual

Implements G4VTwistSurface.

Definition at line 191 of file G4TwistTubsHypeSide.hh.

{
  // approximation with tube surface
 
  return ( fAxisMax[1] - fAxisMin[1] ) * fR0 * fDPhi ;
}

◆ Inside()

EInside G4TwistTubsHypeSide::Inside ( const G4ThreeVector & gp )

virtual

Definition at line 184 of file G4TwistTubsHypeSide.cc.

{
   // Inside returns 
   const G4double halftol
     = 0.5 * G4GeometryTolerance::GetInstance()->GetRadialTolerance();
 
   if (fInside.gp == gp)
   {
      return fInside.inside;
   }
   fInside.gp = gp;
   
   G4ThreeVector p = ComputeLocalPoint(gp);
  
 
   if (p.mag() < DBL_MIN)
   {
      fInside.inside = kOutside;
      return fInside.inside;
   }
   
   G4double rhohype = GetRhoAtPZ(p);
   G4double distanceToOut = fHandedness * (rhohype - p.getRho());
                            // +ve : inside
 
   if (distanceToOut < -halftol)
   {
     fInside.inside = kOutside;
   }
   else
   {
      G4int areacode = GetAreaCode(p);
      if (IsOutside(areacode))
      {
         fInside.inside = kOutside;
      }
      else if (IsBoundary(areacode))
      {
         fInside.inside = kSurface;
      }
      else if (IsInside(areacode))
      {
         if (distanceToOut <= halftol)
         {
            fInside.inside = kSurface;
         }
         else
         {
            fInside.inside = kInside;
         }
      }
      else
      {
         G4cout << "WARNING - G4TwistTubsHypeSide::Inside()" << G4endl
                << "          Invalid option !" << G4endl
                << "          name, areacode, distanceToOut = "
                << GetName() << ", " << std::hex << areacode
                << std::dec << ", " << distanceToOut << G4endl;
      }
   }
 
   return fInside.inside; 
}

◆ SurfacePoint()

G4ThreeVector G4TwistTubsHypeSide::SurfacePoint	(	G4double	phi,
		G4double	z,
		G4bool	isGlobal = `false`
	)

inlinevirtual

Implements G4VTwistSurface.

Definition at line 161 of file G4TwistTubsHypeSide.hh.

{
  G4double rho = std::sqrt(fR02 + z * z * fTan2Stereo) ;
 
  G4ThreeVector SurfPoint (rho*std::cos(phi), rho*std::sin(phi), z) ;
 
  if (isGlobal) { return (fRot * SurfPoint + fTrans); }
  return SurfPoint;
}

Referenced by GetFacets().

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

Public Member Functions

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

◆ G4TwistTubsHypeSide() [1/3]

◆ G4TwistTubsHypeSide() [2/3]

◆ ~G4TwistTubsHypeSide()

◆ G4TwistTubsHypeSide() [3/3]

Member Function Documentation

◆ DistanceToSurface() [1/2]

◆ DistanceToSurface() [2/2]

◆ GetBoundaryMax()

◆ GetBoundaryMin()

◆ GetFacets()

◆ GetNormal()

◆ GetRhoAtPZ()

◆ GetSurfaceArea()

◆ Inside()

◆ SurfacePoint()