Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4TwistTubsHypeSide Class Reference

#include <G4TwistTubsHypeSide.hh>

Inheritance diagram for G4TwistTubsHypeSide:
Collaboration 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)
 
void DebugPrint () const
 
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 *axis0min, G4VTwistSurface *axis1min, G4VTwistSurface *axis0max, G4VTwistSurface *axis1max)
 
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 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 54 of file G4TwistTubsHypeSide.hh.

Constructor & Destructor Documentation

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 51 of file G4TwistTubsHypeSide.cc.

64  : G4VTwistSurface(name, rot, tlate, handedness, axis0, axis1,
65  axis0min, axis1min, axis0max, axis1max),
66  fKappa(kappa), fTanStereo(tanstereo),
67  fTan2Stereo(tanstereo*tanstereo), fR0(r0), fR02(r0*r0), fDPhi(twopi)
68 {
69  if ( (axis0 == kZAxis) && (axis1 == kPhi) )
70  {
71  G4Exception("G4TwistTubsHypeSide::G4TwistTubsHypeSide()",
72  "GeomSolids0002", FatalErrorInArgument,
73  "Should swap axis0 and axis1!");
74  }
75 
76  fInside.gp.set(kInfinity, kInfinity, kInfinity);
77  fInside.inside = kOutside;
78  fIsValidNorm = false;
79 
80  SetCorners();
81  SetBoundaries();
82 
83 }
Definition: geomdefs.hh:54
static const G4double kInfinity
Definition: geomdefs.hh:42
static constexpr double twopi
Definition: G4SIunits.hh:76
G4VTwistSurface(const G4String &name)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41

Here is the call graph for this function:

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 85 of file G4TwistTubsHypeSide.cc.

97  : G4VTwistSurface(name)
98 {
99 
100  fHandedness = handedness; // +z = +ve, -z = -ve
101  fAxis[0] = kPhi;
102  fAxis[1] = kZAxis;
103  fAxisMin[0] = kInfinity; // we cannot fix boundary min of Phi,
104  fAxisMax[0] = kInfinity; // because it depends on z.
105  fAxisMin[1] = EndZ[0];
106  fAxisMax[1] = EndZ[1];
107  fKappa = Kappa;
108  fDPhi = DPhi ;
109 
110  if (handedness < 0) { // inner hyperbolic surface
111  fTanStereo = TanInnerStereo;
112  fR0 = InnerRadius;
113  } else { // outer hyperbolic surface
114  fTanStereo = TanOuterStereo;
115  fR0 = OuterRadius;
116  }
117  fTan2Stereo = fTanStereo * fTanStereo;
118  fR02 = fR0 * fR0;
119 
120  fTrans.set(0, 0, 0);
121  fIsValidNorm = false;
122 
123  fInside.gp.set(kInfinity, kInfinity, kInfinity);
124  fInside.inside = kOutside;
125 
126  SetCorners(EndInnerRadius, EndOuterRadius, DPhi, EndPhi, EndZ) ;
127 
128  SetBoundaries();
129 }
void set(double x, double y, double z)
Definition: geomdefs.hh:54
static const G4double kInfinity
Definition: geomdefs.hh:42
G4double fAxisMax[2]
G4ThreeVector fTrans
G4double fAxisMin[2]
G4VTwistSurface(const G4String &name)

Here is the call graph for this function:

G4TwistTubsHypeSide::~G4TwistTubsHypeSide ( )
virtual

Definition at line 143 of file G4TwistTubsHypeSide.cc.

144 {
145 }
G4TwistTubsHypeSide::G4TwistTubsHypeSide ( __void__ &  a)

Definition at line 134 of file G4TwistTubsHypeSide.cc.

135  : G4VTwistSurface(a), fKappa(0.), fTanStereo(0.), fTan2Stereo(0.),
136  fR0(0.), fR02(0.), fDPhi(0.)
137 {
138 }
G4VTwistSurface(const G4String &name)

Member Function Documentation

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 243 of file G4TwistTubsHypeSide.cc.

