PET projectors `parallelproj.projectors`¶

high-level geometrical forward and back projectors

class parallelproj.projectors.EqualBlockPETProjector(lor_descriptor: EqualBlockPETLORDescriptor, img_shape: tuple[int, int, int], voxel_size: tuple[float, float, float], img_origin: None | Array = None, num_chunks: int = 1)[source]¶

Bases: LinearOperator

geometric non-TOF and TOF sinogram projector for equal block PET scanners

Examples

Non-TOF and TOF projections using a modularized (block) PET scanner geometry

Parameters:

lor_descriptor (EqualBlockPETLORDescriptor) – descriptor of the LOR start / end points
img_shape (tuple[int, int, int]) – shape of the image to be projected
voxel_size (tuple[float, float, float]) – the voxel size of the image to be projected
img_origin (None | Array, optional) – the origin of the image to be projected, by default None means that the center of the image is at world coordinate (0,0,0)
num_chunks (int, optional) – number of chunks to split the block pairs into during projection, by default 1 (all block pairs processed in a single call). Increase this value to reduce peak memory usage at the cost of more projection kernel calls.

Examples

Non-TOF and TOF projections using a modularized (block) PET scanner geometry

property H: AdjointLinearOperator¶: adjoint operator \(A^H\)

__call__(x: Array) → Array¶

alias to apply(x)

Examples

Parameters:: x (Array)
Return type:: Array

adjoint(y: Array) → Array¶

(scaled) adjoint step \(x = \overline{\alpha} A^H y\)

Parameters:: y (Array)
Return type:: Array

Examples

Non-TOF and TOF projections using a modularized (block) PET scanner geometry

adjointness_test(xp: ModuleType, dev: str, verbose: bool = False, iscomplex: bool = False, dtype: type | None = None, **kwargs) → bool¶

test whether the adjoint is correctly implemented

Parameters:

xp (ModuleType) – array module to use
dev (str) – device (cpu or cuda)
verbose (bool, optional) – verbose output
iscomplex (bool, optional) – use complex arrays
dtype (type | None, optional) – data type of the arrays
**kwargs (dict) – passed to np.isclose

Returns:

whether the adjoint is correctly implemented

Return type:

bool

Examples

Non-TOF and TOF projections using a modularized (block) PET scanner geometry

apply(x: Array) → Array¶

(scaled) forward step \(y = \alpha A x\)

Parameters:: x (Array)
Return type:: Array

Examples

property dev: str¶: device

property img_origin: Array¶: image origin - world coordinates of the [0,0,0] voxel

property in_shape: tuple[int, int, int]¶: shape of the input array

property lor_descriptor: EqualBlockPETLORDescriptor¶: LOR descriptor

norm(xp: ModuleType, dev: str, num_iter: int = 30, iscomplex: bool = False, verbose: bool = False) → float¶

estimate norm of the linear operator using power iterations

Parameters:

xp (ModuleType) – array module to use
dev (str) – device (cpu or cuda)
num_iter (int, optional) – number of power iterations
iscomplex (bool, optional) – use complex arrays
verbose (bool, optional) – verbose output

Returns:

the norm of the linear operator

Return type:

float

Examples

property num_chunks: int¶: number of chunks to split the block pairs into during projection

property out_shape: tuple[int, ...]¶: shape of the output array

property scale: float | complex¶: scalar factor applied to the linear operator

show_geometry(ax: Axes3D, color: tuple[float, float, float] = (1.0, 0.0, 0.0), edgecolor: str = 'grey', alpha: float = 0.1) → None[source]¶

show the geometry of the scanner and the FOV of the image

Parameters:

ax (Axes3D) – matplotlib axes object with projection = ‘3d’
color (tuple[float, float, float], optional) – color to use for the FOV cube, by default (1.,0.,0.)
edgecolor (str, optional) – edgecolor to use for the FOV cube, by default ‘grey’
alpha (float, optional) – alpha value of the FOV cube, by default 0.1

Return type:

None

Examples

Non-TOF and TOF projections using a modularized (block) PET scanner geometry

property tof: bool¶: bool indicating whether to use TOF or not

