42 #if !defined(G4GEOM_USE_UPOLYCONE)
54 using namespace CLHEP;
81 for (i=0; i<numZPlanes; i++)
83 if(rInner[i]>rOuter[i])
86 std::ostringstream message;
87 message <<
"Cannot create a Polycone with rInner > rOuter for the same Z"
89 <<
" rInner > rOuter for the same Z !" <<
G4endl
90 <<
" rMin[" << i <<
"] = " << rInner[i]
91 <<
" -- rMax[" << i <<
"] = " << rOuter[i];
92 G4Exception(
"G4Polycone::G4Polycone()",
"GeomSolids0002",
95 if (( i < numZPlanes-1) && ( zPlane[i] == zPlane[i+1] ))
97 if( (rInner[i] > rOuter[i+1])
98 ||(rInner[i+1] > rOuter[i]) )
101 std::ostringstream message;
102 message <<
"Cannot create a Polycone with no contiguous segments."
104 <<
" Segments are not contiguous !" <<
G4endl
105 <<
" rMin[" << i <<
"] = " << rInner[i]
106 <<
" -- rMax[" << i+1 <<
"] = " << rOuter[i+1] <<
G4endl
107 <<
" rMin[" << i+1 <<
"] = " << rInner[i+1]
108 <<
" -- rMax[" << i <<
"] = " << rOuter[i];
109 G4Exception(
"G4Polycone::G4Polycone()",
"GeomSolids0002",
127 Create( phiStart, phiTotal, rz );
147 Create( phiStart, phiTotal, rz );
156 std::ostringstream message;
157 message <<
"Polycone " <<
GetName() <<
"cannot be converted" <<
G4endl
158 <<
"to Polycone with (Rmin,Rmaz,Z) parameters!";
159 G4Exception(
"G4Polycone::G4Polycone()",
"GeomSolids0002",
165 <<
"to optimized polycone with (Rmin,Rmaz,Z) parameters !"
184 if (rz->
Amin() < 0.0)
186 std::ostringstream message;
187 message <<
"Illegal input parameters - " <<
GetName() <<
G4endl
188 <<
" All R values must be >= 0 !";
189 G4Exception(
"G4Polycone::Create()",
"GeomSolids0002",
199 std::ostringstream message;
200 message <<
"Illegal input parameters - " <<
GetName() <<
G4endl
201 <<
" R/Z cross section is zero or near zero: " << rzArea;
202 G4Exception(
"G4Polycone::Create()",
"GeomSolids0002",
209 std::ostringstream message;
210 message <<
"Illegal input parameters - " <<
GetName() <<
G4endl
211 <<
" Too few unique R/Z values !";
212 G4Exception(
"G4Polycone::Create()",
"GeomSolids0002",
218 std::ostringstream message;
219 message <<
"Illegal input parameters - " <<
GetName() <<
G4endl
220 <<
" R/Z segments cross !";
221 G4Exception(
"G4Polycone::Create()",
"GeomSolids0002",
231 if (phiTotal <= 0 || phiTotal >
twopi-1E-10)
247 endPhi = phiStart+phiTotal;
265 next->
r = iterRZ.
GetA();
266 next->
z = iterRZ.
GetB();
267 }
while( ++next, iterRZ.
Next() );
291 if (corner->
r < 1/kInfinity && next->
r < 1/kInfinity)
continue;
299 if (corner->
z > next->
z)
316 }
while( prev=corner, corner=next, corner >
corners );
345 :
G4VCSGfaceted(a), startPhi(0.), endPhi(0.), phiIsOpen(false),
346 numCorner(0), corners(0),original_parameters(0),
378 if (
this == &source)
return *
this;
544 G4int oldprc = os.precision(16);
545 os <<
"-----------------------------------------------------------\n"
546 <<
" *** Dump for solid - " <<
GetName() <<
" ***\n"
547 <<
" ===================================================\n"
548 <<
" Solid type: G4Polycone\n"
551 <<
" ending phi angle : " <<
endPhi/
degree <<
" degrees \n";
555 os <<
" number of Z planes: " << numPlanes <<
"\n"
557 for (i=0; i<numPlanes; i++)
559 os <<
" Z plane " << i <<
": "
562 os <<
" Tangent distances to inner surface (Rmin): \n";
563 for (i=0; i<numPlanes; i++)
565 os <<
" Z plane " << i <<
": "
568 os <<
" Tangent distances to outer surface (Rmax): \n";
569 for (i=0; i<numPlanes; i++)
571 os <<
" Z plane " << i <<
": "
575 os <<
" number of RZ points: " <<
numCorner <<
"\n"
576 <<
" RZ values (corners): \n";
582 os <<
"-----------------------------------------------------------\n";
583 os.precision(oldprc);
601 G4double Aone, Atwo, Afive, phi, zRand, fDPhi, cosu, sinu;
602 G4double rRand1, rmin, rmax, chose, rone, rtwo, qone, qtwo;
603 G4double fDz=(zTwo-zOne)/2., afDz=std::fabs(fDz);
606 rone = (fRmax1-fRmax2)/(2.*fDz);
607 rtwo = (fRmin1-fRmin2)/(2.*fDz);
608 if(fRmax1==fRmax2){qone=0.;}
610 qone = fDz*(fRmax1+fRmax2)/(fRmax1-fRmax2);
612 if(fRmin1==fRmin2){qtwo=0.;}
614 qtwo = fDz*(fRmin1+fRmin2)/(fRmin1-fRmin2);
616 Aone = 0.5*fDPhi*(fRmax2 + fRmax1)*(
sqr(fRmin1-fRmin2)+
sqr(zTwo-zOne));
617 Atwo = 0.5*fDPhi*(fRmin2 + fRmin1)*(
sqr(fRmax1-fRmax2)+
sqr(zTwo-zOne));
618 Afive = fDz*(fRmax1-fRmin1+fRmax2-fRmin2);
619 totArea = Aone+Atwo+2.*Afive;
622 cosu = std::cos(phi);
623 sinu = std::sin(phi);
626 if( (startPhi == 0) && (
endPhi ==
twopi) ) { Afive = 0; }
628 if( (chose >= 0) && (chose < Aone) )
634 rone*sinu*(qone-zRand), zRand);
643 else if(chose >= Aone && chose < Aone + Atwo)
649 rtwo*sinu*(qtwo-zRand), zRand);
658 else if( (chose >= Aone + Atwo + Afive) && (chose < Aone + Atwo + 2.*Afive) )
661 rmin = fRmin2-((zRand-fDz)/(2.*fDz))*(fRmin1-fRmin2);
662 rmax = fRmax2-((zRand-fDz)/(2.*fDz))*(fRmax1-fRmax2);
665 rRand1*std::sin(startPhi), zRand);
670 rmin = fRmin2-((zRand-fDz)/(2.*fDz))*(fRmin1-fRmin2);
671 rmax = fRmax2-((zRand-fDz)/(2.*fDz))*(fRmax1-fRmax2);
674 rRand1*std::sin(
endPhi), zRand);
691 G4double xRand,yRand,zRand,phi,cosphi,sinphi,chose,
692 aOne,aTwo,aFou,rRand,fDz,fSPhi,fDPhi;
693 fDz = std::fabs(0.5*(zTwo-zOne));
697 aOne = 2.*fDz*fDPhi*fRMax;
698 aTwo = 2.*fDz*fDPhi*fRMin;
699 aFou = 2.*fDz*(fRMax-fRMin);
700 totArea = aOne+aTwo+2.*aFou;
702 cosphi = std::cos(phi);
703 sinphi = std::sin(phi);
710 if( (chose >= 0) && (chose < aOne) )
712 xRand = fRMax*cosphi;
713 yRand = fRMax*sinphi;
717 else if( (chose >= aOne) && (chose < aOne + aTwo) )
719 xRand = fRMin*cosphi;
720 yRand = fRMin*sinphi;
724 else if( (chose >= aOne+aTwo) && (chose <aOne+aTwo+aFou) )
726 xRand = rRand*std::cos(fSPhi+fDPhi);
727 yRand = rRand*std::sin(fSPhi+fDPhi);
734 xRand = rRand*std::cos(fSPhi+fDPhi);
735 yRand = rRand*std::sin(fSPhi+fDPhi);
750 G4double xRand,yRand,phi,cosphi,sinphi,rRand1,rRand2,A1,Atot,rCh;
752 cosphi = std::cos(phi);
753 sinphi = std::sin(phi);
757 rRand1 = fRMin1; A1=0.;
762 A1=std::fabs(fRMin2*fRMin2-fRMin1*fRMin1);
766 rRand2=fRMax1; Atot=A1;
771 Atot = A1+std::fabs(fRMax2*fRMax2-fRMax1*fRMax1);
775 if(rCh>A1) { rRand1=rRand2; }
777 xRand = rRand1*cosphi;
778 yRand = rRand1*sinphi;
797 if( (fRMin1 == fRMin2) && (fRMax1 == fRMax2) )
801 return GetPointOnCone(fRMin1,fRMax1,fRMin2,fRMax2,zOne,zTwo,totArea);
810 G4double Area=0,totArea=0,Achose1=0,Achose2=0,phi,cosphi,sinphi,rRand;
815 cosphi = std::cos(phi);
816 sinphi = std::sin(phi);
822 std::vector<G4double> areas;
823 std::vector<G4ThreeVector> points;
828 for(i=0; i<numPlanes-1; i++)
853 areas.push_back(Area);
860 totArea += (areas[0]+areas[numPlanes]);
863 if( (chose>=0.) && (chose<areas[0]) )
869 for (i=0; i<numPlanes-1; i++)
872 Achose2 = (Achose1+areas[i+1]);
873 if(chose>=Achose1 && chose<Achose2)
913 G4bool isConvertible=
true;
919 std::vector<G4double>
Z;
920 std::vector<G4double> Rmin;
921 std::vector<G4double> Rmax;
934 Rmax.push_back (
corners[1].r);icurr=1;
936 else if (Zprev ==
corners[numPlanes-1].z)
938 Rmin.push_back(
corners[numPlanes-1].
r);
950 G4int inextr=0, inextl=0;
951 for (
G4int i=0; i < numPlanes-2; i++)
954 inextl=(icurl <= 0)? numPlanes-1 : icurl-1;
971 Rmin.push_back(
corners[icurl].r);
972 Rmax.push_back(Rmax[countPlanes-2]);
973 Rmax[countPlanes-2]=
corners[icurl].
r;
977 Rmin.push_back(
corners[inextl].r);
978 Rmax.push_back(
corners[icurl].r);
989 isConvertible=
false;
break;
991 icurl=(icurl == 0)? numPlanes-1 : icurl-1;
999 icurl=(icurl == 0)? numPlanes-1 : icurl-1;
1002 Rmax.push_back (
corners[inextr].r);
1006 Z.push_back(Zright);
1015 Rmin.push_back(
corners[icurr].r);
1016 Rmax.push_back(
corners[inextr].r);
1020 Rmin.push_back(
corners[inextr].r);
1021 Rmax.push_back(
corners[icurr].r);
1022 Rmax[countPlanes-2]=
corners[inextr].
r;
1031 Rmin.push_back (
corners[icurr].r);
1043 isConvertible=
false;
break;
1053 inextl=(icurl <= 0)? numPlanes-1 : icurl-1;
1058 Rmin.push_back(
corners[inextl].r);
1063 Rmin.push_back(
corners[inextl].r);
1075 for(
G4int j=0; j < countPlanes; j++)
1088 std::ostringstream message;
1090 <<
"cannot be converted to Polycone with (Rmin,Rmaz,Z) parameters!";
1091 G4Exception(
"G4Polycone::SetOriginalParameters()",
"GeomSolids0002",
1099 for(
G4int j=0; j < numPlanes; j++)
1109 return isConvertible;
G4bool CrossesItself(G4double tolerance)
G4Polyhedron * CreatePolyhedron() const
ThreeVector shoot(const G4int Ap, const G4int Af)
void ComputeDimensions(G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
CLHEP::Hep3Vector G4ThreeVector
G4ThreeVector GetPointOnSurface() const
G4ThreeVector GetPointOnRing(G4double fRMin, G4double fRMax, G4double fRMin2, G4double fRMax2, G4double zOne) const
const G4Polycone & operator=(const G4Polycone &source)
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
G4bool MustBeOutside(const G4ThreeVector &p) const
G4PolyconeSideRZ * corners
std::ostream & StreamInfo(std::ostream &os) const
G4ThreeVector GetPointOnCone(G4double fRmin1, G4double fRmax1, G4double fRmin2, G4double fRmax2, G4double zOne, G4double zTwo, G4double &totArea) const
G4GeometryType GetEntityType() const
G4ThreeVector GetPointOnCut(G4double fRMin1, G4double fRMax1, G4double fRMin2, G4double fRMax2, G4double zOne, G4double zTwo, G4double &totArea) const
G4bool RemoveDuplicateVertices(G4double tolerance)
G4bool RemoveRedundantVertices(G4double tolerance)
G4GLOB_DLL std::ostream G4cout
G4int NumVertices() const
void SetOriginalParameters(G4PolyconeHistorical *pars)
const G4VCSGfaceted & operator=(const G4VCSGfaceted &source)
EInside Inside(const G4ThreeVector &p) const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
virtual void ComputeDimensions(G4Box &, const G4int, const G4VPhysicalVolume *) const
void CopyStuff(const G4Polycone &source)
G4EnclosingCylinder * enclosingCylinder
G4ThreeVector GetPointOnTubs(G4double fRMin, G4double fRMax, G4double zOne, G4double zTwo, G4double &totArea) const
G4bool BisectedBy(G4double a1, G4double b1, G4double a2, G4double b2, G4double tolerance)
G4Polycone(const G4String &name, G4double phiStart, G4double phiTotal, G4int numZPlanes, const G4double zPlane[], const G4double rInner[], const G4double rOuter[])
void Create(G4double phiStart, G4double phiTotal, G4ReduciblePolygon *rz)
G4PolyconeHistorical * original_parameters
virtual EInside Inside(const G4ThreeVector &p) const
G4bool ShouldMiss(const G4ThreeVector &p, const G4ThreeVector &v) const