Tracker

class plasmapy.diagnostics.charged_particle_radiography.synthetic_radiography.Tracker(
grids: AbstractGrid | Iterable[AbstractGrid],
source: Quantity,
detector: Quantity,
dt=None,
dt_range=None,
field_weighting: Literal['volume averaged', 'nearest neighbor'] = 'volume averaged',
detector_hdir=None,
output_directory: Path | None = None,
output_basename: str = 'output',
fraction_exited_threshold: float = 0.999,
verbose: bool = True,
)[source]

Bases: ParticleTracker

Represents a charged particle radiography experiment with simulated or calculated E and B fields given at positions defined by a grid of spatial coordinates. The particle source and detector plane are defined by vectors from the origin of the grid.

Parameters:
  • grids (AbstractGrid or subclass thereof, or list) – of same. A Grid object or list of grid objects containing the required quantities [E_x, E_y, E_z, B_x, B_y, B_z]. If any of these quantities are missing, a warning will be given and that quantity will be assumed to be zero everywhere.

  • source (Quantity, shape (3)) –

    A vector pointing from the origin of the grid to the location of the particle source. This vector will be interpreted as being in either Cartesian, cylindrical, or spherical coordinates based on its units. Valid geometries are:

    • Cartesian (x,y,z) : (meters, meters, meters)

    • Cylindrical (r, theta, z) : (meters, radians, meters)

    • Spherical (r, theta, phi) : (meters, radians, radians)

    In spherical coordinates theta is the polar angle.

  • detector (Quantity, shape (3)) – A vector pointing from the origin of the grid to the center of the detector plane. The vector from the source point to this point defines the normal vector of the detector plane. This vector can also be specified in Cartesian, cylindrical, or spherical coordinates (see the source keyword).

  • dt (Quantity, optional) – An explicitly set time step in units convertible to seconds. Setting this optional keyword overrules the adaptive time step capability and forces the use of this time step throughout.

  • dt_range (tuple of shape (2,) of Quantity, optional) – If specified, the calculated adaptive time step will be clamped between the first and second values.

  • field_weighting (str, default: "volume averaged") –

    String that selects the field weighting algorithm used to determine what fields are felt by the particles. Options are:

    • "nearest neighbor":

      Particles are assigned the fields on the grid vertex closest to them.

    • "volume averaged":

      The fields experienced by a particle are a volume-average of the eight grid points surrounding them.

  • detector_hdir (numpy.ndarray, shape (3), optional) –

    A unit vector (in Cartesian coordinates) defining the horizontal direction on the detector plane. By default, the horizontal axis in the detector plane is defined to be perpendicular to both the source-to-detector vector and the z-axis (unless the source-to-detector axis is parallel to the z axis, in which case the horizontal axis is the x-axis).

    The detector vertical axis is then defined to be orthogonal to both the source-to-detector vector and the detector horizontal axis.

  • output_directory (Path, optional) – Directory for objects that are saved to disk. If a directory is not specified then a memory save routine is used.

  • output_basename (str, optional) – Optional base name for output files.

  • fraction_exited_threshold (float, optional) – The fraction of particles that must leave the grids to terminate the simulation. This does not include particles that have never entered the grids. By default, this is set to 0.999 (corresponding to 99.9%).

  • verbose (bool, optional) – If True, updates on the status of the program will be printed into the standard output while running.

Attributes Summary

fract_entered

The fraction of particles that have entered the grid.

is_adaptive_time_step

Return whether the simulation is calculating an adaptive time step or using the user-provided time step.

is_synchronized_time_step

Return if the simulation is applying the same time step across all particles.

max_deflection

The maximum deflection experienced by one of the particles, determined by comparing their initial and final velocity vectors.

nparticles_removed

Return the number of particles currently being tracked.

nparticles_stopped

Return the number of particles currently being tracked.

nparticles_tracked

Return the number of particles currently being tracked.

num_entered

Count the number of particles that have entered the grids.

num_grids

The number of grids specified at instantiation.

on_any_grid

Binary array for each particle indicating whether it is currently on ANY grid.

particles_on_grid

Returns a boolean mask of shape [ngrids, nparticles] corresponding to whether or not the particle is on the associated grid.

results_dict

A dictionary containing the results of the simulation.

vmax

The maximum velocity of any particle in the simulation.

Methods Summary

add_stopping(method[, materials, I])

Enable particle stopping using experimental stopping powers.

add_wire_mesh(location, extent, nwires, ...)

Add a wire mesh grid between the particle source and the object grid that blocks particles whose paths intersect the wires.

create_particles(nparticles, particle_energy)

