The invention relates generally to techniques for nuclear medical imaging, and in particular to techniques for measuring and correcting for variations in relative efficiencies of Positron Emission Tomography (PET) detectors.
Nuclear medicine is a unique specialty wherein radiation is used to acquire images which show the function and anatomy of organs, bones or tissues of the body. The technique of acquiring nuclear medicine images entails first introducing radiopharmaceuticals into the body, either by injection or ingestion. These radiopharmaceuticals are preferentially taken up by specific organs, bones or tissues of interest (these exemplary organs, bones, or tissues are also more generally referred to herein using the term “objects”.) Upon arriving at their specified area of interest, the radiopharmaceuticals produce gamma photon emissions which emanate from the body and are then captured by a scintillation crystal. The interaction of the gamma photons with the scintillation crystal produces flashes of light which are referred to as “events.” Events are detected by an array of photo detectors (such as photomultiplier tubes) and their spatial locations or positions are then calculated and stored. In this way, an image of the organ or tissue under study is created from detection of the distribution of the radioisotopes in the body.
One particular nuclear medicine imaging technique is known as positron emission tomography, or PET. PET is used to produce images for diagnosing the biochemistry or physiology of a specific organ, tumor or other metabolically active site. The measurement of tissue concentration using a positron emitting radionuclide is based on coincidence detection of the two gamma photons arising from a positron annihilation. When a positron is annihilated by an electron, two 511 keV gamma photons are simultaneously produced and travel in approximately opposite directions. Gamma photons produced by an annihilation event can be detected by a pair of oppositely disposed radiation detectors capable of producing a signal in response to the interaction of the gamma photons with a scintillation crystal. Annihilation events are typically identified by a time coincidence between the detection of the two 511 keV gamma photons in the two oppositely disposed detectors (i.e., the gamma photon emissions are detected virtually simultaneously by each detector). When two oppositely disposed gamma photons each strike an oppositely disposed detector to produce a time coincidence event, they also identify a line(s)-of-response (LOR) along which the annihilation event has occurred. An example of a PET method and apparatus is described in U.S. Pat. No. 6,858,847, which patent is incorporated herein by reference in its entirety.
After being sorted into parallel projections, the LOR defined by the coincidence events are used to reconstruct a three-dimensional distribution of the positron-emitting radionuclide within the patient. In two-dimensional PET, each 2D transverse section or “slice” of the radionuclide distribution is reconstructed independently of adjacent sections. In fully three-dimensional PET, the data are sorted into sets of LOR, where each set is parallel to a particular detector angle, and therefore represents a two dimensional parallel projection p(r, s, φ, θ) of the three dimensional radionuclide distribution within the patient—where “r” and “s” correspond to the radial and axial distances, respectively, of the LOR from the center of the projection view and “φ” and “θ” correspond to the azimuthal and polar angles, respectively, of the projection direction with respect to the z axis in (x, y, z) coordinate space (in other words, φ and θ correspond to a particular LOR direction).
Coincidence events are integrated or collected for each LOR and stored in a sinogram. In this format, a single fixed point in the emitter distribution f(x, y) traces a sinusoid in the sinogram. Each row of a sinogram contains the LOR data for a particular azimuthal angle φ; each element of the row corresponds to a distinct radial offset of the LOR from the center of rotation of the projection. Different sinograms may have corresponded to projections of the tracer distribution at different coordinates along the scanner axis and/or different polar angles with respect to the scanner's axis.
To accurately reconstruct PET data into a usable image, one must know the efficiencies of the detectors that collected them, in order to compensate for their variability. This is often done by placing a thin planar emission source producing negligible scatter in the scanner and comparing the measured responses along lines of response (LOR) normal to the source to the expected uniform responses. Clinical emission data inevitably include scattered radiation along with the unscattered true coincidences (“trues”), however. Scattered radiation and true coincidences are shown in
If M represents a measured response for a detector pair forming an LOR, then neglecting random coincidences, M can be expressed in terms of the incident trues (T) and scattered (S) radiation as:
M=ε
T
T+∫ε(E,φ)S(E,φ)dEdφ
where εT and ε(E,φ) are the detection efficiencies for trues and scatter, respectively, and E and φ are the energy and incidence angle of the scattered radiation. Since there are two photons and two detectors involved for each LOR, the integral actually involves both energies and both incidence angles, but this is left implicit for notational simplicity.
We can write this equation in terms of the integral scatter flux S=∫(E,φ)dEdφ as
M=ε
T
T+ε
S
S
where the average scatter efficiency for the LOR is:
Because S(E,φ) is object dependent. εS may be also.
It is not obvious how one can determine εS directly since, for a source in the field of view (FOV) of the scanner, scatter is nearly always accompanied by true events. Ollinger proposed a partial solution to this problem in “Detector efficiency and compton scatter in full 3D PET,” IEEE Trans. Nuc. Sci., vol. 42, pp, 1168-1173, August 1995. He wrote the scatter detection efficiency (ηSijk in his notation) as the product of three components: the plane efficiency cSK, the detector efficiency εijk and the geometrical efficiency gijk:
εS=ηijkS=ckSεijkgijk.
He proposed to use the same values of εijk for scatter as for trues, since he had no way of independently measuring them. He proposed to estimate the plane efficiencies by comparing measured to computed scatter in the tails of a sinogram of a uniform cylindrical phantom, although the scattered radiation in the tails may have a different energy and angular distribution than the scatter in the LORs passing through the object. He assumed gijk was equal to 1 everywhere, based on heuristic arguments. Finally, Ollinger treated ηSijk as object independent, and thus no components of it were included in the scatter simulation.
The approach to scatter normalization is different in the single scatter simulation (SSS) algorithm, as described in “New, faster, image-based scatter correction for 3D PET,” C.C. Watson, IEEE Trans. Nuc. Sci., vol. 47, pp. 1587-1594, August 2000, the entirety of which is incorporated by reference. From (2) the normalized measured (randoms corrected) data is:
so the trues may be estimated by
Therefore scatter correction does not require independent knowledge of εS and εT, but only of their ratio. This is an advantage since first order effects in the variation of crystal efficiencies cancel out. The SSS algorithm includes an estimate of the detection efficiency of each simulated photon as a function of its energy and incidence angle. It also estimates the trues efficiency and forms their ratio so that the output of the SSS algorithm is an estimate of
for the object-dependent scatter flux.
To the extent that this internal model is accurate, no separate normalization of the scatter is required. This model explicitly accounts for possible variations in the geometrical efficiency for scatter, and allows for its object dependence. There are limitations to this model, however, because it does not include an exact description of the detectors' structure. For example, although most modern PET scanners use pixelated block detectors, gaps that may exist between detector blocks are neglected in the model. These gaps may result in undesirable artifacts in the ultimate image.
It would be useful, therefore, to have a means of measuring and correcting for variations in the εs/εT ratio which is not accounted for in the simulation. Although there may be random residual variations in individual crystal intrinsic efficiency ratios, the focus is on the variations associated with the geometrical structure of the detectors, i.e. those that depend only on the location of the crystals within the block detector, and the location of the block detectors within the scanner. For a uniform ring scanner only one azimuthal projection angle needs to be considered, due to symmetry.
The disclosed method estimates such a correction by making two measurements of a known object, one without significant scatter present like a standard plane source, and a second one that includes scattering media.
A method is disclosed for processing medical imaging emission data, the method comprising: obtaining data representative of a targeted patient tissue region using a positron emission tomography (PET) device; processing the data by applying a scatter correction to the data to obtain processed image data, where the scatter correction comprises applying a measured correction factor for the ratio of scatter detection efficiency to trues detection efficiency of at least one detector of the PET device to a numerical estimate of the scatter. and either: (a) subtracting this normalized estimate of the scatter from the data, or (b) otherwise using the normalized estimate of the scatter in an image reconstruction process; and representing an image of the targeted patient tissue region on a display, wherein the image is obtained using the processed image data.
A system is disclosed for processing medical imaging emission data, the system comprising an imaging system including a PET device having a plurality of detectors, a processing unit connected to and in communication with the imaging system and display, and a machine-readable storage medium encoded with a computer program code such that, when the computer program code is executed by a processor, the processor performs a method comprising: obtaining data representative of a targeted patient tissue region using the PET device: processing the data by applying a scatter correction to the data to obtain processed image data, where the scatter correction comprises applying a measured correction factor for the ratio of scatter detection efficiency to trues detection efficiency of at least one detector of the plurality of detectors to a numerical estimate of the scatter, and either: (a) subtracting this normalized estimate of the scatter from the data, or (b) otherwise using it in the reconstruction process.
The accompanying drawings illustrate preferred embodiments of the disclosure so far devised for the practical application of the principles thereof, and in which:
A system and method are disclosed for normalizing measured data from PET devices in which the normalization takes into consideration block structure and/or gaps inherent in PET devices. The disclosed system and method measure scatter detection efficiency variations due to detector geometry or other reasons, and compensate for such scatter efficiency variations during image processing.
As previously noted, to compensate for the variability between individual detectors of a PET device, the efficiencies of the detectors are needed. Thus, a thin planar emission source 10 producing negligible scatter (
Let M1 be the measurement with bare plane source, and M2 be the measurement with plastic blocks added. If the blocks have a total thickness of dr and a linear attenuation coefficient of μ, then for an LOR normal to the plane source 10,
M1=εTT
M
2=εTe−μΔχT+εSS.
If the attenuation factor e−μΔχ is accurately known, it can be applied to M1 and the product divided into M2 to give
This is just the scatter/trues ratio for the slab phantom scan. Since e−μΔχT should be constant over the normal projection for this phantom, this ratio should be proportional to the SSS estimate of scatter for this projection. Any deviations between the two, aside from an overall scale factor representing e−μΔχT, could be attributed to unmodeled variation in the εS/εT ratio. In fact, there may be other inaccuracies in the simulation as well, but because of the nature of the physical models used in SSS, such errors are expected to be spatially smooth at the scale of the block detectors. Because the primary concern is with those deviations due to block-level variations in detector efficiencies rather than longer-scale modeling errors, it should be adequate to use a simple analytical function such as a polynomial to model the overall trend of the scatter locally. A scatter/trues efficiency correction factor can be determined from the ratio:
To the extent that there are block effects on the scatter/trues efficiency ratio, εC should differ from 1. Once εC is determined in this way, it can be applied as an independent “normalization” of the SSS estimated scatter, so that
Of course, it may be that Cc itself is object dependent to some extent, but it is believed that this is a third-order effect.
Initial measurements were made on a prototype machine having six 38 cm diameter rings of 32 block detectors, each 3×3 cm in dimension, consisting of a 12×12 array of 2.5×2.5×20 mm LSO crystals. There were gaps between the blocks of about 7 mm transaxially and 2 mm axially. Two plane source scans were performed, one with and one without the plastic blocks. The two plastic slabs were PMMA, each 2.5 cm thick. PMMA has a density of 1.19 g/cc and a linear attenuation coefficient at 511 keV of 0.111 cm−1. The total attenuation factor for the slabs is thus 0.574. The energy window was 400-650 keV.
Two approaches were used for analyzing the data. In the first, the data were histogrammed into sinograms with 256 radial bins and 192 projection angles. To reduce noise, the sinograms were summed over all planes and oblique segments. As a consequence, only average transaxial variations in the efficiencies could be estimated. In the second approach the physical LOR data were analyzed before they were combined into the smaller number of virtual LORs in the sinograms. Substantial noise reduction was achieved by exploiting data symmetries. First the data was averaged over all detector block pairs having the same radial and axial offsets.
A further significant noise reduction came from assuming that there are only a few unique LORs within a block-block pair, determined by the location of the two crystals within the block. For the results described below only corner, edge, and central crystals were distinguished, leading to only six unique crystal-crystal, or LOR, types among the 1442=20736 possible pairs for two 12×12 blocks. The scatter/trues ratio of:
for each of these six unique LOR types was averaged over all similar LORs within each unique block pair. Finally, εC was estimated for each LOR by normalizing it to the center-center LOR type, effectively assuming εC=1 for these crystal pairs. These LOR-based εC maps were then sorted into sinograms for application to scatter correction. To compare the sinogram based and LOR based estimations of εC both analyses were performed on a second experimental data set, similar to the one described above, but in this case with a total thickness of PMMA slabs of 10 cm rather than 5 cm for the slab phantom scan.
Although the plane source used was weak and count rates low, they were necessarily different for the two scans, so the possibility of count rate induced variations in detector response which could possible confound the interpretation of the data in terms of scatter/trues efficiency variations cannot be excluded. To address this issue, a Monte Carlo simulation for this scanner geometry, including a detector response model incorporating block gaps, was performed for a 20 cm cylindrical activated water phantom, using Geant4 (for GEometry ANd Tracking). Scattered and unscattered events were tallied separately, so that their detection efficiencies could be compared, similar to the plane source measurements.
The impact of the scatter/trues efficiency compensation on images was evaluated on synthetic data to eliminate effects of noise and other artifacts. A 20 cm uniform 68Ge cylinder centered in the BrainPET was modeled, and the scatter sinogram for this object was computed using the standard SSS algorithm. This scatter sinogram was attenuation corrected and reconstructed using filtered back-projection (FBP) both with and without applying the εC compensation to the scatter first. Gaps in the εC sinogram were filled with the value 1. Because FBP is a linear algorithm, and the trues and scatter are additive, scatter correction may be equivalently performed in either projection or image space. Thus a synthetic image can be estimated containing the uncorrected effects of εC from
where a represents the attenuation correction factors, Iideal is the ideal uniform cylinder image, and ICC is the image that would be expected if block-related scatter/trues efficiency variations were present in the data, but not compensated for during scatter correction. This image was computed using the LOR-based estimate of εC.
The summed sinograms for the sinogram-based estimate are shown in
Only the central portion of this sinogram, where the LORs are nearly normal to the plane source, is of interest.
The correction factor is greater than 1.2 in the vicinity of the edges of the central block-block coincidence diamond. The simulated profile (light grey lines) is derived in a similar way from the Monte Carlo simulation, and is consistent with the results from the measured data. The slight difference in the block gap position at large radius is due to a difference in the way the measured and simulated data were histogrammed. The simulation results are shown in more detail in
In
The model of an uncompensated image, ICC, averaged over all planes is shown in
The sinogram profiles in
The agreement between the measured and simulation derived scatter/trues efficiency correction factors in
From the image example of
Referring now to
At step 210 (
The correction factor may be determined from the formula:
where: εC—scatter/trues efficiency correction factor; M1—measurement with plane source; M2—measurement with scatter media added; Δx—total thickness of the scatter media; μ—linear attenuation coefficient of the scatter media; εS—detection efficiency of scatter; εT detection efficiency of trues; S—incident scatter radiation; and T—incident trues radiation.
The correction factor may comprise a sinogram-based estimate obtained by histogramming the data obtained from the first and second measurements into sinograms and summing the sinograms. In addition, the correction factor may be obtained by replicating a portion of the summed sinogram to obtain a quasi-uniform sinogram.
Alternatively, the correction factor may comprise a line of response (LOR)-based estimate obtained from LOR data obtained from the first and second measurements. The LOR data may be averaged over all detector pairs of the PET device having the same radial and axial offsets.
The correction factor may comprise a sinogram applied to the data as an independent normalization of a single scatter estimation (SSS) scatter sinogram. The correction factor may compensate for scatter detection efficiency variations due to respective gaps between a plurality of detectors of the PET device.
The method described herein may be automated by, for example, tangibly embodying a program of instructions upon a computer readable storage media capable of being read by machine capable of executing the instructions. A general purpose computer is one example of such a machine. A non-limiting exemplary list of appropriate storage media well known in the art would include such devices as a readable or writeable CD, flash memory chips (e.g., thumb drives), various magnetic storage media, and the like.
The features of the method have been disclosed, and further variations will be apparent to persons skilled in the art. Such variations are considered to be within the scope of the appended claims. Reference should be made to the appended claims, rather than the foregoing specification, as indicating the true scope of the disclosed method.
The functions and process steps herein may be performed automatically or wholly or partially in response to user command. An activity (including a step) performed automatically is performed in response to executable instruction or device operation without user direct initiation of the activity.
The systems and processes disclosed herein are not exclusive. Other systems, processes and menus may be derived in accordance with the principles of the invention to accomplish the same objectives. Although this invention has been described with reference to particular embodiments, it is to be understood that the embodiments and variations shown and described herein are for illustration purposes only. Modifications to the current design may be implemented by those skilled in the art, without departing from the scope of the invention. The processes and applications may, in alternative embodiments, be located on one or more (e.g., distributed) processing devices accessing a network linking the elements of
This is a non-provisional of pending U.S. Provisional Patent Application Ser. No. 61/100,921, filed Sep. 29, 2008, by Charles C. Watson, et al., the entirety of which application is incorporated by reference herein.
Number | Date | Country | |
---|---|---|---|
61100921 | Sep 2008 | US |