G4SphericalSurface Class Reference

#include <G4SphericalSurface.hh>

Inheritance diagram for G4SphericalSurface:

G4Surface G4STEPEntity

Public Member Functions

 G4SphericalSurface ()
 G4SphericalSurface (const G4Vector3D &o, const G4Vector3D &xhat, const G4Vector3D &zhat, G4double r, G4double ph1, G4double ph2, G4double th1, G4double th2)
virtual ~G4SphericalSurface ()
G4int operator== (const G4SphericalSurface &s)
G4String GetEntityType () const
virtual const char * NameOf () const
virtual void PrintOn (std::ostream &os=G4cout) const
G4int Intersect (const G4Ray &)
void CalcBBox ()
void Comp (G4Vector3D &v, G4Point3D &min, G4Point3D &max)
virtual G4double HowNear (const G4Vector3D &x) const
virtual G4Vector3D SurfaceNormal (const G4Point3D &p) const
virtual G4int Inside (const G4Vector3D &x) const
virtual G4int WithinBoundary (const G4Vector3D &x) const
virtual G4double Scale () const
virtual G4double Area () const
virtual void resize (G4double r, G4double ph1, G4double ph2, G4double th1, G4double th2)
G4Vector3D GetXAxis () const
G4Vector3D GetZAxis () const
G4double GetRadius () const
G4double GetPhi1 () const
G4double GetPhi2 () const
G4double GetTheta1 () const
G4double GetTheta2 () const
virtual G4Vector3D Normal (const G4Vector3D &p) const

Protected Attributes

G4Vector3D x_axis
G4Vector3D z_axis
G4double radius
G4double phi_1
G4double phi_2
G4double theta_1
G4double theta_2

Detailed Description

Definition at line 49 of file G4SphericalSurface.hh.


Constructor & Destructor Documentation

G4SphericalSurface::G4SphericalSurface (  ) 

Definition at line 39 of file G4SphericalSurface.cc.

References x_axis, and z_axis.

00040   : G4Surface(), radius(1.0), phi_1(0.0), phi_2(2*pi),
00041     theta_1(0.0), theta_2(pi)
00042 {
00043    // default constructor
00044    // default x_axis is ( 1.0, 0.0, 0.0 ), z_axis is ( 0.0, 0.0, 1.0 ),
00045    // default radius is 1.0
00046    // default phi_1 is 0, phi_2 is 2*PI
00047    // default theta_1 is 0, theta_2 is PI
00048 
00049    x_axis = G4Vector3D( 1.0, 0.0, 0.0 );
00050    z_axis = G4Vector3D( 0.0, 0.0, 1.0 );
00051    //   OuterBoundary = new G4BREPPolyline();
00052 }

G4SphericalSurface::G4SphericalSurface ( const G4Vector3D o,
const G4Vector3D xhat,
const G4Vector3D zhat,
G4double  r,
G4double  ph1,
G4double  ph2,
G4double  th1,
G4double  th2 
)

Definition at line 55 of file G4SphericalSurface.cc.

References G4endl, G4Exception(), JustWarning, phi_1, phi_2, G4INCL::Math::pi, radius, theta_1, theta_2, x_axis, and z_axis.

