Quantity System#

Radiotherapy optimization requires more than just physical dose. RBE-weighted dose, biological effect, and LET-weighted dose all depend on combinations of raw dose matrices stored in the Dij object. pyRadPlan’s quantity system manages these dependencies automatically, resolving which matrices need to be combined and caching intermediate results.

Motivation#

Consider RBE-weighted dose in the LQ model:

\[d_\text{RBE} = \frac{-\alpha_x + \sqrt{\alpha_x^2 + 4\beta_x \cdot E}}{2\beta_x}\]

where \(E = \alpha \cdot w + \beta \cdot w^2\) is the biological effect and \(\alpha\), \(\beta\) are dose-influence matrices stored inside dij. Computing this by hand for every optimizer iteration — including correct gradient propagation — is error-prone and backend-specific. The quantity system handles all of this transparently.

Available quantities#

Class	Identifier	Description
`Dose`	`"physical_dose"`	Physical absorbed dose in Gy.
`DoseWeightedLET`	`"let"`	Dose-averaged LET, derived from LET-weighted dose and physical dose.
`LETxDose`	`"let_dose"`	LET-weighted dose matrix.
`AlphaDose`	`"alpha_dose"`	Alpha-weighted dose for LQ model.
`SqrtBetaDose`	`"sqrt_beta_dose"`	Square-root-beta-weighted dose for LQ model.
`Effect`	`"effect"`	Biological effect \(E = \alpha D + \beta D^2\) (LQ model).
`RBExDose`	`"rbe_x_dose"`	RBE-weighted dose in Gy(RBE). Derived from `effect`, `alpha_dose`, `sqrt_beta_dose`.

Selecting a quantity for optimization#

Objectives reference a quantity by identifier. The identifier is resolved against the quantity registry and the available matrices in dij:

from pyRadPlan.optimization.objectives import SquaredDeviation

# Optimize on physical dose (default)
obj = SquaredDeviation(priority=1000, d_ref=60.0, quantity="physical_dose")

# Optimize on RBE-weighted dose (requires alpha_dose and sqrt_beta_dose in dij)
obj = SquaredDeviation(priority=1000, d_ref=60.0, quantity="rbe_x_dose")

Resolution modes#

Each FluenceDependentQuantity checks whether its corresponding matrix is already present in dij ("direct" mode) or must be derived from other quantities ("indirect" mode):

Direct mode: The quantity matrix (e.g. dij.physical_dose) exists and can be multiplied by the fluence vector directly. Fastest path.
Indirect mode: The quantity is not pre-computed but can be derived from other available matrices (e.g. rbe_x_dose is computed from alpha_dose, sqrt_beta_dose, and reference tissue parameters α_x / β_x). The resolver chains quantity computations in the correct dependency order.

The QuantityResolver#

QuantityResolver is the component that binds a set of quantities to a specific dij. During optimization, the planning problem instantiates a resolver and calls it by identifier:

from pyRadPlan.quantities import QuantityResolver

resolver = QuantityResolver(dij)

# Compute physical dose for a given fluence vector
dose = resolver.get("physical_dose")(fluence)

# Compute RBE-weighted dose (resolver handles the dependency chain automatically)
rbe_dose = resolver.get("rbe_x_dose")(fluence)

Gradients and the Array API#

Quantity objects also expose compute_chain_derivative() for propagating derivatives back to the fluence vector. Both forward evaluation and derivative propagation follow the Python Array API standard, so operations can run on NumPy, CuPy, or PyTorch arrays depending on the configured backend and installed optional dependencies.

Adding custom quantities#

New quantities can be added by subclassing FluenceDependentQuantity and registering the class in the pyRadPlan.quantities.QUANTITIES dictionary. A custom quantity must implement the single-scenario forward calculation and chain derivative methods used by the optimization problem; use the existing quantity classes in pyRadPlan/quantities as templates.