250 {
251  //
252  // Decide if and where a line intersects with a hyperbolic
253  // surface (of infinite extent)
254  //
255  // Arguments:
256  // p - (in) Point on trajectory
257  // v - (in) Vector along trajectory
258  // r2 - (in) Square of radius at z = 0
259  // tan2phi - (in) std::tan(stereo)**2
260  // s - (out) Up to two points of intersection, where the
261  // intersection point is p + s*v, and if there are
262  // two intersections, s[0] < s[1]. May be negative.
263  // Returns:
264  // The number of intersections. If 0, the trajectory misses.
265  //
266  //
267  // Equation of a line:
268  //
269  // x = x0 + s*tx y = y0 + s*ty z = z0 + s*tz
270  //
271  // Equation of a hyperbolic surface:
272  //
273  // x**2 + y**2 = r**2 + (z*tanPhi)**2
274  //
275  // Solution is quadratic:
276  //
277  // a*s**2 + b*s + c = 0
278  //
279  // where:
280  //
281  // a = tx**2 + ty**2 - (tz*tanPhi)**2
282  //
283  // b = 2*( x0*tx + y0*ty - z0*tz*tanPhi**2 )
284  //
285  // c = x0**2 + y0**2 - r**2 - (z0*tanPhi)**2
286  //
287 
288  fCurStatWithV.ResetfDone(validate, &gp, &gv);
289 
290  if (fCurStatWithV.IsDone()) {
291  G4int i;
292  for (i=0; i<fCurStatWithV.GetNXX(); i++) {
293  gxx[i] = fCurStatWithV.GetXX(i);
294  distance[i] = fCurStatWithV.GetDistance(i);
295  areacode[i] = fCurStatWithV.GetAreacode(i);
296  isvalid[i] = fCurStatWithV.IsValid(i);
297  }
298  return fCurStatWithV.GetNXX();
299  } else {
300  // initialize
301  G4int i;
302  for (i=0; i<2; i++) {
303  distance[i] = kInfinity;
304  areacode[i] = sOutside;
305  isvalid[i] = false;
306  gxx[i].set(kInfinity, kInfinity, kInfinity);
307  }
308  }
309 
312  G4ThreeVector xx[2];
313 
314  //
315  // special case! p is on origin.
316  //
317 
318  if (p.mag() == 0) {
319  // p is origin.
320  // unique solution of 2-dimension question in r-z plane
321  // Equations:
322  // r^2 = fR02 + z^2*fTan2Stere0
323  // r = beta*z
324  // where
325  // beta = vrho / vz
326  // Solution (z value of intersection point):
327  // xxz = +- std::sqrt (fR02 / (beta^2 - fTan2Stereo))
328  //
329 
330  G4double vz = v.z();
331  G4double absvz = std::fabs(vz);
332  G4double vrho = v.getRho();
333  G4double vslope = vrho/vz;
334  G4double vslope2 = vslope * vslope;
335  if (vrho == 0 || (vrho/absvz) <= (absvz*std::fabs(fTanStereo)/absvz)) {
336  // vz/vrho is bigger than slope of asymptonic line
337  distance[0] = kInfinity;
338  fCurStatWithV.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
339  isvalid[0], 0, validate, &gp, &gv);
340  return 0;
341  }
342 
343  if (vz) {
344  G4double xxz = std::sqrt(fR02 / (vslope2 - fTan2Stereo))
345  * (vz / std::fabs(vz)) ;
346  G4double t = xxz / vz;
347  xx[0].set(t*v.x(), t*v.y(), xxz);
348  } else {
349  // p.z = 0 && v.z =0
350  xx[0].set(v.x()*fR0, v.y()*fR0, 0); // v is a unit vector.
351  }
352  distance[0] = xx[0].mag();
353  gxx[0] = ComputeGlobalPoint(xx[0]);
354 
355  if (validate == kValidateWithTol) {
356  areacode[0] = GetAreaCode(xx[0]);
357  if (!IsOutside(areacode[0])) {
358  if (distance[0] >= 0) isvalid[0] = true;
359  }
360  } else if (validate == kValidateWithoutTol) {
361  areacode[0] = GetAreaCode(xx[0], false);
362  if (IsInside(areacode[0])) {
363  if (distance[0] >= 0) isvalid[0] = true;
364  }
365  } else { // kDontValidate
366  areacode[0] = sInside;
367  if (distance[0] >= 0) isvalid[0] = true;
368  }
369 
370  fCurStatWithV.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
371  isvalid[0], 1, validate, &gp, &gv);
372  return 1;
373  }
374 
375  //
376  // special case end.
377  //
378 
379  G4double a = v.x()*v.x() + v.y()*v.y() - v.z()*v.z()*fTan2Stereo;
380  G4double b = 2.0 * ( p.x() * v.x() + p.y() * v.y() - p.z() * v.z() * fTan2Stereo );
381  G4double c = p.x()*p.x() + p.y()*p.y() - fR02 - p.z()*p.z()*fTan2Stereo;
382  G4double D = b*b - 4*a*c; //discriminant
383  G4int vout = 0;
384 
385  if (std::fabs(a) < DBL_MIN) {
386  if (std::fabs(b) > DBL_MIN) { // single solution
387 
388  distance[0] = -c/b;
389  xx[0] = p + distance[0]*v;
390  gxx[0] = ComputeGlobalPoint(xx[0]);
391 
392  if (validate == kValidateWithTol) {
393  areacode[0] = GetAreaCode(xx[0]);
394  if (!IsOutside(areacode[0])) {
395  if (distance[0] >= 0) isvalid[0] = true;
396  }
397  } else if (validate == kValidateWithoutTol) {
398  areacode[0] = GetAreaCode(xx[0], false);
399  if (IsInside(areacode[0])) {
400  if (distance[0] >= 0) isvalid[0] = true;
401  }
402  } else { // kDontValidate
403  areacode[0] = sInside;
404  if (distance[0] >= 0) isvalid[0] = true;
405  }
406 
407  fCurStatWithV.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
408  isvalid[0], 1, validate, &gp, &gv);
409  vout = 1;
410 
411  } else {
412  // if a=b=0 and c != 0, p is origin and v is parallel to asymptotic line.
413  // if a=b=c=0, p is on surface and v is paralell to stereo wire.
414  // return distance = infinity.
415 
416  fCurStatWithV.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
417  isvalid[0], 0, validate, &gp, &gv);
418 
419  vout = 0;
420  }
421 
422  } else if (D > DBL_MIN) { // double solutions
423 
424  D = std::sqrt(D);
425  G4double factor = 0.5/a;
426  G4double tmpdist[2] = {kInfinity, kInfinity};
427  G4ThreeVector tmpxx[2] ;
428  G4int tmpareacode[2] = {sOutside, sOutside};
429  G4bool tmpisvalid[2] = {false, false};
430  G4int i;
431 
432  for (i=0; i<2; i++) {
433  tmpdist[i] = factor*(-b - D);
434  D = -D;
435  tmpxx[i] = p + tmpdist[i]*v;
436 
437  if (validate == kValidateWithTol) {
438  tmpareacode[i] = GetAreaCode(tmpxx[i]);
439  if (!IsOutside(tmpareacode[i])) {
440  if (tmpdist[i] >= 0) tmpisvalid[i] = true;
441  continue;
442  }
443  } else if (validate == kValidateWithoutTol) {
444  tmpareacode[i] = GetAreaCode(tmpxx[i], false);
445  if (IsInside(tmpareacode[i])) {
446  if (tmpdist[i] >= 0) tmpisvalid[i] = true;
447  continue;
448  }
449  } else { // kDontValidate
450  tmpareacode[i] = sInside;
451  if (tmpdist[i] >= 0) tmpisvalid[i] = true;
452  continue;
453  }
454  }
455 
456  if (tmpdist[0] <= tmpdist[1]) {
457  distance[0] = tmpdist[0];
458  distance[1] = tmpdist[1];
459  xx[0] = tmpxx[0];
460  xx[1] = tmpxx[1];
461  gxx[0] = ComputeGlobalPoint(tmpxx[0]);
462  gxx[1] = ComputeGlobalPoint(tmpxx[1]);
463  areacode[0] = tmpareacode[0];
464  areacode[1] = tmpareacode[1];
465  isvalid[0] = tmpisvalid[0];
466  isvalid[1] = tmpisvalid[1];
467  } else {
468  distance[0] = tmpdist[1];
469  distance[1] = tmpdist[0];
470  xx[0] = tmpxx[1];
471  xx[1] = tmpxx[0];
472  gxx[0] = ComputeGlobalPoint(tmpxx[1]);
473  gxx[1] = ComputeGlobalPoint(tmpxx[0]);
474  areacode[0] = tmpareacode[1];
475  areacode[1] = tmpareacode[0];
476  isvalid[0] = tmpisvalid[1];
477  isvalid[1] = tmpisvalid[0];
478  }
479 
480  fCurStatWithV.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
481  isvalid[0], 2, validate, &gp, &gv);
482  fCurStatWithV.SetCurrentStatus(1, gxx[1], distance[1], areacode[1],
483  isvalid[1], 2, validate, &gp, &gv);
484  vout = 2;
485 
486  } else {
487  // if D<0, no solution
488  // if D=0, just grazing the surfaces, return kInfinity
489 
490  fCurStatWithV.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
491  isvalid[0], 0, validate, &gp, &gv);
492  vout = 0;
493  }
494  return vout;
495 }
void set(double x, double y, double z)
G4int GetAreacode(G4int i) const
static const G4double kInfinity
Definition: geomdefs.hh:42
double x() const
G4ThreeVector ComputeLocalDirection(const G4ThreeVector &gp) const
const char * p
Definition: xmltok.h:285
static const G4int sOutside
double getRho() const
G4ThreeVector ComputeGlobalPoint(const G4ThreeVector &lp) const
int G4int
Definition: G4Types.hh:78
double z() const
G4bool IsOutside(G4int areacode) const
G4double GetDistance(G4int i) const
bool G4bool
Definition: G4Types.hh:79
void SetCurrentStatus(G4int i, G4ThreeVector &xx, G4double &dist, G4int &areacode, G4bool &isvalid, G4int nxx, EValidate validate, const G4ThreeVector *p, const G4ThreeVector *v=0)
G4bool IsValid(G4int i) const
G4bool IsInside(G4int areacode, G4bool testbitmode=false) const
static const G4int sInside
double y() const
G4ThreeVector GetXX(G4int i) const
#define DBL_MIN
Definition: templates.hh:75
double D(double temp)
G4ThreeVector ComputeLocalPoint(const G4ThreeVector &gp) const
double G4double
Definition: G4Types.hh:76
double mag() const
CurrentStatus fCurStatWithV
void ResetfDone(EValidate validate, const G4ThreeVector *p, const G4ThreeVector *v=0)

