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 #include "G4BREPSolidCone.hh"
00036 #include "G4FPlane.hh"
00037 #include "G4FConicalSurface.hh"
00038 #include "G4FCylindricalSurface.hh"
00039 #include "G4CircularCurve.hh"
00040
00041 G4BREPSolidCone::G4BREPSolidCone(const G4String& name,
00042 const G4ThreeVector& origin,
00043 const G4ThreeVector& axis,
00044 const G4ThreeVector& direction,
00045 G4double length,
00046 G4double radius,
00047 G4double large_radius)
00048 : G4BREPSolid(name)
00049 {
00050 nb_of_surfaces = 3;
00051 active=1;
00052
00053
00054 constructorParams.origin = origin;
00055 constructorParams.axis = axis;
00056 constructorParams.direction = direction;
00057 constructorParams.length = length;
00058 constructorParams.radius = radius;
00059 constructorParams.large_radius = large_radius;
00060
00061 InitializeCone();
00062 }
00063
00064 G4BREPSolidCone::G4BREPSolidCone( __void__& a )
00065 : G4BREPSolid(a)
00066 {
00067 }
00068
00069 G4BREPSolidCone::~G4BREPSolidCone()
00070 {
00071 }
00072
00073 G4BREPSolidCone::G4BREPSolidCone(const G4BREPSolidCone& rhs)
00074 : G4BREPSolid(rhs)
00075 {
00076 constructorParams.origin = rhs.constructorParams.origin;
00077 constructorParams.axis = rhs.constructorParams.axis;
00078 constructorParams.direction = rhs.constructorParams.direction;
00079 constructorParams.length = rhs.constructorParams.length;
00080 constructorParams.radius = rhs.constructorParams.radius;
00081 constructorParams.large_radius = rhs.constructorParams.large_radius;
00082
00083 InitializeCone();
00084 }
00085
00086 G4BREPSolidCone& G4BREPSolidCone::operator = (const G4BREPSolidCone& rhs)
00087 {
00088
00089
00090 if (this == &rhs) { return *this; }
00091
00092
00093
00094 G4BREPSolid::operator=(rhs);
00095
00096
00097
00098 constructorParams.origin = rhs.constructorParams.origin;
00099 constructorParams.axis = rhs.constructorParams.axis;
00100 constructorParams.direction = rhs.constructorParams.direction;
00101 constructorParams.length = rhs.constructorParams.length;
00102 constructorParams.radius = rhs.constructorParams.radius;
00103 constructorParams.large_radius = rhs.constructorParams.large_radius;
00104
00105 InitializeCone();
00106
00107 return *this;
00108 }
00109
00110 void G4BREPSolidCone::InitializeCone()
00111 {
00112
00113 SurfaceVec = new G4Surface*[3];
00114 G4Point3D ArcStart1 = G4Point3D(constructorParams.origin
00115 + (constructorParams.radius
00116 * constructorParams.direction));
00117 G4Vector3D tmpaxis(constructorParams.axis);
00118 G4Vector3D tmporigin(constructorParams.origin);
00119 G4Point3D tmppoint;
00120
00121 tmppoint= G4Point3D(constructorParams.origin)
00122 + (constructorParams.length*tmpaxis);
00123 G4Point3D origin2(tmppoint.x(), tmppoint.y(), tmppoint.z());
00124
00125 tmppoint= origin2 + (constructorParams.large_radius*tmpaxis);
00126 G4Point3D ArcStart2(tmppoint.x(), tmppoint.y(), tmppoint.z());
00127
00128 G4Ray::Vcross(tmpaxis, constructorParams.axis, constructorParams.direction);
00129 G4ThreeVector axis2(tmpaxis.x(),tmpaxis.y(), tmpaxis.z());
00130
00131 G4CurveVector CVec;
00132 G4CircularCurve* tmp;
00133
00134 tmp = new G4CircularCurve();
00135 tmp->Init(G4Axis2Placement3D(constructorParams.direction,
00136 axis2, constructorParams.origin),
00137 constructorParams.large_radius);
00138 tmp->SetBounds(ArcStart1, ArcStart1);
00139 CVec.push_back(tmp);
00140
00141 tmp = new G4CircularCurve();
00142 tmp->Init(G4Axis2Placement3D(constructorParams.direction, axis2, origin2),
00143 constructorParams.large_radius);
00144 tmp->SetBounds(ArcStart2, ArcStart2);
00145 CVec.push_back(tmp);
00146
00147 SurfaceVec[0] = new G4FConicalSurface(tmporigin, constructorParams.axis,
00148 constructorParams.length,
00149 constructorParams.radius,
00150 constructorParams.large_radius);
00151 SurfaceVec[0]->SetBoundaries(&CVec);
00152
00153
00154 G4CurveVector CVec2;
00155 tmp = new G4CircularCurve();
00156 tmp->Init(G4Axis2Placement3D(constructorParams.direction,
00157 axis2, constructorParams.origin),
00158 constructorParams.radius);
00159 tmp->SetBounds(ArcStart1, ArcStart1);
00160 CVec2.push_back(tmp);
00161
00162 SurfaceVec[1] = new G4FPlane(tmpaxis, constructorParams.direction, origin2);
00163 SurfaceVec[1]->SetBoundaries(&CVec2);
00164
00165 CVec2[0] = tmp = new G4CircularCurve();
00166 tmp->Init(G4Axis2Placement3D(constructorParams.direction, axis2, origin2),
00167 constructorParams.large_radius);
00168 tmp->SetBounds(ArcStart2, ArcStart2);
00169
00170 SurfaceVec[2] = new G4FPlane(tmpaxis, constructorParams.direction,
00171 constructorParams.origin);
00172 SurfaceVec[2]->SetBoundaries(&CVec2);
00173
00174 Initialize();
00175 }
00176
00177 void G4BREPSolidCone::Initialize()
00178 {
00179
00180
00181
00182 ShortestDistance=1000000;
00183 CheckSurfaceNormals();
00184 if(!Box || !AxisBox) { IsConvex(); }
00185 CalcBBoxes();
00186 }
00187
00188 EInside G4BREPSolidCone::Inside(register const G4ThreeVector& Pt) const
00189 {
00190 G4double dist1 = SurfaceVec[0]->HowNear(Pt);
00191 G4double dist2 = SurfaceVec[1]->ClosestDistanceToPoint(Pt);
00192 G4double dist3 = SurfaceVec[2]->ClosestDistanceToPoint(Pt);
00193 if(dist1 > dist2) dist1 = dist2;
00194 if(dist1 > dist3) dist1 = dist3;
00195 if(dist1 > 0) return kInside;
00196 if(dist1 < 0) return kOutside;
00197 return kSurface;
00198 }
00199
00200 G4ThreeVector G4BREPSolidCone::SurfaceNormal(const G4ThreeVector& Pt) const
00201 {
00202 G4Vector3D n = SurfaceVec[0]->Normal(Pt);
00203 G4ThreeVector norm(n.x(), n.y(), n.z());
00204 return norm;
00205 }
00206
00207 G4double G4BREPSolidCone::DistanceToIn(const G4ThreeVector& Pt) const
00208 {
00209 G4double dist1 = std::fabs(SurfaceVec[0]->HowNear(Pt));
00210 G4double dist2 = std::fabs(SurfaceVec[1]->ClosestDistanceToPoint(Pt));
00211 G4double dist3 = std::fabs(SurfaceVec[2]->ClosestDistanceToPoint(Pt));
00212 if(dist1 > dist2) dist1 = dist2;
00213 if(dist1 > dist3) dist1 = dist3;
00214 return dist1;
00215
00216 }
00217
00218 G4double G4BREPSolidCone::DistanceToIn(register const G4ThreeVector& Pt,
00219 register const G4ThreeVector& V) const
00220 {
00221 Reset();
00222 G4Vector3D Pttmp(Pt);
00223 G4Vector3D Vtmp(V);
00224
00225 G4Ray r(Pttmp, Vtmp);
00226
00227 if(SurfaceVec[0]->Intersect( r ))
00228 {
00229 ShortestDistance = SurfaceVec[0]->GetDistance();
00230 return ShortestDistance;
00231 }
00232 return kInfinity;
00233 }
00234
00235 G4double G4BREPSolidCone::DistanceToOut(register const G4ThreeVector& Pt,
00236 register const G4ThreeVector& V,
00237 const G4bool,
00238 G4bool *validNorm,
00239 G4ThreeVector *) const
00240 {
00241 if(validNorm)
00242 *validNorm = false;
00243 Reset();
00244
00245 G4Vector3D Pttmp(Pt);
00246 G4Vector3D Vtmp(V);
00247
00248
00249 G4Ray r(Pttmp, Vtmp);
00250 if(SurfaceVec[0]->Intersect( r ))
00251 {
00252 ShortestDistance = SurfaceVec[0]->GetDistance();
00253 return ShortestDistance;
00254 }
00255 return kInfinity;
00256 }
00257
00258 G4double G4BREPSolidCone::DistanceToOut(const G4ThreeVector& Pt) const
00259 {
00260 G4double dist1 = std::fabs(SurfaceVec[0]->HowNear(Pt));
00261 G4double dist2 = std::fabs(SurfaceVec[1]->ClosestDistanceToPoint(Pt));
00262 G4double dist3 = std::fabs(SurfaceVec[2]->ClosestDistanceToPoint(Pt));
00263 if(dist1 > dist2) dist1 = dist2;
00264 if(dist1 > dist3) dist1 = dist3;
00265 return dist1;
00266 }
00267
00268 G4VSolid* G4BREPSolidCone::Clone() const
00269 {
00270 return new G4BREPSolidCone(*this);
00271 }
00272
00273 std::ostream& G4BREPSolidCone::StreamInfo(std::ostream& os) const
00274 {
00275
00276
00277 G4BREPSolid::StreamInfo( os )
00278 << "\n origin: " << constructorParams.origin
00279 << "\n axis: " << constructorParams.axis
00280 << "\n direction: " << constructorParams.direction
00281 << "\n length: " << constructorParams.length
00282 << "\n radius: " << constructorParams.radius
00283 << "\n large_radius: " << constructorParams.large_radius
00284 << "\n-----------------------------------------------------------\n";
00285
00286 return os;
00287 }
00288