00061   : G4Surface()
00062 { 
00063   // Require both x_axis and z_axis to be unit vectors
00064 
00065   G4double xhatmag = xhat.mag();
00066   if ( xhatmag != 0.0 )
00067   {
00068     x_axis = xhat * (1/ xhatmag); // this makes the x_axis a unit vector
00069   }
00070   else
00071   {
00072     std::ostringstream message;
00073     message << "x_axis has zero length." << G4endl
00074             << "Default x_axis of (1, 0, 0) is used.";    
00075     G4Exception("G4SphericalSurface::G4SphericalSurface()",
00076                 "GeomSolids1001", JustWarning, message);
00077 
00078     x_axis = G4Vector3D( 1.0, 0.0, 0.0 );
00079   }
00080 
00081   G4double zhatmag = zhat.mag();
00082   
00083   if (zhatmag != 0.0)
00084   {
00085     z_axis = zhat *(1/ zhatmag);  // this makes the z_axis a unit vector
00086   }
00087   else 
00088   {
00089     std::ostringstream message;
00090     message << "z_axis has zero length." << G4endl
00091             << "Default z_axis of (0, 0, 1) is used.";    
00092     G4Exception("G4SphericalSurface::G4SphericalSurface()",
00093                 "GeomSolids1001", JustWarning, message);
00094 
00095     z_axis = G4Vector3D( 0.0, 0.0, 1.0 );
00096   }
00097   
00098   //  Require radius to be non-negative
00099   //
00100   if ( r >= 0.0 )
00101   {
00102     radius = r;
00103   }
00104   else 
00105   {
00106     std::ostringstream message;
00107     message << "Radius cannot be less than zero." << G4endl
00108             << "Default radius of 1.0 is used.";    
00109     G4Exception("G4SphericalSurface::G4SphericalSurface()",
00110                 "GeomSolids1001", JustWarning, message);
00111 
00112     radius = 1.0;
00113   }
00114 
00115   //  Require phi_1 in the range: 0 <= phi_1 < 2*PI
00116   //  and phi_2 in the range: phi_1 < phi_2 <= phi_1 + 2*PI
00117   //
00118   if ( ( ph1 >= 0.0 ) && ( ph1 < 2*pi ) )
00119   {
00120     phi_1 = ph1;
00121   }
00122   else 
00123   {
00124     std::ostringstream message;
00125     message << "Lower azimuthal limit is out of range." << G4endl
00126             << "Default angle of 0 is used.";    
00127     G4Exception("G4SphericalSurface::G4SphericalSurface()",
00128                 "GeomSolids1001", JustWarning, message);
00129 
00130     phi_1 = 0.0;
00131   }
00132 
00133   if ( ( ph2 > phi_1 ) && ( ph2 <=  ( phi_1 + twopi ) ) )
00134   {
00135     phi_2 = ph2;
00136   }
00137   else 
00138   {
00139     std::ostringstream message;
00140     message << "Upper azimuthal limit is out of range." << G4endl
00141             << "Default angle of 2*PI is used.";    
00142     G4Exception("G4SphericalSurface::G4SphericalSurface()",
00143                 "GeomSolids1001", JustWarning, message);
00144 
00145     phi_2 = twopi;
00146   }
00147  
00148   //  Require theta_1 in the range: 0 <= theta_1 < PI 
00149   //  and theta-2 in the range: theta_1 < theta_2 <= theta_1 + PI
00150   //
00151   if ( ( th1 >= 0.0 ) && ( th1 < pi ) )
00152   {
00153     theta_1 = th1;
00154   }
00155   else
00156   {
00157     std::ostringstream message;
00158     message << "Lower polar limit is out of range." << G4endl
00159             << "Default angle of 0 is used.";    
00160     G4Exception("G4SphericalSurface::G4SphericalSurface()",
00161                 "GeomSolids1001", JustWarning, message);
00162 
00163     theta_1 = 0.0;
00164   }  
00165   
00166   if  ( ( th2 > theta_1 ) && ( th2 <= ( theta_1 + pi ) ) )
00167   {
00168     theta_2 =th2;
00169   }
00170   else 
00171   {
00172     std::ostringstream message;
00173     message << "Upper polar limit is out of range." << G4endl
00174             << "Default angle of PI is used.";    
00175     G4Exception("G4SphericalSurface::G4SphericalSurface()",
00176                 "GeomSolids1001", JustWarning, message);
00177 
00178     theta_2 = pi;
00179   } 
00180 }

G4SphericalSurface::~G4SphericalSurface (  )  [virtual]

Definition at line 183 of file G4SphericalSurface.cc.

00184 {
00185 }


Member Function Documentation

G4double G4SphericalSurface::Area (  )  const [virtual]

Definition at line 831 of file G4SphericalSurface.cc.

References phi_1, phi_2, G4INCL::Math::pi, radius, theta_1, and theta_2.

