00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047 #ifndef __G4TWISTTUBSHYPESIDE__
00048 #define __G4TWISTTUBSHYPESIDE__
00049
00050 #include "G4VTwistSurface.hh"
00051 #include "G4Integrator.hh"
00052 #include "G4SimpleIntegration.hh"
00053
00054 class G4TwistTubsHypeSide : public G4VTwistSurface
00055 {
00056 public:
00057
00058 G4TwistTubsHypeSide(const G4String &name,
00059 const G4RotationMatrix &rot,
00060 const G4ThreeVector &tlate,
00061 const G4int handedness,
00062 const G4double kappa,
00063 const G4double tanstereo,
00064 const G4double r0,
00065 const EAxis axis0 = kPhi,
00066 const EAxis axis1 = kZAxis,
00067 G4double axis0min = -kInfinity,
00068 G4double axis1min = -kInfinity,
00069 G4double axis0max = kInfinity,
00070 G4double axis1max = kInfinity);
00071
00072 G4TwistTubsHypeSide(const G4String &name,
00073 G4double EndInnerRadius[2],
00074 G4double EndOuterRadius[2],
00075 G4double DPhi,
00076 G4double EndPhi[2],
00077 G4double EndZ[2],
00078 G4double InnerRadius,
00079 G4double OuterRadius,
00080 G4double Kappa,
00081 G4double TanInnerStereo,
00082 G4double TanOuterStereo,
00083 G4int handedness) ;
00084
00085 virtual ~G4TwistTubsHypeSide();
00086
00087 virtual G4int DistanceToSurface(const G4ThreeVector &gp,
00088 const G4ThreeVector &gv,
00089 G4ThreeVector gxx[],
00090 G4double distance[],
00091 G4int areacode[],
00092 G4bool isvalid[],
00093 EValidate validate = kValidateWithTol);
00094
00095 virtual G4int DistanceToSurface(const G4ThreeVector &gp,
00096 G4ThreeVector gxx[],
00097 G4double distance[],
00098 G4int areacode[]);
00099
00100 virtual G4ThreeVector GetNormal(const G4ThreeVector &xx,
00101 G4bool isGlobal = false) ;
00102 virtual EInside Inside(const G4ThreeVector &gp) ;
00103
00104 virtual G4double GetRhoAtPZ(const G4ThreeVector &p,
00105 G4bool isglobal = false) const ;
00106
00107 virtual G4ThreeVector SurfacePoint(G4double, G4double,
00108 G4bool isGlobal = false) ;
00109 virtual G4double GetBoundaryMin(G4double phi) ;
00110 virtual G4double GetBoundaryMax(G4double phi) ;
00111 virtual G4double GetSurfaceArea() ;
00112 virtual void GetFacets( G4int m, G4int n, G4double xyz[][3],
00113 G4int faces[][4], G4int iside ) ;
00114
00115 public:
00116
00117 G4TwistTubsHypeSide(__void__&);
00118
00119
00120
00121
00122 private:
00123
00124 virtual G4int GetAreaCode(const G4ThreeVector &xx,
00125 G4bool withTol = true);
00126 virtual G4int GetAreaCodeInPhi(const G4ThreeVector &xx,
00127 G4bool withTol = true);
00128 virtual void SetCorners();
00129
00130 virtual void SetCorners(G4double EndInnerRadius[2],
00131 G4double EndOuterRadius[2],
00132 G4double DPhi,
00133 G4double EndPhi[2],
00134 G4double EndZ[2]);
00135 virtual void SetBoundaries();
00136
00137 private:
00138
00139 G4double fKappa;
00140 G4double fTanStereo;
00141 G4double fTan2Stereo;
00142 G4double fR0;
00143 G4double fR02;
00144 G4double fDPhi ;
00145
00146 class Insidetype
00147 {
00148 public:
00149 G4ThreeVector gp;
00150 EInside inside;
00151 };
00152 Insidetype fInside;
00153 };
00154
00155
00156
00157
00158
00159 inline
00160 G4double G4TwistTubsHypeSide::GetRhoAtPZ(const G4ThreeVector &p,
00161 G4bool isglobal) const
00162 {
00163
00164 G4ThreeVector tmpp;
00165 if (isglobal) {
00166 tmpp = fRot.inverse()*p - fTrans;
00167 } else {
00168 tmpp = p;
00169 }
00170 return std::sqrt(fR02 + tmpp.z() * tmpp.z() * fTan2Stereo);
00171 }
00172
00173 inline
00174 G4ThreeVector G4TwistTubsHypeSide::
00175 SurfacePoint(G4double phi , G4double z , G4bool isGlobal)
00176 {
00177 G4double rho = std::sqrt(fR02 + z * z * fTan2Stereo) ;
00178
00179 G4ThreeVector SurfPoint (rho*std::cos(phi), rho*std::sin(phi), z) ;
00180
00181 if (isGlobal) { return (fRot * SurfPoint + fTrans); }
00182 return SurfPoint;
00183 }
00184
00185 inline
00186 G4double G4TwistTubsHypeSide::GetBoundaryMin(G4double z)
00187 {
00188 G4ThreeVector ptmp(0,0,z) ;
00189 G4ThreeVector lowerlimit;
00190 lowerlimit = GetBoundaryAtPZ(sAxis0 & sAxisMin, ptmp);
00191 return std::atan2( lowerlimit.y(), lowerlimit.x() ) ;
00192 }
00193
00194 inline
00195 G4double G4TwistTubsHypeSide::GetBoundaryMax(G4double z )
00196 {
00197 G4ThreeVector ptmp(0,0,z) ;
00198 G4ThreeVector upperlimit;
00199 upperlimit = GetBoundaryAtPZ(sAxis0 & sAxisMax, ptmp);
00200 return std::atan2( upperlimit.y(), upperlimit.x() ) ;
00201 }
00202
00203 inline
00204 G4double G4TwistTubsHypeSide::GetSurfaceArea()
00205 {
00206
00207
00208 return ( fAxisMax[1] - fAxisMin[1] ) * fR0 * fDPhi ;
00209 }
00210
00211 #endif