Here is the call graph for this function:

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

Implements G4VTwistSurface.

Definition at line 501 of file G4TwistTubsHypeSide.cc.

505 {
506  // Find the approximate distance of a point of a hyperbolic surface.
507  // The distance must be an underestimate.
508  // It will also be nice (although not necessary) that the estimate is
509  // always finite no matter how close the point is.
510  //
511  // We arranged G4Hype::ApproxDistOutside and G4Hype::ApproxDistInside
512  // for this function. See these discriptions.
513 
514  const G4double halftol
516 
518 
519  if (fCurStat.IsDone()) {
520  for (G4int i=0; i<fCurStat.GetNXX(); i++) {
521  gxx[i] = fCurStat.GetXX(i);
522  distance[i] = fCurStat.GetDistance(i);
523  areacode[i] = fCurStat.GetAreacode(i);
524  }
525  return fCurStat.GetNXX();
526  } else {
527  // initialize
528  for (G4int i=0; i<2; i++) {
529  distance[i] = kInfinity;
530  areacode[i] = sOutside;
531  gxx[i].set(kInfinity, kInfinity, kInfinity);
532  }
533  }
534 
535 
537  G4ThreeVector xx;
538 
539  //
540  // special case!
541  // If p is on surface, return distance = 0 immediatery .
542  //
543  G4ThreeVector lastgxx[2];
544  for (G4int i=0; i<2; i++) {
545  lastgxx[i] = fCurStatWithV.GetXX(i);
546  }
547 
548  if ((gp - lastgxx[0]).mag() < halftol || (gp - lastgxx[1]).mag() < halftol) {
549  // last winner, or last poststep point is on the surface.
550  xx = p;
551  gxx[0] = gp;
552  distance[0] = 0;
553 
554  G4bool isvalid = true;
555  fCurStat.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
556  isvalid, 1, kDontValidate, &gp);
557 
558  return 1;
559 
560  }
561  //
562  // special case end
563  //
564 
565  G4double prho = p.getRho();
566  G4double pz = std::fabs(p.z()); // use symmetry
567  G4double r1 = std::sqrt(fR02 + pz * pz * fTan2Stereo);
568 
569  G4ThreeVector pabsz(p.x(), p.y(), pz);
570 
571  if (prho > r1 + halftol) { // p is outside of Hyperbolic surface
572 
573  // First point xx1
574  G4double t = r1 / prho;
575  G4ThreeVector xx1(t * pabsz.x(), t * pabsz.y() , pz);
576 
577  // Second point xx2
578  G4double z2 = (prho * fTanStereo + pz) / (1 + fTan2Stereo);
579  G4double r2 = std::sqrt(fR02 + z2 * z2 * fTan2Stereo);
580  t = r2 / prho;
581  G4ThreeVector xx2(t * pabsz.x(), t * pabsz.y() , z2);
582 
583  G4double len = (xx2 - xx1).mag();
584  if (len < DBL_MIN) {
585  // xx2 = xx1?? I guess we
586  // must have really bracketed the normal
587  distance[0] = (pabsz - xx1).mag();
588  xx = xx1;
589  } else {
590  distance[0] = DistanceToLine(pabsz, xx1, (xx2 - xx1) , xx);
591  }
592 
593  } else if (prho < r1 - halftol) { // p is inside of Hyperbolic surface.
594 
595  // First point xx1
596  G4double t;
597  G4ThreeVector xx1;
598  if (prho < DBL_MIN) {
599  xx1.set(r1, 0. , pz);
600  } else {
601  t = r1 / prho;
602  xx1.set(t * pabsz.x(), t * pabsz.y() , pz);
603  }
604 
605  // dr, dz is tangential vector of Hyparbolic surface at xx1
606  // dr = r, dz = z*tan2stereo
607  G4double dr = pz * fTan2Stereo;
608  G4double dz = r1;
609  G4double tanbeta = dr / dz;
610  G4double pztanbeta = pz * tanbeta;
611 
612  // Second point xx2
613  // xx2 is intersection between x-axis and tangential vector
614  G4double r2 = r1 - pztanbeta;
615  G4ThreeVector xx2;
616  if (prho < DBL_MIN) {
617  xx2.set(r2, 0. , 0.);
618  } else {
619  t = r2 / prho;
620  xx2.set(t * pabsz.x(), t * pabsz.y() , 0.);
621  }
622 
623  G4ThreeVector d = xx2 - xx1;
624  distance[0] = DistanceToLine(pabsz, xx1, d, xx);
625 
626  } else { // p is on Hyperbolic surface.
627 
628  distance[0] = 0;
629  xx.set(p.x(), p.y(), pz);
630 
631  }
632 
633  if (p.z() < 0) {
634  G4ThreeVector tmpxx(xx.x(), xx.y(), -xx.z());
635  xx = tmpxx;
636  }
637 
638  gxx[0] = ComputeGlobalPoint(xx);
639  areacode[0] = sInside;
640  G4bool isvalid = true;
641  fCurStat.SetCurrentStatus(0, gxx[0], distance[0], areacode[0],
642  isvalid, 1, kDontValidate, &gp);
643  return 1;
644 }
const XML_Char int len
Definition: expat.h:262
void set(double x, double y, double z)
G4int GetAreacode(G4int i) const
static const G4double kInfinity
Definition: geomdefs.hh:42
double x() const
const char * p
Definition: xmltok.h:285
static const G4int sOutside
double getRho() const
CurrentStatus fCurStat
G4ThreeVector ComputeGlobalPoint(const G4ThreeVector &lp) const
int G4int
Definition: G4Types.hh:78
double z() const
G4double GetDistance(G4int i) const
G4double GetRadialTolerance() const
bool G4bool
Definition: G4Types.hh:79
G4double DistanceToLine(const G4ThreeVector &p, const G4ThreeVector &x0, const G4ThreeVector &d, G4ThreeVector &xx)
void SetCurrentStatus(G4int i, G4ThreeVector &xx, G4double &dist, G4int &areacode, G4bool &isvalid, G4int nxx, EValidate validate, const G4ThreeVector *p, const G4ThreeVector *v=0)
static const G4int sInside
double y() const
G4ThreeVector GetXX(G4int i) const
#define DBL_MIN
Definition: templates.hh:75
G4ThreeVector ComputeLocalPoint(const G4ThreeVector &gp) const
double G4double
Definition: G4Types.hh:76
CurrentStatus fCurStatWithV
static G4GeometryTolerance * GetInstance()
void ResetfDone(EValidate validate, const G4ThreeVector *p, const G4ThreeVector *v=0)

