NuRadioMC.SignalProp.analyticraytracing module

NuRadioMC.SignalProp.analyticraytracing.speed_of_light = 0.299792458

Models in the following list will use the speed-optimized algorithm to calculate the attenuation along the path. Can be overwritten by init.

NuRadioMC.SignalProp.analyticraytracing.get_n_steps(x1, x2, dx)[source]

Returns number of segments necessary for width to be approx dx

NuRadioMC.SignalProp.analyticraytracing.get_segments(x1, x2, dx)[source]

Returns equi.dist. segments (np.linspace). Choose number of segments such that width of segments is approx. dx

NuRadioMC.SignalProp.analyticraytracing.get_z_deep(ice_params)[source]

Calculates the z_deep needed for integral along the homogeneous ice to know the path length or the times. We obtain the depth for which the index of refraction is 0.035% away of that of deep ice. This calculation assumes a monotonically increasing index of refraction with negative depth.

class NuRadioMC.SignalProp.analyticraytracing.ray_tracing_2D(medium, attenuation_model='SP1', log_level=30, n_frequencies_integration=25, use_optimized_start_values=False, overwrite_speedup=None, use_cpp=True)[source]

Bases: ray_tracing_base

initialize 2D analytic ray tracing class

Parameters:
medium: NuRadioMC.utilities.medium class

details of the medium

attenuation_model: string

specifies which attenuation model to use (default ‘SP1’)

log_level: logging.loglevel object

controls verbosity (default WARNING)

n_frequencies_integration: int

specifies for how many frequencies the signal attenuation is being calculated

use_optimized_start_value: bool

if True, the initial C_0 paramter (launch angle) is set to the ray that skims the surface (default: False)

overwrite_speedup: bool

The signal attenuation is calculated using a numerical integration along the ray path. This calculation can be computational inefficient depending on the details of the ice model. An optimization is implemented approximating the integral with a discrete sum with the loss of some accuracy, See PR #507. This optimization is used for all ice models listed in speedup_attenuation_models (i.e., “GL3”). With this argument you can explicitly activate or deactivate (True or False) if you want to use the optimization. (Default: None, i.e., use optimization if ice model is listed in speedup_attenuation_models)

use_cpp: bool

if True, use CPP implementation of minimization routines default: True if CPP version is available

Methods

angular_diff(x_refl, z_refl, pulser_pos, ...)

This is a helper function to find a ray that is subject to specular reflection (no transmission) at the ice-water interface, e.g.

apply_propagation_effects(efield, i_solution)

Apply propagation effects to the electric field Note that the 1/r weakening of the electric field is already accounted for in the signal generation

determine_solution_type(x1, x2, C_0)

returns the type of the solution

ds(t, C_0)

helper to calculate line integral

find_solutions(x1, x2[, plot, reflection, ...])

this function finds all ray tracing solutions

get_C0_from_log(logC0)

transforms the fit parameter C_0 so that the likelihood looks better

get_C_0_from_angle(anglaunch, z_pos[, in_air])

Find parameter C0 corresponding to a given launch angle and z-position of a ray.

get_C_1(x1, C_0)

calculates constant C_1 for a given C_0 and start point x1

get_angle(x, x_start, C_0[, reflection, ...])

calculates the angle with respect to the positive z-axis of the ray path at position x

get_angle_from_logC_0(logC_0, z_pos[, ...])

argument angoff is provided so that the function can be used for minimization in get_C_0_from_angle(), in which case angoff is the angle for which the C_0 is sought and zero is returned when it is found.

get_attenuation(iS, frequency[, ...])

calculates the signal attenuation due to attenuation in the medium (ice)

get_config()

Function that returns the configuration currently used by the raytracer

get_delta_y(C_0, x1, x2[, C0range, ...])

calculates the difference in the y position between the analytic ray tracing path specified by C_0 at the position x2

get_gamma(z)

transforms z coordinate into gamma

get_launch_vector(iS)

calculates the launch vector (in 3D) of solution iS

get_number_of_raytracing_solutions()

Function that returns the maximum number of raytracing solutions that can exist between each given pair of start and end points

get_number_of_solutions()

returns the number of solutions

get_output_parameters()

Returns a list with information about parameters to include in the output data structure that are specific to this raytracer

get_path(x1, x2, C_0[, n_points])

for plotting purposes only, returns the ray tracing path between x1 and x2

get_path_length(x1, x2, C_0[, reflection, ...])

calculates the path length of solution iS

get_path_length_analytic(x1, x2, C_0[, ...])

analytic solution to calculate the distance along the path.

get_path_reflections(x1, x2, C_0[, ...])

calculates the ray path in the presence of reflections at the bottom The full path is constructed by multiple calls to the get_path() function to put together the full path

get_path_segments(x1, x2, C_0[, reflection, ...])

Calculates the different segments of the path that makes up the full ray tracing path One segment per bottom reflection.

get_raytracing_output(i_solution)

Write parameters that are specific to this raytracer into the output data.

get_receive_vector(iS)

calculates the receive vector (in 3D) of solution iS

get_reflection_angle(x1, x2, C_0[, ...])

calculates the angle under which the ray reflects off the surface.

get_reflection_point(C_0, C_1)

calculates the point where the signal gets reflected off the bottom of the ice shelf

