Trace

Trace allows collection of detailed data along the photon packet path. The trace functionality is implemented in xopto.mcbase.mctrace.Trace and also conveniently imported into the xopto.mcml.mctrace and xopto.mcml.mc modules.

The events that are traced by xopto.mcbase.mctrace.Trace: can be customized through the options argument.

  • options=Trace.TRACE_START - only the start state of the photon packets is traced (this is the state of the photon packet after launching and subtraction the weight due to specular reflections at the source-medium or sample boundary).

  • options=Trace.TRACE_END - only the end state of the photon packet is traced (the photon packet is absorbed, leaves the sample through the top or botom surface or escapes the simulation radius).

  • options=Trace.TRACE_START | Trace.TRACE_END - the start and the end states of the photon packet are traced.

  • options=Trace.TRACE_ALL - the entire path of the photon packet is traced (default).

The following events are traced (in addition to the start and end state of the photon packet):

  • refraction or reflection at the layer/voxel boundary.

  • crossing the layer or voxel boundary (even if the materials of the geometrical entities are the same)

  • scattering and absorption events

Each trace entry includes eight floating-point numbers in the following order:

  • x - x coordinate of the photon packet position (m).

  • y - y coordinate of the photon packet position (m).

  • z - z coordinate of the photon packet position (m).

  • px - x component of the propagation direction vector (m).

  • py - y component of the propagation direction vector (m).

  • pz - z component of the propagation direction vector (m).

  • w - weight of the photon packet.

  • pl - total optical path length since the start event. Note that tracing of the optical path length can be tuned off by passing plon=False to Trace().

After running a Monte carlo simulation, the trace data are stored in a structured numpy array of shape (num_packets, max_len), where num_packets is the number of launched photon packets and maxlen is the maximum length of the trace that is set through the constructor parameter maxlen. The actual number of trace entries for each photon packet can be obtained from the n property, which is a numpy vector of length num_packets.

When tracing the entire path of the photon packets, it is important to select an adequate maximum length of the trace. The value will depend on a number of factors including the optical properties of the sample and the configuration of the investigated source and detector. It is advised to iteratively determine the maximum length of the trace by observing the number of trace events xopto.mcbase.mctrace.Trace.n. If the maximum trace length is exceeded during the Monte Carlo simulation, the last trace entry is overwritten with the final state of the photon packet.

The terminal states of the traced photon packets can be easily accessed through the terminal property.

# terminal weights of all the traced photon packets
weights = trace.terminal['w']
# terminal x coordinate of all the traced photon packets
x = trace.terminal['x']

# terminal y coordinate of the first traced photon packets
y = trace.terminal['y'][0]

Trace Filter

Frequently, only the traces of photon packets that end with a particular state are useful. A versatile filter Filter is available for filtering the photon packet traces according to the terminal position, propagation direction, and optical path length of the photon packet. The trace filter Filter is also conveniently imported into the xopto.mcml.mctrace and xopto.mcvox.mctrace modules.

The traces can be filtered by the final state of the photon packet. For a detailed list of available options see Filter. If a filter object is passed to the Trace constructor (parameter filter) the results of the Monte Carlo simulation will be automatically filtered. Alternatively, filters can be easily applied to the existing Trace instances by the xopto.mcbase.mctrace.Filter.__call__() method. By default, the filter modifies the trace (does not create a new trace object for the filtered results) and returns it. This behavior can be changed by passing update=False to the __call__() method. This will create a new Trace instance for the filtered trace.

Note

The filter will also drop all the traces that exceed the maximum trace length, even if they satisfy the filter. The number of dropped traces is returned as the second output of __call__().

In the following example we create a trace with a maximum number of events / length set to 1000 and a filter that only selects traces of photon packets that end at the top sample surface and are escaping within 20 o of the surface normal. Note that the z component of the propagation direction is negative for photon packets that escape the sample through the top surface.

from xopto.mcml import mc
import numpy as np

filt = mc.mctrace.Filter(
    z=(-float('inf'), 0.0)
    pz=(-1.0, -np.cos(np.deg2rad(20.0)))
)
mc.mctrace.Trace(maxlen=1000, filter=filt)

# apply the filter and modify the input trace
some_trace, dropped = filt(some_trace)

# apply the filter and create a new trace instance
filtered_trace, dropped = filt(some_trace, False)

To visualize the paths of photon packets after completing a simulation, use the plot() method. The traces can be plot in several different views:

  • view=’3d’ produces a 3D view of the traces (default).

  • view=’xy’ produces a 2D view of the traces projected onto the x-y plane.

  • view=’xz’ produces a 2D view of the traces projected onto the x-z plane.

  • view=’yz’ produces a 2D view of the traces projected onto the y-z plane.

The show parameter controls if the plot window is shown before returning from plot(). Note that the starting point of the traced path is highlighted with a green marker and the final point of the traced path with a red marker.

trace.plot(show=False) # produces a 3D view
trace.plot(view='xy', show=False) # produces a 2D view in the x-y plane
trace.plot(view='xz', show=False) # produces a 2D view in the x-z plane
trace.plot(view='yz', show=True) # produces a 2D view in the y-z plane