Here is the call graph for this function:

G4double G4TwistTubsHypeSide::GetBoundaryMax ( G4double  phi)
inlinevirtual

Implements G4VTwistSurface.

Definition at line 195 of file G4TwistTubsHypeSide.hh.

196 {
197  G4ThreeVector ptmp(0,0,z) ; // temporary point with z Komponent only
198  G4ThreeVector upperlimit; // upper phi-boundary limit at z = ptmp.z()
199  upperlimit = GetBoundaryAtPZ(sAxis0 & sAxisMax, ptmp);
200  return std::atan2( upperlimit.y(), upperlimit.x() ) ;
201 }
double x() const
virtual G4ThreeVector GetBoundaryAtPZ(G4int areacode, const G4ThreeVector &p) const
static const G4int sAxis0
static const G4int sAxisMax
double y() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4TwistTubsHypeSide::GetBoundaryMin ( G4double  phi)
inlinevirtual

Implements G4VTwistSurface.

Definition at line 186 of file G4TwistTubsHypeSide.hh.

187 {
188  G4ThreeVector ptmp(0,0,z) ; // temporary point with z Komponent only
189  G4ThreeVector lowerlimit; // lower phi-boundary limit at z = ptmp.z()
190  lowerlimit = GetBoundaryAtPZ(sAxis0 & sAxisMin, ptmp);
191  return std::atan2( lowerlimit.y(), lowerlimit.x() ) ;
192 }
double x() const
virtual G4ThreeVector GetBoundaryAtPZ(G4int areacode, const G4ThreeVector &p) const
static const G4int sAxis0
static const G4int sAxisMin
double y() const

