56 using namespace CLHEP;
73 if ( pDx > 0 && pDy > 0 && pDz > 0 )
78 fTalpha = std::tan(pAlpha);
79 fTthetaCphi = std::tan(pTheta)*std::cos(pPhi);
80 fTthetaSphi = std::tan(pTheta)*std::sin(pPhi);
84 std::ostringstream message;
85 message <<
"Invalid Length Parameters for Solid: " << GetName() <<
G4endl
86 <<
" pDx, pDy, pDz = "
87 << pDx <<
", " << pDy <<
", " << pDz;
88 G4Exception(
"G4Para::SetAllParameters()",
"GeomSolids0002",
103 if ((pDx<=0) || (pDy<=0) || (pDz<=0))
105 std::ostringstream message;
106 message <<
"Invalid Length Parameters for Solid: " <<
GetName() <<
G4endl
107 <<
" pDx, pDy, pDz = "
108 << pDx <<
", " << pDy <<
", " << pDz;
125 if (!( pt[0].
z()<0 && pt[0].
z()==pt[1].
z() && pt[0].
z()==pt[2].
z() &&
126 pt[0].
z()==pt[3].
z() && pt[4].
z()>0 && pt[4].
z()==pt[5].
z() &&
127 pt[4].
z()==pt[6].
z() && pt[4].
z()==pt[7].
z() &&
128 (pt[0].
z()+pt[4].
z())==0 &&
129 pt[0].y()==pt[1].y() && pt[2].y()==pt[3].y() &&
130 pt[4].y()==pt[5].y() && pt[6].y()==pt[7].y() &&
131 ( pt[0].y() + pt[2].y() + pt[4].y() + pt[6].y() ) == 0 &&
132 ( pt[0].x() + pt[1].x() + pt[4].x() + pt[5].x() ) == 0) )
134 std::ostringstream message;
135 message <<
"Invalid vertice coordinates for Solid: " <<
GetName();
139 fDx = ((pt[3]).x()-(pt[2]).x())*0.5;
140 fDy = ((pt[2]).y()-(pt[1]).y())*0.5;
142 fTalpha = ((pt[2]).x()+(pt[3]).x()-(pt[1]).x()-(pt[0]).x())*0.25/
fDy ;
156 fTalpha(0.), fTthetaCphi(0.), fTthetaSphi(0.)
172 :
G4CSGSolid(rhs), fDx(rhs.fDx), fDy(rhs.fDy), fDz(rhs.fDz),
173 fTalpha(rhs.fTalpha), fTthetaCphi(rhs.fTthetaCphi),
174 fTthetaSphi(rhs.fTthetaSphi)
186 if (
this == &rhs) {
return *
this; }
280 temp[0] = pt[0].y()+(pt[4].y()-pt[0].y())
281 *(zMin-pt[0].
z())/(pt[4].
z()-pt[0].z()) ;
282 temp[1] = pt[0].y()+(pt[4].y()-pt[0].y())
283 *(zMax-pt[0].
z())/(pt[4].
z()-pt[0].z()) ;
284 temp[2] = pt[2].y()+(pt[6].y()-pt[2].y())
285 *(zMin-pt[2].
z())/(pt[6].
z()-pt[2].z()) ;
286 temp[3] = pt[2].y()+(pt[6].y()-pt[2].y())
287 *(zMax-pt[2].
z())/(pt[6].
z()-pt[2].z()) ;
292 if(temp[i] > yMax) yMax = temp[i] ;
293 if(temp[i] < yMin) yMin = temp[i] ;
316 temp[0] = pt[0].x()+(pt[4].x()-pt[0].x())
317 *(zMin-pt[0].
z())/(pt[4].
z()-pt[0].z()) ;
318 temp[1] = pt[0].x()+(pt[4].x()-pt[0].x())
319 *(zMax-pt[0].
z())/(pt[4].
z()-pt[0].z()) ;
320 temp[2] = pt[2].x()+(pt[6].x()-pt[2].x())
321 *(zMin-pt[2].
z())/(pt[6].
z()-pt[2].z()) ;
322 temp[3] = pt[2].x()+(pt[6].x()-pt[2].x())
323 *(zMax-pt[2].
z())/(pt[6].
z()-pt[2].z()) ;
324 temp[4] = pt[3].x()+(pt[7].x()-pt[3].x())
325 *(zMin-pt[3].
z())/(pt[7].
z()-pt[3].z()) ;
326 temp[5] = pt[3].x()+(pt[7].x()-pt[3].x())
327 *(zMax-pt[3].
z())/(pt[7].
z()-pt[3].z()) ;
328 temp[6] = pt[1].x()+(pt[5].x()-pt[1].x())
329 *(zMin-pt[1].
z())/(pt[5].
z()-pt[1].z()) ;
330 temp[7] = pt[1].x()+(pt[5].x()-pt[1].x())
331 *(zMax-pt[1].
z())/(pt[5].
z()-pt[1].z()) ;
337 if(temp[i] > xMax) xMax = temp[i] ;
338 if(temp[i] < xMin) xMin = temp[i] ;
387 G4bool existsAfterClip=
false;
402 existsAfterClip=
true;
423 existsAfterClip=
true;
429 flag = existsAfterClip ;
444 yt = std::fabs(yt1) ;
480 G4int noSurfaces = 0;
495 xshift = newpx - newpy*
fTalpha;
498 distx = std::fabs(std::fabs(xshift)-
fDx);
499 disty = std::fabs(std::fabs(newpy)-
fDy);
500 distz = std::fabs(std::fabs(p.z())-
fDz);
503 cosntheta = 1/std::sqrt(1+tntheta*tntheta);
515 if ( xshift >= 0.) {sumnorm += nX;}
516 else {sumnorm -= nX;}
521 if ( newpy >= 0.) {sumnorm += nY;}
522 else {sumnorm -= nY;}
527 if ( p.z() >= 0.) {sumnorm += nZ;}
528 else {sumnorm -= nZ;}
530 if ( noSurfaces == 0 )
533 G4Exception(
"G4Para::SurfaceNormal(p)",
"GeomSolids1002",
538 else if ( noSurfaces == 1 ) {norm = sumnorm;}
539 else {norm = sumnorm.unit();}
566 xshift=newpx*calpha+newpy*salpha;
568 distx=std::fabs(std::fabs(xshift)-
fDx*calpha);
569 disty=std::fabs(std::fabs(newpy)-
fDy);
570 distz=std::fabs(std::fabs(p.z())-
fDz);
574 if (distx<distz) {side=
kNX;}
579 if (disty<distz) {side=
kNY;}
589 cosntheta=-1/std::sqrt(1+tntheta*tntheta);
593 cosntheta=1/std::sqrt(1+tntheta*tntheta);
595 norm=
G4ThreeVector(calpha*cosntheta,salpha*cosntheta,-tntheta*cosntheta);
655 smin=(-
fDz-p.z())/v.z();
668 smin=(
fDz-p.z())/v.z();
677 if (std::fabs(p.z())<=
fDz)
723 if (std::fabs(yt)<=
fDy)
736 if (tmin>smin) smin=tmin;
737 if (tmax<smax) smax=tmax;
777 if (std::fabs(xt)<=
fDx)
787 if (tmin>smin) smin=tmin;
788 if (tmax<smax) smax=tmax;
791 if (smax>0&&smin<smax)
817 G4double distz1,distz2,disty1,disty2,distx1,distx2;
838 disty1=(trany-
fDy)*cosy;
839 disty2=(-
fDy-trany)*cosy;
841 if (disty1>safe) safe=disty1;
842 if (disty2>safe) safe=disty2;
846 distx1=(tranx-
fDx)*cosx;
847 distx2=(-
fDx-tranx)*cosx;
849 if (distx1>safe) safe=distx1;
850 if (distx2>safe) safe=distx2;
870 G4double ycomp,calpha,salpha,tntheta,cosntheta;
997 cosntheta=1/std::sqrt(1+tntheta*tntheta);
998 *n=
G4ThreeVector(calpha*cosntheta,salpha*cosntheta,-tntheta*cosntheta);
1023 cosntheta=-1/std::sqrt(1+tntheta*tntheta);
1024 *n=
G4ThreeVector(calpha*cosntheta,salpha*cosntheta,-tntheta*cosntheta);
1046 ycomp=1/std::sqrt(1+fTthetaSphi*fTthetaSphi);
1053 cosntheta=-1/std::sqrt(1+tntheta*tntheta);
1054 *n=
G4ThreeVector(calpha*cosntheta,salpha*cosntheta,-tntheta*cosntheta);
1057 calpha=1/std::sqrt(1+fTalpha*fTalpha);
1060 cosntheta=1/std::sqrt(1+tntheta*tntheta);
1061 *n=
G4ThreeVector(calpha*cosntheta,salpha*cosntheta,-tntheta*cosntheta);
1067 "Undefined side for valid surface normal to solid.");
1082 G4double distz1,distz2,disty1,disty2,distx1,distx2;
1095 G4cout.precision(oldprc) ;
1096 G4Exception(
"G4Para::DistanceToOut(p)",
"GeomSolids1002",
1119 disty1=(
fDy-trany)*cosy;
1120 disty2=(
fDy+trany)*cosy;
1122 if (disty1<safe) safe=disty1;
1123 if (disty2<safe) safe=disty2;
1127 distx1=(
fDx-tranx)*cosx;
1128 distx2=(
fDx+tranx)*cosx;
1130 if (distx1<safe) safe=distx1;
1131 if (distx2<safe) safe=distx2;
1153 vertices->reserve(8);
1185 "Error in allocation of vertices. Out of memory !");
1205 return new G4Para(*
this);
1214 G4int oldprc = os.precision(16);
1215 os <<
"-----------------------------------------------------------\n"
1216 <<
" *** Dump for solid - " <<
GetName() <<
" ***\n"
1217 <<
" ===================================================\n"
1218 <<
" Solid type: G4Para\n"
1219 <<
" Parameters: \n"
1220 <<
" half length X: " <<
fDx/
mm <<
" mm \n"
1221 <<
" half length Y: " <<
fDy/
mm <<
" mm \n"
1222 <<
" half length Z: " <<
fDz/
mm <<
" mm \n"
1228 <<
"-----------------------------------------------------------\n";
1229 os.precision(oldprc);
1244 G4double lambda1, lambda2, chose, aOne, aTwo;
1253 w.z()*v.x() - w.x()*v.z(),
1254 w.x()*v.y() - w.y()*v.x());
1256 aOne = 0.5*Area.mag();
1259 t.z()*u.x() - t.x()*u.z(),
1260 t.x()*u.y() - t.y()*u.x());
1262 aTwo = 0.5*Area.mag();
1268 if( (chose>=0.) && (chose < aOne) )
1272 return (p2+lambda1*v+lambda2*w);
1279 return (p0+lambda1*t+lambda2*u);
1293 G4double chose, aOne, aTwo, aThree, aFour, aFive, aSix;
1323 if( (chose>=0.) && (chose<aOne) )
1325 else if(chose>=aOne && chose<aOne+aTwo)
1327 else if(chose>=aOne+aTwo && chose<aOne+aTwo+aThree)
1329 else if(chose>=aOne+aTwo+aThree && chose<aOne+aTwo+aThree+aFour)
1331 else if(chose>=aOne+aTwo+aThree+aFour && chose<aOne+aTwo+aThree+aFour+aFive)
G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=G4bool(false), G4bool *validNorm=0, G4ThreeVector *n=0) const
void ClipCrossSection(G4ThreeVectorList *pVertices, const G4int pSectionIndex, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
ThreeVector shoot(const G4int Ap, const G4int Af)
G4ThreeVector GetPointOnSurface() const
void ComputeDimensions(G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
EInside Inside(const G4ThreeVector &p) const
void SetAllParameters(G4double pDx, G4double pDy, G4double pDz, G4double pAlpha, G4double pTheta, G4double pPhi)
G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const
static const G4double kInfinity
G4double GetMinYExtent() const
CLHEP::Hep3Vector G4ThreeVector
G4ThreeVector ApproxSurfaceNormal(const G4ThreeVector &p) const
G4bool IsYLimited() const
G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
G4Para(const G4String &pName, G4double pDx, G4double pDy, G4double pDz, G4double pAlpha, G4double pTheta, G4double pPhi)
G4bool IsXLimited() const
std::ostream & StreamInfo(std::ostream &os) const
G4Para & operator=(const G4Para &rhs)
virtual void AddSolid(const G4Box &)=0
G4Polyhedron * CreatePolyhedron() const
static double normal(HepRandomEngine *eptr)
G4double GetMaxXExtent() const
G4double GetMinZExtent() const
G4GLOB_DLL std::ostream G4cout
std::vector< G4ThreeVector > G4ThreeVectorList
G4GeometryType GetEntityType() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
virtual void ComputeDimensions(G4Box &, const G4int, const G4VPhysicalVolume *) const
G4double GetMinXExtent() const
void DescribeYourselfTo(G4VGraphicsScene &scene) const
G4double GetMaxZExtent() const
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4ThreeVectorList * CreateRotatedVertices(const G4AffineTransform &pTransform) const
static const double degree
G4double GetMaxYExtent() const
G4double GetMaxExtent(const EAxis pAxis) const
static const G4double alpha
G4CSGSolid & operator=(const G4CSGSolid &rhs)
G4bool IsZLimited() const
void ClipBetweenSections(G4ThreeVectorList *pVertices, const G4int pSectionIndex, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
G4double GetMinExtent(const EAxis pAxis) const
G4ThreeVector GetPointOnPlane(G4ThreeVector p0, G4ThreeVector p1, G4ThreeVector p2, G4ThreeVector p3, G4double &area) const