Proton Optimisation Tools

A collection of functions useful for proton treatment plan optimisation, implemented in the fredtools.ProtonOptimisation subpackage. The functions are not exported to the top-level fredtools namespace and must be called with the full module path, e.g. fredtools.ProtonOptimisation.preoptimizer.convertCTtoWER.

Pre-optimisation

fredtools.ProtonOptimisation.preoptimizer.convertCTtoWER(img: Image, HU: ArrayLike, WER: ArrayLike, displayInfo: bool = False) → Image

Convert CT map to WER map.

The function converts a 3D Computed Tomography (CT) map with Hounsfield Unit (HU) values, defined as a SimpleITK image object, to an image with Water-Equivalent Ratio values (WER). The two parameters, HU and WER, define the HU to WER conversion, whereas the missing HU values are interpolated linearly.

Parameters:

img (SimpleITK Image) – An object of a SimpleITK image with HU values.
HU (array_like) – An iterable with HU values. It must be of the same size as WER.
WER (array_like) – An iterable with WER values. It must be of the same size as HU.
displayInfo (bool, optional) – Displays a summary of the function results. (def. False)

Returns:

An instance of a SimpleITK image object with WER values.

Return type:

SimpleITK Image

Raises:

ValueError – If the lengths of HU and WER differ, or if the HU values of the image are outside the range of the HU table.

See also

calcWETfromWER: calculate WET image from WER image for point-like source.

fredtools.ProtonOptimisation.preoptimizer.calcWETfromWER(imgWER: Image, SAD: Iterable[int | float | number], imgMask: Image | None = None, displayInfo: bool = False) → Image

Calculate WET image from WER image for point-like source.

The function calculates Water-Equivalent Thickness (WET) for each voxel of an image defined as a SimpleITK image object containing Water-Equivalent Ratio values (WER), inside a mask, defined as a SimpleITK image object describing a binary mask. The WET values are calculated starting from a virtual source, located at [X,Y]=[0,0] position and with Z position defined with a two-element SAD, describing the source point in X and Y. Particularly, the WET is calculated for a virtual source, where rays are deflected in X and Y directions in different distances from the isocenter.

Parameters:

imgWER (SimpleITK Image) – An object of a 3D SimpleITK image with WER values. It must have the same frame of reference as imgMask.
SAD (2-element array_like) – Z coordinates of the virtual point source for deflection in X and Y directions, respectively.
imgMask (SimpleITK Image or None, optional) – An object of a 3D SimpleITK image describing a binary mask, with the same frame of reference as imgWER, or None, then all voxel positions will be calculated (def. None)
displayInfo (bool, optional) – Displays a summary of the function results. (def. False)

Returns:

An instance of a SimpleITK image object with WET values. Voxels outside the mask are set to numpy.nan.

Return type:

SimpleITK Image

Raises:

ValueError – If a parallel beam is requested, which has not been implemented yet.
AttributeError – If imgWER and imgMask do not have the same frame of reference, or if SAD is not a 2-element iterable.

See also

convertCTtoWER: convert CT map to WER map.

fredtools.ProtonOptimisation.preoptimizer.generateIsoLayers(minRange: int | float | number, maxRange: int | float | number, beamParams: DataFrame) → DataFrame

Calculate iso-WET layers and corresponding energies.

The function calculates iso Water-Equivalent Thickness (WET) layers between minimum and maximum range based on predefined parameters of the beam.

Parameters:

minRange (scalar) – Minimum range in [mm] to calculate layers.
maxRange (scalar) – Maximum range in [mm] to calculate layers.
beamParams (pandas.DataFrame) – Parameters of the beam, i.e. dependence of the beam range and width with nominal energies. Must include at least columns: “nomEnergy” (in [MeV]), “rangeProx” and “rangeDist” (in [mm]).

Returns:

An instance of pandas.DataFrame object describing the iso WET layers, with the ranges in [mm] and the nominal energy (“nomEnergy”) in [MeV].

Return type:

pandas.DataFrame

Raises:

ValueError – If beamParams does not contain the required columns, or if minRange or maxRange is outside the range covered by beamParams.

Notes

The layer generation loop is capped at 100 iterations, so at most 100 layers are generated even if the proximal range has not reached minRange.