00832 {
00833   // Returns the Area of a G4SphericalSurface
00834 
00835   return ( 2.0*( theta_2 - theta_1 )*( phi_2 - phi_1)*radius*radius/pi );
00836 }

void G4SphericalSurface::CalcBBox (  )  [virtual]

Reimplemented from G4Surface.

Definition at line 337 of file G4SphericalSurface.cc.

References G4Surface::bbox, G4Surface::origin, and radius.

00338 {
00339   G4double x_min = origin.x() - radius;
00340   G4double y_min = origin.y() - radius;
00341   G4double z_min = origin.z() - radius;
00342   G4double x_max = origin.x() + radius;
00343   G4double y_max = origin.y() + radius;
00344   G4double z_max = origin.z() + radius;  
00345     
00346   G4Point3D Min(x_min, y_min, z_min);
00347   G4Point3D Max(x_max, y_max, z_max);  
00348   bbox = new G4BoundingBox3D( Min, Max); 
00349 }

void G4SphericalSurface::Comp ( G4Vector3D v,
G4Point3D min,
G4Point3D max 
) [inline]

Definition at line 57 of file G4SphericalSurface.icc.

00058 {
00059   if(v.x() > max.x()) max.setX(v.x());
00060   if(v.y() > max.y()) max.setY(v.y());
00061   if(v.z() > max.z()) max.setZ(v.z());
00062 
00063   if(v.x() < min.x()) min.setX(v.x());
00064   if(v.y() < min.y()) min.setY(v.y());
00065   if(v.z() < min.z()) min.setZ(v.z());
00066 }

G4String G4SphericalSurface::GetEntityType (  )  const [inline, virtual]

Reimplemented from G4Surface.

Definition at line 51 of file G4SphericalSurface.icc.

00052 {
00053   return G4String("Spherical_Surface");
00054 }

G4double G4SphericalSurface::GetPhi1 (  )  const [inline]

Definition at line 87 of file G4SphericalSurface.icc.

References phi_1.

00088 {
00089   return phi_1;
00090 }

G4double G4SphericalSurface::GetPhi2 (  )  const [inline]

Definition at line 93 of file G4SphericalSurface.icc.

References phi_2.

00094 {
00095   return phi_2;
00096 }

G4double G4SphericalSurface::GetRadius (  )  const [inline]

Definition at line 81 of file G4SphericalSurface.icc.

References radius.

Referenced by Intersect().

00082 {
00083   return radius;
00084 }

G4double G4SphericalSurface::GetTheta1 (  )  const [inline]

Definition at line 99 of file G4SphericalSurface.icc.

References theta_1.

00100 {
00101   return theta_1;
00102 }

G4double G4SphericalSurface::GetTheta2 (  )  const [inline]

Definition at line 105 of file G4SphericalSurface.icc.

References theta_2.

00106 {
00107   return theta_2;
00108 }

G4Vector3D G4SphericalSurface::GetXAxis (  )  const [inline]

Definition at line 69 of file G4SphericalSurface.icc.

References x_axis.

00070 {
00071   return x_axis;
00072 }

G4Vector3D G4SphericalSurface::GetZAxis (  )  const [inline]

Definition at line 75 of file G4SphericalSurface.icc.

References z_axis.

00076 {
00077   return z_axis;
00078 }

G4double G4SphericalSurface::HowNear ( const G4Vector3D x  )  const [virtual]

Reimplemented from G4Surface.

Definition at line 233 of file G4SphericalSurface.cc.

References G4Surface::origin, and radius.

Referenced by Inside(), and Intersect().

00234 {
00235   //  Distance from the point x to the G4SphericalSurface.
00236   //  The distance will be positive if the point is Inside the 
00237   //  G4SphericalSurface, negative if the point is outside.
00238 
00239   G4Vector3D d = G4Vector3D( x - origin );
00240   G4double rds = d.mag();
00241 
00242   return (radius - rds);
00243 }

G4int G4SphericalSurface::Inside ( const G4Vector3D x  )  const [virtual]

Definition at line 756 of file G4SphericalSurface.cc.

References HowNear().