Generates the angular distributions about the Z-axis, then rotates those distributions to align with the source-to-detector axis.

load_particles(x, v[, particle])

Load arrays of particle positions and velocities.

run()

Runs a particle-tracing simulation.

setup_adaptive_time_step([...])

Set parameters for the adaptive time step candidates.

Attributes Documentation

fract_entered

The fraction of particles that have entered the grid. The denominator of this fraction is based off the number of tracked particles, and therefore does not include stopped or removed particles.

is_adaptive_time_step

Return whether the simulation is calculating an adaptive time step or using the user-provided time step.

is_synchronized_time_step

Return if the simulation is applying the same time step across all particles.

max_deflection

The maximum deflection experienced by one of the particles, determined by comparing their initial and final velocity vectors.

This value can be used to determine the charged particle radiography regime using the dimensionless number defined by Kugland et al. 2012

Returns:

max_deflection – The maximum deflection in radians.

Return type:

float

nparticles_removed

Return the number of particles currently being tracked. That is, they do not have NaN position or velocity.

nparticles_stopped

Return the number of particles currently being tracked. That is, they do not have NaN position or velocity.

nparticles_tracked

Return the number of particles currently being tracked. That is, they do not have NaN position or velocity.

num_entered

Count the number of particles that have entered the grids. This number is calculated by summing the number of non-zero entries in the entered grid array.

num_grids

The number of grids specified at instantiation.

on_any_grid

Binary array for each particle indicating whether it is currently on ANY grid.

particles_on_grid

Returns a boolean mask of shape [ngrids, nparticles] corresponding to whether or not the particle is on the associated grid.

results_dict

A dictionary containing the results of the simulation.

Dictionary keys and descriptions.

Key

Type

Description

"source"

ndarray

The source location vector, in meters.

"detector"

ndarray

The detector location vector, in meters.

"mag"

float

The system magnification.

"nparticles"

int

Number of particles in the simulation.

"max_deflection"

ndarray

The maximum deflection experienced by a particle in the simulation, in radians.

"x"

ndarray, [nparticles,]

The x-coordinate location where each particle hit the detector plane, in meters.

"y"

ndarray, [nparticles,]

The y-coordinate location where each particle hit the detector plane, in meters.

"v"

ndarray, [nparticles, 3]

The velocity of each particle when it hits the detector plane, in meters per second. The velocity is in a coordinate system relative to the detector plane. The components are [normal, horizontal, vertical] relative to the detector plane coordinates.

"x0"

ndarray, [nparticles,]

The x-coordinate location where each particle would have hit the detector plane if the grid fields were zero, in meters. Useful for calculating the source profile.

"y0"

ndarray, [nparticles,]

The y-coordinate location where each particle would have hit the detector plane if the grid fields were zero, in meters. Useful for calculating the source profile.

"v0"

ndarray, [nparticles, 3]

The velocity of each particle when it hit the detector plan if the grid fields were zero, in meters per second. The velocity is in a coordinate system relative to the detector plane. The components are [normal, horizontal, vertical] relative to the detector plane coordinates.

vmax

The maximum velocity of any particle in the simulation.

This quantity is used for determining the grid crossing maximum time step.

Methods Documentation

add_stopping(
method: Literal['NIST', 'Bethe'],
materials: list[str] | None = None,
I: Quantity | None = None,
)

Enable particle stopping using experimental stopping powers.

Interpolation of stopping powers is conducted using data from the NIST PSTAR database. This information is combined with the mass density quantity provided in the grids to calculate the energy loss over the distance travelled during a timestep.

add_wire_mesh(
location,
extent,
nwires,
wire_diameter,
mesh_hdir=None,
mesh_vdir=None,
)[source]

Add a wire mesh grid between the particle source and the object grid that blocks particles whose paths intersect the wires.