get_solution_type(iS)

returns the type of the solution

get_surf_skim_angle(x1)

For a given position x1 = [x,z] and depth profile self.n(), find the angle at which a beam must be emitted to "skim the surface", i.e. arrive horizontally (angle = 90 deg) at the surface; This is used to find the refraction zone.

get_surface_pulse(x1, x2[, infirn, angle, ...])

Calculate the time for a ray to travel from x1 to the surface and arriving at the surface with (a) critical angle or (b) Brewster angle, propagating in a straight line along the surface in air (n=1) in the firn at z=0 (n=self.n(0)) and then reaching the receiver by returning into the firn just at the point to reach the receiver at x2, entering the firn at the surface at the same angle it reached the surface from x1..

get_tof_for_straight_line(x1, x2)

Calculate the time of flight for a hypothatical ray travelling from x1 to x2 in a straight line.

get_travel_time(x1, x2, C_0[, reflection, ...])

calculates the travel time of solution iS

get_travel_time_analytic(x1, x2, C_0[, ...])

analytic solution to calculate the time of flight.

get_turning_point(c)

calculate the turning point, i.e. the maximum of the ray tracing path; parameter is c = self.medium.n_ice ** 2 - C_0 ** -2.

get_y(gamma, C_0, C_1)

analytic form of the ray tracing part given an exponential index of refraction profile

get_y_diff(z_raw, C_0[, in_air])

derivative dy(z)/dz

get_y_turn(C_0, x1)

calculates the y-coordinate of the turning point.

get_y_with_z_mirror(z, C_0[, C_1])

analytic form of the ray tracing part given an exponential index of refraction profile

get_z_from_n(n)

get z from given n - equation from get_n solved for z

get_z_mirrored(x1, x2, C_0)

calculates the mirrored x2 position so that y(z) can be used as a continuous function

get_z_unmirrored(z, C_0)

calculates the unmirrored z position

has_solution()

checks if ray tracing solution exists

is_in_refraction_zone(x1, x2[, C0crit, plot])

Find if receiver at x2 is in the refraction zone of emitter at x1.

n(z)

refractive index as a function of depth

obj_delta_y(logC_0, x1, x2[, reflection, ...])

function to find solution for C0, returns distance in y between function and x2 position result is signed! (important to use a root finder)

obj_delta_y_square(logC_0, x1, x2[, ...])

objective function to find solution for C0

plot_result(x1, x2, C_0, ax)

helper function to visualize results

set_config(config)

Function to change the configuration file used by the raytracer

set_start_and_end_point(x1, x2)

Set the start and end points between which raytracing solutions shall be found It is recommended to also reset the solutions from any previous raytracing to avoid confusing them with the current solution

use_optional_function(function_name, *args, ...)

Use optional function which may be different for each ray tracer.

get_angle_from_C_0

get_attenuation_along_path

get_c

get_launch_angle

get_receive_angle

get_results

reset_solutions

n(z)[source]

refractive index as a function of depth

get_gamma(z)[source]

transforms z coordinate into gamma

get_turning_point(c)[source]

calculate the turning point, i.e. the maximum of the ray tracing path; parameter is c = self.medium.n_ice ** 2 - C_0 ** -2

This is either the point of reflection off the ice surface or the point where the saddle point of the ray (transition from upward to downward going)

Technically, the turning point is set to z=0 if the saddle point is above the surface.

Parameters:
c: float

related to C_0 parameter via c = self.medium.n_ice ** 2 - C_0 ** -2

Returns:
typle (gamma, z coordinate of turning point)
get_y_turn(C_0, x1)[source]

calculates the y-coordinate of the turning point. This is either the point of reflection off the ice surface or the point where the saddle point of the ray (transition from upward to downward going)

Parameters:
C_0: float

C_0 parameter of function

x1: typle

(y, z) start position of ray

get_C_1(x1, C_0)[source]

calculates constant C_1 for a given C_0 and start point x1

get_c(C_0)[source]
get_C0_from_log(logC0)[source]

transforms the fit parameter C_0 so that the likelihood looks better

get_y(gamma, C_0, C_1)[source]

analytic form of the ray tracing part given an exponential index of refraction profile

Parameters:
gamma: (float or array)

gamma is a function of the depth z

C_0: (float)

first parameter

C_1: (float)

second parameter

get_y_diff(z_raw, C_0, in_air=False)[source]

derivative dy(z)/dz

get_y_with_z_mirror(z, C_0, C_1=0)[source]

analytic form of the ray tracing part given an exponential index of refraction profile

this function automatically mirrors z values that are above the turning point, so that this function is defined for all z

Parameters:
z: (float or array)

depth z

C_0: (float)

first parameter

C_1: (float)

second parameter

get_z_mirrored(x1, x2, C_0)[source]

calculates the mirrored x2 position so that y(z) can be used as a continuous function

get_z_unmirrored(z, C_0)[source]

calculates the unmirrored z position

ds(t, C_0)[source]

helper to calculate line integral

get_path_length(x1, x2, C_0, reflection=0, reflection_case=1)[source]

calculates the path length of solution iS

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

analytic: bool

If True the analytic solution is used. If False, a numerical integration is used. (default: True)

Returns:
distance: float

distance from x1 to x2 along the ray path