00757 {
00758   // Return 0 if point x is outside G4SphericalSurface, 1 if Inside.
00759   // Outside means that the distance to the G4SphericalSurface would 
00760   // be negative.
00761   // Use the HowNear function to calculate this distance.
00762 
00763   if ( HowNear( x ) >= 0.0 )  { return 1; }
00764   else                        { return 0; }
00765 }

G4int G4SphericalSurface::Intersect ( const G4Ray  )  [virtual]

Reimplemented from G4Surface.

Definition at line 352 of file G4SphericalSurface.cc.

References G4Surface::closest_hit, G4Surface::distance, G4Ray::GetDir(), G4Surface::GetOrigin(), G4Ray::GetPoint(), GetRadius(), G4Ray::GetStart(), HowNear(), G4Surface::kCarTolerance, Normal(), and WithinBoundary().

00353 {
00354   //  Distance along a Ray (straight line with G4Vector3D) to leave or enter
00355   //  a G4SphericalSurface.  The input variable which_way should be set to +1 
00356   //  to indicate leaving a G4SphericalSurface, -1 to indicate entering a 
00357   //  G4SphericalSurface.
00358   //  p is the point of intersection of the Ray with the G4SphericalSurface.
00359   //  If the G4Vector3D of the Ray is opposite to that of the Normal to
00360   //  the G4SphericalSurface at the intersection point, it will not leave the
00361   //  G4SphericalSurface.
00362   //  Similarly, if the G4Vector3D of the Ray is along that of the Normal 
00363   //  to the G4SphericalSurface at the intersection point, it will not enter 
00364   //  the G4SphericalSurface.
00365   //  This method is called by all finite shapes sub-classed to 
00366   //  G4SphericalSurface.
00367   //  Use the virtual function table to check if the intersection point
00368   //  is within the boundary of the finite shape.
00369   //  A negative result means no intersection.
00370   //  If no valid intersection point is found, set the distance
00371   //  and intersection point to large numbers.
00372 
00373   G4int which_way = (G4int)HowNear(ry.GetStart());
00374   
00375   if(!which_way)  { which_way =-1; }
00376   distance = kInfinity;
00377   
00378   G4Vector3D lv ( kInfinity, kInfinity, kInfinity );
00379 
00380   // p = lv;
00381   //
00382   closest_hit = lv;
00383 
00384   // Origin and G4Vector3D unit vector of Ray.
00385   //
00386   G4Vector3D x= G4Vector3D( ry.GetStart() );
00387   G4Vector3D dhat = ry.GetDir();
00388   G4int isoln = 0, maxsoln = 2;
00389 
00390   // Array of solutions in distance along the Ray
00391   //
00392   G4double ss[2];
00393   ss[0] = -1.0 ; 
00394   ss[1] = -1.0 ;
00395 
00396   // Calculate the two solutions (quadratic equation)
00397   //
00398   G4Vector3D d = G4Vector3D( x - GetOrigin() );
00399   G4double r = GetRadius();
00400 
00401   // Quit with no intersection if the radius of the G4SphericalSurface is zero
00402   //
00403   if ( r <= 0.0 )  { return 0; }
00404   
00405   G4double dsq     = d * d;
00406   G4double rsq     = r * r;
00407   G4double b       = d * dhat;
00408   G4double c       = dsq - rsq;
00409   G4double radical = b * b - c;
00410   
00411   // Quit with no intersection if the radical is negative
00412   //
00413   if ( radical < 0.0 )  { return 0; }
00414   
00415   G4double root = std::sqrt( radical );
00416   ss[0] = -b + root;
00417   ss[1] = -b - root;
00418 
00419   // Order the possible solutions by increasing distance along the Ray
00420   // G4Sort_double( ss, isoln, maxsoln-1 );
00421   //        
00422   if(ss[0] > ss[1])
00423   {
00424     G4double tmp =ss[0];
00425     ss[0] = ss[1];
00426     ss[1] = tmp;
00427   }
00428 
00429   // Now loop over each positive solution, keeping the first one (smallest
00430   // distance along the Ray) which is within the boundary of the sub-shape
00431   // and which also has the correct G4Vector3D with respect to the Normal to
00432   // the G4SphericalSurface at the intersection point
00433   //
00434   for ( isoln = 0; isoln < maxsoln; isoln++ ) 
00435   {
00436     if ( ss[isoln] >= kCarTolerance*0.5 ) 
00437     {
00438       if ( ss[isoln] >= kInfinity )  { return 0; }  // quit if too large
00439       distance = ss[isoln];
00440       closest_hit = ry.GetPoint( distance );
00441       if ( ( ( dhat * Normal( closest_hit ) * which_way ) >= 0.0 ) && 
00442            ( WithinBoundary( closest_hit ) == 1 )                      )
00443       {
00444         distance =  distance*distance;    
00445         return 1;
00446       }
00447     }
00448   }
00449   
00450   // Get here only if there was no solution within the boundary,
00451   // reset distance and intersection point to large numbers
00452   //
00453   distance = kInfinity;
00454   closest_hit = lv;
00455 
00456   return 0;
00457 }          

