|
ATLAS Offline Software
|
#include <LArWheelSolid.h>
|
| 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. More...
|
|
MsgStream & | msg () const |
| The standard message stream. More...
|
|
MsgStream & | msg (const MSG::Level lvl) const |
| The standard message stream. More...
|
|
void | setLevel (MSG::Level lvl) |
| Change the current logging level. More...
|
|
|
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. More...
|
|
Definition at line 90 of file LArWheelSolid.h.
◆ FanBoundExit_t
Enumerator |
---|
NoCross | |
ExitAtInner | |
ExitAtOuter | |
ExitAtFront | |
ExitAtBack | |
ExitAtSide | |
Definition at line 189 of file LArWheelSolid.h.
◆ LArWheelSolid()
Definition at line 32 of file LArWheelSolidInit.cxx.
38 #ifndef PORTABLE_LAR_SHAPE
43 #ifndef PORTABLE_LAR_SHAPE
44 #ifdef LARWHEELSOLID_USE_FANBOUND
47 ATH_MSG_INFO (
"compiled with private find_exit_point" );
126 G4Exception(
"LArWheelSolid",
"UnknownSolidType", FatalException,
127 "Constructor: unknown LArWheelSolid_t");
132 const G4String bs_name =
name +
"-Bounding";
133 #ifdef DEBUG_LARWHEELSOLID
134 const char *venv =
getenv(
"LARWHEELSOLID_VERBOSE");
136 std::cout <<
"The LArWheelSolid build " << __DATE__ <<
" " << __TIME__
138 std::cout <<
"LArWheelSolid verbosity level is " <<
Verbose << std::endl;
175 G4Exception(
"LArWheelSolid",
"UnknownSolidType", FatalException,
176 "Constructor: unknown LArWheelSolid_t");
180 #ifndef PORTABLE_LAR_SHAPE
186 #ifdef DEBUG_LARWHEELSOLID
188 std::cout <<
"Limits: (" <<
m_Zsect.size() <<
")" << std::endl;
190 std::cout <<
i <<
" " <<
m_Zsect[
i] << std::endl;
193 #ifndef PORTABLE_LAR_SHAPE
◆ ~LArWheelSolid()
LArWheelSolid::~LArWheelSolid |
( |
| ) |
|
|
virtual |
◆ CalculateExtent()
G4bool LArWheelSolid::CalculateExtent |
( |
const EAxis |
a, |
|
|
const G4VoxelLimits & |
vl, |
|
|
const G4AffineTransform & |
t, |
|
|
G4double & |
p, |
|
|
G4double & |
q |
|
) |
| const |
◆ 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.
48 const G4double D =
B*
B -
A*
C;
49 LWSDBG(8, std::cout <<
"check D=" << D <<
" out=" <<
out << std::endl);
50 if(D < 0.)
return false;
54 LWSDBG(8, std::cout <<
"Case A=B=0" << std::endl);
58 LWSDBG(8, std::cout <<
"t1=" <<
t1 <<
" - only one solution" << std::endl);
61 const G4double D1 = sqrt(D);
62 const G4double inv_A = 1.0 /
A;
63 t1 = (-
B + D1) * inv_A;
64 t2 = (-
B - D1) * inv_A;
65 LWSDBG(8, std::cout <<
"t1=" <<
t1 <<
" t2=" <<
t2 << std::endl);
◆ clean_tests()
void LArWheelSolid::clean_tests |
( |
void |
| ) |
|
|
private |
◆ CreatePolyhedron()
G4Polyhedron * LArWheelSolid::CreatePolyhedron |
( |
| ) |
const |
◆ DescribeYourselfTo()
void LArWheelSolid::DescribeYourselfTo |
( |
G4VGraphicsScene & |
scene | ) |
const |
◆ 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.
218 LWSDBG(4, std::cout <<
"dti: " << MSG_VECTOR(
p) <<
" "
219 << MSG_VECTOR(
v) << std::endl);
222 #ifdef LARWHEELSOLID_USE_FANBOUND
223 if(FanBound->Inside(
p) == kOutside) {
224 const G4double
d = FanBound->DistanceToIn(
p,
v);
230 if(
v.x() >= 0.)
return kInfinity;
232 const G4double
y2 =
p.y() +
v.y() *
b;
233 const G4double z2 =
p.z() +
v.z() *
b;
237 if(
v.x() <= 0.)
return kInfinity;
239 const G4double
y2 =
p.y() +
v.y() *
b;
240 const G4double z2 =
p.z() +
v.z() *
b;
248 LWSDBG(5, std::cout <<
"dti corrected: " << MSG_VECTOR(
p) << std::endl);
252 LWSDBG(5, std::cout <<
"hit fan dist_p=" << dist_p <<
", m_FHTminusT=" <<
m_FHTminusT << std::endl);
256 #ifdef CHECK_DIRTONORM_ANGLE_ON_SURFACE
258 LWSDBG(5, std::cout <<
"outside fan dist_p=" << dist_p <<
", m_FHTplusT=" <<
m_FHTplusT << std::endl);
260 LWSDBG(5, std::cout <<
"on fan surface dist_p=" << dist_p <<
", m_FHTplusT=" <<
m_FHTplusT <<
", m_FHTminusT=" <<
m_FHTminusT << std::endl);
264 if ( (
p-
d).cosTheta(
v) < -AngularTolerance ) {
273 #ifdef LARWHEELSOLID_USE_FANBOUND
274 q =
p +
v * FanBound->DistanceToOut(
p,
v);
284 <<
" " <<
step << std::endl);
290 const G4double
x1 =
p.x() +
v.x() *
d1,
y1 =
p.y() +
v.y() *
d1;
293 LWSDBG(5, std::cout <<
i <<
">" <<
p <<
" " << dist_p <<
" "
294 <<
p1 <<
" " << dist_p1 << std::endl);
295 G4double
dd = kInfinity;
296 if(dist_p * dist_p1 < 0.){
300 LWSDBG(6, std::cout <<
i <<
"> dd=" <<
dd <<
", d2=" <<
d2 <<
", distance=" <<
distance << std::endl);
303 }
else if(
dd < kInfinity){
312 LWSDBG(5, std::cout <<
"dti exit point: " << MSG_VECTOR(
q) <<
" "
313 << dist_q << std::endl);
314 G4double
dd = kInfinity;
315 if(dist_p * dist_q < 0.){
321 }
else if(
dd < kInfinity){
◆ DistanceToIn() [1/2]
G4double LArWheelSolid::DistanceToIn |
( |
const G4ThreeVector & |
inputP | ) |
const |
Definition at line 17 of file LArWheelSolidDisToIn.cxx.
19 LWSDBG(1, std::cout << TypeStr() <<
" DisToIn" << MSG_VECTOR(inputP) << std::endl);
24 LWSDBG(2, std::cout <<
"Outside BS" << std::endl);
27 G4ThreeVector
p(inputP);
36 const G4double
d = fabs(
GetCalculator()->DistanceToTheNearestFan(
p, p_fan));
39 LWSDBG(2, std::cout <<
"dti result = " <<
result << std::endl);
42 LWSDBG(2, std::cout <<
"already inside, return 0" << MSG_VECTOR(
p) << std::endl);
◆ DistanceToIn() [2/2]
G4double LArWheelSolid::DistanceToIn |
( |
const G4ThreeVector & |
inputP, |
|
|
const G4ThreeVector & |
inputV |
|
) |
| const |
Definition at line 46 of file LArWheelSolidDisToIn.cxx.
49 LWSDBG(1, std::cout << TypeStr() <<
" DisToIn" << MSG_VECTOR(inputP)
50 << MSG_VECTOR(inputV) << std::endl);
54 G4ThreeVector
p(inputP);
55 if(inside_BS == kOutside) {
58 LWSDBG(2, std::cout <<
"Infinity distance to m_BoundingShape" << MSG_VECTOR(inputP) << MSG_VECTOR(inputV) << std::endl);
63 LWSDBG(2, std::cout <<
"shift" << MSG_VECTOR(inputP) << std::endl);
66 const G4double
phi0 =
p.phi();
70 LWSDBG(2, std::cout <<
"already inside fan" << MSG_VECTOR(
p) << std::endl);
72 if(inside_BS == kSurface) {
73 LWSDBG(2, std::cout <<
"On BS surface" << std::endl);
78 G4ThreeVector
v(inputV);
84 #ifdef DEBUG_LARWHEELSOLID
87 std::cout << MSG_VECTOR(inputP)
88 <<
" " << MSG_VECTOR(inputV) << std::endl;
90 std::cout <<
"dti: " <<
d0 <<
", DTI: " <<
distance << std::endl;
96 std::cout <<
"DTI hit at dist. " <<
distance <<
", point "
97 << MSG_VECTOR(
p2) <<
", "
100 std::cout <<
"got infinity from dti" << std::endl;
103 #ifdef LWS_HARD_TEST_DTI
104 if(test_dti(inputP, inputV,
distance)){
106 std::cout << TypeStr() <<
" DisToIn" << MSG_VECTOR(inputP)
107 << MSG_VECTOR(inputV) << std::endl;
111 std::cout << TypeStr() <<
" DisToIn" << MSG_VECTOR(inputP)
112 << MSG_VECTOR(inputV) <<
" " <<
distance << std::endl;
114 #endif // ifdef LWS_HARD_TEST_DTI
116 #endif // ifdef DEBUG_LARWHEELSOLID
◆ DistanceToOut() [1/2]
G4double LArWheelSolid::DistanceToOut |
( |
const G4ThreeVector & |
inputP | ) |
const |
Definition at line 11 of file LArWheelSolidDisToOut.cxx.
13 LWSDBG(1, std::cout << TypeStr() <<
" DisToOut" << MSG_VECTOR(inputP) << std::endl);
15 LWSDBG(2, std::cout <<
"DistanceToOut(p):"
16 <<
" point " << MSG_VECTOR(inputP)
17 <<
" is not inside of the m_BoundingShape."
21 G4ThreeVector
p( inputP );
25 LWSDBG(2, std::cout <<
"already not inside " << MSG_VECTOR(
p) << std::endl);
29 LWSDBG(2, std::cout <<
"dto " <<
d <<
" " <<
d0 << std::endl);
◆ 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.
40 LWSDBG(1, std::cout << TypeStr() <<
" DisToOut" << MSG_VECTOR(inputP)
41 << MSG_VECTOR(inputV) << std::endl);
44 if(inside_BS == kOutside){
45 LWSDBG(2, std::cout <<
"DistanceToOut(p):"
46 <<
" point " << MSG_VECTOR(inputP)
47 <<
" is outside of m_BoundingShape." << std::endl);
48 if(calcNorm) *validNorm =
false;
53 G4ThreeVector
p(inputP);
55 const G4double adtnf_p = fabs(
GetCalculator()->DistanceToTheNearestFan(
p, p_fan));
57 LWSDBG(2, std::cout <<
"DistanceToOut(p, v): point "
59 <<
" is outside of solid." << std::endl);
60 if(calcNorm) *validNorm =
false;
64 G4ThreeVector
v(inputV);
65 const G4double
phi0 =
p.phi() - inputP.phi();
68 #ifdef CHECK_DIRTONORM_ANGLE_ON_SURFACE
69 if(adtnf_p < FHTminusT) {
70 LWSDBG(5, std::cout <<
"inside fan point " << MSG_VECTOR(inputP) <<
", FHTminusT=" << FHTminusT << std::endl);
72 LWSDBG(5, std::cout <<
"on fan surface adtnf_p=" << adtnf_p <<
", m_FHTplusT=" <<
m_FHTplusT <<
", FHTminusT=" << FHTminusT << std::endl);
76 if ( (
p-
d).cosTheta(
v) > AngularTolerance ) {
85 LWSDBG(4, std::cout <<
"dto: " << MSG_VECTOR(
p) <<
" "
86 << MSG_VECTOR(
v) << std::endl);
89 #ifdef LARWHEELSOLID_USE_BS_DTO
91 inputP, inputV, calcNorm, validNorm, sn
95 LWSDBG(5, std::cout <<
"dto exit point too low " << MSG_VECTOR(
q) << std::endl);
97 q.setX(
p.x() +
v.x() *
dy);
99 q.setZ(
p.z() +
v.z() *
dy);
103 LWSDBG(5, std::cout <<
"dto exit " <<
exit << std::endl);
105 LWSDBG(5, std::cout <<
"dto exit point " << MSG_VECTOR(
q) << std::endl);
112 LWSDBG(5, std::cout <<
"dto sections " <<
start <<
" " <<
stop <<
" " <<
step << std::endl);
121 LWSDBG(5, std::cout <<
"at " <<
i <<
" dist to zsec = " <<
d1 << std::endl);
122 const G4double
x1 =
p.x() +
v.x() *
d1,
y1 =
p.y() +
v.y() *
d1;
132 #ifndef LARWHEELSOLID_USE_BS_DTO
139 #ifndef LARWHEELSOLID_USE_BS_DTO
149 LWSDBG(5, std::cout <<
"q=" << MSG_VECTOR(
q) <<
" outside fan cur distance=" <<
distance <<
", m_FHTplusT=" <<
m_FHTplusT << std::endl);
151 #ifndef LARWHEELSOLID_USE_BS_DTO
160 #ifndef LARWHEELSOLID_USE_BS_DTO
167 LWSDBG(5, std::cout <<
"At end_dto distance=" <<
distance << std::endl);
168 #ifdef LARWHEELSOLID_USE_BS_DTO
169 if(calcNorm &&
distance < dto_bs) *validNorm =
false;
172 LWSDBG(5, std::cout <<
"dto calc norm " <<
exit << std::endl);
179 sn->set(0., 0., -1.);
184 sn->set(
q.x(),
q.y(), 0.);
196 #ifdef DEBUG_LARWHEELSOLID
198 std::cout <<
"DTO: " <<
distance <<
" ";
200 std::cout << *validNorm <<
" " << MSG_VECTOR((*sn));
202 std::cout <<
"Norm is not valid";
204 std::cout << std::endl;
206 G4ThreeVector
p2 = inputP + inputV *
distance;
208 std::cout <<
"DTO hit at " << MSG_VECTOR(
p2) <<
", "
212 #ifdef LWS_HARD_TEST_DTO
213 if(test_dto(inputP, inputV,
distance)){
215 std::cout << TypeStr() <<
" DisToOut" << MSG_VECTOR(inputP)
216 << MSG_VECTOR(inputV) << std::endl;
◆ find_exit_point()
Definition at line 216 of file LArFanSection.cxx.
220 LWSDBG(6, std::cout <<
"in fep p" << MSG_VECTOR(
p)<<
", v"<< MSG_VECTOR(
v) <<
", q" << MSG_VECTOR(
q) << std::endl);
228 LWSDBG(6, std::cout <<
"after fs_cross_lower q" << MSG_VECTOR(
q) << std::endl);
230 LWSDBG(6, std::cout <<
"after fs_cross_upper q" << MSG_VECTOR(
q) << std::endl);
243 }
else if(
v.z() < 0.){
254 LWSDBG(7, std::cout <<
"fep side " <<
d <<
" " <<
result <<
" q" << MSG_VECTOR(
q) << std::endl);
256 LWSDBG(7, std::cout <<
"fep out_distlower(q)=" << out_distlower <<
" Tolerance=" <<
s_Tolerance << std::endl);
257 if (out_distlower >= 0.0) {
270 LWSDBG(7, std::cout <<
"fep out_distupper(q)=" << out_distupper <<
" Tolerance=" <<
s_Tolerance << std::endl);
271 if (out_distupper <= 0.0) {
277 assert((
q -
p).
mag() < kInfinity);
◆ 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.
107 LWSDBG(7, std::cout <<
"fcl" << std::endl);
109 const G4double
B =
p.x()*
v.x() +
p.y()*
v.y()
116 LWSDBG(8, std::cout <<
"fcl out_dist(p)=" << out_dist <<
" Tolerance=" <<
s_Tolerance << std::endl);
117 const G4bool
out = out_dist >= 0.0;
119 const G4double zz1 =
p.z() +
v.z() *
t1;
121 LWSDBG(8, std::cout <<
"fcl out on Z " << zz1 << std::endl);
124 const G4double xx1 =
p.x() +
v.x() *
t1;
126 LWSDBG(8, std::cout <<
"fcl out on X " << xx1 << std::endl);
137 const G4double xx2 =
p.x() +
v.x() *
t1 * 0.5;
138 const G4double yy2 =
p.y() +
v.y() *
t1 * 0.5;
139 const G4double dev = fabs(sqrt(xx2 *xx2 + yy2*yy2)
143 LWSDBG(8, std::cout <<
"fcl on the surface" << std::endl);
146 LWSDBG(8, std::cout <<
"fcl out = in" << std::endl);
152 q.setY(
p.y() +
v.y() *
t1);
154 LWSDBG(8, std::cout <<
"fcl got " <<
t1 << std::endl);
157 LWSDBG(8, std::cout <<
"fcl no intersection" << std::endl);
◆ 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.
169 LWSDBG(7, std::cout <<
"fcu" << std::endl);
170 G4double
A =
v.perp2();
171 G4double
B =
p.x()*
v.x() +
p.y()*
v.y();
172 G4double
C =
p.perp2();
175 const G4double &Af =
A, &Bf =
B;
179 if(
check_D(b1, Af, Bf, Cf, Cf >= 0.)){
180 const G4double zz1 =
p.z() +
v.z() * b1;
182 const G4double xx1 =
p.x() +
v.x() * b1;
185 q.setY(
p.y() +
v.y() * b1);
190 LWSDBG(8, std::cout <<
"fcu no cyl intersection" << std::endl);
200 const G4double zz1 =
p.z() +
v.z() *
t1;
201 LWSDBG(8, std::cout <<
"fcu z = " << zz1 <<
", lim: (" <<
m_Zsect.front() <<
", " <<
m_Zmid <<
")" << std::endl);
203 const G4double xx1 =
p.x() +
v.x() *
t1;
204 LWSDBG(8, std::cout <<
"fcu x = " << xx1 <<
", lim: (" <<
m_fs->
xmin <<
", " <<
m_fs->
xmax <<
")" << std::endl);
207 q.setY(
p.y() +
v.y() *
t1);
211 LWSDBG(8, std::cout <<
"fcu no cone intersection" << std::endl);
◆ get_area_at_r()
G4double LArWheelSolid::get_area_at_r |
( |
G4double |
r | ) |
const |
|
private |
Definition at line 511 of file LArWheelSolidTests.cxx.
516 G4ThreeVector(0.,
r,
m_Zmin), G4ThreeVector(0., 0., 1.)
519 G4ThreeVector(0.,
r,
m_Zmax), G4ThreeVector(0., 0., -1.)
◆ get_area_on_face()
G4double LArWheelSolid::get_area_on_face |
( |
void |
| ) |
const |
|
private |
◆ get_area_on_polycone()
G4double LArWheelSolid::get_area_on_polycone |
( |
void |
| ) |
const |
|
private |
◆ get_area_on_side()
G4double LArWheelSolid::get_area_on_side |
( |
void |
| ) |
const |
|
private |
◆ get_length_at_r()
G4double LArWheelSolid::get_length_at_r |
( |
G4double |
r | ) |
const |
|
private |
Definition at line 345 of file LArWheelSolidTests.cxx.
350 G4ThreeVector(0.,
r, -
m_Zmin), G4ThreeVector(0., 0., 1.)
354 G4ThreeVector(0.,
r,
m_Zmax), G4ThreeVector(0., 0., -1.)
◆ get_point_on_accordion_surface()
void LArWheelSolid::get_point_on_accordion_surface |
( |
G4ThreeVector & |
p | ) |
const |
|
private |
Definition at line 99 of file LArWheelSolidTests.cxx.
108 p[1] = rnd->Uniform(rmin, rmax);
110 G4double dphi =
p.phi();
116 if(
d < 0.)
side = -1;
117 if(
d >= 0.)
side = 1;
125 G4ThreeVector D =
p; D[2] = 0.;
137 G4ThreeVector
B =
p + D *
d1;
138 G4double dphi =
B.phi();
154 G4ThreeVector D1 = (
p -
B).
unit();
◆ get_point_on_flat_surface()
void LArWheelSolid::get_point_on_flat_surface |
( |
G4ThreeVector & |
p | ) |
const |
|
private |
◆ get_point_on_polycone_surface()
void LArWheelSolid::get_point_on_polycone_surface |
( |
G4ThreeVector & |
p | ) |
const |
|
private |
Definition at line 165 of file LArWheelSolidTests.cxx.
170 const bool inner = rnd->Uniform() > 0.5?
true:
false;
172 p[0] = 0.;
p[1] = inner? rmin: rmax;
p[2] =
z;
174 G4double dphi =
p.phi();
179 const G4double
r =
p[1];
181 G4ThreeVector A1(0.,
r,
z);
193 G4ThreeVector A2(0.,
r,
z);
211 G4ThreeVector
d = (A2 - A1).
unit();
223 if(inner)
step *= -2;
233 G4ThreeVector B1(0.,
r +
step,
z);
244 G4ThreeVector B2(0.,
r +
step,
z);
256 if(B1i == A1i || B2i == A1i){
268 G4ThreeVector
d1 = (A1 - B1).
unit();
270 G4ThreeVector
d2 = (A2 - B2).
unit();
273 G4ThreeVector
X = X1;
274 G4double phi1 = X1.phi(), phi2 = X2.phi();
277 G4double phiX = rnd->Uniform(phi1, phi2);
◆ GetBoundingShape()
const G4VSolid* LArWheelSolid::GetBoundingShape |
( |
void |
| ) |
const |
|
inline |
◆ GetCalculator()
◆ GetCubicVolume()
G4double LArWheelSolid::GetCubicVolume |
( |
void |
| ) |
|
◆ GetEntityType()
G4GeometryType LArWheelSolid::GetEntityType |
( |
| ) |
const |
Definition at line 65 of file LArWheelSolid.cxx.
69 return G4String(
"LArInnerAbsorberWheel");
72 return G4String(
"LArOuterAbsorberWheel");
75 return G4String(
"LArInnerElecrodWheel");
78 return G4String(
"LArOuterElecrodWheel");
81 return G4String(
"LArInnerAbsorberModule");
84 return G4String(
"LArOuterAbsorberModule");
87 return G4String(
"LArInnerElecrodModule");
90 return G4String(
"LArOuterElecrodModule");
93 return G4String(
"LArInnerGlueWheel");
96 return G4String(
"LArOuterGlueWheel");
99 return G4String(
"LArInnerLeadWheel");
102 return G4String(
"LArOuterLeadWheel");
105 return G4String(
"LArInnerAbsorberCone");
108 return G4String(
"LArInnerElectrodCone");
111 return G4String(
"LArInnerGlueCone");
114 return G4String(
"LArInnerLeadCone");
117 return G4String(
"LArOuterAbsorberFrontCone");
120 return G4String(
"LArOuterElectrodFrontCone");
123 return G4String(
"LArOuterGlueFrontCone");
126 return G4String(
"LArOuterLeadFrontCone");
129 return G4String(
"LArOuterAbsorberBackCone");
132 return G4String(
"LArOuterElectrodBackCone");
135 return G4String(
"LArOuterGlueBackCone");
138 return G4String(
"LArOuterLeadBackCone");
141 G4Exception(
"LArWheelSolid",
"UnknownSolidType", FatalException,
"GetEntityType: Unknown LArWheelType.");
◆ GetExtent()
G4VisExtent LArWheelSolid::GetExtent |
( |
| ) |
const |
◆ GetPointOnSurface()
G4ThreeVector LArWheelSolid::GetPointOnSurface |
( |
void |
| ) |
const |
Definition at line 66 of file LArWheelSolidTests.cxx.
69 rnd =
new TRandom3(0);
72 G4double
r = rnd->Uniform();
74 G4ThreeVector
p(0., 0., 0.);
76 G4double level1 = .980;
77 G4double level2 = .993;
78 const char *
v =
getenv(
"LARWHEELSOLID_TEST_MODE_LEVEL1");
80 const char *v1 =
getenv(
"LARWHEELSOLID_TEST_MODE_LEVEL2");
81 if(v1) level2 =
atof(v1);
84 std::cout <<
"LWS::GPOS " <<
r << std::endl;
89 }
else if(
r <= level2){
94 G4Exception(
"LArWheelSolid",
"Rnd generator error", FatalException,
"GetPointOnSurface: Wrong data from rnd generator");
◆ GetSurfaceArea()
G4double LArWheelSolid::GetSurfaceArea |
( |
void |
| ) |
|
Definition at line 366 of file LArWheelSolidTests.cxx.
373 std::cout <<
"get_area_on_polycone: " << a1/
mm2 << std::endl;
379 std::cout <<
"get_area_on_face: " << a2/
mm2 << std::endl;
385 std::cout <<
"get_area_on_side: " << a3/
mm2 << std::endl;
◆ GetType()
◆ 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.
130 LWSDBG(6, std::cout <<
"iip from " <<
p <<
" to " <<
B
131 <<
" dir " << (
B -
p).
unit()
137 unsigned int niter = 0;
144 if(dist_c * dist_p < 0. || fabs(dist_c) <
m_FHTminusT){
155 LWSDBG(7, std::cout <<
"iip result in " << niter <<
" = " <<
B
156 <<
" " <<
diff.mag() << std::endl);
◆ init_tests()
void LArWheelSolid::init_tests |
( |
void |
| ) |
|
|
private |
Definition at line 610 of file LArWheelSolidTests.cxx.
615 #ifdef DEBUG_LARWHEELSOLID
616 std::cout <<
"LArWheelSolid::init_tests: put " <<
this
639 (GetName() +
"_f_length").c_str(), &fcn_length,
◆ 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.
◆ inner_solid_init()
void LArWheelSolid::inner_solid_init |
( |
const G4String & |
bs_name | ) |
|
|
private |
Definition at line 206 of file LArWheelSolidInit.cxx.
212 std::array<G4double,2> zPlane{}, rInner{}, rOuter{};
214 G4double wheel_thickness = zPlane[1] - zPlane[0];
224 m_Ymid = (rInner[0] + rOuter[1]) * 0.5;
232 bs_name +
"Cone", zPlane[0], zPlane[1],
233 rInner[0], rOuter[0], rInner[1], rOuter[1]
238 2, zPlane.data(), rInner.data(), rOuter.data()
241 #ifdef LARWHEELSOLID_USE_FANBOUND
243 FanBound =
new G4Polycone(bs_name +
"ofFan", phi_min, phi_size,
244 2, zPlane, rInner, rOuter);
246 #ifndef PORTABLE_LAR_SHAPE
253 m_Zsect.push_back(sss + half_wave_length * 0.25);
255 for(G4int
i = 2;
i < num_fs;
i ++){
256 const G4double zi = half_wave_length * (
i - 1) + sss;
257 #if LARWHEELSOLID_ZSECT_MULT > 1
268 rInner[0], rInner[1], rOuter[0], rOuter[1],
◆ Inside()
EInside LArWheelSolid::Inside |
( |
const G4ThreeVector & |
inputP | ) |
const |
Definition at line 15 of file LArWheelSolid.cxx.
17 LWSDBG(10, std::cout << std::setprecision(25));
18 LWSDBG(1, std::cout << TypeStr() <<
" Inside " << MSG_VECTOR(inputP) << std::endl);
20 if(inside_BS == kOutside){
21 LWSDBG(2, std::cout <<
"outside BS" << std::endl);
24 G4ThreeVector
p( inputP );
26 const G4double
d = fabs(
GetCalculator()->DistanceToTheNearestFan(
p, p_fan));
28 LWSDBG(2, std::cout <<
"outside fan d=" <<
d <<
", m_FHTplusT=" <<
m_FHTplusT << std::endl);
32 LWSDBG(2, std::cout <<
"inside fan d=" <<
d <<
", m_FHTminusT=" <<
m_FHTminusT <<
", inside_BS=" <<
inside(inside_BS) << std::endl);
35 LWSDBG(2, std::cout <<
"surface" << std::endl);
◆ Inside_accordion()
EInside LArWheelSolid::Inside_accordion |
( |
const G4ThreeVector & |
p | ) |
const |
|
private |
◆ 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 164 of file AthMessaging.h.
◆ 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 179 of file AthMessaging.h.
180 {
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.
◆ out_iteration_process()
G4double LArWheelSolid::out_iteration_process |
( |
const G4ThreeVector & |
p, |
|
|
G4ThreeVector & |
B, |
|
|
const int |
p_fan |
|
) |
| const |
|
private |
◆ outer_solid_init()
void LArWheelSolid::outer_solid_init |
( |
const G4String & |
bs_name | ) |
|
|
private |
Definition at line 274 of file LArWheelSolidInit.cxx.
280 std::array<G4double,3> zPlane{}, rInner{}, rOuter{};
282 G4double wheel_thickness = zPlane[2] - zPlane[0];
285 const G4double phi_min =
290 m_Ymid = (rInner[0] + rOuter[2]) * 0.5;
292 bool hasFrontSections =
false;
293 bool hasBackSections =
false;
302 bs_name +
"FrontCone", zPlane[0], zPlane[1],
303 rInner[0], rOuter[0], rInner[1], rOuter[1]
305 hasFrontSections =
true;
314 bs_name +
"BackCone", zPlane[1], zPlane[2],
315 rInner[1], rOuter[1], rInner[2], rOuter[2]
317 hasBackSections =
true;
323 3, zPlane.data(), rInner.data(), rOuter.data()
325 hasFrontSections =
true;
326 hasBackSections =
true;
331 #ifdef LARWHEELSOLID_USE_FANBOUND
333 FanBound =
new G4Polycone(bs_name +
"ofFan", phi_min, phi_size,
334 3, zPlane, rInner, rOuter);
336 #ifndef PORTABLE_LAR_SHAPE
342 if(hasFrontSections){
344 m_Zsect.push_back(sss + half_wave_length * 0.25);
350 for(G4int
i = 2;
i < num_fs;
i ++){
351 const G4double zi = half_wave_length * (
i - 1) + sss;
352 #if LARWHEELSOLID_ZSECT_MULT > 1
355 if(hasFrontSections && hasBackSections
359 if((zj <
m_Zmid && hasFrontSections)
360 || (zj >
m_Zmid && hasBackSections)){
365 if(hasFrontSections && hasBackSections
369 if((zi <
m_Zmid && hasFrontSections)
370 || (zi >
m_Zmid && hasBackSections)){
375 m_Zsect.push_back(wheel_thickness - sss - half_wave_length * 0.25);
376 m_Zsect.push_back(wheel_thickness);
382 rInner[0], rInner[1], rOuter[0], rOuter[1],
◆ 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.
263 LWSDBG(6, std::cout <<
"sfmrp " <<
a <<
" " <<
b << std::endl);
264 G4ThreeVector
diff(
b -
a);
272 unsigned int niter = 0;
280 LWSDBG(7, std::cout <<
"sfmrp -> " <<
C <<
" " << fabs(
d1)
281 <<
" " << (
C -
a).
unit() <<
" "
282 << (
C -
a).
mag() << std::endl);
◆ 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.
167 LWSDBG(6, std::cout <<
"sfnp " << MSG_VECTOR(p_in) <<
" "
168 << MSG_VECTOR(p_out) << std::endl);
176 G4double
sign = dist_p_in < 0.? -1. : 1.;
179 unsigned long niter = 0;
184 if(dist_c *
sign <= 0.){
185 LWSDBG(7, std::cout <<
"sfnp0 " << dist_c << std::endl);
192 if(d_prime < 0.)
A =
C;
199 LWSDBG(7, std::cout <<
"sfnp1 " << dist_c << std::endl);
203 if(dist_p_in *
sign < dist_c *
sign){
208 if(dist_p_out *
sign < dist_c *
sign)
C = p_out;
210 LWSDBG(7, std::cout <<
"sfnp2 " << dist_p_out << std::endl);
◆ select_section()
G4int LArWheelSolid::select_section |
( |
const G4double & |
Z | ) |
const |
|
private |
◆ set_failover_point()
void LArWheelSolid::set_failover_point |
( |
G4ThreeVector & |
p, |
|
|
const char * |
m = 0 |
|
) |
| const |
|
private |
◆ set_phi_size()
void LArWheelSolid::set_phi_size |
( |
void |
| ) |
|
|
private |
◆ setLevel()
void AthMessaging::setLevel |
( |
MSG::Level |
lvl | ) |
|
|
inherited |
◆ StreamInfo()
virtual std::ostream& LArWheelSolid::StreamInfo |
( |
std::ostream & |
os | ) |
const |
|
inlinevirtual |
◆ SurfaceNormal()
G4ThreeVector LArWheelSolid::SurfaceNormal |
( |
const G4ThreeVector & |
inputP | ) |
const |
Definition at line 39 of file LArWheelSolid.cxx.
41 LWSDBG(1, std::cout << TypeStr() <<
" SurfaceNormal" << MSG_VECTOR(inputP) << std::endl);
43 if(inside_BS != kInside){
44 LWSDBG(2, std::cout <<
"not inside BS" << std::endl);
47 G4ThreeVector
p( inputP );
51 d.rotateZ(inputP.phi() -
p.phi());
52 LWSDBG(4, std::cout <<
"npnf" << MSG_VECTOR(
d) << std::endl);
54 LWSDBG(2, std::cout <<
"sn " << MSG_VECTOR(
p.unit()) << std::endl);
◆ test()
void LArWheelSolid::test |
( |
void |
| ) |
|
|
private |
Definition at line 396 of file LArWheelSolidTests.cxx.
398 boost::io::ios_all_saver ias(std::cout);
399 const char *on =
getenv(
"LARWHEELSOLID_TEST");
401 std::string test_mode = on;
403 std::cout <<
"============| LArWheelSolid test() routine |=============="
407 std::cout.precision(6);
408 std::cout << std::fixed;
409 const char *
prec =
getenv(
"LARWHEELSOLID_TEST_INTPRECISION");
411 std::cout <<
"Int. precision " << IntPrecision << std::endl;
412 std::cout <<
"test mode " << test_mode << std::endl;
414 std::cout <<
"m_Rmin = " <<
m_Rmin <<
" m_Rmax = " <<
m_Rmax << std::endl
415 <<
"m_Zmin = " <<
m_Zmin <<
" m_Zmax = " <<
m_Zmax << std::endl;
420 if(test_mode.find(
"root") != std::string::npos){
421 F =
new TFile(
"LArWheelSolid_test.root",
"UPDATE");
422 T =
new TNtupleD(GetName(), GetName(),
"x:y:z");
425 const int Nmax(1000000000);
426 char *NN =
getenv(
"LARWHEELSOLID_TEST_NPOINTS");
430 N = strtol(NN, &endptr, 0);
431 if (endptr[0] !=
'\0') {
432 throw std::invalid_argument(
"Could not convert string to int: " + std::string(NN));
436 std::cout <<
"Number of points from LARWHEELSOLID_TEST_NPOINTS environment variable ("<<
N<<
") is too large. Using " << Nmax <<
" instead." << std::endl;
440 std::cout <<
"Number of points from LARWHEELSOLID_TEST_NPOINTS environment variable ("<<
N<<
") is negative!!. Using 0 instead." << std::endl;
443 if(test_mode.find(
"points") == std::string::npos){
446 std::cout <<
N <<
" points" << std::endl;
448 for(
int i = 0;
i <
N; ++
i){
453 std::cout <<
i <<
" "
454 << (ii == kInside?
"inside":
"outside")
458 if(
T)
T->Fill(
p[0],
p[1],
p[2]);
467 if(test_mode.find(
"volume") != std::string::npos){
469 std::cout <<
"GetCubicVolume: " << cv/
CLHEP::mm3 <<
" mm^3" << std::endl;
472 if(test_mode.find(
"area") != std::string::npos){
474 std::cout <<
"GetSurfaceArea: " << sa/
CLHEP::mm2 <<
" mm^2" << std::endl;
477 std::cout <<
"======= end of ArWheelSolid test() routine ============="
481 if(test_mode.find(
"once") != std::string::npos)
exit(0); }
◆ LArWheelSolid_fcn_area
double LArWheelSolid_fcn_area |
( |
double * |
x, |
|
|
double * |
p |
|
) |
| |
|
friend |
Definition at line 530 of file LArWheelSolidTests.cxx.
532 const double &
z =
x[0];
533 const double &
r =
p[0];
534 const double &
index =
p[1];
536 G4ThreeVector
a(0.,
r,
z);
537 double b = solid[
index]->GetCalculator()->AmplitudeOfSurface(
a, -1, 121)
538 - solid[
index]->GetCalculator()->AmplitudeOfSurface(
a, 1, 121);
◆ LArWheelSolid_fcn_area_on_pc
double LArWheelSolid_fcn_area_on_pc |
( |
double * |
x, |
|
|
double * |
p |
|
) |
| |
|
friend |
Definition at line 550 of file LArWheelSolidTests.cxx.
552 const double &
z =
x[0];
553 const double &
index =
p[0];
559 G4ThreeVector
a(0., rmin,
z);
560 result += solid[
index]->GetCalculator()->AmplitudeOfSurface(
a, -1, 232)
561 - solid[
index]->GetCalculator()->AmplitudeOfSurface(
a, 1, 232);
563 result += solid[
index]->GetCalculator()->AmplitudeOfSurface(
a, -1, 343)
564 - solid[
index]->GetCalculator()->AmplitudeOfSurface(
a, 1, 343);
◆ LArWheelSolid_fcn_side_area
double LArWheelSolid_fcn_side_area |
( |
double * |
x, |
|
|
double * |
p |
|
) |
| |
|
friend |
◆ LArWheelSolid_fcn_vol
double LArWheelSolid_fcn_vol |
( |
double * |
x, |
|
|
double * |
p |
|
) |
| |
|
friend |
◆ LArWheelSolid_get_dl
double LArWheelSolid_get_dl |
( |
double * |
x, |
|
|
double * |
par, |
|
|
G4int |
side |
|
) |
| |
|
friend |
Definition at line 570 of file LArWheelSolidTests.cxx.
572 const double &
z =
x[0];
573 const double &
r =
par[0];
576 const double h = 0.001;
579 G4ThreeVector
p(0.,
r,
z +
h);
580 G4double D1 = solid[
index]->GetCalculator()->AmplitudeOfSurface(
p,
side, 5665);
582 D1 -= solid[
index]->GetCalculator()->AmplitudeOfSurface(
p,
side, 5665);
586 G4double D2 = solid[
index]->GetCalculator()->AmplitudeOfSurface(
p,
side, 5665);
588 D2 -= solid[
index]->GetCalculator()->AmplitudeOfSurface(
p,
side, 5665);
591 G4double D = (D2 * 4 - D1) / 3.;
592 G4double
dl = sqrt(1 + D * D);
◆ ATLAS_THREAD_SAFE
std::atomic_flag m_initialized AthMessaging::ATLAS_THREAD_SAFE = ATOMIC_FLAG_INIT |
|
mutableprivateinherited |
◆ m_BoundingShape
G4VSolid* LArWheelSolid::m_BoundingShape |
|
private |
◆ m_Calculator
◆ m_f_area
TF1* LArWheelSolid::m_f_area |
|
protected |
◆ m_f_area_on_pc
TF1 * LArWheelSolid::m_f_area_on_pc |
|
protected |
◆ m_f_length
TF1 * LArWheelSolid::m_f_length |
|
protected |
◆ m_f_side_area
TF1 * LArWheelSolid::m_f_side_area |
|
protected |
◆ m_f_vol
TF1 * LArWheelSolid::m_f_vol |
|
protected |
◆ m_FanHalfThickness
G4double LArWheelSolid::m_FanHalfThickness |
|
private |
◆ m_FanPhiAmplitude
G4double LArWheelSolid::m_FanPhiAmplitude |
|
private |
◆ m_FHTminusT
G4double LArWheelSolid::m_FHTminusT |
|
private |
◆ m_FHTplusT
G4double LArWheelSolid::m_FHTplusT |
|
private |
◆ m_fs
◆ m_imsg
std::atomic<IMessageSvc*> AthMessaging::m_imsg { nullptr } |
|
mutableprivateinherited |
◆ m_IsOuter
G4bool LArWheelSolid::m_IsOuter |
|
private |
◆ m_lvl
std::atomic<MSG::Level> AthMessaging::m_lvl { MSG::NIL } |
|
mutableprivateinherited |
◆ m_MaxPhi
G4double LArWheelSolid::m_MaxPhi |
|
private |
◆ m_MinPhi
G4double LArWheelSolid::m_MinPhi |
|
private |
◆ 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 |
◆ m_PhiPosition
const G4double LArWheelSolid::m_PhiPosition |
|
private |
◆ m_Rmax
G4double LArWheelSolid::m_Rmax |
|
private |
◆ m_Rmin
G4double LArWheelSolid::m_Rmin |
|
private |
◆ m_test_index
double LArWheelSolid::m_test_index |
|
protected |
◆ m_Type
◆ m_Ymid
G4double LArWheelSolid::m_Ymid |
|
private |
◆ m_Ymin
G4double LArWheelSolid::m_Ymin |
|
private |
◆ m_Zmax
G4double LArWheelSolid::m_Zmax |
|
private |
◆ m_Zmid
G4double LArWheelSolid::m_Zmid |
|
private |
◆ m_Zmin
G4double LArWheelSolid::m_Zmin |
|
private |
◆ m_Zsect
std::vector<G4double> LArWheelSolid::m_Zsect |
|
private |
◆ m_Zsect_start_search
G4int LArWheelSolid::m_Zsect_start_search |
|
private |
◆ s_AngularTolerance
const G4double LArWheelSolid::s_AngularTolerance = G4GeometryTolerance::GetInstance()->GetAngularTolerance() / 2 |
|
staticprivate |
◆ s_IterationPrecision
const G4double LArWheelSolid::s_IterationPrecision = 0.001*CLHEP::mm |
|
staticprivate |
◆ s_IterationPrecision2
◆ s_IterationsLimit
const unsigned int LArWheelSolid::s_IterationsLimit = 50 |
|
staticprivate |
◆ s_Tolerance
const G4double LArWheelSolid::s_Tolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance() / 2 |
|
staticprivate |
The documentation for this class was generated from the following files:
std::atomic< MSG::Level > m_lvl
Current logging level.
G4double get_area_on_polycone(void) const
void get_point_on_flat_surface(G4ThreeVector &) const
double GetStraightStartSection() const
G4double search_for_nearest_point(const G4ThreeVector &, const G4double, const G4ThreeVector &, int) const
friend double LArWheelSolid_fcn_side_area(double *, double *)
friend double LArWheelSolid_fcn_area_on_pc(double *, double *)
static const G4double s_Tolerance
double GetFanHalfThickness(LArG4::LArWheelCalculator_t) const
G4ThreeVector GetPointOnSurface(void) const
friend double LArWheelSolid_fcn_area(double *, double *)
G4bool search_for_most_remoted_point(const G4ThreeVector &, const G4ThreeVector &, G4ThreeVector &, const int) const
void get_point_on_accordion_surface(G4ThreeVector &) const
FanBoundExit_t find_exit_point(const G4ThreeVector &p, const G4ThreeVector &v, G4ThreeVector &q) const
std::vector< ALFA_RawDataCollection_p1 > t1
int GetNumberOfHalfWaves() const
static const char * LArWheelCalculatorTypeString(LArG4::LArWheelCalculator_t)
G4bool fs_cross_lower(const G4ThreeVector &p, const G4ThreeVector &v, G4ThreeVector &q) const
const G4double m_PhiPosition
void set_failover_point(G4ThreeVector &p, const char *m=0) const
G4double get_area_on_side(void) const
void get_point_on_polycone_surface(G4ThreeVector &) const
G4double GetSurfaceArea(void)
std::atomic< IMessageSvc * > m_imsg
MessageSvc pointer.
IMessageSvc * getMessageSvc(bool quiet=false)
G4bool fs_cross_upper(const G4ThreeVector &p, const G4ThreeVector &v, G4ThreeVector &q) const
const LArWheelCalculator * GetCalculator(void) const
LArWheelCalculator * m_Calculator
AthMessaging()
Default constructor:
friend double LArWheelSolid_fcn_vol(double *, double *)
double GetWheelInnerRadius(std::array< double, 2 > &rInner) const
int GetNumberOfFans() const
double GetWheelThickness() const
G4bool check_D(G4double &b, G4double A, G4double B, G4double C, G4bool) const
std::vector< G4double > m_Zsect
EInside Inside(const G4ThreeVector &) const
MsgStream & msg() const
The standard message stream.
void GetWheelOuterRadius(std::array< double, 2 > &rOuter) const
double atof(std::string_view str)
Converts a string into a double / float.
G4int m_Zsect_start_search
G4int select_section(const G4double &Z) const
G4double m_FanHalfThickness
static const G4double s_IterationPrecision2
double GetHalfWaveLength() const
double AmplitudeOfSurface(const CLHEP::Hep3Vector &P, int side, int fan_number) const
double DistanceToTheNearestFan(CLHEP::Hep3Vector &p, int &out_fan_number) const
Determines the nearest to the input point fan.
std::vector< ALFA_RawDataContainer_p1 > t2
G4double get_area_on_face(void) const
std::string getenv(const std::string &variableName)
get an environment variable
CLHEP::Hep3Vector NearestPointOnNeutralFibre(const CLHEP::Hep3Vector &p, int fan_number) const
G4double m_FanPhiAmplitude
G4double in_iteration_process(const G4ThreeVector &, G4double, G4ThreeVector &, int) const
const PlainObject unit() const
This is a plugin that makes Eigen look like CLHEP & defines some convenience methods.
std::string m_nm
Message source name.
const LArWheelSolid_t m_Type
void inner_solid_init(const G4String &)
double GetFanStepOnPhi() const
int atoi(std::string_view str)
Helper functions to unpack numbers decoded in string into integers and doubles The strings are requir...
void outer_solid_init(const G4String &)
double DistanceToTheNeutralFibre(const CLHEP::Hep3Vector &p, int fan_number) const
Calculates aproximate, probably underestimate, distance to the neutral fibre of the vertical fan.
void initMessaging() const
Initialize our message level and MessageSvc.
static const G4double s_IterationPrecision
boost::thread_specific_ptr< MsgStream > m_msg_tls
MsgStream instance (a std::cout like with print-out levels)
G4double GetCubicVolume(void)
G4VSolid * m_BoundingShape
static const unsigned int s_IterationsLimit
virtual G4double distance_to_in(G4ThreeVector &, const G4ThreeVector &, int) const
float distance(const Amg::Vector3D &p1, const Amg::Vector3D &p2)
calculates the distance between two point in 3D space
Scalar mag() const
mag method
#define LARWHEELSOLID_ZSECT_MULT
G4double out_iteration_process(const G4ThreeVector &, G4ThreeVector &, const int) const
const char * LArWheelSolidTypeString(LArWheelSolid_t type)