fredtools.ProtonOptimisation.preoptimizer.calcContours(imgMask: Image, level: int | float | number = 0.5, displayInfo: bool = False) → List[Polygon]

Calculate contours from 2D binary mask.

The function calculates list of contours from a 2D image defined as a SimpleITK image object describing a binary or floating mask, along a level value.

Parameters:

imgMask (SimpleITK Image) – Object of a SimpleITK 2D image describing a binary mask.
level (scalar, optional) – Value along which to find contours in the image. (def. 0.5)
displayInfo (bool, optional) – Displays a summary of the function results. (def. False)

Returns:

A list of polygons defined as instances of shapely.polygon objects.

Return type:

list of shapely.polygon

Notes

The imgMask must be a 2D image describing a floating or binary mask. If an image is, for instance, a 3D image describing a slice, it must be squeezed to 2D first. It can be done, for instance by slicing the image, e.g. imgMask[:,:,0] for image where the third axis is single sized.

fredtools.ProtonOptimisation.preoptimizer.convertRayTargetToIsoPlane(rayTarget: ndarray[tuple[Any, ...], dtype[_ScalarT]], SAD: Iterable[int | float | number]) → ndarray[tuple[Any, ...], dtype[_ScalarT]]

Calculate beam position in the isocentre plane.

The function calculates the beam positions in the isocentre plane, based on the target position and distance to the virtual point (SAD).

Parameters:

rayTarget (Nx3 array_like) – Target positions in the format of Nx3 iterable (a single 3-element point is also accepted).
SAD (2-element array_like) – Z coordinates of the virtual point source for deflection in X and Y directions, respectively.

Returns:

An Nx3 array with the ray positions in the isocentre plane.

Return type:

Nx3 numpy array

Notes

The function assumes that the beam goes along +Z direction.

Beam position optimisation

fredtools.ProtonOptimisation.optimiseBeamPositions.optimiseBeamPositions(contourPolygon: Polygon, spotDistance: int | float | number, algorithm: Literal['regular', 'reg', 'hexagonal', 'hex', 'concentric', 'con', 'delaunay', 'del'] = 'regular', **kwargs) → ndarray[tuple[Any, ...], dtype[_ScalarT]]

Calculate and optimise the beam positions in a contour.

The function calculates optimised beam positions in a contour defined as an instance of the shapely.Polygon object. Various optimisation algorithms are implemented. Refer to ‘See Also’ section to read more about each algorithm.

Parameters:

contourPolygon (shapely Polygon) – Object of the shapely.Polygon.
spotDistance (scalar) – The nominal spot distance, in the polygon coordinate unit (usually [mm]), to be used to optimise the beam positions. Depending on the algorithm, the distance between neighbouring spots does not have to be equal to this parameter. Therefore it describes only the nominal distance.
algorithm ({'regular', 'hexagonal', 'concentric', 'delaunay'}, optional) – Algorithm to be used to optimise the beam positions. The short aliases ‘reg’, ‘hex’, ‘con’ and ‘del’ are also accepted. The ‘concentric’ and ‘delaunay’ algorithms are not implemented yet and raise NotImplementedError. (def. ‘regular’)
**kwargs (keyword args, optional) – Additional parameters passed to the selected optimisation algorithm (whichever algorithm is chosen). Refer to the given algorithm routine for more description.

Returns:

Numpy array of size (Nx2) describing the position of N beams, where the first column is X and the second Y directions.

Return type:

numpy array (Nx2)

Raises:

TypeError – If contourPolygon is not an instance of a shapely Polygon.
ValueError – If algorithm cannot be recognised.
NotImplementedError – If the ‘concentric’ or ‘delaunay’ algorithm is requested, as they have not been implemented yet.

See also

optimiseBeamPositionsRegular: Optimise beam positions in a regular grid.
optimiseBeamPositionsHexagonal: Optimise beam positions in a hexagonal grid.

fredtools.ProtonOptimisation.optimiseBeamPositions.optimiseBeamPositionsRegular(contourPolygon: Polygon, spotDistance: int | float | number) → ndarray[tuple[Any, ...], dtype[_ScalarT]]

Calculate the beam positions using regular grid algorithm.