const char * G4SphericalSurface::NameOf (  )  const [virtual]

Definition at line 215 of file G4SphericalSurface.cc.

00216 {
00217   return "G4SphericalSurface";
00218 }

G4Vector3D G4SphericalSurface::Normal ( const G4Vector3D p  )  const [virtual]

Reimplemented from G4Surface.

Definition at line 720 of file G4SphericalSurface.cc.

References CLHEP::detail::n, and G4Surface::origin.

Referenced by Intersect().

00721 { 
00722   // Return the Normal unit vector to the G4SphericalSurface at a point p on
00723   // (or nearly on) the G4SphericalSurface.
00724 
00725   G4Vector3D n = G4Vector3D( p - origin );
00726   G4double nmag = n.mag();
00727   
00728   //  If the point p happens to coincide with the origin (which is possible
00729   //  if the radius is zero), set the Normal to the z-axis unit vector.
00730   //
00731   if ( nmag != 0.0 )  { n = n * (1/ nmag); }
00732   else  { n = G4Vector3D( 0.0, 0.0, 1.0 ); }
00733 
00734   return n;
00735 }

G4int G4SphericalSurface::operator== ( const G4SphericalSurface s  )  [inline]

Definition at line 38 of file G4SphericalSurface.icc.

References G4Surface::origin, phi_1, phi_2, radius, theta_1, theta_2, x_axis, and z_axis.

00039 {
00040   return origin  == surf.origin  &&  
00041   x_axis  == surf.x_axis  &&
00042   z_axis  == surf.z_axis  &&
00043   radius  == surf.radius  && 
00044   phi_1   == surf.phi_1   &&
00045   phi_2   == surf.phi_2   &&
00046   theta_1 == surf.theta_1 &&
00047   theta_2 == surf.theta_2;
00048 }   

void G4SphericalSurface::PrintOn ( std::ostream &  os = G4cout  )  const [virtual]

Definition at line 220 of file G4SphericalSurface.cc.

References G4Surface::origin, phi_1, phi_2, radius, theta_1, theta_2, x_axis, and z_axis.

00221 {  
00222   os << "G4SphericalSurface surface with origin: " << origin << "\t"
00223      << "radius: " << radius << "\n"
00224      << "\t local x_axis: " << x_axis
00225      << "\t local z_axis: " << z_axis << "\n" 
00226      << "\t lower azimuthal limit: " << phi_1 << " radians\n"
00227      << "\t upper azimuthal limit: " << phi_2 << " radians\n"
00228      << "\t lower polar limit    : " << theta_1 << " radians\n"
00229      << "\t upper polar limit    : " << theta_2 << " radians\n";
00230 }

void G4SphericalSurface::resize ( G4double  r,
G4double  ph1,
G4double  ph2,
G4double  th1,
G4double  th2 
) [virtual]

Definition at line 839 of file G4SphericalSurface.cc.

References G4endl, G4Exception(), JustWarning, G4Surface::kAngTolerance, phi_1, phi_2, G4INCL::Math::pi, radius, theta_1, and theta_2.

