CA2327750A1 - Real-time image reconstruction for computed tomography systems - Google Patents

Real-time image reconstruction for computed tomography systems Download PDF

Info

Publication number
CA2327750A1
CA2327750A1 CA002327750A CA2327750A CA2327750A1 CA 2327750 A1 CA2327750 A1 CA 2327750A1 CA 002327750 A CA002327750 A CA 002327750A CA 2327750 A CA2327750 A CA 2327750A CA 2327750 A1 CA2327750 A1 CA 2327750A1
Authority
CA
Canada
Prior art keywords
system matrix
vector
image
reconstruction
matrix
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Abandoned
Application number
CA002327750A
Other languages
French (fr)
Inventor
Vitali Selivanov
Roger Lecomte
Germain Leger
Original Assignee
Universite de Sherbrooke
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Priority claimed from CA002324832A external-priority patent/CA2324832A1/en
Application filed by Universite de Sherbrooke filed Critical Universite de Sherbrooke
Priority to CA002327750A priority Critical patent/CA2327750A1/en
Priority to CA002360656A priority patent/CA2360656A1/en
Publication of CA2327750A1 publication Critical patent/CA2327750A1/en
Abandoned legal-status Critical Current

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/002D [Two Dimensional] image generation
    • G06T11/003Reconstruction from projections, e.g. tomography
    • G06T11/006Inverse problem, transformation from projection-space into object-space, e.g. transform methods, back-projection, algebraic methods
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
    • A61B6/50Clinical applications
    • A61B6/508Clinical applications for non-human patients
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/11Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/12Simultaneous equations, e.g. systems of linear equations
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
    • A61B6/02Devices for diagnosis sequentially in different planes; Stereoscopic radiation diagnosis
    • A61B6/03Computerised tomographs
    • A61B6/037Emission tomography
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2211/00Image generation
    • G06T2211/40Computed tomography
    • G06T2211/428Real-time

Abstract

The present invention discloses a method and a system for real-time reconstruction of tomographic images from projection data-generated by computed tomography scanners. The method is based on the theory of pseudo-inverse matrices applied to tomographic image reconstruction, and it allows immediate, projection-by-projection or event-by-event updating of the reconstructed tomographic image at the moment the measuring instrument supplies additional individual projection data. The data processing involves simple arithmetic, which can be implemented into general-purpose computers in the first case, and digital signal processors, dedicated programmable logic devices or application specific circuits in the second case. Models of the underlying physical processes, such as the photon emission and detection processes or the spatially variant system response, can be easily included into image reconstruction through the system matrix. Corrections for other physical and instrumental factors affecting the accuracy of projection data, such as attenuation in the body or detector efficiency normalization, can be merged with the pseudo-inverted and regularized system matrix beforehand or be performed "on the fly" by simple multiplication of the matrix elements, thus achieving optimal compression. Real-time updating of the reconstructed tomographic image can be performed as soon as measured projection data become available. No storage of the measured raw data in the form of list mode events or sinogram is required, and the user of the tomographic system may choose to store reconstructed images only, thus achieving optimal compression of the tomographic data.

Description