get_path_length_analytic(x1, x2, C_0, reflection=0, reflection_case=1)[source]

analytic solution to calculate the distance along the path. This code is based on the analytic solution found by Ben Hokanson-Fasing and the pyrex implementation.

get_travel_time(x1, x2, C_0, reflection=0, reflection_case=1)[source]

calculates the travel time of solution iS

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

analytic: bool

If True the analytic solution is used. If False, a numerical integration is used. (default: True)

Returns:
time: float

travel time

get_travel_time_analytic(x1, x2, C_0, reflection=0, reflection_case=1)[source]

analytic solution to calculate the time of flight. This code is based on the analytic solution found by Ben Hokanson-Fasing and the pyrex implementation.

get_attenuation_along_path(x1, x2, C_0, frequency, max_detector_freq, reflection=0, reflection_case=1)[source]
get_path_segments(x1, x2, C_0, reflection=0, reflection_case=1)[source]

Calculates the different segments of the path that makes up the full ray tracing path One segment per bottom reflection.

Parameters:
x1: tuple

(y, z) coordinate of start value

x2: tuple

(y, z) coordinate of stop value

C_0: float

C_0 parameter of analytic ray path function

reflection: int (default 0)

the number of bottom reflections to consider

reflection_case: int (default 1)

only relevant if reflection is larger than 0

  • 1: rays start upwards

  • 2: rays start downwards

Returns:
(original x1, x1 of path segment, original x2, x2 of path segment, C_0, C_1 of path segment)
get_angle(x, x_start, C_0, reflection=0, reflection_case=1, in_air=False)[source]

calculates the angle with respect to the positive z-axis of the ray path at position x

Parameters:
x: tuple

(y, z) coordinate to calculate the angle

x_start: tuple

(y, z) start position of the ray

C_0: float

C_0 parameter of analytic ray path function

reflection: int (default 0)

the number of bottom reflections to consider

reflection_case: int (default 1)

only relevant if reflection is larger than 0

  • 1: rays start upwards

  • 2: rays start downwards

get_launch_angle(x1, C_0, reflection=0, reflection_case=1)[source]
get_receive_angle(x1, x2, C_0, reflection=0, reflection_case=1)[source]
get_reflection_angle(x1, x2, C_0, reflection=0, reflection_case=1)[source]

calculates the angle under which the ray reflects off the surface. If not reflection occurs, None is returned

If reflections off the bottom (e.g. Moore’s Bay) are simulated, an array with reflection angles (one for each track segment) is returned

Parameters:
x1: tuple

(y, z) start position of ray

x2: tuple

(y, z) stop position of the ray

C_0: float

C_0 parameter of analytic ray path function

reflection: int (default 0)

the number of bottom reflections to consider

reflection_case: int (default 1)

only relevant if reflection is larger than 0 * 1: rays start upwards * 2: rays start downwards

get_path(x1, x2, C_0, n_points=1000)[source]

for plotting purposes only, returns the ray tracing path between x1 and x2

the result is only valid if C_0 is a solution to the ray tracing problem

Parameters:
x1: array

start position (y, z)

x2: array

stop position (y, z)

C_0: (float)

first parameter

n_points: integer (optional)

the number of coordinates to calculate

Returns:
yy: array

the y coordinates of the ray tracing path

zz: array

the z coordinates of the ray tracing path

get_path_reflections(x1, x2, C_0, n_points=1000, reflection=0, reflection_case=1)[source]

calculates the ray path in the presence of reflections at the bottom The full path is constructed by multiple calls to the get_path() function to put together the full path

Parameters:
x1: tuple

(y, z) coordinate of start value

x2: tuple

(y, z) coordinate of stop value

C_0: float

C_0 parameter of analytic ray path function

n_points: int (default 1000)

the number of points of the numeric path

reflection: int (default 0)

the number of bottom reflections to consider

reflection_case: int (default 1)

only relevant if reflection is larger than 0

  • 1: rays start upwards

  • 2: rays start downwards

Returns:
yy: array

the y coordinates of the ray tracing path

zz: array

the z coordinates of the ray tracing path

get_reflection_point(C_0, C_1)[source]

calculates the point where the signal gets reflected off the bottom of the ice shelf

Returns tuple (y,z)

obj_delta_y_square(logC_0, x1, x2, reflection=0, reflection_case=2)[source]

objective function to find solution for C0

obj_delta_y(logC_0, x1, x2, reflection=0, reflection_case=2)[source]

function to find solution for C0, returns distance in y between function and x2 position result is signed! (important to use a root finder)

get_delta_y(C_0, x1, x2, C0range=None, reflection=0, reflection_case=2)[source]

calculates the difference in the y position between the analytic ray tracing path specified by C_0 at the position x2

determine_solution_type(x1, x2, C_0)[source]

returns the type of the solution

Parameters:
x1: 2dim np.array

start position

x2: 2dim np.array

stop position

C_0: float

C_0 value of ray tracing solution

Returns:
solution_type: int
  • 1: ‘direct’

  • 2: ‘refracted’

  • 3: ‘reflected

find_solutions(x1, x2, plot=False, reflection=0, reflection_case=1)[source]

this function finds all ray tracing solutions

prerequesite is that x2 is above and to the right of x1, this is not a violation of universality because this requirement can be achieved with a simple coordinate transformation