Here is the call graph for this function:

Here is the caller graph for this function:

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

Implements G4VTwistSurface.

Definition at line 929 of file G4TwistTubsHypeSide.cc.

931 {
932 
933  G4double z ; // the two parameters for the surface equation
934  G4double x,xmin,xmax ;
935 
936  G4ThreeVector p ; // a point on the surface, given by (z,u)
937 
938  G4int nnode ;
939  G4int nface ;
940 
941  // calculate the (n-1)*(k-1) vertices
942 
943  G4int i,j ;
944 
945  for ( i = 0 ; i<n ; i++ ) {
946 
947  z = fAxisMin[1] + i*(fAxisMax[1]-fAxisMin[1])/(n-1) ;
948 
949  for ( j = 0 ; j<k ; j++ )
950  {
951  nnode = GetNode(i,j,k,n,iside) ;
952 
953  xmin = GetBoundaryMin(z) ;
954  xmax = GetBoundaryMax(z) ;
955 
956  if (fHandedness < 0) { // inner hyperbolic surface
957  x = xmin + j*(xmax-xmin)/(k-1) ;
958  } else { // outer hyperbolic surface
959  x = xmax - j*(xmax-xmin)/(k-1) ;
960  }
961 
962  p = SurfacePoint(x,z,true) ; // surface point in global coord.system
963 
964  xyz[nnode][0] = p.x() ;
965  xyz[nnode][1] = p.y() ;
966  xyz[nnode][2] = p.z() ;
967 
968  if ( i<n-1 && j<k-1 ) { // clock wise filling
969 
970  nface = GetFace(i,j,k,n,iside) ;
971 
972  faces[nface][0] = GetEdgeVisibility(i,j,k,n,0,1) * ( GetNode(i ,j ,k,n,iside)+1) ;
973  faces[nface][1] = GetEdgeVisibility(i,j,k,n,1,1) * ( GetNode(i+1,j ,k,n,iside)+1) ;
974  faces[nface][2] = GetEdgeVisibility(i,j,k,n,2,1) * ( GetNode(i+1,j+1,k,n,iside)+1) ;
975  faces[nface][3] = GetEdgeVisibility(i,j,k,n,3,1) * ( GetNode(i ,j+1,k,n,iside)+1) ;
976 
977  }
978  }
979  }
980 }
virtual G4double GetBoundaryMin(G4double phi)
double x() const
const char * p
Definition: xmltok.h:285
G4double fAxisMax[2]
virtual G4double GetBoundaryMax(G4double phi)
int G4int
Definition: G4Types.hh:78
double z() const
G4double fAxisMin[2]
G4int GetNode(G4int i, G4int j, G4int m, G4int n, G4int iside)
double y() const
G4int GetFace(G4int i, G4int j, G4int m, G4int n, G4int iside)
double G4double
Definition: G4Types.hh:76
G4int GetEdgeVisibility(G4int i, G4int j, G4int m, G4int n, G4int number, G4int orientation)
virtual G4ThreeVector SurfacePoint(G4double, G4double, G4bool isGlobal=false)