The function calculates beam positions in a contour defined as an instance of the shapely.Polygon object using regular grid algorithm. The algorithm is optimised to place the central beam position at the polygon centroid.

Parameters:

contourPolygon (shapely Polygon) – Object of the shapely.Polygon.
spotDistance (scalar) – The spot distance, in the polygon coordinate unit (usually [mm]), to be used to calculate regular grid beam positions.

Returns:

Numpy array of size (Nx2) describing the position of N beams, where the first column is X and the second Y directions.

Return type:

numpy array (Nx2)

See also

optimiseBeamPositions: Optimise beam positions using various algorithms.

Notes

The regular grid algorithm distributes the beams with the same spacing in X and Y directions. The grid size is calculated to fit the given contour polygon and is moved so that the central beam is at the polygon centroid. All the beam positions which are not inside the polygon, enlarged by a buffer of 0.2 (in the polygon coordinate unit) tolerance, are removed.

fredtools.ProtonOptimisation.optimiseBeamPositions.optimiseBeamPositionsHexagonal(contourPolygon: Polygon, spotDistance: int | float | number, direction: Literal['X', 'Y'] = 'X') → ndarray[tuple[Any, ...], dtype[_ScalarT]]

Calculate the beam positions using the hexagonal grid algorithm.

The function calculates beam positions in a contour defined as an instance of the shapely.Polygon object using hexagonal grid algorithm. The algorithm is optimised to place the central beam position at the polygon centroid.

Parameters:

contourPolygon (shapely Polygon) – Object of the shapely.Polygon.
spotDistance (scalar) – The spot distance, in the polygon coordinate unit (usually [mm]), to be used to calculate hexagonal grid beam positions.
direction ({'X', 'Y'}, optional) – The direction along which the beams should be shifted to create a hexagonal grid. This parameter can be used to align the hexagonal direction to the faster direction of the pencil beam scanning. (def. ‘X’)

Returns:

Numpy array of size (Nx2) describing the position of N beams, where the first column is X and the second Y directions.

Return type:

numpy array (Nx2)

Raises:

ValueError – If direction is not ‘X’ or ‘Y’.

See also

optimiseBeamPositions: Optimise beam positions using various algorithms.

Notes

The hexagonal grid algorithm distributes the beams with the same spacing in X (or in Y) and every second row (or column) of the beam positions is shifted by half of the spotDistance. The grid size is calculated to fit the given contour polygon and is moved so that the central beam is at the polygon centroid. All the beam positions which are not inside the polygon, enlarged by a buffer of 0.2 (in the polygon coordinate unit) tolerance, are removed.

The user can choose in which direction, X or Y, the hexagonal grid should be aligned. This might be important when optimising the beam positions for a given machine where the scanning is faster in one direction than in the other.

fredtools.ProtonOptimisation.optimiseBeamPositions.optimiseBeamPositionsConcentric(contourPolygon: Polygon, spotDistance: int | float | number, **kwargs)

Calculate the beam positions using the concentric rings algorithm.

The function is intended to calculate beam positions in a contour defined as an instance of the shapely.Polygon object using a concentric rings algorithm. The algorithm has not been implemented yet.

Parameters:

contourPolygon (shapely Polygon) – Object of the shapely.Polygon.
spotDistance (scalar) – The spot distance to be used to calculate the beam positions.

Raises:

NotImplementedError – The algorithm has not been implemented yet.

See also

optimiseBeamPositions: Optimise beam positions using various algorithms.

fredtools.ProtonOptimisation.optimiseBeamPositions.optimiseBeamPositionsDelaunay(contourPolygon: Polygon, spotDistance: int | float | number, **kwargs)

Calculate the beam positions using the Delaunay triangulation algorithm.

The function is intended to calculate beam positions in a contour defined as an instance of the shapely.Polygon object using a Delaunay triangulation algorithm. The algorithm has not been implemented yet.

Parameters:

contourPolygon (shapely Polygon) – Object of the shapely.Polygon.
spotDistance (scalar) – The spot distance to be used to calculate the beam positions.

Raises:

NotImplementedError – The algorithm has not been implemented yet.

See also

optimiseBeamPositions: Optimise beam positions using various algorithms.