REAL-TIME IMAGE RECONSTRUCTION FOR
COMPUTED TOMOGRAPHY SYSTEMS
BACKGROUND OF THE INVENTION
1. Field of the invention:
The present invention relates to a method for conducting real-time matrix pseudo-inverse tomographic reconstruction.
2. Brief description of the prior art:
Tomography refers to the cross-sectional imaging of an object from either transmission, emission or reflection data collected from many different directions.
Tomographic imaging deals with reconstructing an image from such data, commonly called projections. From a purely mathematical standpoint, a projection at a given angle is the integral of the image in the direction specified by that angle. The solution to the problem of reconstructing a function from its projections dates back to the paper by Radon in 1917 [1], thus the denomination inverse-Radon transform. However, it is only in the early 1970's, with the invention of the X-ray computed tomographic (CT) scanner [2,3] and the development of reconstruction algorithms [4-7] that tomographic imaging came into widespread practice. Tomographic methods now find applications in a vast number of fields, including radio astronomy, seismology, nondestructive analyses, electron microscopy, and above all medical imaging. In fact, most of the powerful new medical imaging modalities that have been introduced during the last three decades, such as X-ray CT, Single Photon Emission Computed Tomography (SPELT), Positron Emission Tomography (PET), Magnetic Resonance Imaging (MRI) and 3D ultrasound (US), were the result of the application of tomographic principles.
A. Review of current tomographic reconstruction algorithms Comprehensive description of the existing tomographic reconstruction algorithms can be found in textbooks by Herman [8], Dean [9], Natterer [10], Kak and Slaney [11] and Fessler and Ollinger [12], and review articles by Vardi et al.
[13] and Ollinger and Fessler [14]. There are two main classes of tomographic reconstruction methods currently in use today: analytical reconstruction algorithms, which aim at solving the inverse-Radon transform, and iterative algorithms, which attempt to reach a solution by successive estimation of the underlying image and improving image estimates at every iteration.
A. 7 Filtered Backproiection The most common tomographic reconstruction method nowadays is the filtered backprojection (FBP) algorithm, which makes use of the projection-slice theorem [11]. The algorithm proceeds by first transferring each measured projection into Fourier space, filtering with an appropriate filter in frequency space, transferring the result back into the projection space and then backprojecting the result onto the image grid:
x=Back-Project{F,-'~c~f~xF,~b'~j~ (1) where b~ is a j-th projection vector (j=1,...,A), A is a total number of projection angles, F 1 is a 1-D Fourier transform, F 1-1 is a 1-D inverse Fourier transform, c(f~ is a 1-D filter function in the Fourier space, commonly referred to as a "ramp" filter, and z = ~xk : k =1,..., M~ is the image that has to be found, M
being the total number of image pixels (or voxels).
The FBP algorithm produces an image that is a linear combination of the projection data. It can be efficiently implemented in general-purpose computers, its memory requirements are modest and the computation is rather fast. It is the reconstruction method of choice for virtually all X-ray CT and most emission tomography (SPELT, PET) systems currently in operation. Since the FBP
algorithm is linear, it can be reduced to event-by-event reconstruction based on the principle of superposition, by deriving a filter function representing the contribution of an individual event to the image. However, its widespread popularity stems more from historical reasons of computational simplicity than any widely accepted advantage in image quality. In fact, there are several problems with this algorithm:
- data must be assumed to be uniformly distributed on projections (which is generally not the case with the cylindrical scanner geometry), or must be rebinned with equal spacing to accommodate the convolution operation (performed in Fourier space);
- models for the detector response must be space invariant and can only be incorporated into the algorithm as a deconvolution with the attendant noise amplification;
- the ideal ramp filter amplifies high-frequency noise in the projection, and it must often be apodized with a window function to reduce noise in the reconstructed image;
- even though the intensity is known to be non-negative, the algorithm yields negative values, particularly if the data are noisy;
- streak artefact is present; and - event-by-event reconstruction involves some elaborate computation to align, rotate and scale the single event filters before summation to the image matrix.
A.2 Alqebraic Reconstruction There are other tomographic image reconstruction methods called algebraic reconstruction techniques [15,16]. They aim at solving the system of linear equations p11x1 ~' p12x2 +... + plMxM = b1 p21x1 +p22x2 +...+ p2MxM = bz (2) pNlx1+pN2x2-f'...+pN,uxu =bN
where bi : i =1,..., N are the measured projections (sinogram), x = {xk : k =1,..., M~ is the image that has to be found, N and M are the total number of tubes-of-response and image pixels (or voxels), respectively.
These methods are usually iterative in nature and tend to converge to the minimum norm least-squares solution as the iteration number goes to infinity [16].
Some of the drawbacks of the iterative algebraic methods are:
- highly computer intensive, large memory requirements;
- because the algorithms are inherently iterative, real-time reconstruction can hardly be implemented;
- the need to optimize iteration number and the order in which projections are being utilized in the reconstruction;
- non-linearity of the solution due to the iterative nature of reconstruction algorithm, which results in complex quantitative dependence of image estimate on input data; and - the complete measured projection data set is required before starting reconstruction.
A.3 Statistical Reconstruction Statistical image reconstruction methods were introduced in an attempt to overcome some of the drawbacks of the previous methods, in particular taking into account the stochastic nature of the measured data in emission tomography [13,17] and the physical modeling of the detection process [18]. In general, statistical iterative reconstruction methods yield images with much superior visual quality than FBP reconstructed images, especially in the case of low projection statistics. The elimination of negative values and proper modeling of the detector response avoid several of the artifacts present in FBP images. Unfortunately, this superior image quality does not necessarily translate into better quantitation accuracy, as this was recently demonstrated for the pharmacokinetic modeling of dynamic image sequences [19]. The use of statistical reconstruction methods also raises other practical problems that still hamper their utilization outside the research environment. In addition to most of the drawbacks that are common with the iterative algebraic methods, the statistical methods have the following problems:
- in the case of unconstrained iterative estimation, the absence of reliable, operator-independent, criteria to stop iterations: if iteration number is too low, the spatial resolution is sub-optimal; if iteration number is too high, the spatial resolution is artificially enhanced at the expense of unacceptable noise amplification;
- necessity to optimize an extra parameter, called prior, that controls noise or image smoothness in the case of so called Bayesian approach [20]; and - non-linearity of the solution due to the iterative nature of reconstruction algorithm, which results in complex quantitative dependence of image estimate on input data.
SUMMARY OF THE INVENTION
In an attempt to overcome the above discussed drawbacks of the prior art, there is provided, according to the present invention, a method for conducting tomographic image reconstruction utilizing matrix inverse or pseudo-inverse.
As will become evident from the following, the method is also amenable to real time reconstruction by reducing the computation of matrix operations to vector operations corresponding to the reconstruction of individual events or projections by successive update of an image using dedicated or general purpose computing hardware. In accordance with this method, a system matrix is first constructed to relate a tomographic image grid comprised of discrete pixels (voxels) to the measured projections. If the system matrix is sufficiently well conditioned, the image can be reconstructed by inverting the system matrix to find the relationship of the measured projections to the underlying image. To facilitate the matrix inversion with the ill-conditioned matrices typical for the practical imaging geometries, the system matrix is factored by singular value decomposition, wherein this factored system matrix includes a diagonal matrix of singular values.
The singular values are usually presented as a non-increasing sequence called singular value spectrum. The factored system matrix is pseudo-inverted to produce a pseudo-inverted system matrix, in which the vector-columns are representative of the contribution of individual projection data to the reconstructed image. The pseudo-inverted system matrix is regularized by truncating the singular value spectrum at an index Tand removing small singular values ,u;, i = T +1, ... , M, M being the total number of image voxels, or modifying the singular value spectrum to diminish the effect of very small singular values.
Finally, a new reconstructed image is evaluated by summing a previous image estimate and one vector-column of the pseudo-inverted regularized system matrix corresponding to a last event recorded. In operation, reconstruction of the tomographic image involves summation of the previous image estimate and a current update from a new measured event or additional single projection data.
The method of the invention presents, amongst others, the following advantages:
- it allows real-time reconstruction of projection data;
- it allows event-by-event updating of a tomographically reconstructed image;
- it allows individual projection data to be reconstructed independently to update an existing image;
- it allows real-time updating and display of an image;
- it allows to monitor tomographic data acquisition in real-time;
- it avoids data rebinning into projection vectors, since the exact geometry of the measuring instrument can be utilized; and - it may avoid storage of measured projection data (either in list mode, histogram or rebinned sinogram) if further processing of the measured projection data is unnecessary.
The above as well as other advantages and features of the present invention will become more apparent upon reading of the following non restrictive description of a preferred embodiment thereof, given as illustrative example only with reference to the accompanying drawings.
BRIEF DESCRIPTION OF THE DRAWINGS
In the appended drawings:
Figure 1 is an illustration of the various steps involved in the real-time tomographic image reconstruction by the matrix inverse/pseudo-inverse method;
Figure 2 is a schematic diagram of dedicated hardware for performing real-time tomographic image reconstruction based on the matrix pseudo-inverse method;
Figure 3 is a singular value spectrum of the system matrix for an animal PET scanner, 64x64 pixel image;
Figure 4 is a sequence of phantom images demonstrating the use of real time TSVD reconstruction, wherein estimated radioactivity distribution progress is shown with one plane (64x64 pixel images), approximately 23800 additional counts are accumulated with each next image estimate shown, and less than 1/3 of the singular value spectrum is utilized; and Figure 5 is a s equence of phantom images demonstrating the use of real time TSVD reconstruction, wherein estimated radioactivity distribution progress is shown with one plane (64x64 pixel images), projections at eight consecutive angles are included into each next image estimate shown, and less than 1/3 of the singular value spectrum is utilized.
DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENT
An approach that was discussed in the literature [21 ] some time ago, but has not yet gained popularity for tomographic reconstruction is the matrix pseudo-inverse method, whereby the system matrix relating the tomographic image grid to the measured projections must be pseudo-inverted. Its use has been suggested for a limited number of projections. For complete number of projections, i.e. practical problems of tomographic medical imaging, its use has been hindered by the computational burden resulting from the size of the matrix and the inverse problem ill-conditioning which has to be accounted for in a manner preferably independent of an image being reconstructed. Matrix pseudo-inversion can be performed based on singular value decomposition (SVD) of the system matrix with some regularization to diminish the effect of small singular values. The application to the tomographic reconstruction problem is described below, followed by examples of practical hardware implementation of the technique for real-time image reconstruction.
Matrix Pseudo-Inverse Tomogrraphic Reconstruction The image reconstruction problem in tomography can be written in the matrix form Px = b (4) where N
P = p~~ : i =1,..., N; ~ ply =1, j =1,..., M (5) r=i is the system matrix, b = ~b; : i =1,..., N~ is the vector-column of measured projections (sinogram), x = ~xk : k =1,...,M~ is the vector-column (image) that has to be found, N and M are the total number of tubes-of response and image pixels (or voxels), respectively. For example, the matrix P may be obtained in PET
for a given scanner geometry with the approach described in [18].
5 Any matrix P can be factored [22] into P=~T
where U = ~u~~ : i =1,..., N; j =1,...,M~ and V = ~v,~ : i, j =1,..., M~ are orthogonal 10 matrices, W = ~,u,~ : i, j =1,..., M; ,u;~ = 0 if i ~ j ~ is a diagonal matrix containing singular values ,ui ---_ ~ZZ, i =1,...,M . Factored matrix representation (6) is referred to as singular value decomposition (SVD). An algorithm to compute matrix SVD
may be found in [23].
In the case when measured data are noisy and/or (4) is overdetermined/
underdetermined, and the exact solution of (4) does not exist, one may find a minimum norm least-squares solution of (4) using X - p.b (7) where P+ = VW +UT (8) is the pseudo-inverse of P [22]. Singular values are usually presented as a non-increasing sequence called singular value spectrum. Some singular values can be very small (or even zeros if the matrix is singular).
The condition number max,ui c = t (9) min,u~
t of the system matrix can be very high (even infinity for a singular matrix).
Thus, the inverse problem of tomographic image reconstruction is ill-conditioned and the solution exploiting pseudo-inverse (8) will be very sensitive to noise.
One way of solution regularization is to truncate the singular value spectrum at some index T and remove small singular values ,ut , i = T + 1,..., M
from the solution expansion, leaving T N
x = zk : xk = ~vkt,ui ~ ~u jtb j; k =1,...,M (10) i=1 j=1 which is referred to as the truncated SVD (TSVD) solution.
Another regularization approach would be to modify the singular spectrum in order to diminish the effect of the inversion of very small singular values (and zeros if any):
x= xk : zk = ~,Vkl~Ol) 1 ~ ujlbj; k=1,...,M (11) i=1 j=1 where f (~cl ~ is the modified singular value. Such solution is referred to as the modified SVD (MSVD) solution. An example of MSVD can be found in [24].
Real-Time Reconstruction In this section let us concentrate on the TSVD solution. Rearranging the summation order in (10) we get xk = ~ ~ vkif~~ ujtbj = ~ ~vkif~i ujt bj~ k=1,...,M (12) i=lj=1 j=I i=1 Let P+ be the "truncated pseudo-inverse":
T
p+ _ ~Pkj ~ Pkj = ~vktN~~ lu jj a k =1,...,M; j =1,..., N (13) i=1 Then the TSDV solution in the matrix form is x = P+b (14) Let b(t~ be a sinogram containing a total of t counts:
N
b(t~= b;(t~---b; : i=1,...,N; ~bt =t (15) l=I
and z(t) _ ~xk (t) : k =1,..., M~ is the respective solution given for b(t) by (14). Let Ob(S) _ ~~bt : i =1,...,N~ is such that 0, i =1,..., s -1 0b;--_- l,i=s (16) 0, i = s + 1,..., N
and b(S)(t) ---- b(t -1)+ 0 b(S) (17) In other words, b(S)(t)= {b~s)(t): i =1,...,N~ is a sinogram, which differs from b(t-1) by only one count in a bin with index s:
b; (t-1),i=1,...,s-1 b~s)(t)= b; (t-1)+l,i=s (18) b; (t-1),i=s+1,...,N
Assuming that x~°~(t) stands for the image estimate obtained using input data b~°~(t), then the solution given by x(S)(t) = p+b(S)(t) (19) may be found as follows. We have xk~t~= ~pkjb~s~~t~= ~pkjbi~t-1~+pks~bs~t-1~+'1~+
j=1 j=1 N _ N (20) + ~ pkjbi~t-l~= ~pkjbi~t-1~+pkr~ k=1~...,M
j=s+1 j=I
Thus, the new image estimate is simply a sum of the previous image estimate and one column (vector-column) of matrix P+:
x(S)(t)=~ks)(t): xks)(t)=xk~t-1)+p~, k=1,...,M} (21) It is natural to start from the blank image x(o) _ {xk (o): k = i,..., M; xk (o) = o} and use (21 ) afterwards to obtain the current image estimate in real-time based on the previous image and a column of P+ defined by the sinogram bin index where the last event was registered.
MSVD solution can be obtained in real time as well using (21 ) and taking into account that "modified pseudo-inverse" P+ is given in this case by:
M _ P+ _ ~ pkj : pkj ---- ~ vk; f (p.l ) I a jt ; k =1,..., M; j =1,..., N (22a) i=1 The proposed approach is general and can be applied with any number of projection bins and image pixels (or voxels), for 2D as well as for 3D-reconstruction geometry.
by M
P+ - pkj : pkj =- ~ vkt.f (fc; )u j=; k =1,..., M; j =1,..., N (22b) i=I
where f(~; ) is a function of ~; other than f(~; ) _ ~t'~ ~' ' 'T
o, i > T ° which was used with TSVD.

The proposed approach is general and can be applied with any number of projection bins and image pixels (or voxels), for 2D as well as for 3D-reconstruction geometry, with scanners collecting complete or incomplete sets of projections, for non-singular or singular system matrices. A diagram illustrating the process of reaf-5 time tomographic image reconstruction, using positron emission tomography as an example, is shown in Figure 1.
Corrections of aroiection data 10 Physical and instrumental determinants that affect the measurement of projection data must be incorporated into the image reconstruction. Such factors as the spatially variant system response and model of the emission and detection processes can be included into the computation of the system matrix. Other factors related to the specific instrument being used and the subject under study, 15 such as the normalization for detector efficiency or the correction for signal attenuation in tissue (e.g., photon attenuation in the case of emission tomography), can be included by multiplying the elements of regularized pseudo-inverse matrix P+ with the proper factors, either before reconstruction or on the fly during reconstruction. The updated equation (21 ) will change slightly in the case of correction on the fly:
x(S)(t)=~xks)(t): xks)(t)=xk~t-1~+FS x per, k=1,...,M} (23) where FS is the correction factor assigned to the tube-of-response s.
Similarly, random coincidence events in PET could be corrected for by subtracting one column of matrix P+ from the previous image estimate:
x(S)(t)-~kS)(t): xks)(t)=xk~t-1~-p~,k=1,...,M} (24) or, instead, by skipping the next column addition given by (21 ) when the next coincident event is registered. In the case of correction on the fly, (24) would be replaced by:
x(S)(r)=~ks)(r): xks)~t)=xk~t-l~-FS x per, k=1,...,M} (25) Hardware implementation The reconstruction of tomographic images using the proposed method has been shown to be amenable to a simple sum of two vectors representing the previous image estimate and the current update from a new measured event or additional single projection data. The latter vector may be multiplied by a normalization factor accounting for physical or instrumental corrections.
Computation should be performed in floating point arithmetic to satisfy the required precision, but integer arithmetic may be used instead in some special cases. In that case, real values will be scaled to integers that will be used in image update. One significant issue to account for with the integer arithmetic would be the possible value overflow due to the limited number of bits used for the integer representation. In fact, quantitation from regions-of-interest in the reconstructed image should always be carried out in floating point arithmetic to avoid overflow and ensure the required precision.
For real-time reconstruction, the computation is subject to the additional constraint that the duration of the entire calculation must be less than the average time between successive events. This includes fetching the column of P+defined by the bin index of the current event, multiplying by a correction factor if necessary, summing the M elements of the resulting vector to the previous image estimate, and storing the updated image.
Two different acquisition modes must be distinguished here, depending on the way projection data are being generated by the tomographic scanner. In the first mode, the measuring instrument sequentially supplies projection data, which are measured only once in each projection bin. X-ray CT scanners and ultrasound probes fall into this category. In such cases, the average time between events is determined by the scanning speed of the instrument, and is typically in the millisecond range or more. Computation in such case can be performed using general-purpose computers.
In the second mode of acquisition, events are registered randomly in all available projection bins. Instruments operating in counting mode, such as SPECT and PET scanners, fall into this category. Assuming a peak event rate of one million counts/sec, the entire computation would have to be accomplished in less than one microsecond in such case. Even for a small size 2D image (e.g., 64x64 pixels, M=4096), such computational speed can hardly be reached with the current general-purpose computer technology, unless massively parallel processing is employed. The system matrix P+, having dimensions N xM, must also be loaded into memory for fast access. Although extending the reconstruction to larger image sizes or 3D imaging geometry makes such implementation impractical with the current computer technology, this should become feasible in the foreseeable future given the constant progress in memory capacity and processing speed.
Dedicated hardware with highly parallel processing and high storage capacity must be designed to efficiently perform the required computation. The principle of a real-time processing unit to perform tomographic image reconstruction based on the matrix pseudo-inverse method is displayed schematically in Figure 2. Such hardware can easily be implemented using Digital Signal Processors (DSP), programmable logic devices (e.g., Field Programmable Gate Arrays or FPGA) or application specific integrated circuits (ASIC) by persons knowledgeable in the field.
Results To illustrate a typical matrix SVD in PET for a complete set of 2D
projections, the singular value spectrum of the system matrix P for the above mentioned animal PET scanner [25] and image of 64x64 pixels is shown in Figure 3. The matrix condition number is c = 4441,5 in this case.
Sequence of phantom images demonstrating the use of real time TSVD
reconstruction based on (18) is presented in Figures 4 and 5. Data acquired with a phantom of 110 mm diameter made of Lucite and having holes of diameter 2, 3.4, 6.7, 9.7, 13, 15.8, 20.3, 22.7 mm located on a circumference at a distance of 28 mm from the center and filled with '8F-fluorodeoxyglucose were reconstructed. In Figure 4 every image has approximately 23800 additional counts as compared to the previous image estimate. Figure 5 shows image estimates obtained by including projections at consecutive angles. In both cases, less than 1/3 of the singular value spectrum was utilized. An approach for singular value spectrum truncation based on spatial resolution analysis may be found in [26].
Discussion of application in medical imaaina SVD of the system matrix, apart from precise numerical diagnostics of the tomographic reconstruction ill-conditioning with a given detection system geometry, provides a linear and very fast reconstruction means.
Singular value spectrum truncation allows separation of the signal and noise at the reconstruction step. Index T sets the trade-off between noise and resolution. Truncation of the singular spectrum is not the only way of solution regularization. Spectrum modification without truncation may be an appropriate regularizing approach in some situations. It is possible as well that in some special situations when the system matrix corresponding to a particular tomographic system is sufficiently well conditioned no SVD pseudo-inversion is necessary and direct matrix inversion is feasible.
TSVD reconstruction has some drawbacks: negative values in the image estimate, streak artifact with low-count images if index T is lower than a certain threshold value or noise artifact if T is higher than the threshold value. But these features (except for the noise artifact with high T, which is analogous to the artifact developing with high iteration numbers when unconstrained iterative image estimation is performed) are also common to the filtered backprojection image reconstruction, the most popular method used in medical practice today.
TSVD (as well as MSVD) also has very attractive benefits:
1 ) spatially variant system response and model of the signal emission and detection processes can be easily included into image reconstruction through the system matrix;
2) data rebinning is not necessary since the geometry of a given system is utilized;
3) the resolution in reconstructed images may be set based on the spatial resolution analysis in reconstructed images by varying truncation index T
(or modifying singular value spectrum in MSVD case);
4) noise amplification may be controlled by varying truncation index T (or modifying singular value spectrum in MSVD case);
5) image reconstruction can be performed in real time, on an event-by-event or projection-by-projection basis, allowing for instant visualization of the radioactivity distribution while the subject is being scanned; and 6) since measured projection data can be reconstructed "on the fly" as soon as they are acquired by the scanner, storage of the measured data or intermediary calculation results, except for the image estimate being 5 updated, is unnecessary, which allows for optimal compression of the tomographic data.
Although the present invention has been described hereinabove by way of a preferred embodiment thereof, this embodiment can be modified at will, 10 within the scope of the present invention, without departing from the spirit and nature of the subject of the present invention. Moreover, the application of the present invention is not limited to medical imaging only, but may be possible in other imaging techniques utilizing tomographic principles and image reconstruction from projections.
References:
[1] J. Radon, "On the determination of functions from their integrals values along certain manifolds", IEEE Trans. Med. Imag., vol. MI-5, no. 4, pp. 170-176, 1986.
(Translated from the original text that appeared in Berichte Saechsische Akademie des Wissenschaften, vol. 29, pp. 262-277, 1917) [2] G.N. Hounsfield, "Method and apparatus for measuring X- or 0-radiation absorption or transmission at plural angles and analyzing the data", U.S.
Patent 3,778,614, 1973.
[3] G.N. Hounsfield, "Computerized transverse axial scanning (tomography). 1.
Description of system", Brit. J. Radiol., vol. 46(552), pp. 1016-22, 1973.
[4] A.M. Cormack, "Representation of a function by its line integrals, with some radiological applications", J. Appl. Phys., vol. 34, pp. 2722-2727, 1963.
[5] A.M. Cormack, "Representation of a function by its line integrals, with some radiological applications. II", J. Appl. Phys., vol. 35, pp. 2908-2913, 1964.
[6] A.M. Cormack, "Reconstruction of densities from their projections, with applications in radiological physics", Phys. Med. Biol., vol. 18(2), pp. 195-207, 1973.
[7] L.A. Shepp, B.F. Logan, "The Fourier reconstruction of a head section", IEEE
Trans. Nucl. Sci., vol. NS-21, pp.21-43, 1974.
[8] G.T. Herman, Image Reconstruction from Projections, Academic Press, New York 1980.
[9] S.R. Dean, The Radon Transform and Some of Its Applications, John Wiley, New York, 1983.
[10] F. Natterer, The Mathematics of Computerized Tomography, Wiley, New York, 1986.
[11] A.C. Kak, M.Slaney, Principles of Computerized Tomographic Imaging, IEEE
Press, New York, 1988.
[12] J.A. Fessler, J.M. Ollinger, Signal Processing Pitfalls in Positron Emission Tomography, Technical Report 302, Comm. and Sign. Proc. Lab., Dept. of EECS, Univ. of Michigan, 48109-2122, Ann Arbor, MI, 1996.
(http://www.eecs.umich.edu/~fessler).
[13] Y. Vardi, L.A. Shepp, L. Kaufman, "A statistical model for positron emission tomography", J. Am. Stat. Assoc., vol. 80, pp. 8-37, 1985.
[14] J.M. Ollinger and J.A. Fessler, "Positron-Emission Tomography", IEEE
Signal Processing Mag., vol. 14(1 ), pp. 43-55, 1997.
[15] G.T. Herman, L.B. Meyer, "Algebraic reconstruction techniques can be made computationally efficient", IEEE. Trans. Med. Imag., vol. 12, no. 3, pp. 600-609, 1993.
[16] X.-L. Xu, J.-S. Liow, S.C. Strother, "Iterative algebraic reconstruction algorithms for emission computed tomography: A unified framework and its application to positron emission tomography", Med. Phys., vol. 20(6), pp. 1675-1684, 1993.
[17] A.J. Rockmore, A. Macovski, "A maximum likelihood approach to emission image reconstruction from projections", IEEE Trans. Nucl. Sci., vol. 23, pp.

1432, 1976.
[18] V.V. Selivanov, Y. Picard, J. Cadorette, S. Rodrigue, R. Lecomte, "Detector response models for statistical iterative image reconstruction in high resolution PET", IEEE Trans. Nucl. Sci., vol. 47, No. 3, pp. 1168-1175, 2000.

[19] V. Selivanov, D. Lapointe, M. Bentourkia, R. Lecomte, "Cross-validation stopping rule for ML-EM reconstruction of dynamic PET series: Effect on image quality and quantitative accuracy", 1999 IEEE Nuclear Science Symposium &
Medical Imaging Conference Record, CD-ROM, MIC Sessions, M10-18.
[20] E. Levitan, G.T. Herman, "A maximum a posteriori probability expectation maximization algorithm for image reconstruction in emission tomography", IEEE
Trans. Med. Imag., vol. 6, no. 3, pp. 185-192, 1987.
[21] Y.S. Shim, Z.H. Cho, "SVD pseudoinversion image reconstruction", IEEE
Trans. ASSP, vo1.29, no.4, pp.904-909, 1981.
[22] G. Strang, Linear algebra and its applications, Harcourt Brace Jovanovich, 1980, pp. 137-145 [23] W. H. Press, B. P. Flannery, S. A. Teukolsky, W. T. Vetterling, "Numerical recipes in C: The art of scientific computing", Cambridge University Press, 1988, pp. 60-72.
[24] M.F. Smith, "Generalized matrix inverse reconstruction for SPECT using a weighted singular value spectrum", IEEE Trans. Nucl. Sci., vol. 43, no. 3, pp.
2008-2017, 1996.
[25] R. Lecomte, J. Cadorette, P. Richard, S. Rodrigue, D. Rouleau, "Design and engineering aspects of a high resolution positron tomograph for small animal imaging", IEEE Trans. Nucl. Sci., vol. 41, pp. 1446-1452, 1994.
[26] V. Selivanov, R. Lecomte, "Fast PET image reconstruction based on SVD
decomposition of the system matrix", 2000 IEEE Nuclear Science Symposium &
Medical Imaging Conference Record, CD-ROM, MIC Sessions, IEEE664, in press.

Claims (4)

1. A method for conducting real-time matrix-inverse tomographic reconstruction, comprising:
constructing a system matrix relating a tomographic image comprised of discrete pixels (voxels) to the measured projections, said system matrix comprising vector-rows respectively representative of measured projections;

factoring by singular value decomposition the system matrix, wherein the factored system matrix includes a diagonal matrix of singular values;

pseudo-inverting the factored system matrix to produce a pseudo-inverted system matrix, said pseudo-inverted system matrix comprising vector-columns respectively representative of the contribution of individual measured projections to the reconstructed image;

regularizing the pseudo-inverted system matrix by, for instance, truncating the singular value spectrum, presented as a non-increasing sequence of singular values, at an index T and replacing small singular values µ,j = T + 1,...,M with zeros, M being the total number of image pixels (voxels), or modifying the singular value spectrum to diminish the effect of inversion of very small singular values;
measuring an individual projection event from transmission, emission and/or reflection data collected from a subject under study; and estimating a new reconstructed image by summing a previous image estimate and an update vector corresponding to a last event recorded, wherein said update vector represents one vector-column of the pseudo-inverted regularized system matrix corresponding to a last event recorded or one vector-column of the pseudo-inverted regularized system matrix corresponding to a last event recorded scaled according to the correction factor assigned to the respective projection bin;

whereby reconstruction of the tomographic image involves summation of the previous image estimate and a current update from a new measured event.
2. A method for conducting real-time matrix-inverse tomographic reconstruction, comprising:
constructing a system matrix relating a tomographic image comprised of discrete pixels (voxels) to the measured projections, said system matrix comprising vector-rows respectively representative of measured projections;
factoring by singular value decomposition the system matrix, wherein the factored system matrix includes a diagonal matrix of singular values;
pseudo-inverting the factored system matrix to produce a pseudo-inverted system matrix, said pseudo-inverted system matrix comprising vector-columns respectively representative of the contribution of individual measured projections to the reconstructed image;
regularizing the pseudo-inverted system matrix by, for instance, truncating the singular value spectrum, presented as a non-increasing sequence of singular values, at an index T and replacing small singular values µ~,j = T + 1,...,M with zeros, M being the total number of image pixels (voxels), or modifying the singular value spectrum to diminish the effect of inversion of very small singular values;
measuring a single projection data (or datum) from transmission, emission and/or reflection data collected from a subject under study; and estimating a new reconstructed image by summing a previous image estimate and an update vector corresponding to a last projection recorded, wherein said update vector represents one vector-column of the pseudo-inverted regularized system matrix corresponding to a last projection recorded scaled using the value of the data (e.g., number of counts, measured radiation flux, or reflected energy), or one vector-column of the pseudo-inverted regularized system matrix corresponding to a last projection recorded scaled using the value of the data (e.g., number of counts, measured radiation flux, or reflected energy) and using the correction factor assigned to the respective projection bin;

whereby reconstruction of the tomographic image involves summation of the previous image estimate and a current update from a new measured additional single projection data.
3. A method for conducting real-time matrix-inverse tomographic reconstruction, comprising:
constructing a system matrix relating a tomographic image comprised of discrete pixels (voxels) to the measured projections, said system matrix comprising vector-rows respectively representative of measured projections;
inverting the system matrix to produce an inverted system matrix, said inverted system matrix comprising vector-columns respectively representative of the contribution of individual measured projections to the reconstructed image;

measuring an individual projection event from transmission, emission and/or reflection data collected from a subject under study; and estimating a new reconstructed image by summing a previous image estimate and an update vector corresponding to a last event recorded, wherein said update vector represents one vector-column of the inverted system matrix corresponding to a last event recorded or one vector-column of the inverted system matrix corresponding to a last event recorded scaled according to the correction factor assigned to the respective projection bin;
whereby reconstruction of the tomographic image involves summation of the previous image estimate and a current update from a new measured event.
4. A method for conducting real-time matrix-inverse tomographic reconstruction, comprising:

constructing a system matrix relating a tomographic image comprised of discrete pixels (voxels) to the measured projections, said system matrix comprising vector-rows respectively representative of measured projections;
inverting the system matrix to produce an inverted system matrix, said inverted system matrix comprising vector-columns respectively representative of the contribution of individual measured projections to the reconstructed image;
measuring a single projection data (or datum?) from transmission, emission and/or reflection data collected from a subject under study; and estimating a new reconstructed image by summing a previous image estimate and an update vector corresponding to a last projection recorded, wherein said update vector represents one vector-column of the inverted system matrix corresponding to a last projection recorded scaled using the value of the data (e.g., number of counts, measured radiation flux, or reflected energy), or one vector-column of the inverted system matrix corresponding to a last projection recorded scaled using the value of the data (e.g., number of counts, measured radiation flux, or reflected energy) and using the correction factor assigned to the respective projection bin;
whereby reconstruction of the tomographic image involves summation of the previous image estimate and a current update from a new measured additional single projection data.
CA002327750A 2000-10-31 2000-12-07 Real-time image reconstruction for computed tomography systems Abandoned CA2327750A1 (en)

Priority Applications (2)

Application Number Priority Date Filing Date Title
CA002327750A CA2327750A1 (en) 2000-10-31 2000-12-07 Real-time image reconstruction for computed tomography systems
CA002360656A CA2360656A1 (en) 2000-10-31 2001-10-31 Real-time image reconstruction for computed tomography systems

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
CA002324832A CA2324832A1 (en) 2000-10-31 2000-10-31 Real-time image reconstruction for computed tomography systems
CA2,324,832 2000-10-31
CA002327750A CA2327750A1 (en) 2000-10-31 2000-12-07 Real-time image reconstruction for computed tomography systems

Publications (1)

Publication Number Publication Date
CA2327750A1 true CA2327750A1 (en) 2002-04-30

Family

ID=25682199

Family Applications (1)

Application Number Title Priority Date Filing Date
CA002327750A Abandoned CA2327750A1 (en) 2000-10-31 2000-12-07 Real-time image reconstruction for computed tomography systems

Country Status (1)

Country Link
CA (1) CA2327750A1 (en)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110168612A (en) * 2017-01-06 2019-08-23 皇家飞利浦有限公司 The reconstruction control of standardized uptake values (SUV) guidance of result robustness for improving in positron emission tomography (PET) imaging
CN110163821A (en) * 2019-05-09 2019-08-23 安徽大学 A kind of terrible imaging method and system of singular value decomposition compression based on nonnegativity restrictions
CN110176046A (en) * 2019-05-31 2019-08-27 天津大学 The electrical resistance tomography method of constriction coefficient improvement Tikhonov regularization parameter
CN110533599A (en) * 2019-06-13 2019-12-03 中国科学院合肥物质科学研究院 The method for improving polluted gas concentration space distribution two-dimensional ct image reconstruction quality
CN110706298A (en) * 2019-09-10 2020-01-17 天津大学 Regularization weighted least square transmission-reflection dual-mode ultrasonic imaging reconstruction method

Cited By (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110168612A (en) * 2017-01-06 2019-08-23 皇家飞利浦有限公司 The reconstruction control of standardized uptake values (SUV) guidance of result robustness for improving in positron emission tomography (PET) imaging
CN110168612B (en) * 2017-01-06 2023-08-15 皇家飞利浦有限公司 Normalized uptake value directed reconstruction control for improved robustness of results in positron emission tomography imaging
CN110163821A (en) * 2019-05-09 2019-08-23 安徽大学 A kind of terrible imaging method and system of singular value decomposition compression based on nonnegativity restrictions
CN110163821B (en) * 2019-05-09 2022-10-11 安徽大学 Singular value decompression ghost imaging method and system based on non-negative constraint
CN110176046A (en) * 2019-05-31 2019-08-27 天津大学 The electrical resistance tomography method of constriction coefficient improvement Tikhonov regularization parameter
CN110176046B (en) * 2019-05-31 2023-04-07 天津大学 Electrical resistance tomography method with shrinkage coefficient improved Tikhonov regularization parameter
CN110533599A (en) * 2019-06-13 2019-12-03 中国科学院合肥物质科学研究院 The method for improving polluted gas concentration space distribution two-dimensional ct image reconstruction quality
CN110533599B (en) * 2019-06-13 2022-09-13 中国科学院合肥物质科学研究院 Method for improving reconstruction quality of two-dimensional tomographic image of polluted gas concentration spatial distribution
CN110706298A (en) * 2019-09-10 2020-01-17 天津大学 Regularization weighted least square transmission-reflection dual-mode ultrasonic imaging reconstruction method
CN110706298B (en) * 2019-09-10 2023-04-07 天津大学 Regularization weighted least square transmission-reflection dual-mode ultrasonic imaging reconstruction method

Similar Documents

Publication Publication Date Title
EP1938276B1 (en) Distributed iterative image reconstruction
Levkovilz et al. The design and implementation of COSEN, an iterative algorithm for fully 3-D listmode data
US20090123048A1 (en) Image Reconstruction Methods Based on Block Circulant System Matrices
Qi A unified noise analysis for iterative image estimation
Beekman et al. Selection of task-dependent diffusion filters for the post-processing of SPECT images
Van Slambrouck et al. Reconstruction scheme for accelerated maximum likelihood reconstruction: The patchwork structure
Chen et al. Bayesian image reconstruction in positron emission tomography
US20030190065A1 (en) Fast iterative image reconstruction from linograms
US6269323B1 (en) Support vector method for function estimation
CA2327750A1 (en) Real-time image reconstruction for computed tomography systems
CA2324832A1 (en) Real-time image reconstruction for computed tomography systems
Zhang et al. Deep generalized learning model for PET image reconstruction
EP4080459A1 (en) Machine learning-based improvement in iterative image reconstruction
CA2360656A1 (en) Real-time image reconstruction for computed tomography systems
Bouwens et al. Resolution recovery for list-mode reconstruction in SPECT
Valencia Pérez et al. Parallel approach to tomographic reconstruction algorithm using a Nvidia GPU
Koščević et al. Extra-low-dose 2D PET imaging
Houimli et al. Fast and Accurate Thoracic SPECT Image Reconstruction
Ote et al. List-Mode PET Image Reconstruction Using Dykstra-Like Splitting
Lesonen et al. Anatomy-guided multi-resolution image reconstruction in PET
Li et al. ESR-Net: An Efficient Image Super-resolution Network for SPECT Reconstruction
Labbé et al. Determination of a Set of Existing Algorithms For Comparative Evaluation
Yan et al. Multi-GPU cone-beam CT reconstruction
Yu et al. PET image reconstruction: GPU-accelerated particle filter framework
Wang et al. Fast reconstruction for unconstrained cone beam tomosynthesis

Legal Events

Date Code Title Description
FZDE Discontinued