Parameters:
x1: tuple

(y,z) coordinate of start point

x2: tuple

(y,z) coordinate of stop point

reflection: int (default 0)

how many reflections off the reflective layer (bottom of ice shelf) should be simulated

returns an array of the C_0 paramters of the solutions (the array might be empty)
plot_result(x1, x2, C_0, ax)[source]

helper function to visualize results

get_angle_from_C_0(C_0, z_pos, angoff=0, in_air=False)[source]
get_angle_from_logC_0(logC_0, z_pos, angoff=0, in_air=False)[source]

argument angoff is provided so that the function can be used for minimization in get_C_0_from_angle(), in which case angoff is the angle for which the C_0 is sought and zero is returned when it is found.

C_0 has a smallest possible value at 1./self.medium.n_ice . When it approaches this value, very small changes in C_0 correspond to a given change in the angle. In order to prevent the root finding algorithm from crossing into the invalid range of C_0 at values smaller than 1./self.medium.n_ice, the root finding is done with the parameter logC_0 = np.log(C_0 - 1. / self.medium.n_ice), so it is not exactly the log of C_0 as the nome of this method implies. This is the same parameter transformation that is done for find_solutions()

input:

logC_0 = np.log(C_0 - 1. / self.medium.n_ice) angoff = angular offset z_pos = z-position from where ray is emitted

output:

angle corresponding to C_0, minus offset angoff

get_C_0_from_angle(anglaunch, z_pos, in_air=False)[source]

Find parameter C0 corresponding to a given launch angle and z-position of a ray. The parameter is found by means of a root finding procedure

output:

Complete output of optimisation procedure (result.x[0] is the C0 value found by optimisation procedure)

get_z_from_n(n)[source]

get z from given n - equation from get_n solved for z

get_surf_skim_angle(x1)[source]

For a given position x1 = [x,z] and depth profile self.n(), find the angle at which a beam must be emitted to “skim the surface”, i.e. arrive horizontally (angle = 90 deg) at the surface; This is used to find the refraction zone.

Returns:
C0crit: float

C0 of critical angle

thcrit: float

critical angle

is_in_refraction_zone(x1, x2, C0crit=None, plot=False)[source]

Find if receiver at x2 is in the refraction zone of emitter at x1. The refraction zone is the oposite of the shadow zone.

If the C0 of the critical angle, C0crit, is provided, it will not be calculated. This is useful in case find_solutions() is called and C0crit is calculated in the process of determining the initial value for the minimization procedure.

Returns True if x2 is in the refraction zone of x1 - note that the inverse statement is not necessarily true, i.e. when False is returned, it is possible that x2 is in the refraction zone nonetheless

TODO: Why does the reference point not seem to lie exactly on the mirrored path? Instead of returning True/False, it might be useful to return ycheck - x2[0] (in case x2[0]>ycrit), which gives some idea of how close the receiver is to the refraction zone. This could be used to define a “gray zone’ and a ‘far zone’, in which the receiver is most definitely in the shadow zone

get_tof_for_straight_line(x1, x2)[source]

Calculate the time of flight for a hypothatical ray travelling from x1 to x2 in a straight line. Such an array in general is not a solution consistant with Fermat’s principle. It is however useful as a reference time or approximation for signals not explicable with geometric optics.

output:

time of flight for a ray travelling straight from x1 to x2

get_surface_pulse(x1, x2, infirn=False, angle='critical', chdraw=None, label=None)[source]

Calculate the time for a ray to travel from x1 to the surface and arriving at the surface with (a) critical angle or (b) Brewster angle, propagating in a straight line along the surface in air (n=1) in the firn at z=0 (n=self.n(0)) and then reaching the receiver by returning into the firn just at the point to reach the receiver at x2, entering the firn at the surface at the same angle it reached the surface from x1..

Parameters:
x1, x2: arrays

Arrays with x and z positions of emitter x1 and receiver x2

infirn: Boolean.

Set to True if surface ray travels in the firn, set to False (default) if it travels in air.

angle: String

specifying angle at which ray reaches/leaves the surface. Can be ‘Brewster’ or ‘critical’ If neither of these is chosen, a warning is printed and angle is set to ‘critical’

chdraw: string or None

If None, do not draw the path of the ray. If the ray should be drawn, a string consistent with the matplotlib.pyplot library has to be specified, e.g. ‘r:’ to draw a dotted red line. It is assumed that an appropriate figure on which to draw the ray has been created and set as current figure by the user before calling this method.

label: Label for plot
angular_diff(x_refl, z_refl, pulser_pos, receiver_pos, ipulssol, irxsol)[source]

This is a helper function to find a ray that is subject to specular reflection (no transmission) at the ice-water interface, e.g. for Moore’s Bay. For a (virtual) emitter positioned at [x_refl,z_refl], it finds the emission angles such that the rays hit positions pulser_pos and receiver_pos. ipulssol = 0 or 1, respectively means that pulser_pos is hit directly or by means of reflection at the surface, respectively. irxsol is the equivalent parameter for receiver_pos.

angular_diff can be used as the function to find a root of in the following fashion:

result = optimize.root(raytr.angular_diff, x0=x_refl_start, args=(z_refl,pulser_pos,receiver_pos,ipulssol,irxsol))