Here is the call graph for this function:

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

Implements G4VTwistSurface.

Definition at line 150 of file G4TwistTubsHypeSide.cc.

152 {
153  // GetNormal returns a normal vector at a surface (or very close
154  // to surface) point at tmpxx.
155  // If isGlobal=true, it returns the normal in global coordinate.
156  //
157 
158  G4ThreeVector xx;
159  if (isGlobal) {
160  xx = ComputeLocalPoint(tmpxx);
161  if ((xx - fCurrentNormal.p).mag() < 0.5 * kCarTolerance) {
163  }
164  } else {
165  xx = tmpxx;
166  if (xx == fCurrentNormal.p) {
167  return fCurrentNormal.normal;
168  }
169  }
170 
171  fCurrentNormal.p = xx;
172 
173  G4ThreeVector normal( xx.x(), xx.y(), -xx.z() * fTan2Stereo);
174  normal *= fHandedness;
175  normal = normal.unit();
176 
177  if (isGlobal) {
179  } else {
181  }
182  return fCurrentNormal.normal;
183 }
G4ThreeVector ComputeGlobalDirection(const G4ThreeVector &lp) const
double x() const
G4ThreeVector ComputeLocalDirection(const G4ThreeVector &gp) const
G4SurfCurNormal fCurrentNormal
double z() const
static double normal(HepRandomEngine *eptr)
Definition: RandPoisson.cc:77
double y() const
G4ThreeVector ComputeLocalPoint(const G4ThreeVector &gp) const