property tof_parameters: TOFParameters | None¶: TOF parameters

property voxel_size: Array¶: voxel size

property xp: ModuleType¶: array module

class parallelproj.projectors.ListmodePETProjector(event_start_coordinates: Array, event_end_coordinates: Array, img_shape: tuple[int, int, int], voxel_size: tuple[float, float, float], img_origin: None | Array = None)[source]¶

Bases: LinearOperator

non-TOF and TOF listmode projector for regular polygon PET scanners

Examples

PET listmode projector (non-TOF and TOF)

Listmode MLEM, OSEM, and SVRG

Convergence comparison: SGD vs SVRG with regularization (sinogram and listmode)

PDHG and LM-SPHG to optimize the Poisson logL and total variation

Parameters:

event_start_coordinates (Array) – float world coordinates of event LOR start points, shape (num_events, 3)
event_end_coordinates (Array) – float world coordinates of event LOR end points, shape (num_events, 3)
img_shape (tuple[int, int, int]) – shape of the image to be projected
voxel_size (tuple[float, float, float]) – the voxel size of the image to be projected
img_origin (None | Array, optional) – the origin of the image to be projected, by default None means that the center of the image is at world coordinate (0,0,0)

Examples

PET listmode projector (non-TOF and TOF)

Listmode MLEM, OSEM, and SVRG

Convergence comparison: SGD vs SVRG with regularization (sinogram and listmode)

PDHG and LM-SPHG to optimize the Poisson logL and total variation

property H: AdjointLinearOperator¶: adjoint operator \(A^H\)

__call__(x: Array) → Array¶

alias to apply(x)

Examples

Parameters:: x (Array)
Return type:: Array

adjoint(y: Array) → Array¶

(scaled) adjoint step \(x = \overline{\alpha} A^H y\)

Parameters:: y (Array)
Return type:: Array

Examples

PET listmode projector (non-TOF and TOF)

adjointness_test(xp: ModuleType, dev: str, verbose: bool = False, iscomplex: bool = False, dtype: type | None = None, **kwargs) → bool¶

test whether the adjoint is correctly implemented

Parameters:

xp (ModuleType) – array module to use
dev (str) – device (cpu or cuda)
verbose (bool, optional) – verbose output
iscomplex (bool, optional) – use complex arrays
dtype (type | None, optional) – data type of the arrays
**kwargs (dict) – passed to np.isclose

Returns:

whether the adjoint is correctly implemented

Return type:

bool

Examples

apply(x: Array) → Array¶

(scaled) forward step \(y = \alpha A x\)

Parameters:: x (Array)
Return type:: Array

Examples

property event_end_coordinates: Array¶: coordinates of LOR end points

property event_start_coordinates: Array¶: coordinates of LOR start points

property event_tofbins: None | Array¶: TOF bin of each event

property in_shape: tuple[int, int, int]¶: shape of the input array

norm(xp: ModuleType, dev: str, num_iter: int = 30, iscomplex: bool = False, verbose: bool = False) → float¶

estimate norm of the linear operator using power iterations

Parameters:

xp (ModuleType) – array module to use
dev (str) – device (cpu or cuda)
num_iter (int, optional) – number of power iterations
iscomplex (bool, optional) – use complex arrays
verbose (bool, optional) – verbose output

Returns:

the norm of the linear operator

Return type:

float

Examples

property num_events: int¶: number of events

property out_shape: tuple[int, ...]¶: shape of the output array

property scale: float | complex¶: scalar factor applied to the linear operator

property tof: bool¶: bool indicating whether to use TOF projections or not

property tof_parameters: TOFParameters | None¶: TOF parameters

property voxel_size: Array¶: voxel size

property xp: ModuleType¶: array module

class parallelproj.projectors.ParallelViewProjector2D(image_shape: tuple[int, int], radial_positions: Array, view_angles: Array, radius: float, image_origin: tuple[float, float], voxel_size: tuple[float, float])[source]¶

Bases: LinearOperator

2D non-TOF parallel view projector

Examples

2D non-TOF filtered back projection (FBP) of Poisson data