00842 { 
00843   // Resize the G4SphericalSurface to new radius r, new lower and upper 
00844   // azimuthal angle limits ph1 and ph2, and new lower and upper polar 
00845   // angle limits th1 and th2.
00846   
00847   //  Require radius to be non-negative
00848   //
00849   if ( r >= 0.0 )  { radius = r; }
00850   else 
00851   {
00852     std::ostringstream message;
00853     message << "Radius cannot be less than zero." << G4endl
00854             << "Original value of " << radius << " is retained.";    
00855     G4Exception("G4SphericalSurface::resize()",
00856                 "GeomSolids1001", JustWarning, message);
00857   }
00858 
00859   //  Require azimuthal angles to be within bounds
00860   
00861   if ( ( ph1 >= 0.0 ) && ( ph1 < twopi ) )  { phi_1 = ph1; }
00862   else 
00863   {
00864     std::ostringstream message;
00865     message << "Lower azimuthal limit out of range." << G4endl
00866             << "Original value of " << phi_1 << " is retained.";    
00867     G4Exception("G4SphericalSurface::resize()",
00868                 "GeomSolids1001", JustWarning, message);
00869   }
00870   
00871   if ( ( ph2 > phi_1 ) && ( ph2 <= ( phi_1 + twopi ) ) )  { phi_2 = ph2; }
00872   else 
00873   {
00874     ph2 = ( phi_2 <= phi_1 ) ? ( phi_1 + kAngTolerance ) : phi_2;
00875     phi_2 = ph2;
00876     std::ostringstream message;
00877     message << "Upper azimuthal limit out of range." << G4endl
00878             << "Value of " << phi_2 << " is used.";    
00879     G4Exception("G4SphericalSurface::resize()",
00880                 "GeomSolids1001", JustWarning, message);
00881   }
00882 
00883   // Require polar angles to be within bounds
00884   //
00885   if ( ( th1 >= 0.0 ) && ( th1 < pi ) )  { theta_1 = th1; }
00886   else 
00887   {
00888     std::ostringstream message;
00889     message << "Lower polar limit out of range." << G4endl
00890             << "Original value of " << theta_1 << " is retained.";    
00891     G4Exception("G4SphericalSurface::resize()",
00892                 "GeomSolids1001", JustWarning, message);
00893   }
00894   
00895   if ( ( th2 > theta_1 ) && ( th2 <= ( theta_1 + pi ) ) )  { theta_2 = th2; }
00896   else
00897   {
00898     th2 = ( theta_2 <= theta_1 ) ? ( theta_1 + kAngTolerance ) : theta_2;
00899     theta_2 = th2;
00900     std::ostringstream message;
00901     message << "Upper polar limit out of range." << G4endl
00902             << "Value of " << theta_2 << " is used.";    
00903     G4Exception("G4SphericalSurface::resize()",
00904                 "GeomSolids1001", JustWarning, message);
00905   }
00906 }

G4double G4SphericalSurface::Scale (  )  const [virtual]

Definition at line 820 of file G4SphericalSurface.cc.

References radius.

00821 {
00822   // Returns the radius of a G4SphericalSurface unless it is zero, in which
00823   // case returns the arbitrary number 1.0.
00824   // Used for Scale-invariant tests of surface thickness.
00825 
00826   if ( radius == 0.0 )  { return 1.0; }
00827   else                  { return radius; }
00828 }

G4Vector3D G4SphericalSurface::SurfaceNormal ( const G4Point3D p  )  const [virtual]

Implements G4Surface.

Definition at line 738 of file G4SphericalSurface.cc.

References CLHEP::detail::n, and G4Surface::origin.

00739 { 
00740   // Return the Normal unit vector to the G4SphericalSurface at a point p on
00741   // (or nearly on) the G4SphericalSurface.
00742 
00743   G4Vector3D n = G4Vector3D( p - origin );
00744   G4double nmag = n.mag();
00745   
00746   //  If the point p happens to coincide with the origin (which is possible
00747   //  if the radius is zero), set the Normal to the z-axis unit vector.
00748   //
00749   if ( nmag != 0.0 )  { n = n * (1/ nmag); }
00750   else  { n = G4Vector3D( 0.0, 0.0, 1.0 ); }
00751 
00752   return n;
00753 }

