62 axis0min, axis1min, axis0max, axis1max),
66 G4Exception(
"G4TwistTubsSide::G4TwistTubsSide()",
"GeomSolids0002",
101 SetCorners( EndInnerRadius, EndOuterRadius, EndPhi, EndZ) ;
226 for (i=0; i<2; i++) {
251 isvalid[0], 0, validate, &gp, &gv);
271 distance[0] = - c / b;
272 xx[0] = p + distance[0]*v;
278 if (distance[0] >= 0) isvalid[0] =
true;
283 if (distance[0] >= 0) isvalid[0] =
true;
289 if (distance[0] >= 0) isvalid[0] =
true;
293 areacode[0], isvalid[0],
294 0, validate, &gp, &gv);
300 isvalid[0], 1, validate, &gp, &gv);
312 isvalid[0], 0, validate, &gp, &gv);
324 G4bool tmpisvalid[2] = {
false,
false};
327 for (i=0; i<2; i++) {
333 if ( b * D < 0 && std::fabs(bminusD / D) < protection ) {
335 tmpdist[i] = - c/b * ( 1 - acovbb * (1 + 2*acovbb));
337 tmpdist[i] = factor * bminusD;
341 tmpxx[i] = p + tmpdist[i]*v;
346 if (tmpdist[i] >= 0) tmpisvalid[i] =
true;
352 if (tmpdist[i] >= 0) tmpisvalid[i] =
true;
357 if (tmpxx[i].
x() > 0) {
359 if (tmpdist[i] >= 0) tmpisvalid[i] =
true;
368 if (tmpdist[0] <= tmpdist[1]) {
369 distance[0] = tmpdist[0];
370 distance[1] = tmpdist[1];
375 areacode[0] = tmpareacode[0];
376 areacode[1] = tmpareacode[1];
377 isvalid[0] = tmpisvalid[0];
378 isvalid[1] = tmpisvalid[1];
380 distance[0] = tmpdist[1];
381 distance[1] = tmpdist[0];
386 areacode[0] = tmpareacode[1];
387 areacode[1] = tmpareacode[0];
388 isvalid[0] = tmpisvalid[1];
389 isvalid[1] = tmpisvalid[0];
393 isvalid[0], 2, validate, &gp, &gv);
395 isvalid[1], 2, validate, &gp, &gv);
399 for (
G4int k=0; k<2; k++) {
400 if (!isvalid[k])
continue;
403 * xx[k].
z() , xx[k].
z());
404 G4double deltaY = (xx[k] - xxonsurface).mag();
411 for (l=0; l<maxcount; l++) {
414 surfacenormal, xx[k]);
415 deltaY = (xx[k] - xxonsurface).mag();
416 if (deltaY > lastdeltaY) {
425 xxonsurface.
set(xx[k].
x(),
426 fKappa * std::fabs(xx[k].
x()) * xx[k].
z(),
431 message <<
"Exceeded maxloop count!" <<
G4endl
432 <<
" maxloop count " << maxcount;
433 G4Exception(
"G4TwistTubsFlatSide::DistanceToSurface(p,v)",
445 isvalid[0], 0, validate, &gp, &gv);
470 for (i=0; i<2; i++) {
491 for (i=0; i<2; i++) {
495 if ((gp - lastgxx[0]).mag() < halftol
496 || (gp - lastgxx[1]).mag() < halftol) {
548 }
else if (p.
z() < C.
z()) {
554 }
else if (p.
z() < A.
z()) {
563 for (i=0; i<2; i++) {
566 B = x0[i] + ((A.
z() - x0[i].
z()) / d[i].
z()) * d[i];
570 D = x0[i] + ((C.
z() - x0[i].
z()) / d[i].
z()) * d[i];
579 G4double test = (A.
z() - C.
z()) * parity * pside;
604 }
else if (test < 0) {
636 if (std::fabs(distToACB) <= halftol || std::fabs(distToCAD) <= halftol) {
637 xx = (std::fabs(distToACB) < std::fabs(distToCAD) ? xxacb : xxcad);
647 if (distToACB * distToCAD > 0 && distToACB < 0) {
653 if (distToACB * distToCAD > 0) {
656 if (distToACB <= distToCAD) {
657 distance[0] = distToACB;
660 distance[0] = distToCAD;
667 distance[0] = distToACB;
670 distance[0] = distToCAD;
709 if (distToanm * distTocmn > 0 && distToanm < 0) {
712 "Point p is behind the surfaces.");
716 if (std::fabs(distToanm) <= halftol) {
720 }
else if (std::fabs(distTocmn) <= halftol) {
726 if (distToanm <= distTocmn) {
787 if (xx.
z() <=
fAxisMin[zaxis] - ctol) isoutside =
true;
789 }
else if (xx.
z() >
fAxisMax[zaxis] - ctol) {
794 if (xx.
z() >=
fAxisMax[zaxis] + ctol) isoutside =
true;
802 areacode = tmpareacode;
838 "Feature NOT implemented !");
862 x = endInnerRad[zmin]*std::cos(endPhi[zmin]);
863 y = endInnerRad[zmin]*std::sin(endPhi[zmin]);
868 x = endOuterRad[zmin]*std::cos(endPhi[zmin]);
869 y = endOuterRad[zmin]*std::sin(endPhi[zmin]);
874 x = endOuterRad[zmax]*std::cos(endPhi[zmax]);
875 y = endOuterRad[zmax]*std::sin(endPhi[zmax]);
880 x = endInnerRad[zmax]*std::cos(endPhi[zmax]);
881 y = endInnerRad[zmax]*std::sin(endPhi[zmax]);
887 message <<
"Feature NOT implemented !" <<
G4endl
889 <<
" fAxis[1] = " <<
fAxis[1];
902 "Method NOT implemented !");
918 direction = direction.
unit();
924 direction = direction.
unit();
930 direction = direction.
unit();
936 direction = direction.
unit();
942 message <<
"Feature NOT implemented !" <<
G4endl
944 <<
" fAxis[1] = " <<
fAxis[1];
969 for ( i = 0 ; i<
n ; i++ )
974 for ( j = 0 ; j<k ; j++ ) {
976 nnode =
GetNode(i,j,k,n,iside) ;
982 x = xmin + j*(xmax-xmin)/(k-1) ;
984 x = xmax - j*(xmax-xmin)/(k-1) ;
989 xyz[nnode][0] = p.
x() ;
990 xyz[nnode][1] = p.
y() ;
991 xyz[nnode][2] = p.
z() ;
993 if ( i<n-1 && j<k-1 ) {
995 nface =
GetFace(i,j,k,n,iside) ;
G4int GetAreacode(G4int i) const
static const G4int sC0Min1Max
void set(double x, double y, double z)
CurrentStatus fCurStatWithV
static const G4int sInside
virtual G4int GetAreaCode(const G4ThreeVector &xx, G4bool withTol=true)
G4bool IsInside(G4int areacode, G4bool testbitmode=false) const
virtual G4double GetBoundaryMin(G4double phi)
std::vector< ExP01TrackerHit * > a
static const G4double kInfinity
virtual G4ThreeVector GetNormal(const G4ThreeVector &xx, G4bool isGlobal=false)
static const G4int sAxisX
virtual void SetBoundary(const G4int &axiscode, const G4ThreeVector &direction, const G4ThreeVector &x0, const G4int &boundarytype)
void message(RunManager *runmanager)
static const G4int sCorner
virtual void SetBoundaries()
void SetCurrentStatus(G4int i, G4ThreeVector &xx, G4double &dist, G4int &areacode, G4bool &isvalid, G4int nxx, EValidate validate, const G4ThreeVector *p, const G4ThreeVector *v=0)
static double normal(HepRandomEngine *eptr)
virtual G4int DistanceToSurface(const G4ThreeVector &gp, const G4ThreeVector &gv, G4ThreeVector gxx[], G4double distance[], G4int areacode[], G4bool isvalid[], EValidate validate=kValidateWithTol)
virtual void SetCorners()
G4bool IsOutside(G4int areacode) const
G4int GetEdgeVisibility(G4int i, G4int j, G4int m, G4int n, G4int number, G4int orientation)
static const G4int sAxisZ
static const G4int sC0Min1Min
virtual G4double GetBoundaryMax(G4double phi)
virtual G4ThreeVector SurfacePoint(G4double, G4double, G4bool isGlobal=false)
void ResetfDone(EValidate validate, const G4ThreeVector *p, const G4ThreeVector *v=0)
static const G4int sAxis1
G4int GetFace(G4int i, G4int j, G4int m, G4int n, G4int iside)
static const G4int sC0Max1Max
G4double DistanceToLine(const G4ThreeVector &p, const G4ThreeVector &x0, const G4ThreeVector &d, G4ThreeVector &xx)
virtual void GetFacets(G4int m, G4int n, G4double xyz[][3], G4int faces[][4], G4int iside)
G4double DistanceToPlaneWithV(const G4ThreeVector &p, const G4ThreeVector &v, const G4ThreeVector &x0, const G4ThreeVector &n0, G4ThreeVector &xx)
double A(double temperature)
void SetCorner(G4int areacode, G4double x, G4double y, G4double z)
G4ThreeVector GetCorner(G4int areacode) const
virtual ~G4TwistTubsSide()
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
G4ThreeVector ComputeGlobalPoint(const G4ThreeVector &lp) const
virtual G4double DistanceToBoundary(G4int areacode, G4ThreeVector &xx, const G4ThreeVector &p)
G4double GetDistance(G4int i) const
static const G4int sAxis0
G4ThreeVector ComputeLocalPoint(const G4ThreeVector &gp) const
HepRotation & rotateZ(double delta)
virtual void GetBoundaryParameters(const G4int &areacode, G4ThreeVector &d, G4ThreeVector &x0, G4int &boundarytype) const
static const G4int sOutside
G4SurfCurNormal fCurrentNormal
static const G4int sAxisMin
static const G4int sBoundary
G4int GetNode(G4int i, G4int j, G4int m, G4int n, G4int iside)
G4TwistTubsSide(const G4String &name, const G4RotationMatrix &rot, const G4ThreeVector &tlate, G4int handedness, const G4double kappa, const EAxis axis0=kXAxis, const EAxis axis1=kZAxis, G4double axis0min=-kInfinity, G4double axis1min=-kInfinity, G4double axis0max=kInfinity, G4double axis1max=kInfinity)
static const G4int sC0Max1Min
static const G4int sAxisMax
double B(double temperature)
G4double DistanceToPlane(const G4ThreeVector &p, const G4ThreeVector &x0, const G4ThreeVector &n0, G4ThreeVector &xx)
G4ThreeVector GetXX(G4int i) const
virtual G4int AmIOnLeftSide(const G4ThreeVector &me, const G4ThreeVector &vec, G4bool withTol=true)
G4ThreeVector ComputeGlobalDirection(const G4ThreeVector &lp) const
virtual G4ThreeVector GetBoundaryAtPZ(G4int areacode, const G4ThreeVector &p) const
virtual G4double DistanceToPlane(const G4ThreeVector &p, const G4ThreeVector &A, const G4ThreeVector &B, const G4ThreeVector &C, const G4ThreeVector &D, const G4int parity, G4ThreeVector &xx, G4ThreeVector &n)
G4bool IsValid(G4int i) const
G4ThreeVector ComputeLocalDirection(const G4ThreeVector &gp) const