Then get the final value by x_refl = result.x[0] (z_refl is fixed)

The idea is to treat the reflection point as a virtual emitter and then find the x-position at the predefined depth z_refl, for which the emisssion angles to pulser_pos and receiver_pos are the same (i.e. the output of angular_diff is zero). The x-position would be output “result” ofoptimize.root() above.

Returns:
result: float

Is zero if the angles (w.r.t. the vertical) of rays emitted from [x_refl,z_refl] to positions pulser_pos and receiver_pos are the same or greater than zero, if this is not the case. For exact defintion, see “result” in code below

apply_propagation_effects(efield, i_solution)

Apply propagation effects to the electric field Note that the 1/r weakening of the electric field is already accounted for in the signal generation

Parameters:
efield: ElectricField object

The electric field that the effects should be applied to

i_solution: int

Index of the raytracing solution the propagation effects should be based on

Returns:
efield: ElectricField object

The modified ElectricField object

get_attenuation(iS, frequency, max_detector_freq=None)

calculates the signal attenuation due to attenuation in the medium (ice)

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

frequency: array of floats

the frequencies for which the attenuation is calculated

max_detector_freq: float or None

the maximum frequency of the final detector sampling (the simulation is internally run with a higher sampling rate, but the relevant part of the attenuation length calculation is the frequency interval visible by the detector, hence a finer calculation is more important)

Returns:
attenuation: array of floats

the fraction of the signal that reaches the observer (only ice attenuation, the 1/R signal falloff not considered here)

get_config()

Function that returns the configuration currently used by the raytracer

get_launch_vector(iS)

calculates the launch vector (in 3D) of solution iS

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

Returns:
launch_vector: 3dim np.array

the launch vector

get_number_of_raytracing_solutions()

Function that returns the maximum number of raytracing solutions that can exist between each given pair of start and end points

get_number_of_solutions()

returns the number of solutions

get_output_parameters()

Returns a list with information about parameters to include in the output data structure that are specific to this raytracer

! be sure that the first entry is specific to your raytracer !

Returns:
list with entries of form [{‘name’: str, ‘ndim’: int}]

! be sure that the first entry is specific to your raytracer ! ‘name’: Name of the new parameter to include in the data structure ‘ndim’: Dimension of the data structure for the parameter

get_raytracing_output(i_solution)

Write parameters that are specific to this raytracer into the output data.

Parameters:
i_solution: int

The index of the raytracing solution

Returns:
dictionary with the keys matching the parameter names specified in get_output_parameters and the values being
the results from the raytracing
get_receive_vector(iS)

calculates the receive vector (in 3D) of solution iS

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

Returns:
receive_vector: 3dim np.array

the receive vector

get_results()
get_solution_type(iS)

returns the type of the solution

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

Returns:
solution_type: int

integer corresponding to the types in the dictionary solution_types

has_solution()

checks if ray tracing solution exists

reset_solutions()
set_config(config)

Function to change the configuration file used by the raytracer

Parameters:
config: dict

The new configuration settings

set_start_and_end_point(x1, x2)

Set the start and end points between which raytracing solutions shall be found It is recommended to also reset the solutions from any previous raytracing to avoid confusing them with the current solution

Parameters:
x1: np.array of shape (3,), default unit

start point of the ray

x2: np.array of shape (3,), default unit

stop point of the ray

use_optional_function(function_name, *args, **kwargs)

Use optional function which may be different for each ray tracer. If the name of the function is not present for the ray tracer this function does nothing.

Parameters:
function_name: string

name of the function to use

*args: type of the argument required by function

all the neseccary arguments for the function separated by a comma

**kwargs: type of keyword argument of function

all all the neseccary keyword arguments for the function in the form of key=argument and separated by a comma

Examples

use_optional_function('set_shower_axis',np.array([0,0,1]))
use_optional_function('set_iterative_sphere_sizes',sphere_sizes=np.aray([3,1,.5]))
class NuRadioMC.SignalProp.analyticraytracing.ray_tracing(medium, attenuation_model='SP1', log_level=30, n_frequencies_integration=100, n_reflections=0, config=None, detector=None, ray_tracing_2D_kwards={}, use_cpp=True)[source]

Bases: ray_tracing_base

utility class (wrapper around the 2D analytic ray tracing code) to get ray tracing solutions in 3D for two arbitrary points x1 and x2

Methods

apply_propagation_effects(efield, i_solution)

Apply propagation effects to the electric field Note that the 1/r weakening of the electric field is already accounted for in the signal generation

find_solutions()

find all solutions between x1 and x2

get_attenuation(iS, frequency[, ...])

calculates the signal attenuation due to attenuation in the medium (ice)

get_config()

Function that returns the configuration currently used by the raytracer

get_effective_index_birefringence(direction, ...)

Function to find the analytical solutions for the effective refractive indices.

get_focusing(iS[, dz, limit])

calculate the focusing effect in the medium

get_launch_vector(iS)

calculates the launch vector (in 3D) of solution iS

get_number_of_raytracing_solutions()

Function that returns the maximum number of raytracing solutions that can exist between each given pair of start and end points

get_number_of_solutions()

returns the number of solutions

get_output_parameters()

Returns a list with information about parameters to include in the output data structure that are specific to this raytracer