init method

Parameters:

image_shape (tuple[int, int]) – shape of the input image (n1, n2)
radial_positions (Array) – radial positions of the projection views in world coordinates
view_angles (Array) – angles of the projection views in radians
radius (float) – radius of the scanner
image_origin (tuple[float, float]) – world coordinates of the [0,0] voxel
voxel_size (tuple[float, float]) – the voxel size in both directions

Examples

2D non-TOF filtered back projection (FBP) of Poisson data

property H: AdjointLinearOperator¶: adjoint operator \(A^H\)

__call__(x: Array) → Array¶

alias to apply(x)

Examples

Parameters:: x (Array)
Return type:: Array

adjoint(y: Array) → Array¶

(scaled) adjoint step \(x = \overline{\alpha} A^H y\)

Parameters:: y (Array)
Return type:: Array

Examples

2D non-TOF filtered back projection (FBP) of Poisson data

adjointness_test(xp: ModuleType, dev: str, verbose: bool = False, iscomplex: bool = False, dtype: type | None = None, **kwargs) → bool¶

test whether the adjoint is correctly implemented

Parameters:

xp (ModuleType) – array module to use
dev (str) – device (cpu or cuda)
verbose (bool, optional) – verbose output
iscomplex (bool, optional) – use complex arrays
dtype (type | None, optional) – data type of the arrays
**kwargs (dict) – passed to np.isclose

Returns:

whether the adjoint is correctly implemented

Return type:

bool

Examples

apply(x: Array) → Array¶

(scaled) forward step \(y = \alpha A x\)

Parameters:: x (Array)
Return type:: Array

Examples

property dev: str¶: device used for storage of LOR endpoints

property image_origin: Array¶: image origin - world coordinates of the [0,0] voxel

property image_shape: tuple[int, int]¶: image shape

property in_shape: tuple[int, ...]¶: shape of the input array

norm(xp: ModuleType, dev: str, num_iter: int = 30, iscomplex: bool = False, verbose: bool = False) → float¶

estimate norm of the linear operator using power iterations

Parameters:

xp (ModuleType) – array module to use
dev (str) – device (cpu or cuda)
num_iter (int, optional) – number of power iterations
iscomplex (bool, optional) – use complex arrays
verbose (bool, optional) – verbose output

Returns:

the norm of the linear operator

Return type:

float

Examples

property num_rad: int¶: number of radial elements

property num_views: int¶: number of views

property out_shape: tuple[int, ...]¶: shape of the output array

property scale: float | complex¶: scalar factor applied to the linear operator

show_views(views_to_show: None | Array = None, image: None | Array = None, **kwargs: Any) → Figure[source]¶

visualize the geometry of certrain projection views

Parameters:

views_to_show (None | Array) – view numbers to show
image (None | Array) – show an image inside the projector geometry
**kwargs (dict) – passed to matplotlib.pyplot.imshow

Returns:

the matplotlib figure

Return type:

plt.Figure

Examples

property voxel_size: Array¶: voxel size

property xend: Array¶: coordinates of LOR end points

property xp: ModuleType¶: array module

property xstart: Array¶: coordinates of LOR start points

class parallelproj.projectors.ParallelViewProjector3D(image_shape: tuple[int, int, int], radial_positions: Array, view_angles: Array, radius: float, image_origin: tuple[float, float, float], voxel_size: tuple[float, float, float], ring_positions: Array, span: int = 1, max_ring_diff: int | None = None)[source]¶

Bases: LinearOperator

3D non-TOF parallel view projector

Examples

init method

Parameters:

image_shape (tuple[int, int, int]) – shape of the input image (n0, n1, n2) (last direction is axial)
radial_positions (Array) – radial positions of the projection views in world coordinates
view_angles (Array) – angles of the projection views in radians
radius (float) – radius of the scanner
image_origin (tuple[float, float, float]) – world coordinates of the [0,0,0] voxel
voxel_size (tuple[float, float, float]) – the voxel size in all three directions (n0, n1, axial)
ring_positions (Array) – position of the rings in world coordinates
span (int) – span of the sinogram - default is 1
max_ring_diff (int | None) – maximum ring difference - default is None (no limit)

