LArWheelSolid Class Reference

#include <LArWheelSolid.h>

Inheritance diagram for LArWheelSolid:

Collaboration diagram for LArWheelSolid:

Public Member Functions
	LArWheelSolid (const G4String &name, LArWheelSolid_t type, G4int zside=1, LArWheelCalculator *calc=0, const EMECData *emecData=0)
virtual	~LArWheelSolid ()
EInside	Inside (const G4ThreeVector &) const
G4double	DistanceToIn (const G4ThreeVector &, const G4ThreeVector &) const
G4double	DistanceToIn (const G4ThreeVector &) const
G4double	DistanceToOut (const G4ThreeVector &, const G4ThreeVector &, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const
G4double	DistanceToOut (const G4ThreeVector &) const
G4ThreeVector	SurfaceNormal (const G4ThreeVector &) const
G4bool	CalculateExtent (const EAxis, const G4VoxelLimits &, const G4AffineTransform &, G4double &, G4double &) const
G4GeometryType	GetEntityType () const
void	DescribeYourselfTo (G4VGraphicsScene &) const
G4VisExtent	GetExtent () const
G4Polyhedron *	CreatePolyhedron () const
virtual std::ostream &	StreamInfo (std::ostream &os) const
const G4VSolid *	GetBoundingShape (void) const
const LArWheelCalculator *	GetCalculator (void) const
LArWheelSolid_t	GetType (void) const
G4ThreeVector	GetPointOnSurface (void) const
G4double	GetCubicVolume (void)
G4double	GetSurfaceArea (void)
bool	msgLvl (const MSG::Level lvl) const
	Test the output level.
MsgStream &	msg () const
	The standard message stream.
MsgStream &	msg (const MSG::Level lvl) const
	The standard message stream.
void	setLevel (MSG::Level lvl)
	Change the current logging level.

Protected Attributes
TF1 *	m_f_area
TF1 *	m_f_vol
TF1 *	m_f_area_on_pc
TF1 *	m_f_length
TF1 *	m_f_side_area
double	m_test_index

Private Types
enum	FanBoundExit_t {   NoCross , ExitAtInner , ExitAtOuter , ExitAtFront ,   ExitAtBack , ExitAtSide }

Private Member Functions
void	inner_solid_init (const G4String &)
void	outer_solid_init (const G4String &)
void	set_phi_size (void)
virtual G4double	distance_to_in (G4ThreeVector &, const G4ThreeVector &, int) const
G4double	in_iteration_process (const G4ThreeVector &, G4double, G4ThreeVector &, int) const
G4double	search_for_nearest_point (const G4ThreeVector &, const G4double, const G4ThreeVector &, int) const
G4bool	search_for_most_remoted_point (const G4ThreeVector &, const G4ThreeVector &, G4ThreeVector &, const int) const
G4double	out_iteration_process (const G4ThreeVector &, G4ThreeVector &, const int) const
FanBoundExit_t	find_exit_point (const G4ThreeVector &p, const G4ThreeVector &v, G4ThreeVector &q) const
G4bool	fs_cross_lower (const G4ThreeVector &p, const G4ThreeVector &v, G4ThreeVector &q) const
G4bool	fs_cross_upper (const G4ThreeVector &p, const G4ThreeVector &v, G4ThreeVector &q) const
G4bool	check_D (G4double &b, G4double A, G4double B, G4double C, G4bool) const
G4int	select_section (const G4double &Z) const
EInside	Inside_accordion (const G4ThreeVector &) const
void	get_point_on_accordion_surface (G4ThreeVector &) const
void	get_point_on_polycone_surface (G4ThreeVector &) const
void	get_point_on_flat_surface (G4ThreeVector &) const
void	set_failover_point (G4ThreeVector &p, const char *m=0) const
G4double	get_area_on_polycone (void) const
G4double	get_area_on_face (void) const
G4double	get_area_on_side (void) const
G4double	get_area_at_r (G4double r) const
G4double	get_length_at_r (G4double r) const
void	test (void)
void	clean_tests (void)
void	init_tests (void)
void	initMessaging () const
	Initialize our message level and MessageSvc.

Private Attributes
G4bool	m_IsOuter
const LArWheelSolid_t	m_Type
LArWheelCalculator *	m_Calculator
G4double	m_FanHalfThickness
G4double	m_FHTplusT
G4double	m_FHTminusT
G4double	m_FanPhiAmplitude
G4double	m_MinPhi
G4double	m_MaxPhi
const G4double	m_PhiPosition
G4VSolid *	m_BoundingShape
std::vector< G4double >	m_Zsect
G4int	m_Zsect_start_search
LArFanSections *	m_fs
G4double	m_Zmid
G4double	m_Ymin
G4double	m_Zmin
G4double	m_Zmax
G4double	m_Rmin
G4double	m_Rmax
G4double	m_Ymid
std::string	m_nm
	Message source name.
boost::thread_specific_ptr< MsgStream >	m_msg_tls
	MsgStream instance (a std::cout like with print-out levels)
std::atomic< IMessageSvc * >	m_imsg { nullptr }
	MessageSvc pointer.
std::atomic< MSG::Level >	m_lvl { MSG::NIL }
	Current logging level.
std::atomic_flag m_initialized	ATLAS_THREAD_SAFE = ATOMIC_FLAG_INIT
	Messaging initialized (initMessaging)

Static Private Attributes
static const G4double	s_Tolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance() / 2
static const G4double	s_AngularTolerance = G4GeometryTolerance::GetInstance()->GetAngularTolerance() / 2
static const G4double	s_IterationPrecision = 0.001*CLHEP::mm
static const G4double	s_IterationPrecision2 = s_IterationPrecision * s_IterationPrecision
static const unsigned int	s_IterationsLimit = 50

Friends
double	LArWheelSolid_fcn_area (double *x, double *p)
double	LArWheelSolid_fcn_vol (double *x, double *p)
double	LArWheelSolid_fcn_area_on_pc (double *x, double *p)
double	LArWheelSolid_get_dl (double *x, double *par, G4int side)
double	LArWheelSolid_fcn_side_area (double *x, double *p)

Detailed Description

Definition at line 90 of file LArWheelSolid.h.

Member Enumeration Documentation

FanBoundExit_t

enum LArWheelSolid::FanBoundExit_t

private

Enumerator
NoCross
ExitAtInner
ExitAtOuter
ExitAtFront
ExitAtBack
ExitAtSide

Definition at line 189 of file LArWheelSolid.h.

               {
    NoCross, ExitAtInner, ExitAtOuter,
    ExitAtFront, ExitAtBack, ExitAtSide
  } FanBoundExit_t;

Constructor & Destructor Documentation

LArWheelSolid()

LArWheelSolid::LArWheelSolid	(	const G4String &	name,
		LArWheelSolid_t	type,
		G4int	zside = 1,
		LArWheelCalculator *	calc = 0,
		const EMECData *	emecData = 0 )

Definition at line 32 of file LArWheelSolidInit.cxx.

  : G4VSolid(name)
#ifndef PORTABLE_LAR_SHAPE
  , AthMessaging("LArWheelSolid")
#endif
  , m_Type(type), m_Calculator(calc), m_PhiPosition(CLHEP::halfpi), m_fs(0)
{
#ifndef PORTABLE_LAR_SHAPE
#ifdef LARWHEELSOLID_USE_FANBOUND
  ATH_MSG_INFO ( "compiled with G4 FanBound" );
#else
  ATH_MSG_INFO ( "compiled with private find_exit_point" );
#endif
#endif
 
  LArG4::LArWheelCalculator_t calc_type = LArG4::LArWheelCalculator_t(0);
  switch(m_Type){
  case InnerAbsorberWheel:
    calc_type = LArG4::InnerAbsorberWheel;
    break;
  case OuterAbsorberWheel:
    calc_type = LArG4::OuterAbsorberWheel;
    break;
  case InnerElectrodWheel:
    calc_type = LArG4::InnerElectrodWheel;
    break;
  case OuterElectrodWheel:
    calc_type = LArG4::OuterElectrodWheel;
    break;
  case InnerAbsorberModule:
    calc_type = LArG4::InnerAbsorberModule;
    break;
  case OuterAbsorberModule:
    calc_type = LArG4::OuterAbsorberModule;
    break;
  case InnerElectrodModule:
    calc_type = LArG4::InnerElectrodModule;
    break;
  case OuterElectrodModule:
    calc_type = LArG4::OuterElectrodModule;
    break;
  case InnerGlueWheel:
    calc_type = LArG4::InnerGlueWheel;
    break;
  case OuterGlueWheel:
    calc_type = LArG4::OuterGlueWheel;
    break;
  case InnerLeadWheel:
    calc_type = LArG4::InnerLeadWheel;
    break;
  case OuterLeadWheel:
    calc_type = LArG4::OuterLeadWheel;
    break;
  case InnerAbsorberCone:
    calc_type = LArG4::InnerAbsorberWheel;
    break;
  case InnerElectrodCone:
    calc_type = LArG4::InnerElectrodWheel;
    break;
  case InnerGlueCone:
    calc_type = LArG4::InnerGlueWheel;
    break;
  case InnerLeadCone:
    calc_type = LArG4::InnerLeadWheel;
    break;
  case OuterAbsorberFrontCone:
    calc_type = LArG4::OuterAbsorberWheel;
    break;
  case OuterElectrodFrontCone:
    calc_type = LArG4::OuterElectrodWheel;
    break;
  case OuterGlueFrontCone:
    calc_type = LArG4::OuterGlueWheel;
    break;
  case OuterLeadFrontCone:
    calc_type = LArG4::OuterLeadWheel;
    break;
  case OuterAbsorberBackCone:
    calc_type = LArG4::OuterAbsorberWheel;
    break;
  case OuterElectrodBackCone:
    calc_type = LArG4::OuterElectrodWheel;
    break;
  case OuterGlueBackCone:
    calc_type = LArG4::OuterGlueWheel;
    break;
  case OuterLeadBackCone:
    calc_type = LArG4::OuterLeadWheel;
    break;
  default:
    G4Exception("LArWheelSolid", "UnknownSolidType", FatalException,
                "Constructor: unknown LArWheelSolid_t");
  }
 
  if(m_Calculator == 0) m_Calculator = new LArWheelCalculator(*emecData,calc_type, zside);
 
  const G4String bs_name = name + "-Bounding";
#ifdef DEBUG_LARWHEELSOLID
  const char *venv = getenv("LARWHEELSOLID_VERBOSE");
  if(venv) Verbose = atoi(venv);
  std::cout << "The LArWheelSolid build " << __DATE__ << " " << __TIME__
            << std::endl;
  std::cout << "LArWheelSolid verbosity level is " << Verbose << std::endl;
#endif
 
  // Initialize code that depends on wheel type:
  m_FanHalfThickness = GetCalculator()->GetFanHalfThickness();
  m_FHTplusT = m_FanHalfThickness + s_Tolerance;
  m_FHTminusT = m_FanHalfThickness - s_Tolerance;
  switch(m_Type){
  case InnerAbsorberWheel:
  case InnerElectrodWheel:
  case InnerAbsorberModule:
  case InnerElectrodModule:
  case InnerGlueWheel:
  case InnerLeadWheel:
  case InnerAbsorberCone:
  case InnerElectrodCone:
  case InnerGlueCone:
  case InnerLeadCone:
    inner_solid_init(bs_name);
    break;
  case OuterAbsorberWheel:
  case OuterElectrodWheel:
  case OuterAbsorberModule:
  case OuterElectrodModule:
  case OuterGlueWheel:
  case OuterLeadWheel:
  case OuterAbsorberFrontCone:
  case OuterElectrodFrontCone:
  case OuterGlueFrontCone:
  case OuterLeadFrontCone:
  case OuterAbsorberBackCone:
  case OuterElectrodBackCone:
  case OuterGlueBackCone:
  case OuterLeadBackCone:
    outer_solid_init(bs_name);
    break;
  default:
    G4Exception("LArWheelSolid", "UnknownSolidType", FatalException,
                "Constructor: unknown LArWheelSolid_t");
  }
 
  m_Zsect_start_search = (m_Zsect.size() - 1) - 1;
#ifndef PORTABLE_LAR_SHAPE
  init_tests();
  test(); // activated by env. variable
  clean_tests();
#endif
  
#ifdef DEBUG_LARWHEELSOLID
  m_fs->print();
  std::cout << "Limits: (" << m_Zsect.size() << ")" << std::endl;
  for(size_t i = 0; i < m_Zsect.size(); ++ i){
    std::cout << i << " " << m_Zsect[i] << std::endl;
  }
#endif
#ifndef PORTABLE_LAR_SHAPE
  ATH_MSG_DEBUG ( "solid of type "
                  << LArWheelCalculator::LArWheelCalculatorTypeString(calc_type)
                  << " initialized" );
#endif
}

~LArWheelSolid()

LArWheelSolid::~LArWheelSolid ( )

virtual

Definition at line 200 of file LArWheelSolidInit.cxx.

{
  if(m_fs) delete m_fs;
}

Member Function Documentation

CalculateExtent()

G4bool LArWheelSolid::CalculateExtent	(	const EAxis	a,
		const G4VoxelLimits &	vl,
		const G4AffineTransform &	t,
		G4double &	p,
		G4double &	q ) const

Definition at line 58 of file LArWheelSolid.cxx.

{
    return m_BoundingShape->CalculateExtent(a, vl, t, p, q);
}

check_D()

G4bool LArWheelSolid::check_D ( G4double & b, G4double A, G4double B, G4double C, G4bool out ) const

private

Definition at line 41 of file LArFanSection.cxx.

{
    // G4bool out is to be set true if the point is surface-outside
    // then have to discard first intersection
 
    const G4double D = B*B - A*C;
    LWSDBG(8, std::cout << "check D=" << D << " out=" << out << std::endl);
    if(D < 0.) return false;
    G4double t2 = 0.;
    if(A == 0) {
        if(B == 0) {
            LWSDBG(8, std::cout << "Case A=B=0" << std::endl);
            return false;
        }
        t1= -C / (2 * B);
        LWSDBG(8, std::cout << "t1=" << t1 << " - only one solution" << std::endl);
        t2 = t1;
    } else {
        const G4double D1 = sqrt(D);
        const G4double inv_A = 1.0 / A;
        t1 = (-B + D1) * inv_A;
        t2 = (-B - D1) * inv_A;
        LWSDBG(8, std::cout << "t1=" << t1 << " t2=" << t2 << std::endl);     
    }
    
    if(t1 > 0.){
        if(t2 > 0.){
            if(out){
                if(t2 > t1) t1 = t2;
            } else {
                if(t2 < t1) t1 = t2;
            }
        } else if(t2 < 0.){
            if(out) return false;
        } else { // answer is t1
        }
    } else if(t1 < 0.){
        if(t2 > 0.){
            if(out) return false;
            t1 = t2;
        } else if(t2 < 0.){
            return false;
        } else {
            return false;
        }
    } else {
        if(t2 > 0.){
            t1 = t2;
        } else if(t2 < 0.){
            return false;
        } else {
            return false;
        }
    }
    return true;
}

clean_tests()

void LArWheelSolid::clean_tests ( void )

private

Definition at line 486 of file LArWheelSolidTests.cxx.

                                    {
    if(m_f_area) {
        delete m_f_area;
        m_f_area = 0;
    }
    if(m_f_vol) {
        delete m_f_vol;
        m_f_vol = 0;
    }
 
    if(m_f_area_on_pc) {
        delete m_f_area_on_pc;
        m_f_area_on_pc = 0;
    }
 
    if(m_f_length) {
        delete m_f_length;
        m_f_length = 0;
    }
    if(m_f_side_area) {
        delete m_f_side_area;
        m_f_side_area = 0;
    }
}

CreatePolyhedron()

G4Polyhedron * LArWheelSolid::CreatePolyhedron

(

)

const

Definition at line 156 of file LArWheelSolid.cxx.

{
  return m_BoundingShape->CreatePolyhedron();
}

DescribeYourselfTo()

void LArWheelSolid::DescribeYourselfTo

(

G4VGraphicsScene &

scene

)

const

Definition at line 146 of file LArWheelSolid.cxx.

{
  scene.AddSolid(*this);
}

distance_to_in()

G4double LArWheelSolid::distance_to_in ( G4ThreeVector & p, const G4ThreeVector & v, int p_fan ) const

privatevirtual

Definition at line 216 of file LArWheelSolidDisToIn.cxx.

{
    LWSDBG(4, std::cout << "dti:           " << MSG_VECTOR(p) << " "
                        << MSG_VECTOR(v) << std::endl);
 
    G4double distance = 0.;
#ifdef LARWHEELSOLID_USE_FANBOUND
    if(FanBound->Inside(p) == kOutside) {
        const G4double d = FanBound->DistanceToIn(p, v);
        p += v * d;
        distance += d;
    }
#else
    if(p.x() > m_fs->xmax) {
        if(v.x() >= 0.) return kInfinity;
        const G4double b = (m_fs->xmax - p.x()) / v.x();
        const G4double y2 = p.y() + v.y() * b;
        const G4double z2 = p.z() + v.z() * b;
        p.set(m_fs->xmax, y2, z2);
        distance += b;
    } else if(p.x() < m_fs->xmin) {
        if(v.x() <= 0.)    return kInfinity;
        const G4double b = (m_fs->xmin - p.x()) / v.x();
        const G4double y2 = p.y() + v.y() * b;
        const G4double z2 = p.z() + v.z() * b;
        p.set(m_fs->xmin, y2, z2);
        distance += b;
    }
#endif
 
// here p is on surface of or inside the "FanBound",
// distance corrected, misses are accounted for
    LWSDBG(5, std::cout << "dti corrected: " << MSG_VECTOR(p) << std::endl);
 
    G4double dist_p = GetCalculator()->DistanceToTheNeutralFibre(p, p_fan);
    if(fabs(dist_p) < m_FHTminusT) {
        LWSDBG(5, std::cout << "hit fan dist_p=" << dist_p << ", m_FHTminusT=" << m_FHTminusT << std::endl);
        return distance;
    }
 
#ifdef CHECK_DIRTONORM_ANGLE_ON_SURFACE
    if(fabs(dist_p) > m_FHTplusT) {
        LWSDBG(5, std::cout << "outside fan dist_p=" << dist_p << ", m_FHTplusT=" << m_FHTplusT << std::endl);
    } else {
        LWSDBG(5, std::cout << "on fan surface dist_p=" << dist_p << ", m_FHTplusT=" << m_FHTplusT << ", m_FHTminusT=" << m_FHTminusT << std::endl);
 
        const G4ThreeVector d = GetCalculator()->NearestPointOnNeutralFibre(p, p_fan);
        // check dot product between normal and v
        if ( (p-d).cosTheta(v) < -AngularTolerance ) {
            // direction "v" definitely pointing inside
            // return 0.0, it should be in "distance"
            return distance;
        }
    }
#endif
 
    G4ThreeVector q;
#ifdef LARWHEELSOLID_USE_FANBOUND
    q = p + v * FanBound->DistanceToOut(p, v);
#else
    find_exit_point(p, v, q);
#endif
 
    G4int start = select_section(p.z());
    G4int stop = select_section(q.z());
    G4int step = -1;
    if(stop > start) { step = 1;    start ++; stop ++; }
    LWSDBG(5, std::cout << "dti sections " << start << " " << stop
                        << " " << step << std::endl);
    G4ThreeVector p1;
    for(G4int i = start; i != stop; i += step){
// v.z() can't be 0, otherwise start == stop, so the exit point could be only
// at z border of the fan section
        const G4double d1 = (m_Zsect[i] - p.z()) / v.z();
        const G4double x1 = p.x() + v.x() * d1, y1 = p.y() + v.y() * d1;
        p1.set(x1, y1, m_Zsect[i]);
        G4double dist_p1 = GetCalculator()->DistanceToTheNeutralFibre(p1, p_fan);
        LWSDBG(5, std::cout << i << ">" << p << " " << dist_p << " "
                            << p1 << " " << dist_p1 << std::endl);
        G4double dd = kInfinity;
        if(dist_p * dist_p1 < 0.){// it certanly cross current half-wave
            dd = in_iteration_process(p, dist_p, p1, p_fan);
        }
        G4double d2 = search_for_nearest_point(p, dist_p, p1, p_fan);
        LWSDBG(6, std::cout << i << "> dd=" << dd << ", d2=" << d2 << ", distance=" << distance << std::endl);
        if(d2 < kInfinity){
            return distance + d2; // this half-wave is intersected
        } else if(dd < kInfinity){
            return distance + dd;
        }
        distance += d1;
        p.set(x1, y1, m_Zsect[i]);
        dist_p = dist_p1;
    }
 
    G4double dist_q = GetCalculator()->DistanceToTheNeutralFibre(q, p_fan);
    LWSDBG(5, std::cout << "dti exit point: " << MSG_VECTOR(q) << " "
                        << dist_q << std::endl);
    G4double dd = kInfinity;
    if(dist_p * dist_q < 0.){// it certanly cross current half-wave
        dd = in_iteration_process(p, dist_p, q, p_fan);
    }
    G4double d2 = search_for_nearest_point(p, dist_p, q, p_fan);
    if(d2 < kInfinity){
        return distance + d2; // this half-wave is intersected
    } else if(dd < kInfinity){
        return distance + dd;
    }
    return kInfinity;
}

DistanceToIn() [1/2]

G4double LArWheelSolid::DistanceToIn

(

const G4ThreeVector &

inputP

)

const

Definition at line 17 of file LArWheelSolidDisToIn.cxx.

{
    LWSDBG(1, std::cout << TypeStr() << " DisToIn" << MSG_VECTOR(inputP) << std::endl);
    if(m_BoundingShape->Inside(inputP) == kOutside) {
    // here is an approximation - for the point outside m_BoundingShape
    // the solid looks like a m_BoundingShape
    // it's okay since the result could be underestimated
        LWSDBG(2, std::cout << "Outside BS" << std::endl);
        return m_BoundingShape->DistanceToIn(inputP);
    }
    G4ThreeVector p(inputP);
 
    //
    // DistanceToTheNearestFan:
    // rotates point p to the localFan coordinates and returns fan number to out_fan_number parameter
    // returns distance to fan as result
    //
 
    int p_fan = 0;
    const G4double d = fabs(GetCalculator()->DistanceToTheNearestFan(p, p_fan));
    if(d > m_FHTplusT){
        const G4double result = d - m_FanHalfThickness;
        LWSDBG(2, std::cout << "dti result = " << result << std::endl);
        return result;
    }
    LWSDBG(2, std::cout << "already inside, return 0" << MSG_VECTOR(p) << std::endl);
    return 0.;
}

DistanceToIn() [2/2]

G4double LArWheelSolid::DistanceToIn	(	const G4ThreeVector &	inputP,
		const G4ThreeVector &	inputV ) const

Definition at line 46 of file LArWheelSolidDisToIn.cxx.

{
    LWSDBG(1, std::cout << TypeStr() << " DisToIn" << MSG_VECTOR(inputP)
                        << MSG_VECTOR(inputV) << std::endl);
 
    G4double distance = 0.;
    const EInside inside_BS = m_BoundingShape->Inside(inputP);
    G4ThreeVector p(inputP);
    if(inside_BS == kOutside) {
        distance = m_BoundingShape->DistanceToIn(inputP, inputV);
        if(distance == kInfinity) {
            LWSDBG(2, std::cout << "Infinity distance to m_BoundingShape" << MSG_VECTOR(inputP) << MSG_VECTOR(inputV) << std::endl);
            return kInfinity;
        }
        p += inputV * distance;
        assert(m_BoundingShape->Inside(p) != kOutside);
        LWSDBG(2, std::cout << "shift" << MSG_VECTOR(inputP) << std::endl);
    }
 
    const G4double phi0 = p.phi();
    int p_fan = 0;
    const G4double d = GetCalculator()->DistanceToTheNearestFan(p, p_fan);
    if(fabs(d) < m_FHTminusT){
        LWSDBG(2, std::cout << "already inside fan" << MSG_VECTOR(p) << std::endl);
        // if initial point is on BS surface and inside fan volume it is a real surface
        if(inside_BS == kSurface) {
            LWSDBG(2, std::cout << "On BS surface" << std::endl);
            return m_BoundingShape->DistanceToIn(inputP, inputV);
        }
        return distance;
    }
    G4ThreeVector v(inputV);
    v.rotateZ(p.phi() - phi0);
 
    const G4double d0 = distance_to_in(p, v, p_fan);
    distance += d0;
 
#ifdef DEBUG_LARWHEELSOLID
    if(Verbose > 2){
        if(Verbose > 3){
            std::cout << MSG_VECTOR(inputP)
                      << " " << MSG_VECTOR(inputV) << std::endl;
        }
        std::cout << "dti: " << d0 << ", DTI: " << distance << std::endl;
    }
    if(Verbose > 3){
        if(d0 < kInfinity){
            G4ThreeVector p2 = inputP + inputV*distance;
            EInside i = Inside(p2);
            std::cout << "DTI hit at dist. " << distance << ", point "
                      << MSG_VECTOR(p2) << ", "
                      << inside(i) << std::endl;
        } else {
            std::cout << "got infinity from dti" << std::endl;
        }
    }
#ifdef LWS_HARD_TEST_DTI
    if(test_dti(inputP, inputV, distance)){
        if(Verbose == 1){
            std::cout << TypeStr() << " DisToIn" << MSG_VECTOR(inputP)
                      << MSG_VECTOR(inputV) << std::endl;
        }
    }
    if(Verbose == 1){
        std::cout << TypeStr() << " DisToIn" << MSG_VECTOR(inputP)
                  << MSG_VECTOR(inputV) << " " << distance << std::endl;
    }
#endif // ifdef LWS_HARD_TEST_DTI
 
#endif // ifdef DEBUG_LARWHEELSOLID
 
  return distance;
}

DistanceToOut() [1/2]

G4double LArWheelSolid::DistanceToOut

(

const G4ThreeVector &

inputP

)

const

Definition at line 11 of file LArWheelSolidDisToOut.cxx.

{
    LWSDBG(1, std::cout << TypeStr() << " DisToOut" << MSG_VECTOR(inputP) << std::endl);
    if(m_BoundingShape->Inside(inputP) != kInside){
        LWSDBG(2, std::cout << "DistanceToOut(p):"
                            << " point " << MSG_VECTOR(inputP)
                            << " is not inside of the m_BoundingShape."
                            << std::endl);
        return 0.;
    }
    G4ThreeVector p( inputP );
    int p_fan = 0;
    const G4double d = m_FanHalfThickness - fabs(GetCalculator()->DistanceToTheNearestFan(p, p_fan));
    if(d < s_Tolerance){
        LWSDBG(2, std::cout << "already not inside " << MSG_VECTOR(p) << std::endl);
        return 0.;
    }
    const G4double d0 = m_BoundingShape->DistanceToOut(inputP);
    LWSDBG(2, std::cout << "dto " << d << " " << d0 << std::endl);
    if(d > d0) return d0;
    else return d;
}

DistanceToOut() [2/2]

G4double LArWheelSolid::DistanceToOut	(	const G4ThreeVector &	inputP,
		const G4ThreeVector &	inputV,
		const G4bool	calcNorm = false,
		G4bool *	validNorm = 0,
		G4ThreeVector *	n = 0 ) const

Definition at line 34 of file LArWheelSolidDisToOut.cxx.

{
    LWSDBG(1, std::cout << TypeStr() << " DisToOut" << MSG_VECTOR(inputP)
                        << MSG_VECTOR(inputV) << std::endl);
 
    const EInside inside_BS = m_BoundingShape->Inside(inputP);
    if(inside_BS == kOutside){
        LWSDBG(2, std::cout << "DistanceToOut(p):"
                            << " point " << MSG_VECTOR(inputP)
                            << " is outside of m_BoundingShape." << std::endl);
        if(calcNorm) *validNorm = false;
        return 0.;
    }
 
    // If here inside or on surface of BS
    G4ThreeVector p(inputP);
    int p_fan = 0;
    const G4double adtnf_p = fabs(GetCalculator()->DistanceToTheNearestFan(p, p_fan));
    if(adtnf_p >= m_FHTplusT) {
        LWSDBG(2, std::cout << "DistanceToOut(p, v): point "
                            << MSG_VECTOR(inputP)
                            << " is outside of solid." << std::endl);
        if(calcNorm) *validNorm = false;
        return 0.;
    }
 
    G4ThreeVector v(inputV);
    const G4double phi0 = p.phi() - inputP.phi();
    v.rotateZ(phi0);
 
#ifdef CHECK_DIRTONORM_ANGLE_ON_SURFACE
    if(adtnf_p < FHTminusT) {
        LWSDBG(5, std::cout << "inside fan point " << MSG_VECTOR(inputP) << ", FHTminusT=" << FHTminusT << std::endl);
    } else {
        LWSDBG(5, std::cout << "on fan surface adtnf_p=" << adtnf_p << ", m_FHTplusT=" << m_FHTplusT << ", FHTminusT=" << FHTminusT << std::endl);
 
        const G4ThreeVector d = GetCalculator()->NearestPointOnNeutralFibre(p, p_fan);
        // check dot product between normal and v
        if ( (p-d).cosTheta(v) > AngularTolerance ) {
            // direction "v" definitely pointing outside
            // return 0.0
            return 0.;
        }
    }
#endif
 
// former distance_to_out starts here
    LWSDBG(4, std::cout << "dto: " << MSG_VECTOR(p) << " "
                        << MSG_VECTOR(v) << std::endl);
 
    G4ThreeVector q(p);
#ifdef LARWHEELSOLID_USE_BS_DTO
    const G4double dto_bs = m_BoundingShape->DistanceToOut(
        inputP, inputV, calcNorm, validNorm, sn
    );
    q = p + v * dto_bs;
    if(q.y() < m_Ymin){
        LWSDBG(5, std::cout << "dto exit point too low " << MSG_VECTOR(q) << std::endl);
        const G4double dy = (m_Ymin - p.y()) / v.y();
        q.setX(p.x() + v.x() * dy);
        q.setY(m_Ymin);
        q.setZ(p.z() + v.z() * dy);
    }
#else
    FanBoundExit_t exit = find_exit_point(p, v, q);
    LWSDBG(5, std::cout << "dto exit " << exit << std::endl);
#endif
    LWSDBG(5, std::cout << "dto exit point " << MSG_VECTOR(q) << std::endl);
 
    G4double distance = 0.;
    G4int start = select_section(p.z());
    G4int stop = select_section(q.z());
    G4int step = -1;
    if(stop > start){ step = 1; start ++; stop ++; }
    LWSDBG(5, std::cout << "dto sections " << start << " " << stop << " " << step << std::endl);
 
    G4double tmp;
    G4ThreeVector p1, C;
 
    for(G4int i = start; i != stop; i += step){
        const G4double d1 = (m_Zsect[i] - p.z()) / v.z();
// v.z() can't be 0, otherwise start == stop, so the exit point could be only
// at z border of the fan section
        LWSDBG(5, std::cout << "at " << i << " dist to zsec = " << d1 << std::endl);
        const G4double x1 = p.x() + v.x() * d1, y1 = p.y() + v.y() * d1;
        p1.set(x1, y1, m_Zsect[i]);
        const G4double dd = fabs(GetCalculator()->DistanceToTheNeutralFibre(p1, p_fan));
        if(dd > m_FHTplusT){
            tmp = out_iteration_process(p, p1, p_fan);
            //while(search_for_most_remoted_point(p, out_section, C)){
            if(search_for_most_remoted_point(p, p1, C, p_fan)){
                tmp = out_iteration_process(p, C, p_fan);
            }
            distance += tmp;
#ifndef LARWHEELSOLID_USE_BS_DTO
            exit = NoCross;
#endif
            goto end_dto;
        }
        if(search_for_most_remoted_point(p, p1, C, p_fan)){
            distance += out_iteration_process(p, C, p_fan);
#ifndef LARWHEELSOLID_USE_BS_DTO
            exit = NoCross;
#endif
            goto end_dto;
        }
        distance += d1;
        p.set(x1, y1, m_Zsect[i]);
    }
 
    if(fabs(GetCalculator()->DistanceToTheNeutralFibre(q, p_fan)) > m_FHTplusT){
        LWSDBG(5, std::cout << "q=" << MSG_VECTOR(q) << " outside fan cur distance=" << distance << ", m_FHTplusT=" << m_FHTplusT << std::endl);
        tmp = out_iteration_process(p, q, p_fan);
#ifndef LARWHEELSOLID_USE_BS_DTO
        exit = NoCross;
#endif
    } else {
        tmp = (q - p).mag();
    }
    //while(search_for_most_remoted_point(out, out1, C, p_fan)){
    if(search_for_most_remoted_point(p, q, C, p_fan)){
        tmp = out_iteration_process(p, C, p_fan);
#ifndef LARWHEELSOLID_USE_BS_DTO
        exit = NoCross;
#endif
    }
    distance += tmp;
// former distance_to_out ends here
  end_dto:
    LWSDBG(5, std::cout << "At end_dto distance=" << distance  << std::endl);
#ifdef LARWHEELSOLID_USE_BS_DTO
    if(calcNorm && distance < dto_bs) *validNorm = false;
#else
    if(calcNorm){
        LWSDBG(5, std::cout << "dto calc norm " << exit << std::endl);
        switch(exit){
            case ExitAtBack:
                sn->set(0., 0., 1.);
                *validNorm = true;
                break;
            case ExitAtFront:
                sn->set(0., 0., -1.);
                *validNorm = true;
                break;
            case ExitAtOuter:
                q.rotateZ(-phi0);
                sn->set(q.x(), q.y(), 0.);
                if(q.z() <= m_Zmid) sn->setZ(- q.perp() * m_fs->Amax);
                sn->setMag(1.0);
                *validNorm = true;
                break;
            default:
                *validNorm = false;
                break;
        }
    }
#endif
 
#ifdef DEBUG_LARWHEELSOLID
    if(Verbose > 2){
        std::cout << "DTO: " << distance << " ";
        if (*validNorm) {
          std::cout << *validNorm << " " << MSG_VECTOR((*sn));
        } else {
          std::cout << "Norm is not valid";
        }
        std::cout << std::endl;
        if(Verbose > 3){
            G4ThreeVector p2 = inputP + inputV * distance;
            EInside i = Inside(p2);
            std::cout << "DTO hit at " << MSG_VECTOR(p2) << ", "
                      << inside(i) << std::endl;
        }
    }
#ifdef LWS_HARD_TEST_DTO
    if(test_dto(inputP, inputV, distance)){
        if(Verbose == 1){
            std::cout << TypeStr() << " DisToOut" << MSG_VECTOR(inputP)
                      << MSG_VECTOR(inputV) << std::endl;
        }
    }
#endif
#endif
    return distance;
}

find_exit_point()

LArWheelSolid::FanBoundExit_t LArWheelSolid::find_exit_point ( const G4ThreeVector & p, const G4ThreeVector & v, G4ThreeVector & q ) const

private

Definition at line 216 of file LArFanSection.cxx.

{
    LWSDBG(6, std::cout << "in fep p" << MSG_VECTOR(p)<< ", v"<< MSG_VECTOR(v) << ", q" << MSG_VECTOR(q) << std::endl);
 
/* by construction, cannot have true from both upper and lower */
/* the only problem is the points on surface but "slightly outside" */
/* fs_cross_* account for (x, z) range */
// lower has to be checked first, since outer might find more distant
// intersection in the acceptable (x, z) range
    if(fs_cross_lower(p, v, q))   return ExitAtInner;
    LWSDBG(6, std::cout << "after fs_cross_lower q" << MSG_VECTOR(q) << std::endl);
    if(fs_cross_upper(p, v, q))   return ExitAtOuter;
    LWSDBG(6, std::cout << "after fs_cross_upper q" << MSG_VECTOR(q) << std::endl);
 
    FanBoundExit_t result = ExitAtSide;
    G4double d;
    if(v.x() > 0.) d = (m_fs->xmax - p.x()) / v.x();
    else if(v.x() < 0.) d = (m_fs->xmin - p.x()) / v.x();
    else d = kInfinity;
 
    G4double dz;
    FanBoundExit_t resultz = NoCross;
    if(v.z() > 0.){
        dz = (m_Zsect.back() - p.z()) / v.z();
        resultz = ExitAtBack;
    } else if(v.z() < 0.){
        dz = (m_Zsect.front() - p.z()) / v.z();
        resultz = ExitAtFront;
    } else {
        dz = kInfinity;
    }
    if(d > dz){
        d = dz;
        result = resultz;
    }
    q = p + v * d;
    LWSDBG(7, std::cout << "fep side " << d << " " << result << " q" << MSG_VECTOR(q) << std::endl);
    const G4double out_distlower = m_fs->Amin*q.z() + m_fs->Bmin - q.perp(); // > 0 - below lower cone
    LWSDBG(7, std::cout << "fep out_distlower(q)=" << out_distlower << " Tolerance=" << s_Tolerance << std::endl);
    if (out_distlower >= 0.0) {
        // side intersection point is below lower cone
        // initial point p was at exit boundary
        q = p;
        return NoCross;
    }
 
    if (m_IsOuter && q.z() >= m_Zmid && q.z() <= m_Zsect.back()+s_Tolerance && q.perp2() >= m_fs->Cflat2) {
        // outside of upper cylinder
        q = p;
        return NoCross;
    }
    const G4double out_distupper = m_fs->Amax*q.z() + m_fs->Bmax - q.perp(); // < 0 - above upper cone
    LWSDBG(7, std::cout << "fep out_distupper(q)=" << out_distupper << " Tolerance=" << s_Tolerance << std::endl);
    if (out_distupper <= 0.0) {
        // side intersection point is above upper cone
        // initial point p was at exit boundary
        q = p;
        return NoCross;
    }
    assert((q - p).mag() < kInfinity);
    return result;
}

fs_cross_lower()

bool LArWheelSolid::fs_cross_lower ( const G4ThreeVector & p, const G4ThreeVector & v, G4ThreeVector & q ) const

private

Definition at line 103 of file LArFanSection.cxx.

{
    LWSDBG(7, std::cout << "fcl" << std::endl);
    const G4double A = v.perp2() - m_fs->Amin2*v.z()*v.z();
    const G4double B = p.x()*v.x() + p.y()*v.y()
                     - m_fs->Amin2*p.z()*v.z() - m_fs->ABmin*v.z();
    const G4double C = p.perp2() - m_fs->Amin2*p.z()*p.z()
                     - 2.*m_fs->ABmin*p.z() - m_fs->Bmin2; 
                     
    G4double t1(0.0);
    const G4double out_dist = m_fs->Amin*p.z() + m_fs->Bmin - p.perp();
    LWSDBG(8, std::cout << "fcl out_dist(p)=" << out_dist << " Tolerance=" << s_Tolerance << std::endl);
    const G4bool out = out_dist >= 0.0;
    if(check_D(t1, A, B, C, out)){
        const G4double zz1 = p.z() + v.z() * t1;
        if(zz1 < m_Zsect.front() || zz1 > m_Zsect.back()){
            LWSDBG(8, std::cout << "fcl out on Z " << zz1 << std::endl);
            return false;
        }
        const G4double xx1 = p.x() + v.x() * t1;
        if(xx1 < m_fs->xmin || xx1 > m_fs->xmax){
            LWSDBG(8, std::cout << "fcl out on X " << xx1 << std::endl);
            return false;
        }
        if(out_dist == 0.){ // entry point is exactly on the cone
        // here we got t1 > 0 from check_D, founded point seems to be in x and z ranges
        // if the track leaves the surface, then the entry is the intersection,
        // and the distance is 0
        // if the track is on the surface, then there is no lower cone intersection
 
        // estimate deviation of the track from the surface
        // (exact calculations are too complicated)
            const G4double xx2 = p.x() + v.x() * t1 * 0.5;
            const G4double yy2 = p.y() + v.y() * t1 * 0.5;
            const G4double dev = fabs(sqrt(xx2 *xx2 + yy2*yy2)
                               - m_fs->Amin*(p.z() + zz1)*0.5
                               - m_fs->Bmin);
            if(dev < s_Tolerance){
                LWSDBG(8, std::cout << "fcl on the surface" << std::endl);
                return false;
            } else {
                LWSDBG(8, std::cout << "fcl out = in" << std::endl);
                q = p;
                return true;
            }
        }
        q.setX(xx1);
        q.setY(p.y() + v.y() * t1);
        q.setZ(zz1);
        LWSDBG(8, std::cout << "fcl got " << t1 << std::endl);
        return true;
    }
    LWSDBG(8, std::cout << "fcl no intersection" << std::endl);
    return false;
}

fs_cross_upper()

bool LArWheelSolid::fs_cross_upper ( const G4ThreeVector & p, const G4ThreeVector & v, G4ThreeVector & q ) const

private

Definition at line 165 of file LArFanSection.cxx.

{
    LWSDBG(7, std::cout << "fcu" << std::endl);
    G4double A = v.perp2();
    G4double B = p.x()*v.x() + p.y()*v.y();
    G4double C = p.perp2();
 
    if(m_IsOuter){
        const G4double &Af = A, &Bf = B;
        const G4double Cf = C - m_fs->Cflat2;
        
        G4double b1;
        if(check_D(b1, Af, Bf, Cf, Cf >= 0.)){
            const G4double zz1 = p.z() + v.z() * b1;
            if(zz1 >= m_Zmid && zz1 <= m_Zsect.back()){
                const G4double xx1 = p.x() + v.x() * b1;
                if(xx1 < m_fs->xmin || xx1 > m_fs->xmax) return false;
                q.setX(xx1);
                q.setY(p.y() + v.y() * b1);
                q.setZ(zz1);
                return true;
            }
        }
        LWSDBG(8, std::cout << "fcu no cyl intersection" << std::endl);
    }
 
    A -= m_fs->Amax2*v.z()*v.z();
    B -= m_fs->Amax2*p.z()*v.z() + m_fs->ABmax*v.z();
    C -= m_fs->Amax2*p.z()*p.z() + 2.*m_fs->ABmax*p.z() + m_fs->Bmax2;
 
    G4double t1;
    const G4bool out = m_fs->Amax*p.z() + m_fs->Bmax <= p.perp();
    if(check_D(t1, A, B, C, out)){
        const G4double zz1 = p.z() + v.z() * t1;
        LWSDBG(8, std::cout << "fcu z = " << zz1 << ", lim: (" << m_Zsect.front() << ", " << m_Zmid << ")" << std::endl);
        if(zz1 < m_Zsect.front() || zz1 > m_Zmid) return false;
        const G4double xx1 = p.x() + v.x() * t1;
        LWSDBG(8, std::cout << "fcu x = " << xx1 << ", lim: (" << m_fs->xmin << ", " << m_fs->xmax << ")" << std::endl);
        if(xx1 < m_fs->xmin || xx1 > m_fs->xmax) return false;
        q.setX(xx1);
        q.setY(p.y() + v.y() * t1);
        q.setZ(zz1);
        return true;
    }
    LWSDBG(8, std::cout << "fcu no cone intersection" << std::endl);
    return false;
}

get_area_at_r()

G4double LArWheelSolid::get_area_at_r ( G4double r ) const

private

Definition at line 511 of file LArWheelSolidTests.cxx.

{
  m_f_area->SetParameter(0, r);
 
  double zmin = m_BoundingShape->DistanceToIn(
                                               G4ThreeVector(0., r, m_Zmin), G4ThreeVector(0., 0., 1.)
                                               );
  double zmax = m_BoundingShape->DistanceToIn(
                                               G4ThreeVector(0., r, m_Zmax), G4ThreeVector(0., 0., -1.)
                                               );
  zmax = m_Zmax - zmax;
 
  double result = m_f_area->Integral(zmin, zmax);
 
  return result;
}

get_area_on_face()

G4double LArWheelSolid::get_area_on_face ( void ) const

private

Definition at line 333 of file LArWheelSolidTests.cxx.

{
  G4double result = 0.;
  G4double rmin, rmax;
  get_r(m_BoundingShape, m_Zmin, rmin, rmax);
  result += rmax - rmin;
  get_r(m_BoundingShape, m_Zmax, rmin, rmax);
  result += rmax - rmin;
  result *= GetCalculator()->GetFanHalfThickness() * 2.;
  return result;
}

get_area_on_polycone()

G4double LArWheelSolid::get_area_on_polycone ( void ) const

private

Definition at line 328 of file LArWheelSolidTests.cxx.

{
  return m_f_area_on_pc->Integral(m_Zmin, m_Zmax);
}

get_area_on_side()

G4double LArWheelSolid::get_area_on_side ( void ) const

private

Definition at line 361 of file LArWheelSolidTests.cxx.

{
  return m_f_side_area->Integral(m_Rmin, m_Rmax, IntPrecision);
}

get_length_at_r()

G4double LArWheelSolid::get_length_at_r ( G4double r ) const

private

Definition at line 345 of file LArWheelSolidTests.cxx.

{
  m_f_length->SetParameter(0, r);
 
  double zmin = m_BoundingShape->DistanceToIn(
    G4ThreeVector(0., r, -m_Zmin), G4ThreeVector(0., 0., 1.)
  );
  zmin -= m_Zmin * 2.;
  double zmax = m_BoundingShape->DistanceToIn(
                                               G4ThreeVector(0., r, m_Zmax), G4ThreeVector(0., 0., -1.)
                                               );
  zmax = m_Zmax - zmax;
  double result = m_f_length->Integral(zmin, zmax);
  return result;
}

get_point_on_accordion_surface()

void LArWheelSolid::get_point_on_accordion_surface ( G4ThreeVector & p ) const

private

Definition at line 99 of file LArWheelSolidTests.cxx.

{
  p[0] = 0.;
  p[1] = 0.;
  p[2] = rnd->Uniform(m_Zmin, m_Zmax);
 
  G4double rmin, rmax;
  get_r(m_BoundingShape, p[2], rmin, rmax);
 
  p[1] = rnd->Uniform(rmin, rmax);
  p.setPhi(rnd->Uniform(m_MinPhi, m_MaxPhi));
  G4double dphi = p.phi();
  int p_fan = 0;
  G4double d = GetCalculator()->DistanceToTheNearestFan(p, p_fan);
  dphi -= p.phi();
 
  G4int side = 0;
  if(d < 0.) side = -1;
  if(d >= 0.) side = 1;
 
  G4double a = GetCalculator()->AmplitudeOfSurface(p, side, p_fan);
  p[0] = a;
 
  p.rotateZ(dphi);
 
  if(m_BoundingShape->Inside(p) == kOutside){
    G4ThreeVector D = p; D[2] = 0.;
    G4double d1 = m_BoundingShape->DistanceToIn(p, D);
    if(d1 > 10.*CLHEP::m){
      D *= -1;
      d1 = m_BoundingShape->DistanceToIn(p, D);
    }
    if(d1 > 10.*CLHEP::m){
      set_failover_point(p, "acc fail0");
      return;
    }
    d1 *= 2.;
 
    G4ThreeVector B = p + D * d1;
    G4double dphi = B.phi();
    int B_fan = 0;
    GetCalculator()->DistanceToTheNearestFan(B, B_fan);
    dphi -= B.phi();
 
    B[0] = GetCalculator()->AmplitudeOfSurface(B, side, B_fan);
    B.rotateZ(dphi);
    EInside Bi = m_BoundingShape->Inside(B);
    if(Bi == kSurface){
      p = B;
      return;
    }
    if(Bi == kOutside){
      set_failover_point(p, "acc fail1");
      return;
    }
    G4ThreeVector D1 = (p - B).unit();
    G4ThreeVector X = B + D1 * m_BoundingShape->DistanceToOut(B, D1);
    if(Inside(X) == kSurface){
      p = X;
    } else { // failed
      set_failover_point(p, "acc fail2");
      return;
    }
  }
}

get_point_on_flat_surface()

void LArWheelSolid::get_point_on_flat_surface ( G4ThreeVector & p ) const

private

Definition at line 288 of file LArWheelSolidTests.cxx.

{
  p[0] = 0.;
  p[1] = 0.;
  p[2] = rnd->Uniform() > 0.5? m_Zmin: m_Zmax;
 
  G4double rmin, rmax;
  get_r(m_BoundingShape, p[2], rmin, rmax);
 
  p[1] = rnd->Uniform(rmin, rmax);
  p.setPhi(rnd->Uniform(m_MinPhi, m_MaxPhi));
  G4double dphi = p.phi();
  int p_fan = 0;
  GetCalculator()->DistanceToTheNearestFan(p, p_fan);
  dphi -= p.phi();
 
  p[0] = rnd->Uniform(
                      GetCalculator()->AmplitudeOfSurface(p, -1, p_fan),
                      GetCalculator()->AmplitudeOfSurface(p, 1, p_fan)
                      );
 
  p.rotateZ(dphi);
 
  if(m_BoundingShape->Inside(p) != kSurface){
    set_failover_point(p, "flat fail");
  }
}

get_point_on_polycone_surface()

void LArWheelSolid::get_point_on_polycone_surface ( G4ThreeVector & p ) const

private

Definition at line 165 of file LArWheelSolidTests.cxx.

{
  const G4double z = rnd->Uniform(m_Zmin, m_Zmax);
  G4double rmin, rmax;
  get_r(m_BoundingShape, z, rmin, rmax);
  const bool inner = rnd->Uniform() > 0.5? true: false;
 
  p[0] = 0.; p[1] = inner? rmin: rmax; p[2] = z;
  p.setPhi(rnd->Uniform(m_MinPhi, m_MaxPhi));
  G4double dphi = p.phi();
  int p_fan = 0;
  GetCalculator()->DistanceToTheNearestFan(p, p_fan);
  dphi -= p.phi();
 
  const G4double r = p[1];
 
  G4ThreeVector A1(0., r, z);
  A1[0] = GetCalculator()->AmplitudeOfSurface(A1, -1, p_fan);
  A1.rotateZ(dphi);
  EInside A1i = m_BoundingShape->Inside(A1);
  //    EInside A1a = Inside_accordion(A1);
  //std::cout << "A1: " << A1i << " " << A1a << std::endl;
  if(A1i == kSurface){
    //std::cout << "got A1" << std::endl;
    p = A1;
    return;
  }
 
  G4ThreeVector A2(0., r, z);
  A2[0] = GetCalculator()->AmplitudeOfSurface(A2, 1, p_fan);
  A2.rotateZ(dphi);
  EInside A2i = m_BoundingShape->Inside(A2);
  //    EInside A2a = Inside_accordion(A2);
  //std::cout << "A2: " << A2i << " " << A2a << std::endl;
  if(A2i == kSurface){
    //std::cout << "got A2" << std::endl;
    p = A2;
    return;
  }
 
  if(A1i != A2i){
    if(A2i == kOutside){
      std::swap(A1, A2);
      std::swap(A1i, A2i);
    }
    // here A1 is outside BP, A2 is inside BP
    G4ThreeVector d = (A2 - A1).unit();
    p = A1 + d * m_BoundingShape->DistanceToIn(A1, d);
    //std::cout << "got A1<->A2" << std::endl;
    return;
  }
  // here A1i == A2i
 
  G4double step;
  if(A1i == kInside){
    G4double d1 = m_BoundingShape->DistanceToOut(A1);
    G4double d2 = m_BoundingShape->DistanceToOut(A2);
    step = d1 > d2? d1 : d2;
    if(inner) step *= -2;
    else step *= 2;
  } else {
    G4double d1 = m_BoundingShape->DistanceToIn(A1);
    G4double d2 = m_BoundingShape->DistanceToIn(A2);
    step = d1 > d2? d1 : d2;
    if(inner) step *= 2;
    else step *= -2;
  }
 
  G4ThreeVector B1(0., r + step, z);
  B1[0] = GetCalculator()->AmplitudeOfSurface(B1, -1, p_fan);
  B1.rotateZ(dphi);
  EInside B1i = m_BoundingShape->Inside(B1);
  //    EInside B1a = Inside_accordion(B1);
  //std::cout << "B1: " << B1i << " " << B1a << std::endl;
  if(B1i == kSurface){
    //std::cout << "got B1" << std::endl;
    p = B1;
    return;
  }
  G4ThreeVector B2(0., r + step, z);
  B2[0] = GetCalculator()->AmplitudeOfSurface(B2, 1, p_fan);
  B2.rotateZ(dphi);
  EInside B2i = m_BoundingShape->Inside(B2);
  //    EInside B2a = Inside_accordion(B2);
  //std::cout << "B2: " << B2i << " " << B2a << std::endl;
  if(B2i == kSurface){
    //std::cout << "got B2" << std::endl;
    p = B2;
    return;
  }
 
  if(B1i == A1i || B2i == A1i){ // failed
    set_failover_point(p, "pol fail1");
    return;
  }
  if(A1i == kInside){
    std::swap(A1, B1);
    std::swap(A2, B2);
    std::swap(A1i, B1i);
    std::swap(A2i, B2i);
  }
  // here A* outside, B* inside, all on accordion surface
 
  G4ThreeVector d1 = (A1 - B1).unit();
  G4ThreeVector X1 = B1 + d1 * m_BoundingShape->DistanceToOut(B1, d1);
  G4ThreeVector d2 = (A2 - B2).unit();
  G4ThreeVector X2 = B2 + d2 * m_BoundingShape->DistanceToOut(B2, d2);
 
  G4ThreeVector X = X1;
  G4double phi1 = X1.phi(), phi2 = X2.phi();
  // X1 corresponds to side = -, X2 to side = +
  if(phi1 > 0. && phi2 < 0.) phi2 += CLHEP::twopi;
  G4double phiX = rnd->Uniform(phi1, phi2);
  if(phiX > CLHEP::pi) phiX -= CLHEP::twopi;
  X.setPhi(phiX);
 
  if(Inside(X) == kSurface){
    p = X;
  } else { // failed
    set_failover_point(p, "pol fail2");
  }
}

GetBoundingShape()

const G4VSolid * LArWheelSolid::GetBoundingShape ( void ) const

inline

Definition at line 129 of file LArWheelSolid.h.

{ return m_BoundingShape; }

GetCalculator()

const LArWheelCalculator * LArWheelSolid::GetCalculator ( void ) const

inline

Definition at line 130 of file LArWheelSolid.h.

{ return m_Calculator; }

GetCubicVolume()

G4double LArWheelSolid::GetCubicVolume

(

void

)

Definition at line 316 of file LArWheelSolidTests.cxx.

{
  // sagging ignored, effect should be negligible
    double result = m_f_vol->Integral(m_Rmin, m_Rmax, IntPrecision)
 
#ifndef LOCAL_DEBUG
    * GetCalculator()->GetNumberOfFans()
#endif
    ;
  return result;
}

GetEntityType()

G4GeometryType LArWheelSolid::GetEntityType

(

)

const

Definition at line 65 of file LArWheelSolid.cxx.

{
  switch(m_Type){
  case InnerAbsorberWheel:
    return G4String("LArInnerAbsorberWheel");
    break;
  case OuterAbsorberWheel:
    return G4String("LArOuterAbsorberWheel");
    break;
  case InnerElectrodWheel:
    return G4String("LArInnerElecrodWheel");
    break;
  case OuterElectrodWheel:
    return G4String("LArOuterElecrodWheel");
    break;
  case InnerAbsorberModule:
    return G4String("LArInnerAbsorberModule");
    break;
  case OuterAbsorberModule:
    return G4String("LArOuterAbsorberModule");
    break;
  case InnerElectrodModule:
    return G4String("LArInnerElecrodModule");
    break;
  case OuterElectrodModule:
    return G4String("LArOuterElecrodModule");
    break;
  case InnerGlueWheel:
    return G4String("LArInnerGlueWheel");
    break;
  case OuterGlueWheel:
    return G4String("LArOuterGlueWheel");
    break;
  case InnerLeadWheel:
    return G4String("LArInnerLeadWheel");
    break;
  case OuterLeadWheel:
    return G4String("LArOuterLeadWheel");
    break;
  case InnerAbsorberCone:
    return G4String("LArInnerAbsorberCone");
    break;
  case InnerElectrodCone:
    return G4String("LArInnerElectrodCone");
    break;
  case InnerGlueCone:
    return G4String("LArInnerGlueCone");
    break;
  case InnerLeadCone:
    return G4String("LArInnerLeadCone");
    break;
  case OuterAbsorberFrontCone:
    return G4String("LArOuterAbsorberFrontCone");
    break;
  case OuterElectrodFrontCone:
    return G4String("LArOuterElectrodFrontCone");
    break;
  case OuterGlueFrontCone:
    return G4String("LArOuterGlueFrontCone");
    break;
  case OuterLeadFrontCone:
    return G4String("LArOuterLeadFrontCone");
    break;
  case OuterAbsorberBackCone:
    return G4String("LArOuterAbsorberBackCone");
    break;
  case OuterElectrodBackCone:
    return G4String("LArOuterElectrodBackCone");
    break;
  case OuterGlueBackCone:
    return G4String("LArOuterGlueBackCone");
    break;
  case OuterLeadBackCone:
    return G4String("LArOuterLeadBackCone");
    break;
  default:
    G4Exception("LArWheelSolid", "UnknownSolidType", FatalException,"GetEntityType: Unknown LArWheelType.");
  }
  return G4String("");
}

GetExtent()

G4VisExtent LArWheelSolid::GetExtent

(

)

const

Definition at line 151 of file LArWheelSolid.cxx.

{
  return m_BoundingShape->GetExtent();
}

GetPointOnSurface()

G4ThreeVector LArWheelSolid::GetPointOnSurface

(

void

)

const

Definition at line 66 of file LArWheelSolidTests.cxx.

{
  if(rnd == 0){
    rnd = new TRandom3(0);
  }
 
  G4double r = rnd->Uniform();
 
  G4ThreeVector p(0., 0., 0.);
 
  G4double level1 = .980;
  G4double level2 = .993;
  const char *v = getenv("LARWHEELSOLID_TEST_MODE_LEVEL1");
  if(v) level1 = atof(v);
  const char *v1 = getenv("LARWHEELSOLID_TEST_MODE_LEVEL2");
  if(v1) level2 = atof(v1);
 
#if LOCAL_DEBUG > 1
  std::cout << "LWS::GPOS " << r << std::endl;
#endif
 
  if(r <= level1){
    get_point_on_accordion_surface(p);
  } else if(r <= level2){
    get_point_on_polycone_surface(p);
  } else if(r <= 1.){
    get_point_on_flat_surface(p);
  } else {
    G4Exception("LArWheelSolid", "Rnd generator error", FatalException, "GetPointOnSurface: Wrong data from rnd generator");
  }
  return p;
}

GetSurfaceArea()

G4double LArWheelSolid::GetSurfaceArea

(

void

)

Definition at line 366 of file LArWheelSolidTests.cxx.

{
  double result = 0.;
 
  double a1 = get_area_on_polycone();
  result += a1;
#ifdef LOCAL_DEBUG
  std::cout << "get_area_on_polycone: " << a1/mm2 << std::endl;
#endif
 
  double a2 = get_area_on_face();
  result += a2;
#ifdef LOCAL_DEBUG
  std::cout << "get_area_on_face: " << a2/mm2 << std::endl;
#endif
 
  double a3 = get_area_on_side();
  result += a3;
#ifdef LOCAL_DEBUG
  std::cout << "get_area_on_side: " << a3/mm2 << std::endl;
#endif
 
  // sagging ignored, effect should be negligible
  return result
#ifndef LOCAL_DEBUG
    * GetCalculator()->GetNumberOfFans()
#endif
    ;
}

GetType()

LArWheelSolid_t LArWheelSolid::GetType ( void ) const

inline

Definition at line 131 of file LArWheelSolid.h.

{ return m_Type; }

in_iteration_process()

G4double LArWheelSolid::in_iteration_process ( const G4ThreeVector & p, G4double dist_p, G4ThreeVector & B, int p_fan ) const

private

Definition at line 126 of file LArWheelSolidDisToIn.cxx.

{
    LWSDBG(6, std::cout << "iip from " << p << " to " << B
                        << " dir " << (B - p).unit()
                        << std::endl);
 
    G4ThreeVector A, C, diff;
    A = p;
    G4double dist_c;
    unsigned int niter = 0;
  //  assert(fabs(GetCalculator()->DistanceToTheNeutralFibre(A)) > m_FHTplusT);
  //  assert(GetCalculator()->DistanceToTheNeutralFibre(A) == dist_p);
    do {
        C = A + B;
        C *= 0.5;
        dist_c = GetCalculator()->DistanceToTheNeutralFibre(C, p_fan);
        if(dist_c * dist_p < 0. || fabs(dist_c) < m_FHTminusT){
            B = C;
        } else {
            A = C;
        }
        niter ++;
        diff = A - B;
    } while(diff.mag2() > s_IterationPrecision2 && niter < s_IterationsLimit);
    assert(niter < s_IterationsLimit);
    assert(fabs(GetCalculator()->DistanceToTheNeutralFibre(B, p_fan)) < m_FHTplusT);
    diff = p - B;
    LWSDBG(7, std::cout << "iip result in " << niter << " = " << B
                        << " " << diff.mag() << std::endl);
    return diff.mag();
}

init_tests()

void LArWheelSolid::init_tests ( void )

private

Definition at line 610 of file LArWheelSolidTests.cxx.

{
  m_test_index = double(solid.size());
  solid[m_test_index] = this;
 
#ifdef DEBUG_LARWHEELSOLID
    std::cout << "LArWheelSolid::init_tests: put " << this
              << " with index " << m_test_index << std::endl;
#endif
 
  m_f_area = new TF1(
                   (GetName() + "_f_area").c_str(), &LArWheelSolid_fcn_area,
                   m_Zmin, m_Zmax, 2
                   );
  m_f_area->FixParameter(1, m_test_index);
 
  m_f_vol = new TF1(
                  (GetName() + "_f_vol").c_str(), &LArWheelSolid_fcn_vol,
                  m_Rmin, m_Rmax, 1
                  );
  m_f_vol->FixParameter(0, m_test_index);
 
  m_f_area_on_pc = new TF1(
                         (GetName() + "_f_area_pc").c_str(), &LArWheelSolid_fcn_area_on_pc,
                         m_Zmin, m_Zmax, 1
                         );
  m_f_area_on_pc->FixParameter(0, m_test_index);
 
  m_f_length = new TF1(
                     (GetName() + "_f_length").c_str(), &fcn_length,
                     m_Zmin, m_Zmax, 2
                     );
  m_f_length->FixParameter(1, m_test_index);
 
  m_f_side_area = new TF1(
                        (GetName() + "_f_side_area").c_str(), &LArWheelSolid_fcn_side_area,
                        m_Rmin, m_Rmax, 1
                        );
  m_f_side_area->FixParameter(0, m_test_index);
}

initMessaging()

void AthMessaging::initMessaging ( ) const

privateinherited

Initialize our message level and MessageSvc.

This method should only be called once.

Definition at line 39 of file AthMessaging.cxx.

{
  m_imsg = Athena::getMessageSvc();
  // If user did not set an explicit level, set a default
  if (m_lvl == MSG::NIL) {
    m_lvl = m_imsg ?
      static_cast<MSG::Level>( m_imsg.load()->outputLevel(m_nm) ) :
      MSG::INFO;
  }
}

inner_solid_init()

void LArWheelSolid::inner_solid_init ( const G4String & bs_name )

private

Definition at line 206 of file LArWheelSolidInit.cxx.

{
    m_IsOuter = false;
    m_FanPhiAmplitude = 0.065; // internal technical constant, should not go in DB
    set_phi_size();
 
    std::array<G4double,2> zPlane{}, rInner{}, rOuter{};
    zPlane[1] = GetCalculator()->GetWheelThickness();
    G4double wheel_thickness = zPlane[1] - zPlane[0];
    GetCalculator()->GetWheelInnerRadius(rInner);
    GetCalculator()->GetWheelOuterRadius(rOuter);
    const G4double phi_min = m_PhiPosition - m_FanPhiAmplitude
                           - GetCalculator()->GetFanStepOnPhi() * 2;
 
    m_Zmin = zPlane[0]; m_Zmax = zPlane[1];
    m_Rmin = rInner[0]; m_Rmax = rOuter[1];
    m_Ymin = m_Rmin * 0.9;
    m_Zmid = zPlane[1];
    m_Ymid = (rInner[0] + rOuter[1]) * 0.5;
 
    if(m_Type == InnerAbsorberCone
    || m_Type == InnerElectrodCone
    || m_Type == InnerGlueCone
    || m_Type == InnerLeadCone
    ){
        m_BoundingShape = new G4ShiftedCone(
            bs_name + "Cone", zPlane[0], zPlane[1],
            rInner[0], rOuter[0], rInner[1], rOuter[1]
        );
    } else {
        m_BoundingShape = new G4Polycone(
            bs_name + "Polycone", m_MinPhi, m_MaxPhi - m_MinPhi,
            2, zPlane.data(), rInner.data(), rOuter.data()
        );
    }
#ifdef LARWHEELSOLID_USE_FANBOUND
    const G4double phi_size = (m_FanPhiAmplitude + GetCalculator()->GetFanStepOnPhi() * 2) * 2;
    FanBound = new G4Polycone(bs_name + "ofFan", phi_min, phi_size,
                              2, zPlane, rInner, rOuter);
#endif
#ifndef PORTABLE_LAR_SHAPE
    ATH_MSG_INFO(m_BoundingShape->GetName() + " is the m_BoundingShape");
#endif
    
    const G4double half_wave_length = GetCalculator()->GetHalfWaveLength();
    const G4double sss = GetCalculator()->GetStraightStartSection();
    m_Zsect.push_back(0.);
    m_Zsect.push_back(sss + half_wave_length * 0.25);
    const G4int num_fs = GetCalculator()->GetNumberOfHalfWaves() + 1;
    for(G4int i = 2; i < num_fs; i ++){
          const G4double zi = half_wave_length * (i - 1) + sss;
#if LARWHEELSOLID_ZSECT_MULT > 1
          for(G4int j = LARWHEELSOLID_ZSECT_MULT - 1; j > 0; -- j){
            m_Zsect.push_back(zi - half_wave_length * j / LARWHEELSOLID_ZSECT_MULT);
          }
#endif
          m_Zsect.push_back(zi);
    }
    m_Zsect.push_back(wheel_thickness - m_Zsect[1]);
    m_Zsect.push_back(wheel_thickness - m_Zsect[0]);
 
    m_fs = new LArFanSections(
        rInner[0], rInner[1], rOuter[0], rOuter[1],
        m_Rmax*cos(phi_min), m_Zsect.front(), m_Zsect.back()
    );
}

Inside()

EInside LArWheelSolid::Inside

(

const G4ThreeVector &

inputP

)

const

Definition at line 15 of file LArWheelSolid.cxx.

{
    LWSDBG(10, std::cout << std::setprecision(25));
    LWSDBG(1, std::cout << TypeStr() << " Inside " << MSG_VECTOR(inputP) << std::endl);
    const EInside inside_BS = m_BoundingShape->Inside(inputP);
    if(inside_BS == kOutside){
        LWSDBG(2, std::cout << "outside BS" << std::endl);
        return kOutside;
    }
    G4ThreeVector p( inputP );
    int p_fan = 0;
    const G4double d = fabs(GetCalculator()->DistanceToTheNearestFan(p, p_fan));
    if(d > m_FHTplusT){
        LWSDBG(2, std::cout << "outside fan d=" << d << ", m_FHTplusT=" << m_FHTplusT << std::endl);
        return kOutside;
    }
    if(d < m_FHTminusT){
        LWSDBG(2, std::cout << "inside fan d=" << d << ", m_FHTminusT=" << m_FHTminusT << ", inside_BS=" << inside(inside_BS) << std::endl);
        return inside_BS;
    }
    LWSDBG(2, std::cout << "surface" << std::endl);
    return kSurface;
}

Inside_accordion()

EInside LArWheelSolid::Inside_accordion ( const G4ThreeVector & p ) const

private

Definition at line 29 of file LArWheelSolidTests.cxx.

{
  G4ThreeVector p1 = p;
  int p1_fan = 0;
  G4double d = m_FanHalfThickness - fabs(GetCalculator()->DistanceToTheNearestFan(p1, p1_fan));
  if(d > s_Tolerance) return kInside;
  if(d > -s_Tolerance) return kSurface;
  return kOutside;
}

msg() [1/2]

MsgStream & AthMessaging::msg ( ) const

inlineinherited

The standard message stream.

Returns a reference to the default message stream May not be invoked before sysInitialize() has been invoked.

Definition at line 163 of file AthMessaging.h.

{
  MsgStream* ms = m_msg_tls.get();
  if (!ms) {
    if (!m_initialized.test_and_set()) initMessaging();
    ms = new MsgStream(m_imsg,m_nm);
    m_msg_tls.reset( ms );
  }
 
  ms->setLevel (m_lvl);
  return *ms;
}

msg() [2/2]

MsgStream & AthMessaging::msg ( const MSG::Level lvl ) const

inlineinherited

The standard message stream.

Returns a reference to the default message stream May not be invoked before sysInitialize() has been invoked.

Definition at line 178 of file AthMessaging.h.

{ return msg() << lvl; }

msgLvl()

bool AthMessaging::msgLvl ( const MSG::Level lvl ) const

inlineinherited

Test the output level.

Parameters

lvl	The message level to test against

Returns

boolean Indicating if messages at given level will be printed

Return values

true	Messages at level "lvl" will be printed

Definition at line 151 of file AthMessaging.h.

{
  if (m_lvl <= lvl) {
    msg() << lvl;
    return true;
  } else {
    return false;
  }
}

out_iteration_process()

G4double LArWheelSolid::out_iteration_process ( const G4ThreeVector & p, G4ThreeVector & B, const int p_fan ) const

private

Definition at line 229 of file LArWheelSolidDisToOut.cxx.

{
    LWSDBG(6, std::cout << "oip: " << p << " " << B);
    assert(fabs(GetCalculator()->DistanceToTheNeutralFibre(p, p_fan)) < m_FHTplusT);
    assert(fabs(GetCalculator()->DistanceToTheNeutralFibre(B, p_fan)) > m_FHTminusT);
    G4ThreeVector A(p), C, diff;
    unsigned int niter = 0;
    do {
        C = A + B;
        C *= 0.5;
        if(fabs(GetCalculator()->DistanceToTheNeutralFibre(C, p_fan)) < m_FHTplusT){
            A = C;
        } else {
            B = C;
        }
        niter ++;
        diff = A - B;
    } while(diff.mag2() > s_IterationPrecision2 && niter < s_IterationsLimit);
    assert(fabs(GetCalculator()->DistanceToTheNeutralFibre(B, p_fan)) > m_FHTplusT);
    assert(niter < s_IterationsLimit);
    diff = p - B;
    LWSDBG(7, std::cout << " -> " << B << " " << diff.mag());
    LWSDBG(6, std::cout << std::endl);
    return diff.mag();
}

outer_solid_init()

void LArWheelSolid::outer_solid_init ( const G4String & bs_name )

private

Definition at line 274 of file LArWheelSolidInit.cxx.

{
    m_IsOuter = true;
    m_FanPhiAmplitude = 0.02; // internal technical constant, should not go in DB
    set_phi_size();
 
    std::array<G4double,3> zPlane{}, rInner{}, rOuter{};
    zPlane[2] = GetCalculator()->GetWheelThickness();
    G4double wheel_thickness = zPlane[2] - zPlane[0];
    zPlane[1] = GetCalculator()->GetWheelInnerRadius(rInner);
    GetCalculator()->GetWheelOuterRadius(rOuter);
    const G4double phi_min =
        m_PhiPosition - m_FanPhiAmplitude -
        GetCalculator()->GetFanStepOnPhi() * 2;
 
    m_Zmid = zPlane[1];
    m_Ymid = (rInner[0] + rOuter[2]) * 0.5;
 
    bool hasFrontSections = false;
    bool hasBackSections = false;
    if(m_Type == OuterAbsorberFrontCone
    || m_Type == OuterElectrodFrontCone
    || m_Type == OuterGlueFrontCone
    || m_Type == OuterLeadFrontCone
    ){
        m_Zmin = zPlane[0]; m_Zmax = zPlane[1];
        m_Rmin = rInner[0]; m_Rmax = rOuter[1];
        m_BoundingShape = new G4ShiftedCone(
            bs_name + "FrontCone", zPlane[0], zPlane[1],
            rInner[0], rOuter[0], rInner[1], rOuter[1]
        );
        hasFrontSections = true;
    } else if(m_Type == OuterAbsorberBackCone
    || m_Type == OuterElectrodBackCone
    || m_Type == OuterGlueBackCone
    || m_Type == OuterLeadBackCone
    ){
        m_Zmin = zPlane[1]; m_Zmax = zPlane[2];
        m_Rmin = rInner[1]; m_Rmax = rOuter[2];
        m_BoundingShape = new G4ShiftedCone(
            bs_name + "BackCone", zPlane[1], zPlane[2],
            rInner[1], rOuter[1], rInner[2], rOuter[2]
        );
        hasBackSections = true;
    } else {
        m_Zmin = zPlane[0]; m_Zmax = zPlane[2];
        m_Rmin = rInner[0]; m_Rmax = rOuter[2];
        m_BoundingShape = new G4Polycone(
            bs_name + "Polycone", m_MinPhi, m_MaxPhi - m_MinPhi,
            3, zPlane.data(), rInner.data(), rOuter.data()
        );
        hasFrontSections = true;
        hasBackSections = true;
    }
 
    m_Ymin = m_Rmin * 0.9;
 
#ifdef LARWHEELSOLID_USE_FANBOUND
  const G4double phi_size = (m_FanPhiAmplitude + GetCalculator()->GetFanStepOnPhi() * 2) * 2;
  FanBound = new G4Polycone(bs_name + "ofFan", phi_min, phi_size,
                            3, zPlane, rInner, rOuter);
#endif
#ifndef PORTABLE_LAR_SHAPE
    ATH_MSG_INFO(m_BoundingShape->GetName() + " is the m_BoundingShape");
#endif
    const G4double half_wave_length = GetCalculator()->GetHalfWaveLength();
    const G4double sss = GetCalculator()->GetStraightStartSection();
 
    if(hasFrontSections){
        m_Zsect.push_back(0.);
        m_Zsect.push_back(sss + half_wave_length * 0.25);
    } else {
        m_Zsect.push_back(m_Zmid);
    }
    const G4int num_fs = GetCalculator()->GetNumberOfHalfWaves() + 1;
 
    for(G4int i = 2; i < num_fs; i ++){
        const G4double zi = half_wave_length * (i - 1) + sss;
#if LARWHEELSOLID_ZSECT_MULT > 1
        for(G4int j = LARWHEELSOLID_ZSECT_MULT - 1; j > 0; -- j){
            G4double zj = zi - half_wave_length * j / LARWHEELSOLID_ZSECT_MULT;
            if(hasFrontSections && hasBackSections
            && m_Zsect.back() < m_Zmid && zj >= m_Zmid){
                m_Zsect.push_back(m_Zmid);
            }
            if((zj < m_Zmid && hasFrontSections)
            || (zj > m_Zmid && hasBackSections)){
                m_Zsect.push_back(zj);
            }
        }
#endif
        if(hasFrontSections && hasBackSections
        && m_Zsect.back() < m_Zmid && zi >= m_Zmid){
            m_Zsect.push_back(m_Zmid);
        }
        if((zi < m_Zmid && hasFrontSections)
        || (zi > m_Zmid && hasBackSections)){
            m_Zsect.push_back(zi);
        }
    }
    if(hasBackSections){
        m_Zsect.push_back(wheel_thickness - sss - half_wave_length * 0.25);
        m_Zsect.push_back(wheel_thickness);
    } else {
        m_Zsect.push_back(m_Zmid);
    }
 
    m_fs = new LArFanSections(
        rInner[0], rInner[1], rOuter[0], rOuter[1],
        m_Rmax*cos(phi_min), m_Zsect.front(), m_Zmid
    );
}

search_for_most_remoted_point()

G4bool LArWheelSolid::search_for_most_remoted_point ( const G4ThreeVector & a, const G4ThreeVector & b, G4ThreeVector & C, const int p_fan ) const

private

Definition at line 259 of file LArWheelSolidDisToOut.cxx.

{
    LWSDBG(6, std::cout << "sfmrp " << a << " " << b << std::endl);
    G4ThreeVector diff(b - a);
 
    if(diff.mag2() <= s_IterationPrecision2) return false;
    G4ThreeVector A(a), B(b), l(diff.unit() * s_IterationPrecision);
  // find the most remoted point on the line AB
  // and check if it is outside vertical fan
  // small vector along the segment AB
    G4double d1, d2;
    unsigned int niter = 0;
  // searching for maximum of (GetCalculator()->DistanceToTheNeutralFibre)^2 along AB
    do {
        C = A + B;
        C *= 0.5;
        d1 = GetCalculator()->DistanceToTheNeutralFibre(C, p_fan);
        if(fabs(d1) > m_FHTplusT){
      // here out_iteration_process gives the answer
            LWSDBG(7, std::cout << "sfmrp -> " << C << " " << fabs(d1)
                                << "  " << (C - a).unit() << " "
                                << (C - a).mag() << std::endl);
            return true;
        }
    // sign of derivative
    //d1 = GetCalculator()->DistanceToTheNeutralFibre(C + l);
        d2 = GetCalculator()->DistanceToTheNeutralFibre(C - l, p_fan);
        if(d1 * d1 - d2 * d2 > 0.) A = C;
        else B = C;
        niter ++;
        diff = A - B;
    } while(diff.mag2() > s_IterationPrecision2 && niter < s_IterationsLimit);
  // the line entirely lies inside fan
    assert(niter < s_IterationsLimit);
    return false;
}

search_for_nearest_point()

G4double LArWheelSolid::search_for_nearest_point ( const G4ThreeVector & p_in, const G4double dist_p_in, const G4ThreeVector & p_out, int p_fan ) const

private

Definition at line 162 of file LArWheelSolidDisToIn.cxx.

{
    LWSDBG(6, std::cout << "sfnp " << MSG_VECTOR(p_in) << " "
                        << MSG_VECTOR(p_out) << std::endl);
 
    G4ThreeVector A, B, C, l, diff;
    A = p_in;
    B = p_out;
    diff = B - A;
    l = diff.unit() * s_IterationPrecision;
  // this is to correctly take the sign of the distance into account
    G4double sign = dist_p_in < 0.? -1. : 1.;
    G4double d_prime;
    G4double dist_c;
    unsigned long niter = 0;
    do {
        C = A + B;
        C *= 0.5;
        dist_c = GetCalculator()->DistanceToTheNeutralFibre(C, p_fan);
        if(dist_c * sign <= 0.){ // we are in coditions for in_iteration_process
            LWSDBG(7, std::cout << "sfnp0 " << dist_c << std::endl);
            return in_iteration_process(p_in, dist_p_in, C, p_fan);
        }
    // calculate sign of derivative of distance to the neutral fibre
    // hope this substitution is acceptable
        diff = C - l;
        d_prime = (dist_c - GetCalculator()->DistanceToTheNeutralFibre(diff, p_fan)) * sign;
        if(d_prime < 0.) A = C;
        else B = C;
        niter ++;
        diff = A - B;
    } while(diff.mag2() > s_IterationPrecision2 && niter < s_IterationsLimit);
    assert(niter < s_IterationsLimit);
    if(fabs(dist_c) < m_FHTminusT){
        LWSDBG(7, std::cout << "sfnp1 " << dist_c << std::endl);
        return in_iteration_process(p_in, dist_p_in, C, p_fan);
    }
  //  let's check at p_in and p_out
    if(dist_p_in * sign < dist_c * sign){
        C = p_in;
        dist_c = dist_p_in;
    }
    G4double dist_p_out = GetCalculator()->DistanceToTheNeutralFibre(p_out, p_fan);
    if(dist_p_out *sign < dist_c * sign) C = p_out;
    if(fabs(dist_p_out) < m_FHTminusT){
        LWSDBG(7, std::cout << "sfnp2 " << dist_p_out << std::endl);
        return in_iteration_process(p_in, dist_p_in, C, p_fan);
    }
    return kInfinity;
}

select_section()

G4int LArWheelSolid::select_section ( const G4double & Z ) const

private

Definition at line 164 of file LArWheelSolid.cxx.

{
    for(G4int i = m_Zsect_start_search; i > 0; -- i){
        if(Z > m_Zsect[i]) return i;
    }
    return 0;
}

set_failover_point()

void LArWheelSolid::set_failover_point ( G4ThreeVector & p, const char * m = 0 ) const

private

Definition at line 39 of file LArWheelSolidTests.cxx.

  {
    p[0] = 0.; p[1] = m_Rmin; p[2] = m_Zmin;
    int p_fan = 0;
    GetCalculator()->DistanceToTheNearestFan(p, p_fan);
 
#ifdef LOCAL_DEBUG
    if(m) std::cout << m << std::endl;
#endif
  }

set_phi_size()

void LArWheelSolid::set_phi_size ( void )

private

Definition at line 389 of file LArWheelSolidInit.cxx.

{
    if(GetCalculator()->GetisModule()){
        m_MinPhi = m_PhiPosition - CLHEP::pi * (1. / 8.) - m_FanPhiAmplitude;
        m_MaxPhi = m_PhiPosition + CLHEP::pi * (1. / 8.) + m_FanPhiAmplitude;
    } else {
        m_MinPhi = 0.;
        m_MaxPhi = CLHEP::twopi;
    }
}

setLevel()

void AthMessaging::setLevel ( MSG::Level lvl )

inherited

Change the current logging level.

Use this rather than msg().setLevel() for proper operation with MT.

Definition at line 28 of file AthMessaging.cxx.

{
  m_lvl = lvl;
}

StreamInfo()

virtual std::ostream & LArWheelSolid::StreamInfo ( std::ostream & os ) const

inlinevirtual

Definition at line 127 of file LArWheelSolid.h.

{ return os; }

SurfaceNormal()

G4ThreeVector LArWheelSolid::SurfaceNormal

(

const G4ThreeVector &

inputP

)

const

Definition at line 39 of file LArWheelSolid.cxx.

{
    LWSDBG(1, std::cout << TypeStr() << " SurfaceNormal" << MSG_VECTOR(inputP) << std::endl);
    EInside inside_BS = m_BoundingShape->Inside(inputP);
    if(inside_BS != kInside){
        LWSDBG(2, std::cout << "not inside BS" << std::endl);
        return m_BoundingShape->SurfaceNormal(inputP);
    }
    G4ThreeVector p( inputP );
    int p_fan = 0;
    GetCalculator()->DistanceToTheNearestFan(p, p_fan);
    G4ThreeVector d = GetCalculator()->NearestPointOnNeutralFibre(p, p_fan);
    d.rotateZ(inputP.phi() -  p.phi()); // rotate back to initial position
    LWSDBG(4, std::cout << "npnf" << MSG_VECTOR(d) << std::endl);
    p = inputP - d;
    LWSDBG(2, std::cout << "sn " << MSG_VECTOR(p.unit()) << std::endl);
    return(p.unit());
}

test()

void LArWheelSolid::test ( void )

private

Definition at line 396 of file LArWheelSolidTests.cxx.

{
  boost::io::ios_all_saver ias(std::cout);
  const char *on = getenv("LARWHEELSOLID_TEST");
  if(on == 0) return;
  std::string test_mode = on;
 
  std::cout << "============| LArWheelSolid test() routine |=============="
            << std::endl;
  std::cout << "Solid of type " << LArWheelSolidTypeString(m_Type)
            << std::endl;
  std::cout.precision(6);
  std::cout << std::fixed;
  const char *prec = getenv("LARWHEELSOLID_TEST_INTPRECISION");
  if(prec) IntPrecision = atof(prec);
  std::cout << "Int. precision " << IntPrecision << std::endl;
  std::cout << "test mode " << test_mode << std::endl;
 
  std::cout << "m_Rmin = " << m_Rmin << " m_Rmax = " << m_Rmax << std::endl
            << "m_Zmin = " << m_Zmin << " m_Zmax = " << m_Zmax << std::endl;
 
  //    TFile *F = new TFile("LArWheelSolid_test.root", "RECREATE");
  TFile *F = 0;
  TNtupleD *T = 0;
  if(test_mode.find("root") != std::string::npos){
    F = new TFile("LArWheelSolid_test.root", "UPDATE");
    T = new TNtupleD(GetName(), GetName(), "x:y:z");
  }
  int N = 1000000;
  const int Nmax(1000000000);
  char *NN = getenv("LARWHEELSOLID_TEST_NPOINTS");
 
  if(NN) {
    char *endptr;
    N = strtol(NN, &endptr, 0);
    if (endptr[0] != '\0') {
      throw std::invalid_argument("Could not convert string to int: " + std::string(NN));
    }
  }
  if (Nmax<N) {
    std::cout << "Number of points from LARWHEELSOLID_TEST_NPOINTS environment variable ("<<N<<") is too large. Using " << Nmax << " instead." << std::endl;
    N=Nmax;
  }
  if (N<0) {
    std::cout << "Number of points from LARWHEELSOLID_TEST_NPOINTS environment variable ("<<N<<") is negative!!. Using 0 instead." << std::endl;
    N=0;
  }
  if(test_mode.find("points") == std::string::npos){
    N = 0;
  } else {
    std::cout << N << " points" << std::endl;
  }
  for(int i = 0; i < N; ++ i){
    G4ThreeVector p = GetPointOnSurface();
#ifdef LOCAL_DEBUG
    EInside ii = Inside(p);
    if(ii != kSurface){
      std::cout << i << " "
                << (ii == kInside? "inside": "outside")
                << std::endl;
    }
#endif
    if(T) T->Fill(p[0], p[1], p[2]);
  }
  if(F){
    T->Write();
    F->Write();
    F->Close();
    delete F;
  }
 
  if(test_mode.find("volume") != std::string::npos){
    double cv = GetCubicVolume();
    std::cout << "GetCubicVolume: " << cv/CLHEP::mm3 << " mm^3" << std::endl;
  }
 
  if(test_mode.find("area") != std::string::npos){
    double sa = GetSurfaceArea();
    std::cout << "GetSurfaceArea: " << sa/CLHEP::mm2 << " mm^2" << std::endl;
  }
 
  std::cout << "======= end of ArWheelSolid test() routine ============="
            << std::endl;
 
  if(m_Type == OuterAbsorberWheel) {
    if(test_mode.find("once") != std::string::npos) exit(0); }
 
  ias.restore();
}

Friends And Related Symbol Documentation

LArWheelSolid_fcn_area

double LArWheelSolid_fcn_area ( double * x, double * p )

friend

Definition at line 530 of file LArWheelSolidTests.cxx.

{
  const double &z = x[0];
  const double &r = p[0];
  const double &index = p[1];
 
  G4ThreeVector a(0., r, z);
  double b = solid[index]->GetCalculator()->AmplitudeOfSurface(a, -1, 121) // sagging ignored, effect should be negligible, use arbitrary fan number
    - solid[index]->GetCalculator()->AmplitudeOfSurface(a, 1, 121);
  return b;
}

LArWheelSolid_fcn_area_on_pc

double LArWheelSolid_fcn_area_on_pc ( double * x, double * p )

friend

Definition at line 550 of file LArWheelSolidTests.cxx.

{
  const double &z = x[0];
  const double &index = p[0];
 
  G4double rmin, rmax;
  get_r(solid[index]->m_BoundingShape, z, rmin, rmax);
 
  double result = 0.;
  G4ThreeVector a(0., rmin, z);
  result += solid[index]->GetCalculator()->AmplitudeOfSurface(a, -1, 232) // sagging ignored, effect should be negligible, use arbitrary fan number
    - solid[index]->GetCalculator()->AmplitudeOfSurface(a, 1, 232);
  a[1] = rmax;
  result += solid[index]->GetCalculator()->AmplitudeOfSurface(a, -1, 343)
    - solid[index]->GetCalculator()->AmplitudeOfSurface(a, 1, 343);
 
  return result;
}

LArWheelSolid_fcn_side_area

double LArWheelSolid_fcn_side_area ( double * x, double * p )

friend

Definition at line 602 of file LArWheelSolidTests.cxx.

{
  const double &r = x[0];
  const double &index = p[0];
 
  return solid[index]->get_length_at_r(r);
}

LArWheelSolid_fcn_vol

double LArWheelSolid_fcn_vol ( double * x, double * p )

friend

Definition at line 542 of file LArWheelSolidTests.cxx.

{
  const double &r = x[0];
  const double &index = p[0];
 
  return solid[index]->get_area_at_r(r);
}

LArWheelSolid_get_dl

double LArWheelSolid_get_dl ( double * x, double * par, G4int side )

friend

Definition at line 570 of file LArWheelSolidTests.cxx.

{
  const double &z = x[0];
  const double &r = par[0];
  const double &index = par[1];
 
  const double h = 0.001;
 
  //check what happens if z+h > m_Zmax etc
  G4ThreeVector p(0., r, z + h);
  G4double D1 = solid[index]->GetCalculator()->AmplitudeOfSurface(p, side, 5665);  // sagging ignored, effect should be negligible, use arbitrary fan number
  p[2] = z - h;
  D1 -= solid[index]->GetCalculator()->AmplitudeOfSurface(p, side, 5665);
  D1 /= 2 * h;
 
  p[2] = z + h / 2;
  G4double D2 = solid[index]->GetCalculator()->AmplitudeOfSurface(p, side, 5665);
  p[2] = z - h / 2;
  D2 -= solid[index]->GetCalculator()->AmplitudeOfSurface(p, side, 5665);
  D2 /= h;
 
  G4double D = (D2 * 4 - D1) / 3.;
  G4double dl = sqrt(1 + D * D);
 
  return dl;
}

Member Data Documentation

ATLAS_THREAD_SAFE

std::atomic_flag m_initialized AthMessaging::ATLAS_THREAD_SAFE = ATOMIC_FLAG_INIT

mutableprivateinherited

Messaging initialized (initMessaging)

Definition at line 141 of file AthMessaging.h.

m_BoundingShape

G4VSolid* LArWheelSolid::m_BoundingShape

private

Definition at line 154 of file LArWheelSolid.h.

m_Calculator

LArWheelCalculator* LArWheelSolid::m_Calculator

private

Definition at line 148 of file LArWheelSolid.h.

m_f_area

TF1* LArWheelSolid::m_f_area

protected

Definition at line 227 of file LArWheelSolid.h.

m_f_area_on_pc

TF1 * LArWheelSolid::m_f_area_on_pc

protected

Definition at line 227 of file LArWheelSolid.h.

m_f_length

TF1 * LArWheelSolid::m_f_length

protected

Definition at line 227 of file LArWheelSolid.h.

m_f_side_area

TF1 * LArWheelSolid::m_f_side_area

protected

Definition at line 227 of file LArWheelSolid.h.

m_f_vol

TF1 * LArWheelSolid::m_f_vol

protected

Definition at line 227 of file LArWheelSolid.h.

m_FanHalfThickness

G4double LArWheelSolid::m_FanHalfThickness

private

Definition at line 149 of file LArWheelSolid.h.

m_FanPhiAmplitude

G4double LArWheelSolid::m_FanPhiAmplitude

private

Definition at line 150 of file LArWheelSolid.h.

m_FHTminusT

G4double LArWheelSolid::m_FHTminusT

private

Definition at line 149 of file LArWheelSolid.h.

m_FHTplusT

G4double LArWheelSolid::m_FHTplusT

private

Definition at line 149 of file LArWheelSolid.h.

m_fs

LArFanSections* LArWheelSolid::m_fs

private

Definition at line 162 of file LArWheelSolid.h.

m_imsg

std::atomic<IMessageSvc*> AthMessaging::m_imsg { nullptr }

mutableprivateinherited

MessageSvc pointer.

Definition at line 135 of file AthMessaging.h.

{ nullptr };

m_IsOuter

G4bool LArWheelSolid::m_IsOuter

private

Definition at line 146 of file LArWheelSolid.h.

m_lvl

std::atomic<MSG::Level> AthMessaging::m_lvl { MSG::NIL }

mutableprivateinherited

Current logging level.

Definition at line 138 of file AthMessaging.h.

{ MSG::NIL };

m_MaxPhi

G4double LArWheelSolid::m_MaxPhi

private

Definition at line 152 of file LArWheelSolid.h.

m_MinPhi

G4double LArWheelSolid::m_MinPhi

private

Definition at line 151 of file LArWheelSolid.h.

m_msg_tls

boost::thread_specific_ptr<MsgStream> AthMessaging::m_msg_tls

mutableprivateinherited

MsgStream instance (a std::cout like with print-out levels)

Definition at line 132 of file AthMessaging.h.

m_nm

std::string AthMessaging::m_nm

privateinherited

Message source name.

Definition at line 129 of file AthMessaging.h.

m_PhiPosition

const G4double LArWheelSolid::m_PhiPosition

private

Definition at line 153 of file LArWheelSolid.h.

m_Rmax

G4double LArWheelSolid::m_Rmax

private

Definition at line 169 of file LArWheelSolid.h.

m_Rmin

G4double LArWheelSolid::m_Rmin

private

Definition at line 169 of file LArWheelSolid.h.

m_test_index

double LArWheelSolid::m_test_index

protected

Definition at line 229 of file LArWheelSolid.h.

m_Type

const LArWheelSolid_t LArWheelSolid::m_Type

private

Definition at line 147 of file LArWheelSolid.h.

m_Ymid

G4double LArWheelSolid::m_Ymid

private

Definition at line 171 of file LArWheelSolid.h.

m_Ymin

G4double LArWheelSolid::m_Ymin

private

Definition at line 167 of file LArWheelSolid.h.

m_Zmax

G4double LArWheelSolid::m_Zmax

private

Definition at line 169 of file LArWheelSolid.h.

m_Zmid

G4double LArWheelSolid::m_Zmid

private

Definition at line 165 of file LArWheelSolid.h.

m_Zmin

G4double LArWheelSolid::m_Zmin

private

Definition at line 169 of file LArWheelSolid.h.

m_Zsect

std::vector<G4double> LArWheelSolid::m_Zsect

private

Definition at line 159 of file LArWheelSolid.h.

m_Zsect_start_search

G4int LArWheelSolid::m_Zsect_start_search

private

Definition at line 160 of file LArWheelSolid.h.

s_AngularTolerance

const G4double LArWheelSolid::s_AngularTolerance = G4GeometryTolerance::GetInstance()->GetAngularTolerance() / 2

staticprivate

Definition at line 141 of file LArWheelSolid.h.

s_IterationPrecision

const G4double LArWheelSolid::s_IterationPrecision = 0.001*CLHEP::mm

staticprivate

Definition at line 142 of file LArWheelSolid.h.

s_IterationPrecision2

const G4double LArWheelSolid::s_IterationPrecision2 = s_IterationPrecision * s_IterationPrecision

staticprivate

Definition at line 143 of file LArWheelSolid.h.

s_IterationsLimit

const unsigned int LArWheelSolid::s_IterationsLimit = 50

staticprivate

Definition at line 144 of file LArWheelSolid.h.

s_Tolerance

const G4double LArWheelSolid::s_Tolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance() / 2

staticprivate

Definition at line 140 of file LArWheelSolid.h.

---

The documentation for this class was generated from the following files:

• LArWheelSolid.h
• LArFanSection.cxx
• LArWheelSolid.cxx
• LArWheelSolidDisToIn.cxx
• LArWheelSolidDisToOut.cxx
• LArWheelSolidInit.cxx
• LArWheelSolidTests.cxx