#include <G4PolyconeSide.hh>
Inheritance diagram for G4PolyconeSide:
Definition at line 67 of file G4PolyconeSide.hh.
G4PolyconeSide::G4PolyconeSide | ( | const G4PolyconeSideRZ * | prevRZ, | |
const G4PolyconeSideRZ * | tail, | |||
const G4PolyconeSideRZ * | head, | |||
const G4PolyconeSideRZ * | nextRZ, | |||
G4double | phiStart, | |||
G4double | deltaPhi, | |||
G4bool | phiIsOpen, | |||
G4bool | isAllBehind = false | |||
) |
Definition at line 57 of file G4PolyconeSide.cc.
References allBehind, cone, corners, deltaPhi, G4GeometryTolerance::GetInstance(), G4GeometryTolerance::GetSurfaceTolerance(), length, ncorners, nextRS, nextZS, phiIsOpen, prevRS, prevZS, G4PolyconeSideRZ::r, r, rNorm, rNormEdge, rS, startPhi, G4PolyconeSideRZ::z, z, zNorm, zNormEdge, and zS.
Referenced by Clone().
00065 : ncorners(0), corners(0) 00066 { 00067 kCarTolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance(); 00068 fSurfaceArea = 0.0; 00069 fPhi.first = G4ThreeVector(0,0,0); 00070 fPhi.second= 0.0; 00071 00072 // 00073 // Record values 00074 // 00075 r[0] = tail->r; z[0] = tail->z; 00076 r[1] = head->r; z[1] = head->z; 00077 00078 phiIsOpen = thePhiIsOpen; 00079 if (phiIsOpen) 00080 { 00081 deltaPhi = theDeltaPhi; 00082 startPhi = thePhiStart; 00083 00084 // 00085 // Set phi values to our conventions 00086 // 00087 while (deltaPhi < 0.0) deltaPhi += twopi; 00088 while (startPhi < 0.0) startPhi += twopi; 00089 00090 // 00091 // Calculate corner coordinates 00092 // 00093 ncorners = 4; 00094 corners = new G4ThreeVector[ncorners]; 00095 00096 corners[0] = G4ThreeVector( tail->r*std::cos(startPhi), 00097 tail->r*std::sin(startPhi), tail->z ); 00098 corners[1] = G4ThreeVector( head->r*std::cos(startPhi), 00099 head->r*std::sin(startPhi), head->z ); 00100 corners[2] = G4ThreeVector( tail->r*std::cos(startPhi+deltaPhi), 00101 tail->r*std::sin(startPhi+deltaPhi), tail->z ); 00102 corners[3] = G4ThreeVector( head->r*std::cos(startPhi+deltaPhi), 00103 head->r*std::sin(startPhi+deltaPhi), head->z ); 00104 } 00105 else 00106 { 00107 deltaPhi = twopi; 00108 startPhi = 0.0; 00109 } 00110 00111 allBehind = isAllBehind; 00112 00113 // 00114 // Make our intersecting cone 00115 // 00116 cone = new G4IntersectingCone( r, z ); 00117 00118 // 00119 // Calculate vectors in r,z space 00120 // 00121 rS = r[1]-r[0]; zS = z[1]-z[0]; 00122 length = std::sqrt( rS*rS + zS*zS); 00123 rS /= length; zS /= length; 00124 00125 rNorm = +zS; 00126 zNorm = -rS; 00127 00128 G4double lAdj; 00129 00130 prevRS = r[0]-prevRZ->r; 00131 prevZS = z[0]-prevRZ->z; 00132 lAdj = std::sqrt( prevRS*prevRS + prevZS*prevZS ); 00133 prevRS /= lAdj; 00134 prevZS /= lAdj; 00135 00136 rNormEdge[0] = rNorm + prevZS; 00137 zNormEdge[0] = zNorm - prevRS; 00138 lAdj = std::sqrt( rNormEdge[0]*rNormEdge[0] + zNormEdge[0]*zNormEdge[0] ); 00139 rNormEdge[0] /= lAdj; 00140 zNormEdge[0] /= lAdj; 00141 00142 nextRS = nextRZ->r-r[1]; 00143 nextZS = nextRZ->z-z[1]; 00144 lAdj = std::sqrt( nextRS*nextRS + nextZS*nextZS ); 00145 nextRS /= lAdj; 00146 nextZS /= lAdj; 00147 00148 rNormEdge[1] = rNorm + nextZS; 00149 zNormEdge[1] = zNorm - nextRS; 00150 lAdj = std::sqrt( rNormEdge[1]*rNormEdge[1] + zNormEdge[1]*zNormEdge[1] ); 00151 rNormEdge[1] /= lAdj; 00152 zNormEdge[1] /= lAdj; 00153 }
G4PolyconeSide::~G4PolyconeSide | ( | ) | [virtual] |
G4PolyconeSide::G4PolyconeSide | ( | const G4PolyconeSide & | source | ) |
Definition at line 186 of file G4PolyconeSide.cc.
References CopyStuff().
00187 : G4VCSGface(), ncorners(0), corners(0) 00188 { 00189 CopyStuff( source ); 00190 }
G4PolyconeSide::G4PolyconeSide | ( | __void__ & | ) |
Definition at line 160 of file G4PolyconeSide.cc.
References r, rNormEdge, z, and zNormEdge.
00161 : startPhi(0.), deltaPhi(0.), phiIsOpen(false), allBehind(false), 00162 cone(0), rNorm(0.), zNorm(0.), rS(0.), zS(0.), length(0.), 00163 prevRS(0.), prevZS(0.), nextRS(0.), nextZS(0.), ncorners(0), corners(0), 00164 kCarTolerance(0.), fSurfaceArea(0.) 00165 { 00166 r[0] = r[1] = 0.; 00167 z[0] = z[1] = 0.; 00168 rNormEdge[0]= rNormEdge[1] = 0.; 00169 zNormEdge[0]= zNormEdge[1] = 0.; 00170 }
void G4PolyconeSide::CalculateExtent | ( | const EAxis | axis, | |
const G4VoxelLimits & | voxelLimit, | |||
const G4AffineTransform & | tranform, | |||
G4SolidExtentList & | extentList | |||
) | [virtual] |
Implements G4VCSGface.
Definition at line 552 of file G4PolyconeSide.cc.
References G4SolidExtentList::AddSurface(), G4ClippablePolygon::AddVertexInOrder(), G4AffineTransform::ApplyPointTransform(), G4ClippablePolygon::ClearAllVertices(), DBL_MIN, deltaPhi, FindLineIntersect(), kMaxMeshSections, kMeshAngleDefault, kMinMeshSections, nextRS, nextZS, G4ClippablePolygon::PartialClip(), phiIsOpen, prevRS, prevZS, r, rNorm, rS, G4ClippablePolygon::SetNormal(), startPhi, G4AffineTransform::TransformAxis(), z, and zS.
00556 { 00557 G4ClippablePolygon polygon; 00558 00559 // 00560 // Here we will approximate (ala G4Cons) and divide our conical section 00561 // into segments, like G4Polyhedra. When doing so, the radius 00562 // is extented far enough such that the segments always lie 00563 // just outside the surface of the conical section we are 00564 // approximating. 00565 // 00566 00567 // 00568 // Choose phi size of our segment(s) based on constants as 00569 // defined in meshdefs.hh 00570 // 00571 G4int numPhi = (G4int)(deltaPhi/kMeshAngleDefault) + 1; 00572 if (numPhi < kMinMeshSections) 00573 numPhi = kMinMeshSections; 00574 else if (numPhi > kMaxMeshSections) 00575 numPhi = kMaxMeshSections; 00576 00577 G4double sigPhi = deltaPhi/numPhi; 00578 00579 // 00580 // Determine radius factor to keep segments outside 00581 // 00582 G4double rFudge = 1.0/std::cos(0.5*sigPhi); 00583 00584 // 00585 // Decide which radius to use on each end of the side, 00586 // and whether a transition mesh is required 00587 // 00588 // {r0,z0} - Beginning of this side 00589 // {r1,z1} - Ending of this side 00590 // {r2,z0} - Beginning of transition piece connecting previous 00591 // side (and ends at beginning of this side) 00592 // 00593 // So, order is 2 --> 0 --> 1. 00594 // ------- 00595 // 00596 // r2 < 0 indicates that no transition piece is required 00597 // 00598 G4double r0, r1, r2, z0, z1; 00599 00600 r2 = -1; // By default: no transition piece 00601 00602 if (rNorm < -DBL_MIN) 00603 { 00604 // 00605 // This side faces *inward*, and so our mesh has 00606 // the same radius 00607 // 00608 r1 = r[1]; 00609 z1 = z[1]; 00610 z0 = z[0]; 00611 r0 = r[0]; 00612 00613 r2 = -1; 00614 00615 if (prevZS > DBL_MIN) 00616 { 00617 // 00618 // The previous side is facing outwards 00619 // 00620 if ( prevRS*zS - prevZS*rS > 0 ) 00621 { 00622 // 00623 // Transition was convex: build transition piece 00624 // 00625 if (r[0] > DBL_MIN) r2 = r[0]*rFudge; 00626 } 00627 else 00628 { 00629 // 00630 // Transition was concave: short this side 00631 // 00632 FindLineIntersect( z0, r0, zS, rS, 00633 z0, r0*rFudge, prevZS, prevRS*rFudge, z0, r0 ); 00634 } 00635 } 00636 00637 if ( nextZS > DBL_MIN && (rS*nextZS - zS*nextRS < 0) ) 00638 { 00639 // 00640 // The next side is facing outwards, forming a 00641 // concave transition: short this side 00642 // 00643 FindLineIntersect( z1, r1, zS, rS, 00644 z1, r1*rFudge, nextZS, nextRS*rFudge, z1, r1 ); 00645 } 00646 } 00647 else if (rNorm > DBL_MIN) 00648 { 00649 // 00650 // This side faces *outward* and is given a boost to 00651 // it radius 00652 // 00653 r0 = r[0]*rFudge; 00654 z0 = z[0]; 00655 r1 = r[1]*rFudge; 00656 z1 = z[1]; 00657 00658 if (prevZS < -DBL_MIN) 00659 { 00660 // 00661 // The previous side is facing inwards 00662 // 00663 if ( prevRS*zS - prevZS*rS > 0 ) 00664 { 00665 // 00666 // Transition was convex: build transition piece 00667 // 00668 if (r[0] > DBL_MIN) r2 = r[0]; 00669 } 00670 else 00671 { 00672 // 00673 // Transition was concave: short this side 00674 // 00675 FindLineIntersect( z0, r0, zS, rS*rFudge, 00676 z0, r[0], prevZS, prevRS, z0, r0 ); 00677 } 00678 } 00679 00680 if ( nextZS < -DBL_MIN && (rS*nextZS - zS*nextRS < 0) ) 00681 { 00682 // 00683 // The next side is facing inwards, forming a 00684 // concave transition: short this side 00685 // 00686 FindLineIntersect( z1, r1, zS, rS*rFudge, 00687 z1, r[1], nextZS, nextRS, z1, r1 ); 00688 } 00689 } 00690 else 00691 { 00692 // 00693 // This side is perpendicular to the z axis (is a disk) 00694 // 00695 // Whether or not r0 needs a rFudge factor depends 00696 // on the normal of the previous edge. Similar with r1 00697 // and the next edge. No transition piece is required. 00698 // 00699 r0 = r[0]; 00700 r1 = r[1]; 00701 z0 = z[0]; 00702 z1 = z[1]; 00703 00704 if (prevZS > DBL_MIN) r0 *= rFudge; 00705 if (nextZS > DBL_MIN) r1 *= rFudge; 00706 } 00707 00708 // 00709 // Loop 00710 // 00711 G4double phi = startPhi, 00712 cosPhi = std::cos(phi), 00713 sinPhi = std::sin(phi); 00714 00715 G4ThreeVector v0( r0*cosPhi, r0*sinPhi, z0 ), 00716 v1( r1*cosPhi, r1*sinPhi, z1 ), 00717 v2, w0, w1, w2; 00718 transform.ApplyPointTransform( v0 ); 00719 transform.ApplyPointTransform( v1 ); 00720 00721 if (r2 >= 0) 00722 { 00723 v2 = G4ThreeVector( r2*cosPhi, r2*sinPhi, z0 ); 00724 transform.ApplyPointTransform( v2 ); 00725 } 00726 00727 do 00728 { 00729 phi += sigPhi; 00730 if (numPhi == 1) phi = startPhi+deltaPhi; // Try to avoid roundoff 00731 cosPhi = std::cos(phi), 00732 sinPhi = std::sin(phi); 00733 00734 w0 = G4ThreeVector( r0*cosPhi, r0*sinPhi, z0 ); 00735 w1 = G4ThreeVector( r1*cosPhi, r1*sinPhi, z1 ); 00736 transform.ApplyPointTransform( w0 ); 00737 transform.ApplyPointTransform( w1 ); 00738 00739 G4ThreeVector deltaV = r0 > r1 ? w0-v0 : w1-v1; 00740 00741 // 00742 // Build polygon, taking special care to keep the vertices 00743 // in order 00744 // 00745 polygon.ClearAllVertices(); 00746 00747 polygon.AddVertexInOrder( v0 ); 00748 polygon.AddVertexInOrder( v1 ); 00749 polygon.AddVertexInOrder( w1 ); 00750 polygon.AddVertexInOrder( w0 ); 00751 00752 // 00753 // Get extent 00754 // 00755 if (polygon.PartialClip( voxelLimit, axis )) 00756 { 00757 // 00758 // Get dot product of normal with target axis 00759 // 00760 polygon.SetNormal( deltaV.cross(v1-v0).unit() ); 00761 00762 extentList.AddSurface( polygon ); 00763 } 00764 00765 if (r2 >= 0) 00766 { 00767 // 00768 // Repeat, for transition piece 00769 // 00770 w2 = G4ThreeVector( r2*cosPhi, r2*sinPhi, z0 ); 00771 transform.ApplyPointTransform( w2 ); 00772 00773 polygon.ClearAllVertices(); 00774 00775 polygon.AddVertexInOrder( v2 ); 00776 polygon.AddVertexInOrder( v0 ); 00777 polygon.AddVertexInOrder( w0 ); 00778 polygon.AddVertexInOrder( w2 ); 00779 00780 if (polygon.PartialClip( voxelLimit, axis )) 00781 { 00782 polygon.SetNormal( deltaV.cross(v0-v2).unit() ); 00783 00784 extentList.AddSurface( polygon ); 00785 } 00786 00787 v2 = w2; 00788 } 00789 00790 // 00791 // Next vertex 00792 // 00793 v0 = w0; 00794 v1 = w1; 00795 } while( --numPhi > 0 ); 00796 00797 // 00798 // We are almost done. But, it is important that we leave no 00799 // gaps in the surface of our solid. By using rFudge, however, 00800 // we've done exactly that, if we have a phi segment. 00801 // Add two additional faces if necessary 00802 // 00803 if (phiIsOpen && rNorm > DBL_MIN) 00804 { 00805 cosPhi = std::cos(startPhi); 00806 sinPhi = std::sin(startPhi); 00807 00808 G4ThreeVector a0( r[0]*cosPhi, r[0]*sinPhi, z[0] ), 00809 a1( r[1]*cosPhi, r[1]*sinPhi, z[1] ), 00810 b0( r0*cosPhi, r0*sinPhi, z[0] ), 00811 b1( r1*cosPhi, r1*sinPhi, z[1] ); 00812 00813 transform.ApplyPointTransform( a0 ); 00814 transform.ApplyPointTransform( a1 ); 00815 transform.ApplyPointTransform( b0 ); 00816 transform.ApplyPointTransform( b1 ); 00817 00818 polygon.ClearAllVertices(); 00819 00820 polygon.AddVertexInOrder( a0 ); 00821 polygon.AddVertexInOrder( a1 ); 00822 polygon.AddVertexInOrder( b0 ); 00823 polygon.AddVertexInOrder( b1 ); 00824 00825 if (polygon.PartialClip( voxelLimit , axis)) 00826 { 00827 G4ThreeVector normal( sinPhi, -cosPhi, 0 ); 00828 polygon.SetNormal( transform.TransformAxis( normal ) ); 00829 00830 extentList.AddSurface( polygon ); 00831 } 00832 00833 cosPhi = std::cos(startPhi+deltaPhi); 00834 sinPhi = std::sin(startPhi+deltaPhi); 00835 00836 a0 = G4ThreeVector( r[0]*cosPhi, r[0]*sinPhi, z[0] ), 00837 a1 = G4ThreeVector( r[1]*cosPhi, r[1]*sinPhi, z[1] ), 00838 b0 = G4ThreeVector( r0*cosPhi, r0*sinPhi, z[0] ), 00839 b1 = G4ThreeVector( r1*cosPhi, r1*sinPhi, z[1] ); 00840 transform.ApplyPointTransform( a0 ); 00841 transform.ApplyPointTransform( a1 ); 00842 transform.ApplyPointTransform( b0 ); 00843 transform.ApplyPointTransform( b1 ); 00844 00845 polygon.ClearAllVertices(); 00846 00847 polygon.AddVertexInOrder( a0 ); 00848 polygon.AddVertexInOrder( a1 ); 00849 polygon.AddVertexInOrder( b0 ); 00850 polygon.AddVertexInOrder( b1 ); 00851 00852 if (polygon.PartialClip( voxelLimit, axis )) 00853 { 00854 G4ThreeVector normal( -sinPhi, cosPhi, 0 ); 00855 polygon.SetNormal( transform.TransformAxis( normal ) ); 00856 00857 extentList.AddSurface( polygon ); 00858 } 00859 } 00860 00861 return; 00862 }
G4VCSGface* G4PolyconeSide::Clone | ( | ) | [inline, virtual] |
Implements G4VCSGface.
Definition at line 101 of file G4PolyconeSide.hh.
References G4PolyconeSide().
00101 { return new G4PolyconeSide( *this ); }
void G4PolyconeSide::CopyStuff | ( | const G4PolyconeSide & | source | ) | [protected] |
Definition at line 212 of file G4PolyconeSide.cc.
References allBehind, cone, corners, deltaPhi, fSurfaceArea, kCarTolerance, length, ncorners, nextRS, nextZS, phiIsOpen, prevRS, prevZS, r, rNorm, rNormEdge, rS, startPhi, z, zNorm, zNormEdge, and zS.
Referenced by G4PolyconeSide(), and operator=().
00213 { 00214 r[0] = source.r[0]; 00215 r[1] = source.r[1]; 00216 z[0] = source.z[0]; 00217 z[1] = source.z[1]; 00218 00219 startPhi = source.startPhi; 00220 deltaPhi = source.deltaPhi; 00221 phiIsOpen = source.phiIsOpen; 00222 allBehind = source.allBehind; 00223 00224 kCarTolerance = source.kCarTolerance; 00225 fSurfaceArea = source.fSurfaceArea; 00226 00227 cone = new G4IntersectingCone( *source.cone ); 00228 00229 rNorm = source.rNorm; 00230 zNorm = source.zNorm; 00231 rS = source.rS; 00232 zS = source.zS; 00233 length = source.length; 00234 prevRS = source.prevRS; 00235 prevZS = source.prevZS; 00236 nextRS = source.nextRS; 00237 nextZS = source.nextZS; 00238 00239 rNormEdge[0] = source.rNormEdge[0]; 00240 rNormEdge[1] = source.rNormEdge[1]; 00241 zNormEdge[0] = source.zNormEdge[0]; 00242 zNormEdge[1] = source.zNormEdge[1]; 00243 00244 if (phiIsOpen) 00245 { 00246 ncorners = 4; 00247 corners = new G4ThreeVector[ncorners]; 00248 00249 corners[0] = source.corners[0]; 00250 corners[1] = source.corners[1]; 00251 corners[2] = source.corners[2]; 00252 corners[3] = source.corners[3]; 00253 } 00254 }
G4double G4PolyconeSide::Distance | ( | const G4ThreeVector & | p, | |
G4bool | outgoing | |||
) | [virtual] |
Implements G4VCSGface.
Definition at line 389 of file G4PolyconeSide.cc.
References DistanceAway().
00390 { 00391 G4double normSign = outgoing ? -1 : +1; 00392 G4double distFrom, distOut2; 00393 00394 // 00395 // We have two tries for each hemisphere. Try the closest first. 00396 // 00397 distFrom = normSign*DistanceAway( p, false, distOut2 ); 00398 if (distFrom > -0.5*kCarTolerance ) 00399 { 00400 // 00401 // Good answer 00402 // 00403 if (distOut2 > 0) 00404 return std::sqrt( distFrom*distFrom + distOut2 ); 00405 else 00406 return std::fabs(distFrom); 00407 } 00408 00409 // 00410 // Try second side. 00411 // 00412 distFrom = normSign*DistanceAway( p, true, distOut2 ); 00413 if (distFrom > -0.5*kCarTolerance) 00414 { 00415 00416 if (distOut2 > 0) 00417 return std::sqrt( distFrom*distFrom + distOut2 ); 00418 else 00419 return std::fabs(distFrom); 00420 } 00421 00422 return kInfinity; 00423 }
G4double G4PolyconeSide::DistanceAway | ( | const G4ThreeVector & | p, | |
G4bool | opposite, | |||
G4double & | distOutside2, | |||
G4double * | rzNorm = 0 | |||
) | [protected] |
Definition at line 906 of file G4PolyconeSide.cc.
References deltaPhi, GetPhi(), length, phiIsOpen, r, rNorm, rNormEdge, rS, sqr(), startPhi, z, zNorm, zNormEdge, and zS.
Referenced by Distance(), Inside(), Intersect(), and Normal().
00910 { 00911 // 00912 // Convert our point to r and z 00913 // 00914 G4double rx = p.perp(), zx = p.z(); 00915 00916 // 00917 // Change sign of r if opposite says we should 00918 // 00919 if (opposite) rx = -rx; 00920 00921 // 00922 // Calculate return value 00923 // 00924 G4double deltaR = rx - r[0], deltaZ = zx - z[0]; 00925 G4double answer = deltaR*rNorm + deltaZ*zNorm; 00926 00927 // 00928 // Are we off the surface in r,z space? 00929 // 00930 G4double q = deltaR*rS + deltaZ*zS; 00931 if (q < 0) 00932 { 00933 distOutside2 = q*q; 00934 if (edgeRZnorm) *edgeRZnorm = deltaR*rNormEdge[0] + deltaZ*zNormEdge[0]; 00935 } 00936 else if (q > length) 00937 { 00938 distOutside2 = sqr( q-length ); 00939 if (edgeRZnorm) 00940 { 00941 deltaR = rx - r[1]; 00942 deltaZ = zx - z[1]; 00943 *edgeRZnorm = deltaR*rNormEdge[1] + deltaZ*zNormEdge[1]; 00944 } 00945 } 00946 else 00947 { 00948 distOutside2 = 0; 00949 if (edgeRZnorm) *edgeRZnorm = answer; 00950 } 00951 00952 if (phiIsOpen) 00953 { 00954 // 00955 // Finally, check phi 00956 // 00957 G4double phi = GetPhi(p); 00958 while( phi < startPhi ) phi += twopi; 00959 00960 if (phi > startPhi+deltaPhi) 00961 { 00962 // 00963 // Oops. Are we closer to the start phi or end phi? 00964 // 00965 G4double d1 = phi-startPhi-deltaPhi; 00966 while( phi > startPhi ) phi -= twopi; 00967 G4double d2 = startPhi-phi; 00968 00969 if (d2 < d1) d1 = d2; 00970 00971 // 00972 // Add result to our distance 00973 // 00974 G4double dist = d1*rx; 00975 00976 distOutside2 += dist*dist; 00977 if (edgeRZnorm) 00978 { 00979 *edgeRZnorm = std::max(std::fabs(*edgeRZnorm),std::fabs(dist)); 00980 } 00981 } 00982 } 00983 00984 return answer; 00985 }
G4double G4PolyconeSide::Extent | ( | const G4ThreeVector | axis | ) | [virtual] |
Implements G4VCSGface.
Definition at line 488 of file G4PolyconeSide.cc.
References cone, DBL_MIN, deltaPhi, GetPhi(), phiIsOpen, r, startPhi, z, G4IntersectingCone::ZHi(), and G4IntersectingCone::ZLo().
00489 { 00490 if (axis.perp2() < DBL_MIN) 00491 { 00492 // 00493 // Special case 00494 // 00495 return axis.z() < 0 ? -cone->ZLo() : cone->ZHi(); 00496 } 00497 00498 // 00499 // Is the axis pointing inside our phi gap? 00500 // 00501 if (phiIsOpen) 00502 { 00503 G4double phi = GetPhi(axis); 00504 while( phi < startPhi ) phi += twopi; 00505 00506 if (phi > deltaPhi+startPhi) 00507 { 00508 // 00509 // Yeah, looks so. Make four three vectors defining the phi 00510 // opening 00511 // 00512 G4double cosP = std::cos(startPhi), sinP = std::sin(startPhi); 00513 G4ThreeVector a( r[0]*cosP, r[0]*sinP, z[0] ); 00514 G4ThreeVector b( r[1]*cosP, r[1]*sinP, z[1] ); 00515 cosP = std::cos(startPhi+deltaPhi); sinP = std::sin(startPhi+deltaPhi); 00516 G4ThreeVector c( r[0]*cosP, r[0]*sinP, z[0] ); 00517 G4ThreeVector d( r[1]*cosP, r[1]*sinP, z[1] ); 00518 00519 G4double ad = axis.dot(a), 00520 bd = axis.dot(b), 00521 cd = axis.dot(c), 00522 dd = axis.dot(d); 00523 00524 if (bd > ad) ad = bd; 00525 if (cd > ad) ad = cd; 00526 if (dd > ad) ad = dd; 00527 00528 return ad; 00529 } 00530 } 00531 00532 // 00533 // Check either end 00534 // 00535 G4double aPerp = axis.perp(); 00536 00537 G4double a = aPerp*r[0] + axis.z()*z[0]; 00538 G4double b = aPerp*r[1] + axis.z()*z[1]; 00539 00540 if (b > a) a = b; 00541 00542 return a; 00543 }
void G4PolyconeSide::FindLineIntersect | ( | G4double | x1, | |
G4double | y1, | |||
G4double | tx1, | |||
G4double | ty1, | |||
G4double | x2, | |||
G4double | y2, | |||
G4double | tx2, | |||
G4double | ty2, | |||
G4double & | x, | |||
G4double & | y | |||
) | [static, protected] |
Definition at line 1062 of file G4PolyconeSide.cc.
Referenced by CalculateExtent().
01067 { 01068 // 01069 // The solution is a simple linear equation 01070 // 01071 G4double deter = tx1*ty2 - tx2*ty1; 01072 01073 G4double s1 = ((x2-x1)*ty2 - tx2*(y2-y1))/deter; 01074 G4double s2 = ((x2-x1)*ty1 - tx1*(y2-y1))/deter; 01075 01076 // 01077 // We want the answer to not depend on which order the 01078 // lines were specified. Take average. 01079 // 01080 x = 0.5*( x1+s1*tx1 + x2+s2*tx2 ); 01081 y = 0.5*( y1+s1*ty1 + y2+s2*ty2 ); 01082 }
G4double G4PolyconeSide::GetPhi | ( | const G4ThreeVector & | p | ) | [protected] |
Definition at line 871 of file G4PolyconeSide.cc.
Referenced by DistanceAway(), Extent(), and PointOnCone().
00872 { 00873 G4double val=0.; 00874 00875 if (fPhi.first != p) 00876 { 00877 val = p.phi(); 00878 fPhi.first = p; 00879 fPhi.second = val; 00880 } 00881 else 00882 { 00883 val = fPhi.second; 00884 } 00885 return val; 00886 }
G4ThreeVector G4PolyconeSide::GetPointOnFace | ( | ) | [virtual] |
Implements G4VCSGface.
Definition at line 1100 of file G4PolyconeSide.cc.
References deltaPhi, G4UniformRand, r, startPhi, and z.
01101 { 01102 G4double x,y,zz; 01103 G4double rr,phi,dz,dr; 01104 dr=r[1]-r[0];dz=z[1]-z[0]; 01105 phi=startPhi+deltaPhi*G4UniformRand(); 01106 rr=r[0]+dr*G4UniformRand(); 01107 01108 x=rr*std::cos(phi); 01109 y=rr*std::sin(phi); 01110 01111 // PolyconeSide has a Ring Form 01112 // 01113 if (dz==0.) 01114 { 01115 zz=z[0]; 01116 } 01117 else 01118 { 01119 if(dr==0.) // PolyconeSide has a Tube Form 01120 { 01121 zz = z[0]+dz*G4UniformRand(); 01122 } 01123 else 01124 { 01125 zz = z[0]+(rr-r[0])*dz/dr; 01126 } 01127 } 01128 01129 return G4ThreeVector(x,y,zz); 01130 }
EInside G4PolyconeSide::Inside | ( | const G4ThreeVector & | p, | |
G4double | tolerance, | |||
G4double * | bestDistance | |||
) | [virtual] |
Implements G4VCSGface.
Definition at line 429 of file G4PolyconeSide.cc.
References DistanceAway(), kInside, kOutside, and kSurface.
00432 { 00433 // 00434 // Check both sides 00435 // 00436 G4double distFrom[2], distOut2[2], dist2[2]; 00437 G4double edgeRZnorm[2]; 00438 00439 distFrom[0] = DistanceAway( p, false, distOut2[0], edgeRZnorm ); 00440 distFrom[1] = DistanceAway( p, true, distOut2[1], edgeRZnorm+1 ); 00441 00442 dist2[0] = distFrom[0]*distFrom[0] + distOut2[0]; 00443 dist2[1] = distFrom[1]*distFrom[1] + distOut2[1]; 00444 00445 // 00446 // Who's closest? 00447 // 00448 G4int i = std::fabs(dist2[0]) < std::fabs(dist2[1]) ? 0 : 1; 00449 00450 *bestDistance = std::sqrt( dist2[i] ); 00451 00452 // 00453 // Okay then, inside or out? 00454 // 00455 if ( (std::fabs(edgeRZnorm[i]) < tolerance) 00456 && (distOut2[i] < tolerance*tolerance) ) 00457 return kSurface; 00458 else if (edgeRZnorm[i] < 0) 00459 return kInside; 00460 else 00461 return kOutside; 00462 }
G4bool G4PolyconeSide::Intersect | ( | const G4ThreeVector & | p, | |
const G4ThreeVector & | v, | |||
G4bool | outgoing, | |||
G4double | surfTolerance, | |||
G4double & | distance, | |||
G4double & | distFromSurface, | |||
G4ThreeVector & | normal, | |||
G4bool & | isAllBehind | |||
) | [virtual] |
Implements G4VCSGface.
Definition at line 260 of file G4PolyconeSide.cc.
References allBehind, cone, DBL_MIN, DistanceAway(), G4IntersectingCone::LineHitsCone(), PointOnCone(), rNorm, and zNorm.
00268 { 00269 G4double s1, s2; 00270 G4double normSign = outgoing ? +1 : -1; 00271 00272 isAllBehind = allBehind; 00273 00274 // 00275 // Check for two possible intersections 00276 // 00277 G4int nside = cone->LineHitsCone( p, v, &s1, &s2 ); 00278 if (nside == 0) return false; 00279 00280 // 00281 // Check the first side first, since it is (supposed to be) closest 00282 // 00283 G4ThreeVector hit = p + s1*v; 00284 00285 if (PointOnCone( hit, normSign, p, v, normal )) 00286 { 00287 // 00288 // Good intersection! What about the normal? 00289 // 00290 if (normSign*v.dot(normal) > 0) 00291 { 00292 // 00293 // We have a valid intersection, but it could very easily 00294 // be behind the point. To decide if we tolerate this, 00295 // we have to see if the point p is on the surface near 00296 // the intersecting point. 00297 // 00298 // What does it mean exactly for the point p to be "near" 00299 // the intersection? It means that if we draw a line from 00300 // p to the hit, the line remains entirely within the 00301 // tolerance bounds of the cone. To test this, we can 00302 // ask if the normal is correct near p. 00303 // 00304 G4double pr = p.perp(); 00305 if (pr < DBL_MIN) pr = DBL_MIN; 00306 G4ThreeVector pNormal( rNorm*p.x()/pr, rNorm*p.y()/pr, zNorm ); 00307 if (normSign*v.dot(pNormal) > 0) 00308 { 00309 // 00310 // p and intersection in same hemisphere 00311 // 00312 G4double distOutside2; 00313 distFromSurface = -normSign*DistanceAway( p, false, distOutside2 ); 00314 if (distOutside2 < surfTolerance*surfTolerance) 00315 { 00316 if (distFromSurface > -surfTolerance) 00317 { 00318 // 00319 // We are just inside or away from the 00320 // surface. Accept *any* value of distance. 00321 // 00322 distance = s1; 00323 return true; 00324 } 00325 } 00326 } 00327 else 00328 distFromSurface = s1; 00329 00330 // 00331 // Accept positive distances 00332 // 00333 if (s1 > 0) 00334 { 00335 distance = s1; 00336 return true; 00337 } 00338 } 00339 } 00340 00341 if (nside==1) return false; 00342 00343 // 00344 // Well, try the second hit 00345 // 00346 hit = p + s2*v; 00347 00348 if (PointOnCone( hit, normSign, p, v, normal )) 00349 { 00350 // 00351 // Good intersection! What about the normal? 00352 // 00353 if (normSign*v.dot(normal) > 0) 00354 { 00355 G4double pr = p.perp(); 00356 if (pr < DBL_MIN) pr = DBL_MIN; 00357 G4ThreeVector pNormal( rNorm*p.x()/pr, rNorm*p.y()/pr, zNorm ); 00358 if (normSign*v.dot(pNormal) > 0) 00359 { 00360 G4double distOutside2; 00361 distFromSurface = -normSign*DistanceAway( p, false, distOutside2 ); 00362 if (distOutside2 < surfTolerance*surfTolerance) 00363 { 00364 if (distFromSurface > -surfTolerance) 00365 { 00366 distance = s2; 00367 return true; 00368 } 00369 } 00370 } 00371 else 00372 distFromSurface = s2; 00373 00374 if (s2 > 0) 00375 { 00376 distance = s2; 00377 return true; 00378 } 00379 } 00380 } 00381 00382 // 00383 // Better luck next time 00384 // 00385 return false; 00386 }
G4ThreeVector G4PolyconeSide::Normal | ( | const G4ThreeVector & | p, | |
G4double * | bestDistance | |||
) | [virtual] |
Implements G4VCSGface.
Definition at line 468 of file G4PolyconeSide.cc.
References DistanceAway(), rNorm, and zNorm.
00470 { 00471 if (p == G4ThreeVector(0.,0.,0.)) { return p; } 00472 00473 G4double dFrom, dOut2; 00474 00475 dFrom = DistanceAway( p, false, dOut2 ); 00476 00477 *bestDistance = std::sqrt( dFrom*dFrom + dOut2 ); 00478 00479 G4double rds = p.perp(); 00480 if (rds!=0.) { return G4ThreeVector(rNorm*p.x()/rds,rNorm*p.y()/rds,zNorm); } 00481 return G4ThreeVector( 0.,0., zNorm ).unit(); 00482 }
G4PolyconeSide & G4PolyconeSide::operator= | ( | const G4PolyconeSide & | source | ) |
Definition at line 196 of file G4PolyconeSide.cc.
References cone, CopyStuff(), corners, and phiIsOpen.
00197 { 00198 if (this == &source) { return *this; } 00199 00200 delete cone; 00201 if (phiIsOpen) { delete [] corners; } 00202 00203 CopyStuff( source ); 00204 00205 return *this; 00206 }
G4bool G4PolyconeSide::PointOnCone | ( | const G4ThreeVector & | hit, | |
G4double | normSign, | |||
const G4ThreeVector & | p, | |||
const G4ThreeVector & | v, | |||
G4ThreeVector & | normal | |||
) | [protected] |
Definition at line 993 of file G4PolyconeSide.cc.
References cone, corners, DBL_MIN, deltaPhi, GetPhi(), G4IntersectingCone::HitOn(), phiIsOpen, rNorm, startPhi, and zNorm.
Referenced by Intersect().
00998 { 00999 G4double rx = hit.perp(); 01000 // 01001 // Check radial/z extent, as appropriate 01002 // 01003 if (!cone->HitOn( rx, hit.z() )) return false; 01004 01005 if (phiIsOpen) 01006 { 01007 G4double phiTolerant = 2.0*kCarTolerance/(rx+kCarTolerance); 01008 // 01009 // Check phi segment. Here we have to be careful 01010 // to use the standard method consistent with 01011 // PolyPhiFace. See PolyPhiFace::InsideEdgesExact 01012 // 01013 G4double phi = GetPhi(hit); 01014 while( phi < startPhi-phiTolerant ) phi += twopi; 01015 01016 if (phi > startPhi+deltaPhi+phiTolerant) return false; 01017 01018 if (phi > startPhi+deltaPhi-phiTolerant) 01019 { 01020 // 01021 // Exact treatment 01022 // 01023 G4ThreeVector qx = p + v; 01024 G4ThreeVector qa = qx - corners[2], 01025 qb = qx - corners[3]; 01026 G4ThreeVector qacb = qa.cross(qb); 01027 01028 if (normSign*qacb.dot(v) < 0) return false; 01029 } 01030 else if (phi < phiTolerant) 01031 { 01032 G4ThreeVector qx = p + v; 01033 G4ThreeVector qa = qx - corners[1], 01034 qb = qx - corners[0]; 01035 G4ThreeVector qacb = qa.cross(qb); 01036 01037 if (normSign*qacb.dot(v) < 0) return false; 01038 } 01039 } 01040 01041 // 01042 // We have a good hit! Calculate normal 01043 // 01044 if (rx < DBL_MIN) 01045 normal = G4ThreeVector( 0, 0, zNorm < 0 ? -1 : 1 ); 01046 else 01047 normal = G4ThreeVector( rNorm*hit.x()/rx, rNorm*hit.y()/rx, zNorm ); 01048 return true; 01049 }
G4double G4PolyconeSide::SurfaceArea | ( | ) | [virtual] |
Implements G4VCSGface.
Definition at line 1087 of file G4PolyconeSide.cc.
References deltaPhi, r, sqr(), and z.
01088 { 01089 if(fSurfaceArea==0) 01090 { 01091 fSurfaceArea = (r[0]+r[1])* std::sqrt(sqr(r[0]-r[1])+sqr(z[0]-z[1])); 01092 fSurfaceArea *= 0.5*(deltaPhi); 01093 } 01094 return fSurfaceArea; 01095 }
G4bool G4PolyconeSide::allBehind [protected] |
Definition at line 138 of file G4PolyconeSide.hh.
Referenced by CopyStuff(), G4PolyconeSide(), and Intersect().
G4IntersectingCone* G4PolyconeSide::cone [protected] |
Definition at line 140 of file G4PolyconeSide.hh.
Referenced by CopyStuff(), Extent(), G4PolyconeSide(), Intersect(), operator=(), PointOnCone(), and ~G4PolyconeSide().
G4ThreeVector* G4PolyconeSide::corners [protected] |
Definition at line 154 of file G4PolyconeSide.hh.
Referenced by CopyStuff(), G4PolyconeSide(), operator=(), PointOnCone(), and ~G4PolyconeSide().
G4double G4PolyconeSide::deltaPhi [protected] |
Definition at line 135 of file G4PolyconeSide.hh.
Referenced by CalculateExtent(), CopyStuff(), DistanceAway(), Extent(), G4PolyconeSide(), GetPointOnFace(), PointOnCone(), and SurfaceArea().
G4double G4PolyconeSide::length [protected] |
Definition at line 144 of file G4PolyconeSide.hh.
Referenced by CopyStuff(), DistanceAway(), and G4PolyconeSide().
G4int G4PolyconeSide::ncorners [protected] |
G4double G4PolyconeSide::nextRS [protected] |
Definition at line 147 of file G4PolyconeSide.hh.
Referenced by CalculateExtent(), CopyStuff(), and G4PolyconeSide().
G4double G4PolyconeSide::nextZS [protected] |
Definition at line 147 of file G4PolyconeSide.hh.
Referenced by CalculateExtent(), CopyStuff(), and G4PolyconeSide().
G4bool G4PolyconeSide::phiIsOpen [protected] |
Definition at line 137 of file G4PolyconeSide.hh.
Referenced by CalculateExtent(), CopyStuff(), DistanceAway(), Extent(), G4PolyconeSide(), operator=(), PointOnCone(), and ~G4PolyconeSide().
G4double G4PolyconeSide::prevRS [protected] |
Definition at line 145 of file G4PolyconeSide.hh.
Referenced by CalculateExtent(), CopyStuff(), and G4PolyconeSide().
G4double G4PolyconeSide::prevZS [protected] |
Definition at line 145 of file G4PolyconeSide.hh.
Referenced by CalculateExtent(), CopyStuff(), and G4PolyconeSide().
G4double G4PolyconeSide::r[2] [protected] |
Definition at line 134 of file G4PolyconeSide.hh.
Referenced by CalculateExtent(), CopyStuff(), DistanceAway(), Extent(), G4PolyconeSide(), GetPointOnFace(), and SurfaceArea().
G4double G4PolyconeSide::rNorm [protected] |
Definition at line 142 of file G4PolyconeSide.hh.
Referenced by CalculateExtent(), CopyStuff(), DistanceAway(), G4PolyconeSide(), Intersect(), Normal(), and PointOnCone().
G4double G4PolyconeSide::rNormEdge[2] [protected] |
Definition at line 150 of file G4PolyconeSide.hh.
Referenced by CopyStuff(), DistanceAway(), and G4PolyconeSide().
G4double G4PolyconeSide::rS [protected] |
Definition at line 143 of file G4PolyconeSide.hh.
Referenced by CalculateExtent(), CopyStuff(), DistanceAway(), and G4PolyconeSide().
G4double G4PolyconeSide::startPhi [protected] |
Definition at line 135 of file G4PolyconeSide.hh.
Referenced by CalculateExtent(), CopyStuff(), DistanceAway(), Extent(), G4PolyconeSide(), GetPointOnFace(), and PointOnCone().
G4double G4PolyconeSide::z[2] [protected] |
Definition at line 134 of file G4PolyconeSide.hh.
Referenced by CalculateExtent(), CopyStuff(), DistanceAway(), Extent(), G4PolyconeSide(), GetPointOnFace(), and SurfaceArea().
G4double G4PolyconeSide::zNorm [protected] |
Definition at line 142 of file G4PolyconeSide.hh.
Referenced by CopyStuff(), DistanceAway(), G4PolyconeSide(), Intersect(), Normal(), and PointOnCone().
G4double G4PolyconeSide::zNormEdge[2] [protected] |
Definition at line 150 of file G4PolyconeSide.hh.
Referenced by CopyStuff(), DistanceAway(), and G4PolyconeSide().
G4double G4PolyconeSide::zS [protected] |
Definition at line 143 of file G4PolyconeSide.hh.
Referenced by CalculateExtent(), CopyStuff(), DistanceAway(), and G4PolyconeSide().