pymkm.physics.specific_energy

Computation of specific energy quantities for MKM and SMK models.

This module defines the SpecificEnergy to compute: - Single-event specific energy z₁(b) - Saturation-corrected z′₁(b) using z₀ (square-root or quadratic) - Dose-averaged values (z̄, z̄′)

The sensitive region is modeled as a water cylinder perpendicular to the ion path. Radial dose profiles are integrated using impact parameter b to match MKM/SMK formalisms.

Example usage:

from pymkm.physics.particle_track import ParticleTrack
from pymkm.dosimetry.specific_energy import SpecificEnergy

# Load or generate a ParticleTrack instance
track = ParticleTrack(...)  # contains D(r) and penumbra radius

# Define geometry
region_radius = 0.5  # in micrometers
sz = SpecificEnergy(track, region_radius)

# Compute single-event specific energy z1(b)
z1_array, b_array = sz.single_event_specific_energy()

# Compute saturation parameter z0 from beta0
z0 = sz.compute_saturation_parameter(domain_radius=0.5, nucleus_radius=3.0, beta0=0.05)

# Apply saturation correction
z1_prime_array = sz.saturation_corrected_single_event_specific_energy(z0, z1_array)

# Compute dose-averaged specific energy (corrected)
z_prime_bar = sz.dose_averaged_specific_energy(z1_array, b_array, z1_prime_array, model="square_root")
class pymkm.physics.specific_energy.SpecificEnergy(particle_track: ParticleTrack, region_radius: float)[source]

Bases: object

Compute microdosimetric specific energy quantities from a single ion track.

Supports MKM/SMK calculations of z₁(b), z′₁(b), z₀, and dose-averaged values. The sensitive region is modeled as a cylinder perpendicular to the particle track.

static compute_saturation_parameter(domain_radius: float, nucleus_radius: float, beta0: float) → float[source]

Compute the saturation parameter z₀ for MKM overkill correction.

Parameters:

  • domain_radius (float) – Radius of the sensitive domain (e.g. sub-nuclear volume), in μm.

  • nucleus_radius (float) – Radius of the cell nucleus, in μm.

  • beta0 (float) – β₀ coefficient of the LQ model at low LET (in Gy⁻²).

Returns:

Saturation parameter z₀ in Gy.

Return type:

float

Raises:

ValueError – If any input is not strictly positive.

dose_averaged_specific_energy(z_array: ndarray, b_array: ndarray, z_corrected: ndarray | None = None, model: str = 'square_root', integration_method: str = 'trapz') → float[source]

Compute dose-averaged specific energy: z̄ or z̄′.

If no saturation correction is applied, computes:

z̄ = [ ∫ z₁(b)² * b db ] / [ ∫ z₁(b) * b db ]

If corrected values z′₁(b) are provided:

  • Square-root model: uses [z′₁(b)]² in the numerator

  • Quadratic model: uses z₁(b) or z₁_sat(b) in the numerator

Parameters:

  • z_array (np.ndarray) – Uncorrected single-event specific energy values z₁(b) [Gy].

  • b_array (np.ndarray) – Impact parameter values [μm], must be sorted.

  • z_corrected (Optional[np.ndarray]) – Optional corrected values (e.g. z′₁(b) or z₁_sat(b)).

  • model (str) – Saturation model used if z_corrected is given (‘square_root’ or ‘quadratic’).

  • integration_method (str) – Integration rule to use: ‘trapz’, ‘simps’, or ‘quad’.

Returns:

Dose-averaged specific energy [Gy].

Return type:

float

Raises:

ValueError – If model or integration method is invalid.

saturation_corrected_single_event_specific_energy(z0: float, z_array: ndarray, model: str = 'square_root') → ndarray[source]

Compute z′₁(b): saturation-corrected specific energy from z₁(b) using z₀.

Applies overkill correction to single-event specific energy using either:

  • Square-root model: z′₁(b) = z₀ · sqrt(1 - exp(-[z₁(b)/z₀]²))

  • Quadratic model: z₁_sat(b) = [z′₁(b)]² / z₁(b)

Parameters:

  • z0 (float) – Saturation parameter z₀ (Gy).

  • z_array (np.ndarray) – Array of uncorrected z₁(b) values (Gy).

  • model (str) – Saturation correction model to apply: ‘square_root’ or ‘quadratic’.

Returns:

Array of corrected specific energy values: z′₁(b) or z₁_sat(b) (Gy).

Return type:

np.ndarray

Raises:

ValueError – If the model is not supported.

single_event_specific_energy(impact_parameters: ndarray | None = None, base_points_b: int | None = None, base_points_r: int | None = None, parallel: bool = False, return_time: bool = False) → Tuple[ndarray, ndarray] | Tuple[ndarray, ndarray, float][source]

Compute z₁(b): single-event specific energy as a function of impact parameter.

Parameters:

  • impact_parameters (Optional[np.ndarray]) – Optional array of impact parameters b (in μm).

  • base_points_b (Optional[int]) – Number of sampling points for the b grid.

  • base_points_r (Optional[int]) – Number of radial integration points per b.

  • parallel (bool) – If True, evaluates all b values in parallel using threads.

  • return_time (bool) – If True, returns the computation time in seconds.

Returns:

  • If return_time is False: tuple (z_array, b_array)

  • If return_time is True: tuple (z_array, b_array, elapsed_time)

where: - z_array: specific energy per event at each b [Gy] - b_array: impact parameter values [μm] - elapsed_time: wall-clock time in seconds [s]

Return type:

Union[Tuple[np.ndarray, np.ndarray], Tuple[np.ndarray, np.ndarray, float]]

Raises:

ValueError – If impact_parameters are invalid.