get_path(iS[, n_points])

helper function that returns the 3D ray tracing path of solution iS

get_path_length(iS[, analytic])

calculates the path length of solution iS

get_path_properties_birefringence(i_solution)

Function to extract important information about the birefringent propagation along the path.

get_polarization_birefringence(N1, N2, ...)

Function for the normalized e-field vector of a wave for the direction of propagation in spherical coordinates with special cases.

get_polarization_birefringence_simple(n, ...)

Function for the normalized e-field vector of a wave for the direction of propagation in cartesian coordinates without special cases.

get_pulse_propagation_birefringence(pulse, ...)

Function for the time trace propagation according to the polarization change due to birefringence.

get_raytracing_output(i_solution)

Write parameters that are specific to this raytracer into the output data.

get_receive_vector(iS)

calculates the receive vector (in 3D) of solution iS

get_reflection_angle(iS)

calculates the angle of reflection at the surface (in case of a reflected ray)

get_solution_type(iS)

returns the type of the solution

get_travel_time(iS[, analytic])

calculates the travel time of solution iS

has_solution()

checks if ray tracing solution exists

on_sky_birefringence(theta, phi, polarization)

Function for the normalized e-field vector from cartesian to spherical coordinates.

reset_solutions()

Resets the raytracing solutions back to None.

set_config(config)

Change the configuration file used by the raytracer

set_solution(raytracing_results)

Read an already calculated raytracing solution from the input array

set_start_and_end_point(x1, x2)

Set the start and end points of the raytracing

use_optional_function(function_name, *args, ...)

Use optional function which may be different for each ray tracer.

get_ray_path

get_results

class initilization

Parameters:
medium: medium class

class describing the index-of-refraction profile

attenuation_model: string

signal attenuation model

log_name: string

name under which things should be logged

log_level: logging object

specify the log level of the ray tracing class

  • logging.ERROR

  • logging.WARNING

  • logging.INFO

  • logging.DEBUG

default is WARNING

n_frequencies_integration: int

the number of frequencies for which the frequency dependent attenuation length is being calculated. The attenuation length for all other frequencies is obtained via linear interpolation.

n_reflections: int (default 0)

in case of a medium with a reflective layer at the bottom, how many reflections should be considered

config: dict

a dictionary with the optional config settings. If None, the config is intialized with default values, which is needed to avoid any “key not available” errors. The default settings are

  • self._config = {‘propagation’: {}}

  • self._config[‘propagation’][‘attenuate_ice’] = True

  • self._config[‘propagation’][‘focusing_limit’] = 2

  • self._config[‘propagation’][‘focusing’] = False

  • self._config[‘propagation’][‘birefringence’] = False

detector: detector object
ray_tracing_2D_kwards: dict

Additional arguments which are passed to ray_tracing_2D

use_cpp: bool

if True, use CPP implementation of minimization routines default: True if CPP version is available

Methods

apply_propagation_effects(efield, i_solution)

Apply propagation effects to the electric field Note that the 1/r weakening of the electric field is already accounted for in the signal generation

find_solutions()

find all solutions between x1 and x2

get_attenuation(iS, frequency[, ...])

calculates the signal attenuation due to attenuation in the medium (ice)

get_config()

Function that returns the configuration currently used by the raytracer

get_effective_index_birefringence(direction, ...)

Function to find the analytical solutions for the effective refractive indices.

get_focusing(iS[, dz, limit])

calculate the focusing effect in the medium

get_launch_vector(iS)

calculates the launch vector (in 3D) of solution iS

get_number_of_raytracing_solutions()

Function that returns the maximum number of raytracing solutions that can exist between each given pair of start and end points

get_number_of_solutions()

returns the number of solutions

get_output_parameters()

Returns a list with information about parameters to include in the output data structure that are specific to this raytracer

get_path(iS[, n_points])

helper function that returns the 3D ray tracing path of solution iS

get_path_length(iS[, analytic])

calculates the path length of solution iS

get_path_properties_birefringence(i_solution)

Function to extract important information about the birefringent propagation along the path.

get_polarization_birefringence(N1, N2, ...)

Function for the normalized e-field vector of a wave for the direction of propagation in spherical coordinates with special cases.

get_polarization_birefringence_simple(n, ...)

Function for the normalized e-field vector of a wave for the direction of propagation in cartesian coordinates without special cases.

get_pulse_propagation_birefringence(pulse, ...)

Function for the time trace propagation according to the polarization change due to birefringence.

get_raytracing_output(i_solution)

Write parameters that are specific to this raytracer into the output data.

get_receive_vector(iS)

calculates the receive vector (in 3D) of solution iS

get_reflection_angle(iS)

calculates the angle of reflection at the surface (in case of a reflected ray)

get_solution_type(iS)

returns the type of the solution

get_travel_time(iS[, analytic])

calculates the travel time of solution iS

has_solution()

checks if ray tracing solution exists

on_sky_birefringence(theta, phi, polarization)

Function for the normalized e-field vector from cartesian to spherical coordinates.

reset_solutions()

Resets the raytracing solutions back to None.

set_config(config)

Change the configuration file used by the raytracer

set_solution(raytracing_results)

Read an already calculated raytracing solution from the input array

set_start_and_end_point(x1, x2)