G4int G4SphericalSurface::WithinBoundary ( const G4Vector3D x  )  const [virtual]

Definition at line 768 of file G4SphericalSurface.cc.

References phi_1, phi_2, G4INCL::Math::pi, theta_1, theta_2, x_axis, and z_axis.

Referenced by Intersect().

00769 { 
00770   // Return 1 if point x is on the G4SphericalSurface, otherwise return zero
00771   // (x is assumed to lie on the surface of the G4SphericalSurface, so one 
00772   // only checks the angular limits)
00773 
00774   G4Vector3D y_axis = G4Vector3D( z_axis.cross( x_axis ) );
00775 
00776   // Components of x in the local coordinate system of the G4SphericalSurface
00777   //
00778   G4double px = x * x_axis;
00779   G4double py = x * y_axis;
00780   G4double pz = x * z_axis;
00781 
00782   // Check if within polar angle limits
00783   //
00784   G4double theta = std::acos( pz / x.mag() );  // acos in range 0 to PI
00785   
00786   // Normal case
00787   //
00788   if ( theta_2 <= pi ) 
00789   {
00790     if ( ( theta < theta_1 ) || ( theta > theta_2 ) )  { return 0; }
00791   }
00792   else 
00793   {
00794     theta += pi;
00795     if ( ( theta < theta_1 ) || ( theta > theta_2 ) )  { return 0; }
00796   }
00797 
00798   // Now check if within azimuthal angle limits
00799   //
00800   G4double phi = std::atan2( py, px );  // atan2 in range -PI to PI
00801   
00802   if ( phi < 0.0 )  { phi += twopi; }
00803   
00804   // Normal case
00805   //
00806   if ( ( phi >= phi_1 )  &&  ( phi <= phi_2 ) )  { return 1; }
00807   
00808   // This is for the case that phi_2 is greater than 2*PI
00809   //
00810   phi += twopi;
00811   
00812   if ( ( phi >= phi_1 )  &&  ( phi <= phi_2 ) )  { return 1; }
00813 
00814   // Get here if not within azimuthal limits
00815  
00816   return 0;
00817 }


Field Documentation

G4double G4SphericalSurface::phi_1 [protected]

Definition at line 213 of file G4SphericalSurface.hh.

Referenced by Area(), G4SphericalSurface(), GetPhi1(), operator==(), PrintOn(), resize(), and WithinBoundary().

G4double G4SphericalSurface::phi_2 [protected]

Definition at line 217 of file G4SphericalSurface.hh.

Referenced by Area(), G4SphericalSurface(), GetPhi2(), operator==(), PrintOn(), resize(), and WithinBoundary().

G4double G4SphericalSurface::radius [protected]

Definition at line 210 of file G4SphericalSurface.hh.

Referenced by Area(), CalcBBox(), G4SphericalSurface(), GetRadius(), HowNear(), operator==(), PrintOn(), resize(), and Scale().

G4double G4SphericalSurface::theta_1 [protected]

Definition at line 221 of file G4SphericalSurface.hh.

Referenced by Area(), G4SphericalSurface(), GetTheta1(), operator==(), PrintOn(), resize(), and WithinBoundary().

G4double G4SphericalSurface::theta_2 [protected]

Definition at line 225 of file G4SphericalSurface.hh.

Referenced by Area(), G4SphericalSurface(), GetTheta2(), operator==(), PrintOn(), resize(), and WithinBoundary().

G4Vector3D G4SphericalSurface::x_axis [protected]

Definition at line 202 of file G4SphericalSurface.hh.

Referenced by G4SphericalSurface(), GetXAxis(), operator==(), PrintOn(), and WithinBoundary().

G4Vector3D G4SphericalSurface::z_axis [protected]

Definition at line 206 of file G4SphericalSurface.hh.

Referenced by G4SphericalSurface(), GetZAxis(), operator==(), PrintOn(), and WithinBoundary().


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:53:24 2013 for Geant4 by  doxygen 1.4.7