Selecting the optimal projector =============================== | This page outlines which projector is more optimal for specific data types. You can find references to the projector types below: | Type 1: http://cit.fer.hr/index.php/CIT/article/view/2969 | Type 2: https://doi.org/10.1118/1.3501884 | Type 3: https://doi.org/10.1088/0031-9155/59/3/561 | Type 4: https://doi.org/10.1088/0031-9155/56/13/004 | Type 5: https://doi.org/10.1109/TCI.2017.2675705 | Type 6: https://doi.org/10.1109/23.552756 PET data -------- For PET data, projector_type = 3 is the most optimal one, although also the slowest to compute. This is a volume-of-intersection-based projector, where the voxels are approximated as spheres and the tube of intersection as a cylinder. The user can determine the size of the sphere with ``options.voxel_radius``, where a value of 1 means that the sphere is just large enough to fully contain the cubic voxel. Larger values have larger spheres, while smaller values lead to smaller spheres. The radius (mm) of the tube can be set with ``options.tube_radius``. An alternative to the above is projector_type = 2, which is an orthogonal distance-based ray-tracer. This is practically equal in speed to the above. It computes the orthogonal distance from the voxel where the ray currently is to the nearby voxels. If the distance is lower than the threshold value, then the distance is included, otherwise it is omitted. You can adjust this distance with ``options.tube_width_xy`` in the XY-plane, i.e. 2D, and with ``options.tube_width_z`` in 3D. If ``options.tube_width_z`` is set, then ``options.tube_width_xy`` is ignored. ``options.tube_width_z`` always operates in all three dimensions, while ``options.tube_width_xy`` only operates in 2D and is thus a faster alternative. The fastest methods are projector_type = 1 or projector_type = 4. The first one computes the exact line of intersection length in each voxel using an improved version of the Siddon's algorithm. The latter is an interpolation-based projector, similar to Joseph's method, i.e. it linearly interpolates values after a pre-determined step size. This step size can be adjusted with ``options.dL`` and is the relative size of one voxel. I.e. ``options.dL = 1`` would use the length of one voxel as the interpolation length. A multi-ray version is also available, though recommended only for projector_type = 1. The number of transaxial rays can be adjusted with ``options.n_rays_transaxial`` and axial rays with ``options.n_rays_axial``. The total number of rays is thus the multiplication of these two. While projector_types 1 and 4 are the fastest, they produce lower quality images than projector types 2 and 3. However, you can compensate for this somewhat by using point spread function (PSF) blurring. This can be enabled by setting ``options.use_psf`` to true. The FWHM can be adjusted with ``options.FWHM`` and must include values for all three dimensions. Note that projector type 4 is OpenCL/CUDA only! This means that only implementations 2, 3, and 5 support it (note that Python uses only implementation 2!). Note that you can also use hybrid methods if you wish, such as ``projector_type = 13``, but these are not recommended with PET data. The first value corresponds to the forward projection (1 in this case) and the second to the backprojection (3 in this case). Projector types 5 and 6 and their hybrids are NOT supported with PET data! PSF modeling ^^^^^^^^^^^^ As already mentioned above, PSF blurring/modeling can be enabled by setting ``options.use_psf`` to true. This uses a shift-invariant Gaussian blurring with the FWHM values taken from ``options.FWHM``. Each dimension should have its own FWHM value (in mm), but the values can be identical. The PSF is thus identical along the ray and only potentially varies along a specific dimension, although even then the variation is always constant. The PSF is applied in both forward and backward projection cases. It is also applied to custom reconstruction cases calling the projector operators, if selected. CT data ------- All optimal CT projectors are OpenCL/CUDA only! This means that only implementations 2, 3 and 5 support them and the CPU version of implementation 2 does not! For CT data, the most optimal choice is generally the branchless distance-driven (BDD) projector, i.e. projector_type = 5. However, this method is slow and can be practically useless in certain cases (depending on the scanner geometry). It has reduced usefulness when using regularization. BDD is the branchless version of the distance-driven projector. In general, the DD methods compute the area of the intersection. For forward projection, we project rays from the source to the four corners of a detector pixel. The area inside these four points in the image volume is then computed for each slice. Backprojection works similarly, but we project the lines to the corners of each voxel and then compute the area on the detector. In general, backprojection tends to have a greater beneficial effect compared to the other projectors than forward projection. See below for hybrid projectors. Note that BDD is not yet supported with helical CT data with a curved panel. Projector type 4 is a method that should work well in most cases. It is a ray-based interpolation-based projector in the forward projection and a voxel-based projector in the backprojection. The forward projection is identical to the PET version, i.e. it linearly interpolates values after a pre-determined step size. This step size can be adjusted with ``options.dL`` and is the relative size of one voxel. I.e. ``options.dL = 1`` would use the length of one voxel as the interpolation length. The backprojection projects rays from the source through each voxel to the detector and then interpolates this value on the detector using linear interpolation. This projector is not exactly adjoint, but the difference is generally less than 1%. Using a longer interpolation length leads to faster computation, but too large values lead to loss of accuracy. Generally, values of 0.5 or 1 are a good choice, but values up to 2 should work fine. An alternative method to projector type 4 is type 1, which computes the exact intersection length of each ray with each voxel. Note, however, that this should always be combined with another backprojector! See the hybrid methods below. This is a faster forward projector than type 4 if ``options.dL = 0.5`` or smaller, but slower if ``options.dL = 1`` or higher. Type 4 can have faster convergence than type 1 (i.e. the same number of iterations can lead to a noisier image with type 4), but in general there shouldn't be significant quality differences (note that this applies only to hybrid cases!). Use of type 1 on its own, i.e. as ``options.projector_type = 1`` is NOT recommended! Projector types 2 and 3 are not recommended for CT data, but they work with implementation 4, i.e. they support CPU reconstruction outside of OpenCL (as does projector type 1). Projector type 1 also works with implementation 1, if you need the actual system matrix, see :doc:`reconstruction`. Hybrid projectors can be a good choice for CT data. For example, ``projector_type = 14`` is a method where the forward projection uses an improved version of Siddon's ray-tracer as outlined above and the voxel-based backprojector. Alternatively, BDD can be combined with either, e.g. ``projector_type = 45``. Note that generally it is not recommended to use hybrid projectors where BDD is the forward projector, such as ``projector_type = 54``. A good combination of quality and speed is to use ``projector_type = 45``, or ``projector_type = 15``. SPECT data ---------- For SPECT, currently only projector types 1, 2, and 6 are supported. Projector type 1 is a ray-based projector algorithmically identical to the PET version, i.e. the exact intersection length between a voxel and the ray is computed and converted to probability. This projector can use either one or more rays. However, at the moment the number of rays has to be such that the square root of the total number of rays is an integer. This means that ``options.nRays = 9`` is accepted, but ``options.nRays = 7`` is not. The pattern of the rays traced is similar for all detector pixels and can be selected either arbitrarily or from a cone determined by the collimator geometry. This serves as an approximation for the collimator cone of response in the case that the collimator holes do not align with the detector pixels. However, practical tests have shown that this approximation has a minimal effect on image quality. If the number of rays is 1, no collimator correction is applied. Projector type 6, on the other hand, is a rotation-based (rotate-and-sum) projector where the image is rotated and then reconstructed as a parallel beam case with a point spread function. This works with any "traditional" SPECT scanner with rotating panels. This projector always applies the resolution recovery/collimator response function to the reconstruction. The intrinsic resolution of the scanner is required for this projector. Alternatively, projector type 2 is an orthogonal distance-based projector. The basic principle is the same as in the PET version and in projector type 1: each ray is traced from the "source" to the detector. This projector, however, also "spreads" the ray along the cone of the SPECT collimator. This means that at each distance from the point to the collimator, the length of this spread is computed and the voxel weighted accordingly. This means that the spread is larger the further we are from the collimator and more voxels are taken into account. The weighting is applied as a Gaussian weighting, where the standard deviation is determined by the cone (collimator). The Siddon ray tracing projector is usable for any parallel hole or pinhole collimator. Any other type of collimator such as converging or diverging hole requires additional detector coordinate and/or sinogram manipulation. Orthogonal distance ray tracer and the rotate-and-sum projector support only parallel-hole collimators. If you use traditional SPECT scanners with rotating flat-panels, projector type 6 is both high-quality and very fast. If you want to reconstruct data without collimator effect, use projector type 1 with one ray. Otherwise, use projector type 2. Other data ---------- Projector type 1 is recommended. It is the most robust and flexible method and should work in all voxel-based ray-tracing cases. Projector type 3 might also work, but this depends on the data. See PET data above for details on projector type 3. Type 4 should also be applicable to all cases, but, as mentioned above, only works in OpenCL/CUDA environments. If your data is similar to that of CT data (i.e. individual projections on a flat panel), then using CT projectors should be fine. In such a case, see the CT data above.