Set the start and end points of the raytracing

use_optional_function(function_name, *args, ...)

Use optional function which may be different for each ray tracer.

get_ray_path

get_results

reset_solutions()[source]

Resets the raytracing solutions back to None. This is useful to do when changing the start and end points in order to not accidentally use results from previous raytracings.

set_start_and_end_point(x1, x2)[source]

Set the start and end points of the raytracing

Parameters:
x1: 3dim np.array

start point of the ray

x2: 3dim np.array

stop point of the ray

set_solution(raytracing_results)[source]

Read an already calculated raytracing solution from the input array

Parameters:
raytracing_results: dict

The dictionary containing the raytracing solution.

find_solutions()[source]

find all solutions between x1 and x2

get_solution_type(iS)[source]

returns the type of the solution

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

Returns:
solution_type: int

integer corresponding to the types in the dictionary solution_types

get_path(iS, n_points=1000)[source]

helper function that returns the 3D ray tracing path of solution iS

Parameters:
iS: int

ray tracing solution

n_points: int

number of points of path

get_effective_index_birefringence(direction, nx, ny, nz)[source]

Function to find the analytical solutions for the effective refractive indices. The calculations are described here: https://link.springer.com/article/10.1140/epjc/s10052-023-11238-y

Parameters:
direction: numpy.array

propagation direction of the wave

nx: float

the index of refraction in the x-direction

ny: float

the index of refraction in the y-direction

nz: float

the index of refraction in the z-direction

Returns:
output format: np.array([n1, n2])
meaning: effective refractive indices
get_polarization_birefringence_simple(n, direction, nx, ny, nz)[source]

Function for the normalized e-field vector of a wave for the direction of propagation in cartesian coordinates without special cases. For the function with special cases see get_polarization_birefringence. For a birefringent medium, the e-field vector is calculated from the diagonalized dielectric tensor and the propagation direction. The calculations are described here: https://link.springer.com/article/10.1140/epjc/s10052-023-11238-y

Parameters:
n: float

the effective index of refraction in the propagation direction calculated by get_effective_index_birefringence

direction: numpy.array

propagation direction of the wave

nx: float

the index of refraction in the x-direction

ny: float

the index of refraction in the y-direction

nz: float

the index of refraction in the z-direction

Returns:
efieldnp.ndarray of shape (3,)

normalized e-field vector in cartesian coordinates

get_polarization_birefringence(N1, N2, direction, nx, ny, nz)[source]

Function for the normalized e-field vector of a wave for the direction of propagation in spherical coordinates with special cases. For a birefringent medium, the e-field vector is calculated from the diagonalized dielectric tensor and the propagation direction. The calculations are described here: https://link.springer.com/article/10.1140/epjc/s10052-023-11238-y

Parameters:
N1: float

the first effective index of refraction in the propagation direction calculated by get_effective_index_birefringence

N2: float

the second effective index of refraction in the propagation direction calculated by get_effective_index_birefringence

direction: numpy.array

propagation direction of the wave

nx: float

the index of refraction in the x-direction

ny: float

the index of refraction in the y-direction

nz: float

the index of refraction in the z-direction

Returns:
efieldnp.ndarray of shape (2, 3)

normalized e-field vector in spherical coordinates for both birefringence solutions

on_sky_birefringence(theta, phi, polarization)[source]

Function for the normalized e-field vector from cartesian to spherical coordinates. The function does the same as the following radiotool functions, only faster: from radiotools import coordinatesystems cs = coordinatesystems.cstrafo(theta, phi) sky = cs.transform_from_ground_to_onsky(p)

Parameters:
theta: float

zenith angle of the propagation direction

phi: float

azimuth angle of the propagation direction

polarization: np.array([px, py, pz])

normalized e-field vector in cartesian coordinates

Returns:
efieldnp.ndarray of shape (3,)

normalized e-field vector in spherical coordinates

get_pulse_propagation_birefringence(pulse, samp_rate, i_solution, bire_model='southpole_A')[source]

Function for the time trace propagation according to the polarization change due to birefringence. The trace propagation is explained in this paper: https://link.springer.com/article/10.1140/epjc/s10052-023-11238-y

Parameters:
pulse: np.ndarray

3d array with the frequency spectrum of np.array([eR, eTheta, ePhi]), usually provided by the apply_propagation_effects function

samp_rate: float

sampling rate of the time traces

i_solution: int

choose which ray-tracing solution should be propagated

bire_model: string

choose the interpolation to fit the measured refractive index data options include (A, B, C, D, E) description can be found under: NuRadioMC/NuRadioMC/utilities/birefringence_models/model_description

Returns:
final pulse: numpy.array([eR, eTheta, ePhi])

[0] - eR - final frequency spectrum of the radial component - not altered by the function [1] - eTheta - final frequency spectrum of the theta component [2] - ePhi - final frequency spectrum of the phi component

get_path_properties_birefringence(i_solution, bire_model='southpole_A')[source]

Function to extract important information about the birefringent propagation along the path. The important properties include effective refractive indices, polarization eigenvectors, and incremental time delays

Parameters:
i_solution: int

choose which ray-tracing solution should be propagated

bire_model: string