Examples

property H: AdjointLinearOperator¶: adjoint operator \(A^H\)

__call__(x: Array) → Array¶

alias to apply(x)

Examples

Parameters:: x (Array)
Return type:: Array

adjoint(y: Array) → Array¶

(scaled) adjoint step \(x = \overline{\alpha} A^H y\)

Parameters:: y (Array)
Return type:: Array

Examples

adjointness_test(xp: ModuleType, dev: str, verbose: bool = False, iscomplex: bool = False, dtype: type | None = None, **kwargs) → bool¶

test whether the adjoint is correctly implemented

Parameters:

xp (ModuleType) – array module to use
dev (str) – device (cpu or cuda)
verbose (bool, optional) – verbose output
iscomplex (bool, optional) – use complex arrays
dtype (type | None, optional) – data type of the arrays
**kwargs (dict) – passed to np.isclose

Returns:

whether the adjoint is correctly implemented

Return type:

bool

Examples

apply(x: Array) → Array¶

(scaled) forward step \(y = \alpha A x\)

Parameters:: x (Array)
Return type:: Array

Examples

property image_origin: Array¶: image origin - world coordinates of the [0,0,0] voxel

property image_shape: tuple[int, int, int]¶: image shape

property in_shape: tuple[int, int, int]¶: shape of the input array

property max_ring_diff: int¶: maximum ring difference

norm(xp: ModuleType, dev: str, num_iter: int = 30, iscomplex: bool = False, verbose: bool = False) → float¶

estimate norm of the linear operator using power iterations

Parameters:

xp (ModuleType) – array module to use
dev (str) – device (cpu or cuda)
num_iter (int, optional) – number of power iterations
iscomplex (bool, optional) – use complex arrays
verbose (bool, optional) – verbose output

Returns:

the norm of the linear operator

Return type:

float

Examples

property out_shape: tuple[int, int, int]¶: shape of the output array

property scale: float | complex¶: scalar factor applied to the linear operator

property voxel_size: Array¶: the voxel size in all directions

property xend: Array¶: coordinates of LOR end points

property xp: ModuleType¶: array module

property xstart: Array¶: coordinates of LOR start points

class parallelproj.projectors.RegularPolygonPETProjector(lor_descriptor: RegularPolygonPETLORDescriptor, img_shape: tuple[int, int, int], voxel_size: tuple[float, float, float], img_origin: None | Array = None, views: None | Array = None, cache_lor_endpoints: bool = True)[source]¶

Bases: LinearOperator

geometric non-TOF and TOF sinogram projector for regular polygon PET scanners

Examples

Michelograms and axial sinogram compression

PET non-TOF sinogram projector

PET TOF sinogram projector

Convergence comparison: MLEM vs OSEM vs SVRG

Convergence comparison: SGD vs SVRG with quadratic regularization

PDHG and SPDHG for PET reconstruction with a directional TV prior

Listmode MLEM, OSEM, and SVRG

Convergence comparison: SGD vs SVRG with regularization (sinogram and listmode)

PDHG and LM-SPHG to optimize the Poisson logL and total variation

pytorch parallelproj projection layer

Parameters:

lor_descriptor (RegularPolygonPETLORDescriptor) – descriptor of the LOR start / end points
img_shape (tuple[int, int, int]) – shape of the image to be projected
voxel_size (tuple[float, float, float]) – the voxel size of the image to be projected
img_origin (None | Array, optional) – the origin of the image to be projected, by default None means that the center of the image is at world coordinate (0,0,0)
views (None | Array, optional) – sinogram views to be projected, by default None means that all views are being projected
cache_lor_endpoints (bool, optional) – whether to cache the LOR endpoints, by default True setting it to False will save memory but will slow down computations

Examples

Michelograms and axial sinogram compression

PET non-TOF sinogram projector

PET TOF sinogram projector

Convergence comparison: MLEM vs OSEM vs SVRG

Convergence comparison: SGD vs SVRG with quadratic regularization

PDHG and SPDHG for PET reconstruction with a directional TV prior

