@@ -111,9 +111,16 @@ double mass_attenuation_coef( float atomic_number, float energy );
111111double transmition_coefficient_generic ( float atomic_number, float areal_density,
112112 float energy );
113113
114+ // When debugging we will grab a static mutex so we dont get jumbled stdout
115+ #define DEBUG_RAYTRACE_CALCS 0
116+
114117
115118void example_integration ();
116119
120+ /* * Runs some simple test cases for #cylinder_exit_position, and cause assert(0) on error. */
121+ void test_cylinder_exit_position ();
122+
123+
117124int DistributedSrcCalc_integrand ( const int *ndim, const double xx[],
118125 const int *ncomp, double ff[], void *userdata );
119126
@@ -137,6 +144,31 @@ double exit_point_of_sphere_z( const double source_point[3],
137144 double observation_dist,
138145 bool postiveSolution = true );
139146
147+ /* * Starting from a 'source_point' within the volume of the cylinder, and heading towards the 'detector_point' (think center
148+ of the detector face), returns the total attenuation coefficient along the path, including recursing into any sub-tubes, as well as
149+ sets the point where the ray leaves the cylinder.
150+
151+ Note: the cylinder is always oriented along the z-axis, and centered at {0,0,0}.
152+
153+ @param[in] radius The outer radius of the cylinder
154+ @param[in] half_length The half-length of the cylinder
155+ @param[in] source The {x, y, z} source location; must be within cylinder volume, and transformed so that {0,0,0} is center of cylinder.
156+ @param[in] detector The {x, y, z} point on the detector face we care about (so center of detector, unless you are integrating over
157+ the detector face), in the coordinate system where cylinder is centered at {0,0,0}.
158+ @param[out] exit_point The final exit point from the cylinder, where the path will no longer go through the volume.
159+ @returns the attenuation coefficient from the path through the cylinder and its sub-cylinders. E.g.,
160+ \code{.cpp}
161+ double exit_point[3];
162+ const double distance_in_m = cylinder_exit_position( 0.5*m, 100*cm, {0,0.1*m,20*cm}, {0,0,10*m}, exit_point );
163+ const double trans_fraction = exp( -trans_coef ); //trans_fraction will be between 0 and 1.
164+ \endcode
165+ */
166+ double cylinder_exit_position ( const double radius, const double half_length,
167+ const double source[3 ],
168+ const double detector[3 ],
169+ double exit_point[3 ] );
170+
171+
140172// distance(...): returns distance between two points specified in terms of x, y,
141173// and z.
142174double distance ( const double point_a[3 ], const double point_b[3 ] );
@@ -170,7 +202,7 @@ struct DistributedSrcCalc
170202 void eval_cyl_end_on ( const double xx[], const int *ndimptr,
171203 double ff[], const int *ncompptr ) const ;
172204
173- void eval_cyl_side_on ( const double xx[], const int *ndimptr,
205+ void eval_cylinder ( const double xx[], const int *ndimptr,
174206 double ff[], const int *ncompptr ) const ;
175207
176208 void eval_rect ( const double xx[], const int *ndimptr,
0 commit comments