choose the interpolation to fit the measured refractive index data options include (A, B, C, D, E) description can be found under: NuRadioMC/NuRadioMC/utilities/birefringence_models/model_description

Returns:
path_properties: dict

a dictionary containing the following keys:

  • ‘path’: np.ndarray - propagation path in x, y, z with the same granularity as the nirefringent propagation

  • ‘nominal_refractive_index’: np.ndarray - nominal refractive index if only density effects were taken into account

  • ‘refractive_index_x’: np.ndarray - refractive index for the x-direction

  • ‘refractive_index_y’: np.ndarray - refractive index for the y-direction

  • ‘refractive_index_z’: np.ndarray - refractive index for the z-direction

  • ‘first_refractive_index’: np.ndarray - effective refractive index of the first birefringent state along the full path

  • ‘second_refractive_index’: np.ndarray - effective refractive index of the second birefringent state along the full path

  • ‘first_polarization_vector’: np.ndarray - polarization vector of the first birefringent state in spherical coordinates along the full path

  • ‘second_polarization_vector’: np.ndarray - polarization vector of the second birefringent state in spherical coordinates along the full path

  • ‘first_time_delay’: np.ndarray - incremental time delays of the first birefringent state along the full path

  • ‘second_time_delay’: np.ndarray - incremental time delays of the second birefringent state along the full path

get_launch_vector(iS)[source]

calculates the launch vector (in 3D) of solution iS

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

Returns:
launch_vector: 3dim np.array

the launch vector

get_receive_vector(iS)[source]

calculates the receive vector (in 3D) of solution iS

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

Returns:
receive_vector: 3dim np.array

the receive vector

get_reflection_angle(iS)[source]

calculates the angle of reflection at the surface (in case of a reflected ray)

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

Returns:
reflection_angle: float or None

the reflection angle (for reflected rays) or None for direct and refracted rays

get_path_length(iS, analytic=True)[source]

calculates the path length of solution iS

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

analytic: bool

If True the analytic solution is used. If False, a numerical integration is used. (default: True)

Returns:
distance: float

distance from x1 to x2 along the ray path

get_travel_time(iS, analytic=True)[source]

calculates the travel time of solution iS

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

analytic: bool

If True the analytic solution is used. If False, a numerical integration is used. (default: True)

Returns:
time: float

travel time

get_attenuation(iS, frequency, max_detector_freq=None)[source]

calculates the signal attenuation due to attenuation in the medium (ice)

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

frequency: array of floats

the frequencies for which the attenuation is calculated

max_detector_freq: float or None

the maximum frequency of the final detector sampling (the simulation is internally run with a higher sampling rate, but the relevant part of the attenuation length calculation is the frequency interval visible by the detector, hence a finer calculation is more important)

Returns:
attenuation: array of floats

the fraction of the signal that reaches the observer (only ice attenuation, the 1/R signal falloff not considered here)

get_focusing(iS, dz=-0.01, limit=2.0)[source]

calculate the focusing effect in the medium

Parameters:
iS: int

choose for which solution to compute the launch vector, counting starts at zero

dz: float

the infinitesimal change of the depth of the receiver, 1cm by default

limit: float

The maximum signal focusing.

Returns:
focusing: float

gain of the signal at the receiver due to the focusing effect

get_ray_path(iS)[source]
get_output_parameters()[source]

Returns a list with information about parameters to include in the output data structure that are specific to this raytracer

! be sure that the first entry is specific to your raytracer !

Returns:
list with entries of form [{‘name’: str, ‘ndim’: int}]

! be sure that the first entry is specific to your raytracer ! ‘name’: Name of the new parameter to include in the data structure ‘ndim’: Dimension of the data structure for the parameter

get_config()

Function that returns the configuration currently used by the raytracer

get_number_of_raytracing_solutions()

Function that returns the maximum number of raytracing solutions that can exist between each given pair of start and end points

get_number_of_solutions()

returns the number of solutions

get_raytracing_output(i_solution)[source]

Write parameters that are specific to this raytracer into the output data.

Parameters:
i_solution: int

The index of the raytracing solution

Returns:
dictionary with the keys matching the parameter names specified in get_output_parameters and the values being
the results from the raytracing
get_results()
has_solution()

checks if ray tracing solution exists

use_optional_function(function_name, *args, **kwargs)

Use optional function which may be different for each ray tracer. If the name of the function is not present for the ray tracer this function does nothing.

Parameters:
function_name: string

name of the function to use

*args: type of the argument required by function

all the neseccary arguments for the function separated by a comma

**kwargs: type of keyword argument of function

all all the neseccary keyword arguments for the function in the form of key=argument and separated by a comma

Examples

use_optional_function('set_shower_axis',np.array([0,0,1]))
use_optional_function('set_iterative_sphere_sizes',sphere_sizes=np.aray([3,1,.5]))
apply_propagation_effects(efield, i_solution)[source]

Apply propagation effects to the electric field Note that the 1/r weakening of the electric field is already accounted for in the signal generation

Parameters:
efield: ElectricField object

The electric field that the effects should be applied to

i_solution: int

Index of the raytracing solution the propagation effects should be based on

Returns:
efield: ElectricField object

The modified ElectricField object

set_config(config)[source]

Change the configuration file used by the raytracer

Parameters:
config: dict or None

The new configuration settings If None, the default config settings will be applied