Here is the call graph for this function:

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

Definition at line 160 of file G4TwistTubsHypeSide.hh.

162 {
163  // Get Rho at p.z() on Hyperbolic Surface.
164  G4ThreeVector tmpp;
165  if (isglobal) {
166  tmpp = fRot.inverse()*p - fTrans;
167  } else {
168  tmpp = p;
169  }
170  return std::sqrt(fR02 + tmpp.z() * tmpp.z() * fTan2Stereo);
171 }
const char * p
Definition: xmltok.h:285
G4RotationMatrix fRot
double z() const
HepRotation inverse() const
G4ThreeVector fTrans

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4TwistTubsHypeSide::GetSurfaceArea ( )
inlinevirtual

Implements G4VTwistSurface.

Definition at line 204 of file G4TwistTubsHypeSide.hh.

205 {
206  // approximation with tube surface
207 
208  return ( fAxisMax[1] - fAxisMin[1] ) * fR0 * fDPhi ;
209 }
G4double fAxisMax[2]
G4double fAxisMin[2]
EInside G4TwistTubsHypeSide::Inside ( const G4ThreeVector gp)
virtual

Definition at line 188 of file G4TwistTubsHypeSide.cc.

189 {
190  // Inside returns
191  const G4double halftol
193 
194  if (fInside.gp == gp) {
195  return fInside.inside;
196  }
197  fInside.gp = gp;
198 
200 
201 
202  if (p.mag() < DBL_MIN) {
203  fInside.inside = kOutside;
204  return fInside.inside;
205  }
206 
207  G4double rhohype = GetRhoAtPZ(p);
208  G4double distanceToOut = fHandedness * (rhohype - p.getRho());
209  // +ve : inside
210 
211  if (distanceToOut < -halftol) {
212 
213  fInside.inside = kOutside;
214 
215  } else {
216 
217  G4int areacode = GetAreaCode(p);
218  if (IsOutside(areacode)) {
219  fInside.inside = kOutside;
220  } else if (IsBoundary(areacode)) {
221  fInside.inside = kSurface;
222  } else if (IsInside(areacode)) {
223  if (distanceToOut <= halftol) {
224  fInside.inside = kSurface;
225  } else {
226  fInside.inside = kInside;
227  }
228  } else {
229  G4cout << "WARNING - G4TwistTubsHypeSide::Inside()" << G4endl
230  << " Invalid option !" << G4endl
231  << " name, areacode, distanceToOut = "
232  << GetName() << ", " << std::hex << areacode << std::dec << ", "
233  << distanceToOut << G4endl;
234  }
235  }
236 
237  return fInside.inside;
238 }
const char * p
Definition: xmltok.h:285
G4bool IsBoundary(G4int areacode, G4bool testbitmode=false) const
double getRho() const
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
G4bool IsOutside(G4int areacode) const
G4double GetRadialTolerance() const
G4bool IsInside(G4int areacode, G4bool testbitmode=false) const
virtual G4double GetRhoAtPZ(const G4ThreeVector &p, G4bool isglobal=false) const
virtual G4String GetName() const
#define DBL_MIN
Definition: templates.hh:75
#define G4endl
Definition: G4ios.hh:61
G4ThreeVector ComputeLocalPoint(const G4ThreeVector &gp) const
double G4double
Definition: G4Types.hh:76
double mag() const
static G4GeometryTolerance * GetInstance()

Here is the call graph for this function:

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

Implements G4VTwistSurface.

Definition at line 175 of file G4TwistTubsHypeSide.hh.

176 {
177  G4double rho = std::sqrt(fR02 + z * z * fTan2Stereo) ;
178 
179  G4ThreeVector SurfPoint (rho*std::cos(phi), rho*std::sin(phi), z) ;
180 
181  if (isGlobal) { return (fRot * SurfPoint + fTrans); }
182  return SurfPoint;
183 }
G4RotationMatrix fRot
G4ThreeVector fTrans
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:


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