38 #if !(defined(G4GEOM_USE_UELLIPTICALCONE) && defined(G4GEOM_USE_SYS_USOLIDS))
66 using namespace CLHEP;
83 if ( (pxSemiAxis <= 0.) || (pySemiAxis <= 0.) || (pzMax <= 0.) )
86 message <<
"Invalid semi-axis or height for solid: " <<
GetName()
87 <<
"\n X semi-axis, Y semi-axis, height = "
88 << pxSemiAxis <<
", " << pySemiAxis <<
", " << pzMax;
89 G4Exception(
"G4EllipticalCone::G4EllipticalCone()",
"GeomSolids0002",
96 message <<
"Invalid z-coordinate for cutting plane for solid: " <<
GetName()
97 <<
"\n Z top cut = " << pzTopCut;
98 G4Exception(
"G4EllipticalCone::G4EllipticalCone()",
"GeomSolids0002",
113 xSemiAxis(0.), ySemiAxis(0.), zheight(0.), zTopCut(0.),
114 cosAxisMin(0.), invXX(0.), invYY(0.)
132 :
G4VSolid(rhs), halfCarTol(rhs.halfCarTol),
133 fCubicVolume(rhs.fCubicVolume), fSurfaceArea(rhs.fSurfaceArea),
134 xSemiAxis(rhs.xSemiAxis), ySemiAxis(rhs.ySemiAxis),
135 zheight(rhs.zheight), zTopCut(rhs.zTopCut),
136 cosAxisMin(rhs.cosAxisMin), invXX(rhs.invXX), invYY(rhs.invYY)
148 if (
this == &rhs) {
return *
this; }
179 pMin.
set(-xmax,-ymax,-zcut);
180 pMax.
set( xmax, ymax, zcut);
184 if (pMin.
x() >= pMax.
x() || pMin.
y() >= pMax.
y() || pMin.
z() >= pMax.
z())
187 message <<
"Bad bounding box (min >= max) for solid: "
189 <<
"\npMin = " << pMin
190 <<
"\npMax = " <<
pMax;
191 G4Exception(
"G4EllipticalCone::BoundingLimits()",
"GeomMgt0001",
219 return exist = (pMin <
pMax) ?
true :
false;
224 static const G4int NSTEPS = 48;
226 static const G4double sinHalf = std::sin(0.5*ang);
227 static const G4double cosHalf = std::cos(0.5*ang);
228 static const G4double sinStep = 2.*sinHalf*cosHalf;
229 static const G4double cosStep = 1. - 2.*sinHalf*sinHalf;
242 baseA[
k].set(sxmax*cosCur,symax*sinCur,-zcut);
243 baseB[
k].set(sxmin*cosCur,symin*sinCur, zcut);
246 sinCur = sinCur*cosStep + cosCur*sinStep;
247 cosCur = cosCur*cosStep - sinTmp*sinStep;
250 std::vector<const G4ThreeVectorList *> polygons(2);
251 polygons[0] = &baseA;
252 polygons[1] = &baseB;
299 if (nsurf == 1)
return norm;
300 else if (nsurf > 1)
return norm.
unit();
307 G4int oldprc = message.precision(16);
308 message <<
"Point p is not on surface (!?) of solid: "
310 message <<
"Position:\n";
311 message <<
" p.x() = " << p.
x()/
mm <<
" mm\n";
312 message <<
" p.y() = " << p.
y()/
mm <<
" mm\n";
313 message <<
" p.z() = " << p.
z()/
mm <<
" mm";
315 G4Exception(
"G4EllipticalCone::SurfaceNormal(p)",
"GeomSolids1002",
394 yi = p.
y() + q*v.
y();
439 yi = p.
y() + q*v.
y();
469 return distMin = std::fabs(lambda);
484 return distMin = std::fabs(lambda);
528 return distMin = std::fabs(-B/(2.*A));
531 G4double plus = (-B+std::sqrt(discr))/(2.*A);
532 G4double minus = (-B-std::sqrt(discr))/(2.*A);
541 if ( truenorm*v >= 0)
601 return (dist > 0) ? dist : 0.;
616 enum surface_e {kPlaneSurf, kCurvedSurf, kNoSurf}
surface;
629 distMin = std::fabs(lambda);
631 if (!calcNorm) {
return distMin; }
633 distMin = std::fabs(lambda);
645 distMin = std::fabs(lambda);
646 if (!calcNorm) {
return distMin; }
648 distMin = std::fabs(lambda);
665 if(!calcNorm) {
return distMin = std::fabs(-B/(2.*A)); }
670 G4double plus = (-B+std::sqrt(discr))/(2.*A);
671 G4double minus = (-B-std::sqrt(discr))/(2.*A);
677 lambda = std::fabs(plus) < std::fabs(minus) ? plus : minus;
687 if ( std::fabs(lambda) < distMin )
691 distMin = std::fabs(lambda);
699 if( truenorm.
dot(v) > 0 )
733 truenorm /= truenorm.
mag();
741 G4int oldprc = message.precision(16);
742 message <<
"Undefined side for valid surface normal to solid."
745 <<
" p.x() = " << p.
x()/
mm <<
" mm" <<
G4endl
746 <<
" p.y() = " << p.
y()/
mm <<
" mm" <<
G4endl
747 <<
" p.z() = " << p.
z()/
mm <<
" mm" <<
G4endl
749 <<
" v.x() = " << v.
x() <<
G4endl
750 <<
" v.y() = " << v.
y() <<
G4endl
751 <<
" v.z() = " << v.
z() <<
G4endl
752 <<
"Proposed distance :" <<
G4endl
753 <<
" distMin = " << distMin/
mm <<
" mm";
754 message.precision(oldprc);
755 G4Exception(
"G4EllipticalCone::DistanceToOut(p,v,..)",
777 G4int oldprc = message.precision(16);
778 message <<
"Point p is outside (!?) of solid: " <<
GetName() <<
"\n"
780 <<
" p.x() = " << p.
x()/
mm <<
" mm\n"
781 <<
" p.y() = " << p.
y()/
mm <<
" mm\n"
782 <<
" p.z() = " << p.
z()/
mm <<
" mm";
783 message.precision(oldprc) ;
784 G4Exception(
"G4Ellipsoid::DistanceToOut(p)",
"GeomSolids1002",
793 return (dist > 0) ? dist : 0.;
802 return G4String(
"G4EllipticalCone");
820 G4int oldprc = os.precision(16);
821 os <<
"-----------------------------------------------------------\n"
822 <<
" *** Dump for solid - " <<
GetName() <<
" ***\n"
823 <<
" ===================================================\n"
824 <<
" Solid type: G4EllipticalCone\n"
829 <<
" height z: " <<
zheight/
mm <<
" mm \n"
830 <<
" half length in z: " <<
zTopCut/
mm <<
" mm \n"
831 <<
"-----------------------------------------------------------\n";
832 os.precision(oldprc);
853 G4double sside = s0*(kmax*kmax - kmin*kmin);
854 G4double ssurf[3] = { szmin, sside, szmax };
855 for (
auto i=1; i<3; ++i) { ssurf[i] += ssurf[i-1]; }
861 if (select <= ssurf[1]) k = 1;
862 if (select <= ssurf[0]) k = 0;
886 G4double mu_max = R*std::sqrt(hh + R*R);
889 for (
auto i=0; i<1000; ++i)
929 fCubicVolume = (kmax - kmin)*(kmax*kmax + kmax*kmin + kmin*kmin)*v0;
948 +
CLHEP::pi*x0*y0*(kmin*kmin + kmax*kmax);
994 #endif // !defined(G4GEOM_USE_UELLIPTICALCONE) || !defined(G4GEOM_USE_SYS_USOLIDS)