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
00048
00049
00050
00051 #ifndef G4TET_HH
00052 #define G4TET_HH
00053
00054 #include "G4VSolid.hh"
00055
00056 class G4Tet : public G4VSolid
00057 {
00058
00059 public:
00060
00061 G4Tet(const G4String& pName,
00062 G4ThreeVector anchor,
00063 G4ThreeVector p2,
00064 G4ThreeVector p3,
00065 G4ThreeVector p4,
00066 G4bool *degeneracyFlag=0);
00067
00068 virtual ~G4Tet();
00069
00070 void ComputeDimensions(G4VPVParameterisation* p,
00071 const G4int n,
00072 const G4VPhysicalVolume* pRep);
00073
00074 G4bool CalculateExtent(const EAxis pAxis,
00075 const G4VoxelLimits& pVoxelLimit,
00076 const G4AffineTransform& pTransform,
00077 G4double& pmin, G4double& pmax) const;
00078
00079
00080 EInside Inside(const G4ThreeVector& p) const;
00081
00082 G4ThreeVector SurfaceNormal( const G4ThreeVector& p) const;
00083
00084 G4double DistanceToIn(const G4ThreeVector& p, const G4ThreeVector& v) const;
00085
00086 G4double DistanceToIn(const G4ThreeVector& p) const;
00087
00088 G4double DistanceToOut(const G4ThreeVector& p, const G4ThreeVector& v,
00089 const G4bool calcNorm=false,
00090 G4bool *validNorm=0, G4ThreeVector *n=0) const;
00091
00092 G4double DistanceToOut(const G4ThreeVector& p) const;
00093
00094 G4double GetCubicVolume();
00095 G4double GetSurfaceArea();
00096
00097 G4GeometryType GetEntityType() const;
00098
00099 G4VSolid* Clone() const;
00100
00101 std::ostream& StreamInfo(std::ostream& os) const;
00102
00103 G4ThreeVector GetPointOnSurface() const;
00104
00105
00106
00107 void DescribeYourselfTo (G4VGraphicsScene& scene) const;
00108 G4VisExtent GetExtent () const;
00109 G4Polyhedron* CreatePolyhedron () const;
00110 G4NURBS* CreateNURBS () const;
00111 G4Polyhedron* GetPolyhedron () const;
00112
00113 public:
00114
00115 G4Tet(__void__&);
00116
00117
00118
00119
00120 G4Tet(const G4Tet& rhs);
00121 G4Tet& operator=(const G4Tet& rhs);
00122
00123
00124 const char* CVSHeaderVers()
00125 { return "$Id: G4Tet.hh 67011 2013-01-29 16:17:41Z gcosmo $"; }
00126 const char* CVSFileVers()
00127 { return CVSVers; }
00128 void PrintWarnings(G4bool flag)
00129 { warningFlag=flag; }
00130 static G4bool CheckDegeneracy(G4ThreeVector anchor,
00131 G4ThreeVector p2,
00132 G4ThreeVector p3,
00133 G4ThreeVector p4);
00134 std::vector<G4ThreeVector> GetVertices() const;
00135
00136
00137 protected:
00138
00139 G4ThreeVectorList*
00140 CreateRotatedVertices(const G4AffineTransform& pTransform) const;
00141
00142
00143
00144 private:
00145
00146 G4double fCubicVolume, fSurfaceArea;
00147
00148 mutable G4Polyhedron* fpPolyhedron;
00149
00150 G4ThreeVector GetPointOnFace(G4ThreeVector p1, G4ThreeVector p2,
00151 G4ThreeVector p3, G4double& area) const;
00152 static const char CVSVers[];
00153
00154 private:
00155
00156 G4ThreeVector fAnchor, fP2, fP3, fP4, fMiddle;
00157 G4ThreeVector fNormal123, fNormal142, fNormal134, fNormal234;
00158
00159 G4bool warningFlag;
00160
00161 G4double fCdotN123, fCdotN142, fCdotN134, fCdotN234;
00162 G4double fXMin, fXMax, fYMin, fYMax, fZMin, fZMax;
00163 G4double fDx, fDy, fDz, fTol, fMaxSize;
00164 };
00165
00166 #endif