Parameters:
  • location (Quantity, shape (3)) –

    A vector pointing from the origin of the grid to the center of the mesh grid. This location must be between the source and the object grid.

    This vector will be interpreted as being in either cartesian, cylindrical, or spherical coordinates based on its units. Valid geometries are:

    • Cartesian (x,y,z) : (meters, meters, meters)

    • Cylindrical (r, theta, z) : (meters, radians, meters)

    • Spherical (r, theta, phi) : (meters, radians, radians)

    In spherical coordinates theta is the polar angle.

  • extent (Tuple of 1 or 2 Quantity) – The size of the mesh grid (in the mesh plane). If one value is provided, the mesh is circular and the value provided is interpreted as the diameter. If two values are provided, the mesh is rectangular and the values are interpreted as the width and height respectively.

  • nwires (Tuple of 1 or 2 ints, or a single int) – The number of wires in the horizontal and vertical directions. If only one value is provided, the number in the two directions is assumed to be equal. Note that a wire will cross the center of the mesh only when nwires is odd.

  • wire_diameter (Quantity) – The diameter of the wires.

  • mesh_hdir (numpy.ndarray, shape (3), optional) – A unit vector (in Cartesian coordinates) defining the horizontal direction on the mesh plane. Modifying this vector can rotate the mesh in the plane or tilt the mesh plane relative to the source-detector axis. By default, mesh_hdir is set equal to detector_hdir (see detector_hdir keyword in __init__).

  • mesh_vdir (numpy.ndarray, shape (3), optional) – A unit vector (in Cartesian coordinates) defining the vertical direction on the mesh plane. Modifying this vector can tilt the mesh relative to the source-detector axis. By default, mesh_vdir is defined to be perpendicular to mesh_hdir and the detector plane normal (such that the mesh is parallel to the detector plane).

Raises:

ValueError – Raises a ValueError if the provided mesh location is not between the source and the object grid.

create_particles(
nparticles,
particle_energy,
max_theta=None,
particle: Particle = Particle('p+'),
distribution: Literal['monte-carlo', 'uniform'] = 'monte-carlo',
random_seed=None,
) None[source]

Generates the angular distributions about the Z-axis, then rotates those distributions to align with the source-to-detector axis.

By default, particles are generated over almost the entire pi/2. However, if the detector is far from the source, many of these particles will never be observed. The max_theta keyword allows these extraneous particles to be neglected to focus computational resources on the particles who will actually hit the detector.

Parameters:
  • nparticles (integer) – The number of particles to include in the simulation. The default is 1e5.

  • particle_energy (Quantity) – The energy of the particle, in units convertible to eV. All particles are given the same energy.

  • max_theta (Quantity, optional) – The largest velocity vector angle (measured from the source-to-detector axis) for which particles should be generated. Decreasing this angle can eliminate particles that would never reach the detector region of interest. If no value is given, a guess will be made based on the size of the grid. Units must be convertible to radians.

  • particle (Particle or string representation of same, optional) – Representation of the particle species as either a Particle object or a string representation. The default particle is protons.

  • distribution (str) –

    A keyword which determines how particles will be distributed in velocity space. Options are:

    • ’monte-carlo’: velocities will be chosen randomly,

      such that the flux per solid angle is uniform.

    • ’uniform’: velocities will be distributed such that,

      left unperturbed,they will form a uniform pattern on the detection plane. This method requires that nparticles be a perfect square. If it is not, nparticles will be set as the largest perfect square smaller than the provided nparticles.

    Simulations run in the 'uniform' mode will imprint a grid pattern on the image, but will well-sample the field grid with a smaller number of particles. The default is 'monte-carlo'.

  • random_seed (int, optional) – A random seed to be used when generating random particle distributions, e.g. with the monte-carlo distribution.

load_particles(
x,
v,
particle: Particle = Particle('p+'),
) None[source]

Load arrays of particle positions and velocities.

Parameters:
  • x (Quantity, shape (N,3)) – Positions for N particles

  • v (Quantity, shape (N,3)) – Velocities for N particles

  • particle (particle-like, optional) – Representation of the particle species as either a Particle object or a string representation. The default particle is protons.

run() None[source]

Runs a particle-tracing simulation.

Timesteps are adaptively calculated based on the local grid resolution of the particles and the electric and magnetic fields they are experiencing. After all particles have left the grid, they are advanced to the detector plane where they can be used to construct a synthetic diagnostic image.

Return type:

None

setup_adaptive_time_step(
time_steps_per_gyroperiod: int | None = 12,
Courant_parameter: float | None = 0.5,
) None

Set parameters for the adaptive time step candidates.

Parameters:
  • time_steps_per_gyroperiod (int, optional) – The minimum ratio of the particle gyroperiod to the timestep. Higher numbers correspond to higher temporal resolution. The default is twelve.

  • Courant_parameter (float, optional) – The Courant parameter is the minimum ratio of the timestep to the grid crossing time, grid cell length / particle velocity. Lower Courant numbers correspond to higher temporal resolution.

Notes

Two candidates are calculated for the adaptive time step: a time step based on the gyroradius of the particle and a time step based on the resolution of the grid. The candidate associated with the gyroradius of the particle takes a time_steps_per_gyroperiod parameter that specifies how many times the orbit of a gyrating particles will be subdivided. The other candidate, associated with the spatial resolution of the grid object, calculates a time step using the time it would take the fastest particle to cross some fraction of a grid cell length. This fraction is the Courant number.