Aspects of the present disclosure relate in general to medical imaging, and more particularly to improved single photon emission computed tomography (SPECT) reconstruction of 3D images.
Medical imaging of metabolic and biochemical activity within a patient is known as functional imaging. Functional imaging techniques include, for example, nuclear imaging such as Positron Emission Tomography (PET), Single Photon Computed Tomography (SPECT), functional magnetic resonance imaging (fMRI), and functional computed tomography (fCT). An overview of SPECT, PET systems, their combination with computer tomography (CT) systems as well as iterative image reconstruction for emission tomography is given in chapter 7, chapter 11, and chapter 21 of M. Wernick and J. Aarsvold, “Emission tomography: the fundamentals of PET and SPECT,” Elsevier Academic Press, 2004, the contents of which are herein incorporated by reference.
In general, SPECT imaging is performed by using a gamma camera to acquire multiple two-dimensional (2D) projections and then using a computer to perform tomographic image reconstruction to obtain a three-dimensional (3D) image. For example, a gamma photon-emitting radioisotope may be introduced into a patient's body, and any of various techniques can be used to bind the radioisotope to a location of interest in the body. The patient lies on a bed, and one or more gamma cameras are attached to the gantry which rotates, causing the gamma camera(s) to rotate around the patient. Detectors of the gamma camera(s) acquire projection data at each orientation by detecting gamma photons emitted by the radioisotope.
SPECT imaging devices are generally equipped with a collimator such as a parallel-hole, fan-beam, or cone-beam, or coded aperture (CA) collimator. In the case of parallel hole collimation the point response function focuses on a small area, so that each detector pixel mainly receives photons from the incoming gamma photon beam in the normal (i.e., perpendicular) direction. In the case of CA collimation, a thin shield with numerous holes (apertures) drilled in a certain pattern is placed parallel to the surface of a detector. With CA collimation, a detector pixel receives flux from various locations, so the projections from different regions of the source overlap, introducing complexities to the task of tomographic image reconstruction.
In some embodiments of the present disclosure, a method of generating a three dimensional (3D) image includes acquiring projection data for a portion of the body of a patient at multiple views using one or more detectors, the projection data including multiple two dimensional (2D) projections. A 3D image is initialized. For each view among the plurality of views, the 3D image is transformed using a view transformation corresponding to said view to generate an initial transformed image corresponding to said view, and multiple iterations of an MLEM process are performed based on at least the initial transformed image and the projection data. The MLEM process is initialized with the initial transformed image. The 3D image is updated based on an output of the MLEM process.
In some embodiments of the present disclosure, a method of image processing includes transforming a first 3D image using a view transformation corresponding to a view, to generate a transformed image. A second 3D image is initialized with the transformed image. For each iteration of a plurality of N iterations, wherein N is an integer greater than one, maximum likelihood expectation maximization (MLEM) processing is performed. The iterative MLEM processing refines the second 3D image by: forward projecting the second 3D image to generate a data model; comparing the data model with a two dimensional (2D) projection data set corresponding to said view, to generate a correction factor corresponding to said view; backprojecting the correction factor corresponding to said view and applying a normalization factor in image space, to generate an image update factor corresponding to said view and to said iteration, and updating the second 3D image using the image update factor corresponding to said view and to said iteration. The processing flow then proceeds to the next iteration.
In some embodiments, a machine-readable storage medium tangibly embodies a program of instructions executable by a processor to cause the processor to perform operations of the image processing/image generation methods described above.
The following will be apparent from elements of the figures, which are provided for illustrative purposes and are not necessarily to scale.
This description of the exemplary embodiments is intended to be read in connection with the accompanying drawings, which are to be considered part of the entire written description.
Referring to
In image reconstruction, object space and data space are related to each other through a system matrix. Thus, for any projection operation, one can use the appropriate system matrix and its transpose to transform objects between object space and data space. In general, a forward projection is an application of the appropriate system matrix an object in object space. The result of a forward projection is a “projected object” in data space. As an example in nuclear imaging, a forward projection is the linear operation that transforms the functional activity density into the total data model of predicted detection events. Corresponding to the forward projection, the backward projection from the data space into object space can be described as an application of the transpose of the appropriate system matrix.
In some embodiments, projection data are acquired for a portion of the body of patient 145 at multiple views, e.g., M views, using detector(s) 130. Each view corresponds to a particular rotation and/or translation of a detector connected to gantry 120. The projection data include multiple 2D projections.
Jv0=TvI (1)
In equation (1), Jv0 is the initial transformed image, and Tv is the view transformation corresponding to the first view. Because the subscript v for the terms Jv0 and Tv denotes the view index, and because the first view has been selected, v equals 1 at this stage.
Multiple iterations, e.g., N iterations, of a maximum likelihood expectation maximization (MLEM) process 430 are performed. MLEM is described at, e.g., L. A. Shepp and Y. Vardi, “Maximum Likelihood Reconstruction for Emission Tomography,” IEEE Transactions on Medical Imaging, vol. 1, pp. 113-122 (1982), the entire contents of which are hereby incorporated by reference herein. Within MLEM process 430, a 3D image J is initialized (block 432) with the initial transformed image, i.e.:
J−Jv0 (2)
At the first iteration out of the N iterations of the MLEM loop shown in
qv=HvJ (3)
In equation (3), qv represents the data model corresponding to the current view, and Hv is a system matrix corresponding to the current view.
At block 436, the data model corresponding to the current view (i.e., v=1) is compared with projection data corresponding to the current view, to generate a correction factor corresponding to the current view. This comparison may occur in data space on a pixel-by-pixel basis, e.g., by computing a ratio of the data model and projection data as follows:
rv=pv/qv (4)
In equation (4), pv represents the projection data corresponding to the current view, and rv represents the correction factor corresponding to the current view.
At block 438, the correction factor for the current view is backprojected to 3D image space, and a normalization is performed, to generate an image update factor corresponding to said view and to said iteration, e.g., as follows:
Fv=HvTrv/Nv (5)
In equation (5), Fv is the image update factor, HvT is the transpose of Hv and Nv is a normalization factor, e.g., a constant dependent on the current view. For convenience, dependence on the iteration is not denoted in the symbol Fv.
Then, the processing flow proceeds to block 432, where image J is updated using the image update factor corresponding to the current view and the current iteration on a per-voxel basis, e.g., as follows:
Jnew=FvJ (6)
Then, the processing flow proceeds to the next iteration of MLEM loop 430. In this manner, N iterations of MLEM loop 430 are performed to refine image J.
After N iterations of MLEM loop 430, a final image update factor for the current view is generated, e.g., as follows:
F_finalv−Jnew/Jv0 (7)
In other words, a final image update factor for each view may be computed as a ratio, on a per-voxel basis, of the final version of image J (i.e., after N iteration of MLEM process 430) to the initial version of J.
Then, a next view is selected, and flow proceeds to block 420. In this manner, similar processing as described as above is performed for each view among the M views. After block 440 is performed for the Mth view, the final image update factors corresponding to all the views are combined, e.g., as follows:
Fall=(ΣvTTvF_finalv)/Nall (8)
In equation (8), Nall is a constant that is used as a normalization factor. This is one possible implementation of normalization, and other techniques for normalization may be used as well.
Flow proceeds to block 410, where image I is updated, e.g., on a voxel-by-voxel basis as follows:
Inew=FallI (9)
Thus, image I is refined by the outermost loop in
In some embodiments, storage medium 154 tangibly embodies a program of instructions executable by processor 152 to cause processor 152 to perform operations in processes 500 and 600 as well as various other processing described herein.
It is understood by those familiar with the art that techniques described herein may be implemented in hardware, firmware, or software encoded (e.g., as instructions executable by a processor) on a non-transitory computer-readable storage medium.
The apparatuses and processes are not limited to the specific embodiments described herein. In addition, components of each apparatus and each process can be practiced independent and separate from other components and processes described herein.
The previous description of embodiments is provided to enable any person skilled in the art to practice the disclosure. The various modifications to these embodiments will be readily apparent to those skilled in the art, and the generic principles defined herein may be applied to other embodiments without the use of inventive faculty. The present disclosure is not intended to be limited to the embodiments shown herein, but is to be accorded the widest scope consistent with the principles and novel features disclosed herein.
This application claims priority under 35 U.S.C. § 119(e) from U.S. Provisional Application Ser. No. 62/012,445 filed Jun. 16, 2014, the entirety of which is hereby incorporated by reference herein.
Filing Document | Filing Date | Country | Kind |
---|---|---|---|
PCT/IB2015/054533 | 6/16/2015 | WO | 00 |
Publishing Document | Publishing Date | Country | Kind |
---|---|---|---|
WO2015/193802 | 12/23/2015 | WO | A |
Number | Name | Date | Kind |
---|---|---|---|
20030190065 | Hamill | Oct 2003 | A1 |
20050152590 | Thieret | Jul 2005 | A1 |
20090110255 | Vija | Apr 2009 | A1 |
20100014730 | Hahn | Jan 2010 | A1 |
20140177785 | Funk | Jun 2014 | A1 |
20160166218 | Williams | Jun 2016 | A1 |
Number | Date | Country |
---|---|---|
101031938 | Sep 2007 | CN |
102364524 | Feb 2012 | CN |
102525501 | Jul 2012 | CN |
103180883 | Jun 2013 | CN |
P0400894 | Aug 2006 | HU |
226835 | Nov 2009 | HU |
P1200013 | Jul 2013 | HU |
Entry |
---|
Ing-Tsung, Hsiao et al. “An accelerated convergent ordered subsets algorithm for emission tomography; An accelerated convergent ordered subsets algorithm for emission tomography”, Physics in Medicine and Biology, Institute of Physics Publishing, Bristol GB, vol. 49, No. 11, Jun. 7, 2004, pp. 2145-2156. |
PCT International Search Report and Written Opinion dated Sep. 4, 2015 (10 pages). |
Search Report for Corresponding Hungarian Patent Application No. P1700011, dated Mar. 24, 2017. |
L.A. Shepp and Y. Vardi, “Maximum Likelihood Reconstruction for Emission Tomography,” IEEE Transations on Medical Imaging, vol. 1, pp. 113-122 (1982). |
Number | Date | Country | |
---|---|---|---|
20170098316 A1 | Apr 2017 | US |
Number | Date | Country | |
---|---|---|---|
62012445 | Jun 2014 | US |