Listmode MLEM, OSEM, and SVRG

Convergence comparison: SGD vs SVRG with regularization (sinogram and listmode)

PDHG and LM-SPHG to optimize the Poisson logL and total variation

pytorch parallelproj projection layer

property H: AdjointLinearOperator¶: adjoint operator \(A^H\)

__call__(x: Array) → Array¶

alias to apply(x)

Examples

Parameters:: x (Array)
Return type:: Array

adjoint(y: Array) → Array¶

(scaled) adjoint step \(x = \overline{\alpha} A^H y\)

Parameters:: y (Array)
Return type:: Array

Examples

PET non-TOF sinogram projector

adjointness_test(xp: ModuleType, dev: str, verbose: bool = False, iscomplex: bool = False, dtype: type | None = None, **kwargs) → bool¶

test whether the adjoint is correctly implemented

Parameters:

xp (ModuleType) – array module to use
dev (str) – device (cpu or cuda)
verbose (bool, optional) – verbose output
iscomplex (bool, optional) – use complex arrays
dtype (type | None, optional) – data type of the arrays
**kwargs (dict) – passed to np.isclose

Returns:

whether the adjoint is correctly implemented

Return type:

bool

Examples

apply(x: Array) → Array¶

(scaled) forward step \(y = \alpha A x\)

Parameters:: x (Array)
Return type:: Array

Examples

clear_cached_lor_endpoints() → None[source]¶

clear cached LOR endpoints

Examples

Convergence comparison: MLEM vs OSEM vs SVRG

Convergence comparison: SGD vs SVRG with quadratic regularization

PDHG and SPDHG for PET reconstruction with a directional TV prior

Listmode MLEM, OSEM, and SVRG

Convergence comparison: SGD vs SVRG with regularization (sinogram and listmode)

Return type:: None

convert_sinogram_to_listmode(sinogram: Array, shuffle: bool = False) → tuple[Array, Array, Array | None][source]¶

convert a non-TOF or TOF emission sinogram to listmode events

Parameters:

sinogram (Array) – an integer (TOF or non-TOF) emission sinogram
shuffle (bool, optional) – if True, randomly shuffle the order of the output events, by default False. Shuffling is implemented via numpy.random.permutation(num_events), which draws from numpy’s global random state. Use numpy.random.seed() before calling this method for reproducible results.

Returns:

event_start_coordinates, event_end_coordinates, event_tofbin in case of non-TOF, event_tofbin is None

Return type:

tuple[Array, Array, Array | None]

Examples

Listmode MLEM, OSEM, and SVRG

Convergence comparison: SGD vs SVRG with regularization (sinogram and listmode)

property img_origin: Array¶: image origin - world coordinates of the [0,0,0] voxel

property in_shape: tuple[int, int, int]¶: shape of the input array

property lor_descriptor: RegularPolygonPETLORDescriptor¶: LOR descriptor

norm(xp: ModuleType, dev: str, num_iter: int = 30, iscomplex: bool = False, verbose: bool = False) → float¶

estimate norm of the linear operator using power iterations

Parameters:

xp (ModuleType) – array module to use
dev (str) – device (cpu or cuda)
num_iter (int, optional) – number of power iterations
iscomplex (bool, optional) – use complex arrays
verbose (bool, optional) – verbose output

Returns:

the norm of the linear operator

Return type:

float

Examples

property out_shape: tuple[int, ...]¶: shape of the output array

property scale: float | complex¶: scalar factor applied to the linear operator

show_geometry(ax: Axes3D, color: tuple[float, float, float] = (1.0, 0.0, 0.0), edgecolor: str = 'grey', alpha: float = 0.1) → None[source]¶

show the geometry of the scanner and the FOV of the image

Parameters:

ax (Axes3D) – matplotlib axes object with projection = ‘3d’
color (tuple[float, float, float], optional) – color to use for the FOV cube, by default (1.,0.,0.)
edgecolor (str, optional) – edgecolor to use for the FOV cube, by default ‘grey’
alpha (float, optional) – alpha value of the FOV cube, by default 0.1

Return type: