Systems and methods for attenuation compensation in nuclear medicine imaging based on emission data

Systems and methods for attenuation compensation in nuclear medicine imaging based on emission data are provided. One method includes acquiring emission data at a plurality of energy windows for a person having administered thereto a radiopharmaceutical comprising at least one radioactive isotope. The method also includes performing a preliminary reconstruction of the acquired emission data to create one or more preliminary images of a peak energy window and a scatter energy window and determining a body outline of the person from at least one of the reconstructed preliminary image of the peak energy window or of the scatter energy window. The method further includes identifying a heart contour and segmenting at least the left lung. The method additionally includes defining an attenuation map based on the body outline and segmented left lung and reconstructing an image of a region of interest of the person using an iterative joint estimation reconstruction.

BACKGROUND OF THE INVENTION

The subject matter disclosed herein relates generally to nuclear medicine imaging systems, and more particularly to single photon emission computed tomography (SPECT) imaging systems and compensating for emission attenuation in SPECT systems, especially in cardiac imaging, using emission data.

Different types of imaging techniques are known and used for medical diagnostic imaging. For example, diagnostic nuclear imaging, such as SPECT imaging, is used to study radionuclide distribution in a subject, such as a patient. Typically, one or more radiopharmaceuticals or radioisotopes are injected into the patient. Gamma camera detector heads, typically including a collimator, are placed adjacent to a surface of the patient to capture and record emitted radiation to thereby acquire image data. Different configurations are known wherein the gamma cameras may remain in a fixed location/orientation (e.g., focused detector modules) relative to an object of interest during a scan or may be rotated about the patient. Image reconstruction techniques, such as backprojection, may then be used to construct images of radiotracer uptake distribution within internal structures of the subject based upon the acquired image or acquired data, such as list data.

While such conventional systems may provide quality images with good diagnostic value, photon attenuation is a major physical factor affecting the quality of reconstructed images in SPECT systems. Such attenuation may occur, for example, due to tissues between the sources of emissions and the system detectors. However, in SPECT imaging, and specifically in cardiology, it is important to obtain an accurate emission image (a three-dimensional 3D map of the radioisotope distribution within the imaged patient) in the presence of attenuation (in large part due to Compton scattered radiation) caused by the patient's body.

In cardiac imaging, photon attenuation accounts for up to 85% loss of emitted photons from the myocardium area. Moreover, data inconsistencies with models used in image reconstruction from a quantitative point of view are also spatially variant (e.g., 70-85% error within myocardium only in some cases). Thus, known reconstruction methods require knowledge of the attenuation map, for example, the 3D model of the patient tissue in areas affecting the radiation arriving at the detector. These methods currently mostly rely on direction transmission measurements that may include a radioactive source that is often ineffective or measurements from an x-ray computed-tomography (CT) system that are costly, as well as can add radiation dose to the patient, additional imaging time, geometrical mis-registration and resolution differences. Models may be used to characterize the attenuation, although actual attenuation may differ substantially. Moreover, because of the high variability of patient sizes and shapes, a “patient standard” can yield a poor reconstruction result.

BRIEF DESCRIPTION OF THE INVENTION

In accordance with an embodiment, a method for image reconstruction is provided. The method includes acquiring emission data at a plurality of energy windows for a person having administered thereto a radiopharmaceutical comprising at least one radioactive isotope, wherein the energy windows comprise (i) at least a peak energy window centered around a peak emission of the isotope and (ii) at least one scatter energy window at an energy range lower than the peak energy window. The method also includes performing a preliminary reconstruction of the acquired emission data to create one or more preliminary images of the peak energy window and the scatter energy window and determining a body outline of the person from at least one of the reconstructed preliminary image of the peak energy window or the reconstructed preliminary image of the scatter energy window. The method farther includes identifying a heart contour of the person from the reconstructed preliminary image of the peak energy window and segmenting at least the left lung of the person from the reconstructed preliminary image of the scatter energy window using the identified heart contour as a landmark. The method additionally includes defining an attenuation map based on at least the determined body outline and the segmented left lung and reconstructing an image of a region of interest of the person using an iterative joint estimation reconstruction including updating the attenuation map and the image of the peak energy window, wherein the joint estimation reconstruction comprises using data acquired in the plurality of energy windows.

In accordance with another embodiment, a nuclear medicine (NM) imaging system is provided that includes a gantry and a plurality of nuclear medicine (NM) cameras coupled to the gantry and configured to acquire emission data at a plurality of energy windows for a person having administered thereto a radiopharmaceutical comprising at least one radioactive isotope, wherein the energy windows comprise (i) at least a peak energy window centered around a peak emission of the isotope and (ii) at least one scatter energy window at an energy range lower than the peak energy window. The NM imaging system also includes an image reconstruction module configured to (i) perform a preliminary reconstruction of the acquired emission data to create one or more preliminary images of the peak energy window and the scatter energy window, (ii) determine a body outline of the person from at least one of the reconstructed preliminary image of the peak energy window or the reconstructed preliminary image of the scatter energy window, (iii) identify a heart contour of the person from the reconstructed preliminary image of the peak energy window, (iv) segment at least the left lung of the person from the reconstructed preliminary image of the scatter energy window using the identified heart contour as a landmark, (v) define an attenuation map based on at least the determined body outline and the segmented left lung and (vi) reconstruct an image of a region of interest of the person using an iterative joint estimation reconstruction including updating the attenuation map and the image of the peak energy window, wherein the joint estimation reconstruction comprises using data acquired in the plurality of energy windows.

In accordance with yet another embodiment, a non-transitory computer readable storage medium for performing image reconstruction using a processor is provided. The non-transitory computer readable storage medium includes instructions to command the processor to acquire emission data at a plurality of energy windows for a person having administered thereto a radiopharmaceutical comprising at least one radioactive isotope, wherein the energy windows comprise (i) at least a peak energy window centered around a peak emission of the isotope and (ii) at least one scatter energy window at an energy range lower than the peak energy window. The non-transitory computer readable storage medium also includes instructions to command the processor to perform a preliminary reconstruction of the acquired emission data to create one or more preliminary images of the peak energy window and the scatter energy window and determine a body outline of the person from at least one of the reconstructed preliminary image of the peak energy window or the reconstructed preliminary image of the scatter energy window. The non-transitory computer readable storage medium also includes instructions to command the processor to identify a heart contour of the person from the reconstructed preliminary image of the peak energy window, segment at least the left lung of the person from the reconstructed preliminary image of the scatter energy window using the identified heart contour as a landmark and define an attenuation map based on at least the determined body outline and the segmented left lung. The non-transitory computer readable storage medium further includes instructions to command the processor to reconstruct an image of a region of interest of the person using an iterative joint estimation reconstruction including updating the attenuation map and the image of the peak energy window, wherein the joint estimation reconstruction comprises using data acquired in the plurality of energy windows.

DETAILED DESCRIPTION OF THE INVENTION

Described herein are systems and methods to determine and compensate for attenuation within nuclear medicine imaging systems, in particular, single photon emission computer tomography (SPECT) imaging systems. Various embodiments use only emission data acquired by the SPECT system to estimate and compensate for attenuation.

FIG. 1is a block diagram of an exemplary nuclear medicine imaging system20constructed in accordance with various embodiments, which in this embodiment is SPECT imaging system. The system20in one embodiment includes an integrated gantry22that further includes a rotor24oriented about a gantry central bore26. The rotor24is configured to support one or more nuclear medicine (NM) cameras28and30(two are shown for illustration). In various embodiments the NM cameras28and30may be, for example, general purpose gamma cameras or non-general purpose gamma cameras, such as focused pinhole gamma camera modules configured for cardiac imaging. The NM cameras28and30may be formed from different types of suitable material, which may be direct conversion materials or indirect conversion materials, which may be pixelated detectors or cameras. For example, in indirect conversion material, the scintillator, which is typically made of a crystalline material, such as sodium iodide (NaI), converts the received gamma radiation to lower energy light energy (e.g., in an ultraviolet range). In these systems, photomultiplier tubes then receive this light and generate image data corresponding to photons impacting specific discrete picture element (pixel) regions. In direct conversion material, such as cadmium zinc telluride (CZT), the impinging photons are converted directly into electrical signals.

The rotor24is further configured to rotate axially about an examination axis including a patient table34that may include a bed that is slidingly coupled to a bed support system to support a patient36, which may be coupled directly to a floor or may be coupled to the gantry22through a base coupled to the gantry22. The bed may include a stretcher slidingly coupled to an upper surface of the bed. The patient table34is configured to facilitate ingress and egress of the patient36into an examination position that is substantially aligned with the examination axis. During an imaging scan, the patient table34may be controlled to move the bed and/or stretcher axially into and out of the bore26. The operation and control of the imaging system20may be performed in any manner known in the art. It should be noted that the various embodiments may be implemented in connection with imaging systems that include rotating gantries or stationary gantries.

A collimator38may be provided in combination with the NM cameras28and30. For example, a collimator38may be coupled to front detecting faces of each of the NM cameras28and30. The collimators38may be any suitable type of collimator known in the art.

The outputs from the NM cameras28and30are communicated to a processing unit40, which may be any suitable computer or computing device. As used herein, the term “computer” or “module” may include any processor-based or microprocessor-based system including systems using microcontrollers, reduced instruction set computers (RISC), ASICs, logic circuits, and any other circuit or processor capable of executing the functions described herein. The above examples are exemplary only, and are thus not intended to limit in any way the definition and/or meaning of the term “computer”.

The processing unit40may include an attenuation compensation module50to perform attenuation compensation as described in more detail herein. The attenuation compensation module50may be implemented in hardware, software, or a combination of hardware and software.

It should be noted that the imaging system20may also be a multi-modality imaging system, such as an NM/MR imaging system. During an imaging scan, the patient table34may be controlled by a table controller unit44that is part of a controller42. The table controller unit44may control the patient table34to move the patient table34axially into and out of the bore26. The NM cameras28and30may be located at multiple positions (e.g., in an L-mode configuration) with respect to the patient36. It should be noted that although the NM cameras28and30are configured for movable operation along (or about) the gantry22, the NM cameras28and30may be fixed thereto. The controller42also includes a gantry motor controller46that controls movement of the gantry22, for example, rotational movement about the patient or movement of the NM cameras28and30, such as pivoting movement or movement towards/away from the patient36.

Thus, the controller42may control the movement and positioning of the patient table34with respect to the gamma cameras28and30and the movement and positioning of the NM cameras28and30with respect to the patient36to position the desired anatomy (e.g., organ(s)) of the patient36within the fields of view (FOVs) of the NM cameras28and30, which may be performed prior to acquiring an image of the organ of interest. The table controller44and the gantry motor controller46each may be automatically commanded by the processing unit40, manually controlled by an operator, or a combination thereof. The imaging data may be combined and reconstructed into an image as described in more detail below, which may comprise 2D images, a 3D volume or a 3D volume over time (4D).

A Data Acquisition System (DAS)48receives analog and/or digital electrical signal data produced by the NM cameras28and30and decodes the data for subsequent processing as described in more detail herein. An image reconstruction processor52receives the data from the DAS48and reconstructs an image using any reconstruction process known in the art with attenuation compensation as described herein. A data storage device54may be provided to store data from the DAS48or reconstructed image data. An input device56also may be provided to receive user inputs and a display58may be provided to display reconstructed images.

In operation, prior to data collection, a radioisotope, such as a radiopharmaceutical substance (sometimes referred to as a radiotracer), is administered to the patient36, and may be bound or taken up by particular tissues or organs. Typical radioisotopes include various radioactive forms of elements, although many in SPECT imaging are based upon an isotope of technetium (99Tc) that emits gamma radiation during decay. Various additional substances may be selectively combined with such radioisotopes to target specific areas or tissues of the body.

Gamma radiation emitted by the radioisotope, temporarily present at a location within the patient is detected by the NM cameras28and30. Although the NM cameras28and30are illustrated inFIG. 1as planar devices positioned above the patient36, the NM cameras28and30may be positioned below the patient36, both above and below the patient36, next to the patient36and may wrap at least partially around the patient36.

The imaging system20in some embodiments may be coupled to one of more networks to allow for the transfer of system data to and from the imaging system20, as well as to permit transmission and storage of image data and processed images. For example, a local area networks, wide area networks, wireless networks, and so forth may allow for storage of image data on radiology department information systems or on hospital information systems. Such network connections further allow for transmission of image data to remote post-processing systems, physician offices, and so forth.

The various embodiments described herein may be used, for example, in conjunction with dedicated SPECT systems for imaging particular organs of interest, such as for cardiac imaging and evaluation. Some of these systems are characterized by a limited field of view (FOV) aimed to contain the organ of interest and/or non-parallel collimation. Such systems are sometime referred to as “shift variant” imaging systems. Here, shift variance means that system response to an object (e.g., point source) differs depending on location of the object in the FOV. Among the differences are differences in geometrical shape of system response, system sensitivity, and attenuation path from an emitting object to the system detectors.

The various embodiments may be used in connection with different SPECT imaging configurations, such as shown inFIGS. 2 and 3. In the embodiment ofFIG. 2, an imaging system scanner60comprises the gantry22that supports a dual head camera (shown in an L-configuration). The camera comprises two camera sections, illustrated as the NM cameras28and30disposed in the gantry22to acquire data over approximately 180 degrees around the patient36. In the embodiment illustrated, the imaging system scanner60is configured for cardiac imaging, and the embodiments described herein allow for characterization and correction of scatter and attenuation of emissions62originating at locations in and around the heart64. In general, such emissions will traverse at least some regions of the heart64, as well as soft tissues66of the body, and particularly the left lung68. It should be noted that the gamma cameras NM and30may be mounted to the gantry22with mounting mechanisms70that allow for movement in addition to about the patient36, such as pivoting movement or translation towards or away from the patient36.

FIG. 3illustrates another configuration of an imaging scanner80allowing that defines a multi-pinhole acquisition system that include a plurality of modules82, which in this embodiment are pinhole gamma camera modules. The modules82are positioned and oriented around the patient volume to collect emissions62that traverse similar tissues of the patient36. It should be noted that in the case of pinhole acquisition systems, the pinholes of the modules82may be adjusted such that the pinholes are focused on the volume on interest and may be stationary during image acquisition, for example, the Discovery NM 530 c available from GE Healthcare.

It should also be noted that other types and configurations of cameras may be employed, such as a camera of the type disclosed in U.S. Pat. No. 6,242,743.

The modules82may take different forms as illustrated inFIGS. 4 and 5. For example, the pinhole configuration as illustrated inFIG. 4includes a detector90having pinhole collimator92in combination therewith, for example, coupled to a detecting face94of the detector90. The module82may pivot or rotate.

The module82may also include different types of collimation, such as a parallel hole collimator96as shown inFIG. 5. However, other types of collimation may be used including diverging and converging types of collimation as known in the art. In the embodiment ofFIG. 5, a pivot98is provided.

FIG. 6illustrates a more detailed illustration of the tissues traversed by the radiation emissions in these scenarios. The body of the patient36will extend to a skin-air boundary and have a general contour100in the imaging volume from which SPECT data is acquired. Within the body and in the case of cardiac imaging, the heart64will have a contour102that represents the boundary or transition between the tissues of the heart64and those of surrounding anatomies. The lung68(which is the left lung) will have a further contour104representing the extent and the transition between the lung tissues and those surrounding tissues. During SPECT imaging data acquisition, emissions62may radiate in all directions and traverse some or all of these tissues and be scattered and attenuated differently by each. For example, certain lines of radiation106may traverse cardiac and soft tissues only, while other lines of radiation108may traverse cardiac tissues and lung tissues with little soft tissue therebetween, followed by soft tissue up to the body contour. Still further lines of radiation110(direct radiation) may traverse cardiac tissues, soft tissues, and then further traverse the lung and more soft tissue before exiting the body. Some radiation112may scatter as well, wherein an electron114in the tissue scatters the gamma. The various embodiments use these contours for characterizing the scatter and/or attenuation of the emissions for image data processing and image reconstruction.

Various embodiments provide attenuation correction in SPECT using only emission data. A process flow120is shown inFIG. 7that generally illustrates attenuation compensation performed in some embodiments. A detailed description will then follow. The process flow120includes performing a preliminary reconstruction122based on acquired emission data. The acquired emission data that is used in the process flow120includes a full spectrum of energies in various embodiments as shown inFIG. 8to create at least preliminary images of the main emission or peak energy window and the one or more scatter windows. Thus, emission data is acquired for a patient at a plurality of energy windows or levels (e.g., list mode data) such that multiple energy windows can be retrospectively defined. For example, a peak energy window, illustrated as the main energy window140corresponding to the dominant energy peak in the energy response144and one or more scatter energy windows142corresponding to lesser or no energy peaks in the energy response144are acquired. In various embodiments, the main energy window is generally centered around the peak emission of the isotope and the scatter energy windows are at energy ranges lower than the peak energy window. It should be noted the energy response144defines a profile that may identify photons with small angle scatter and large angle scatter. Thus, scatter information, such as in lower energy windows that detect scatter with different scatter angles may be used in various embodiments. It should be noted that scatter with a smaller scatter angle has a smaller deviation and smaller energy loss, while scatter with a larger scatter angle has a larger deviations and larger energy loss.

Referring again toFIG. 7, this preliminary reconstruction122is performed without correction for physical effects, namely no attenuation correction. This preliminary reconstruction122generally defines a boundary of interest, for example, the outer boundary of the patient. It should be noted that the preliminary reconstruction122may be performed using only the main energy window or optionally include scatter data from one or more of the scatter energy windows (to improve outer boundary detection). Thus, a rough estimation124of body contours is determined, for example, by segmenting the body outline using a reconstructed preliminary image of the peak energy window and optionally one or more scatter energy windows.

An identification126of the heart contour of the patient is then determined using from the reconstructed preliminary image. This identification126may be performed using any method known in the art. A segmentation128of at least the left lung is then performed, which may include using the identified heart contour to assist the lung segmentation. For example, the boundary between the left ventricle and the lung may be identified and then a seeding and growing process may be used to identify the boundary of the left lung. Thus, the left lung may be segmented from the reconstructed preliminary image of the scatter energy window using the identified heart contour as a landmark.

Binary maps generated from the rough estimation of the body contour and the segmented left lung filled with linear attenuation coefficients are then used as an input to a joint estimation reconstruction130. The inputs define an initial approximation or guess of the attenuation map, namely a preliminary attenuation map. The joint estimation reconstruction130is an iterative process wherein two updates are performed at each iteration. First, an estimate of the attenuation map is used to perform attenuation correction, which is then used to update the emission data. Thus, at each step, the emission estimate is updated based on the attenuation map from the previous iteration step, which is then used to update the attenuation map in the current step. The joint estimation reconstruction130is accordingly performed with attenuation and scatter compensation to generate a reconstructed image132.

Thus, in various embodiments, the body outline is identified using the peak energy window, the scatter energy window or a combination thereof, for example, a summation of the peak energy window and the scatter energy window. The heart contour is identified from the reconstruction “peak”. The lung(s) are identified from scatter data, such as using one or more scatter energy windows. Various operations or steps to identify the different landmarks and compensate for attenuation may be performed, for example, as described in more detail below.

More particularly, various embodiments provide a method150as shown inFIG. 9for attenuation compensation, particularly in SPECT imaging, especially cardiac SPECT imaging. The method150includes acquiring emission data for multiple energy windows at152, which is used to segment a body outline with a preliminary reconstruction at154. Thus, using emission data for multiple energy windows (peak energy window+scatter energy window(s)), a preliminary reconstruction without attenuation correction is performed. In one embodiment, a main emission or peak energy window reconstruction using any suitable SPECT reconstruction method may be used to determine a rough estimation of the body contour. As described in more detail herein, scatter energy window reconstruction may be used to supplement the main emission or peak energy window reconstruction. Thus, no x-ray CT data is used in the method150.

In some embodiments, additional projection views are acquired at152. For example, additional projection views are acquired from the supine direction to a standard 180 degree acquisition arc as shown inFIG. 10, such as using rotational SPECT. Thus, additional views may be acquired to resolve the body outline by rotating the NM cameras or detectors (e.g., the NM cameras28and30) additional gantry steps, for example, which in one embodiment is a distance about equal to the size of the NM camera or detector. For example,FIG. 10illustrates three gantry positioned for the NM cameras28and30. It should be noted that the NM cameras28and30are rotated through a plurality of gantry steps and only three are shown for illustration. As can be seen, locations1and2are part of the standard 180 degree acquisition while location3acquires additional views.

In some embodiments, for example in a focused collimation system170(having a focused geometry) as shown inFIG. 11, persistence data172may be used. The persistence data172is acquired, for example, during positioning of the patient and is not discarded in these embodiments. Thus, this persistence data172provides data similar to that of a scout data.

Referring again toFIG. 9, thereafter the heart contour is segmented at156from the reconstructed preliminary image of the peak energy window and optionally one or more scatter energy windows. The identification of the heart may be determined using any suitable known heart segmentation method. The heart contour optionally may be used as a landmark to determine an interface between the heart and the left lung to assist in lung segmentation at158. For example, the heart contour may be enclosed in an ellipsoid (e.g., graphical overlay) with the left lung identified at an interface using a known direction from the left ventricle of the heart.

Thus, at158, the left lung is segmented using a scatter window reconstruction, namely from the reconstructed preliminary image of the scatter energy window(s). For example, using scatter data, a scatter window reconstruction may be performed using a regular reconstruction, such as a main window reconstruction with straight line projections from the emission to the detector such that geometry changes are ignored. In some embodiments, a special projector for scatter reconstruction may be used such as a model, for example, a Monte-Carlo based method to model the scatter geometry to improve lung contrast. It should be noted that all or a subset of the voxels may be updated at a time.

Thus, a segmentation based determination of the lung may be used to obtain a binary map by using the voxel values and predetermined threshold values, such as to identify tissue. In some embodiments, the segmentation may be assisted by a knowledge set such as an a priori constructed lung model.

Then, an iterative joint estimation reconstruction using data from the determined contours is performed at160wherein each iteration includes two updates. In various embodiments, a preliminary attenuation map is defined based on the determined body outline and the segmented left lung. In particular, in each iteration the emission data is reconstructed with an attenuation correction estimate and the attenuation map is updated. For example, a maximum likelihood process may be used for the emission update and a conjugate gradient-like process may be used for the attenuation map update. It should be noted that in some embodiments specifically constructed priors (e.g., adding regularizations), such as joint entropy or other intermediate filters (based on neighbor voxels), or bi-normal distribution may be used to provide smoothness to the resulting images and form the images to develop desired properties.

Thus, various embodiments provide a reconstruction process in two main stages. First, an initial estimate of the attenuation map is created from a series of reconstruction and segmentation steps. Second, this estimate, along with SPECT emission projections, are used as an input (initial approximation) for an iterative joint estimation process, when SPECT data reconstruction with attenuation and, optionally, scatter compensation, and attenuation map estimate are interchangeably updated and refined until a pre-defined criterion is met. It should be noted that the various steps of the method150may be achieved in a single step, by reconstructing an optimized scatter window.

Reconstruction with attenuation and scatter compensation may assign a value to the scatter/attenuation for different trajectories through the tissues traversed by each trajectory. For example,FIG. 12is a diagrammatical illustration of such mapping in the case of cardiac imaging. In this illustration, the mapping180is compiled for anatomies of interest and shown disposed in a three-dimensional segment consisting of discrete volume elements or voxels182. Based upon the density and position of the various tissues, the voxels may indicate more or less scatter/attenuation. The mapping may be determined from the body and tissue contour and volume determinations as described in more detail herein, and used in the reconstruction of images from the acquired SPECT data.

In operation, and for example, the rough estimation of the body contour may be performed in various ways, depending on geometry of acquisition. For conventional rotational SPECT acquisitions, reconstruction of counts in the main emission window, and/or counts detected in the scatter window is used. For alternative geometries, such as characterized by stationary acquisition/limited acquisition arc/small field of view (FOV), additional data, such as scout or persistence data from remote detector positions may be added. In this case, several auxiliary views may be appended to the projection data corresponding to standard acquisition orbit/geometry and reconstructed together, to ensure full visibility of the body outline as described herein. This initial reconstruction will be then segmented into “body” and “outside air” classifications.

Following this body contour estimation phase, and in the case of cardiac imaging, a rough estimation of the left lung volume is performed. This estimation may be based upon “seeding” from the edge of the reconstructed left ventricular surface to provide an additional landmark for lung identification and segmentation. The data resulting from this phase defines an initial estimation of the attenuation maps. The attenuation map is reconstructed on the same voxel grid and volume as the emission data, and from the same data. Thus, the attenuation maps obtained from this process are intrinsically registered to the emission data.

For example, as illustrated inFIG. 13, a preliminary reconstruction of the acquired emission data may be performed to create at least preliminary images190and192of the peak energy window and scatter energy window, respectively. The heart may be identified as shown in the image194, including providing an overlay195such that the left ventricle contour of the heart is segmented. The arrows in image194represent lung search “seeds” to segment the left lung198as illustrated in the image196.

Thus, with the attenuation map determined, image reconstruction with compensation for scatter/attenuation of SPECT emissions may be performed. In various embodiments, the determination of the contours as shown, for example, inFIG. 13, provides the initial estimate of the attenuation map, which will be enhanced further in the process of joint estimation as described herein.

In operation, the second stage of the reconstruction process results in final SPECT images reconstructed with compensation for effects of attenuation and optionally scatter. In accordance with one embodiment, the reconstruction process is iterative, and may be provided as described below. As an initialization step, activity uptake distribution is assumed, in a standard manner, to be uniform in accordance with the relationship:
X(0)=M(0)=H′1pEq. 1
wherein SPECT emission projections, initial estimate of attenuation map and, optionally, scatter estimate (same volume as the volume of emission reconstruction) are the inputs.

Following initialization, joint estimation is performed. For each iteration of the joint estimation process, two subsequent updates are performed. First, activity concentration estimate x (SPECT reconstruction) is advanced following, for example, a conventional penalized likelihood scheme. In this update, a current estimate of the attenuation map xkis used. In the second update, the attenuation map estimate is refined using the just obtained activity concentration estimate. The update of the attenuation map, which does not obey Poisson statistics, is not driven by likelihood maximization, but by a general optimization scheme such as a coordinate descent. So, a single iteration of a joint estimation algorithm may be described by the relationships:

x(k+1)=arg⁢⁢minx⁢(L⁡(g,x)+β⁢⁢P⁡(x,m(k))Eq.⁢2m(k+1)=arg⁢⁢minm⁢(F(x(k+1),m))Eq.⁢3
where x represents the activity distribution and in the attenuation map.

From Equation 2, the value of x is updated in accordance with the relationship:

The values of m(k)are then thresholded in accordance with the relationship single threshold:

mj(k+1)={mj(k+1)=c⋒;mj(k+1)≥t(k+1)mj(k+1)=c⋓;elseEq.⁢6
(or using other segmentation techniques described known in the art) into air and soft tissue compartments and filled with linear attenuation coefficients from a look-up table, in accordance with the radiopharmaceutical used in the process of data formation.

Combining the stage of initial independent reconstruction of attenuation map with joint estimation, “cross-talk” is reduced or eliminated. “Cross-talk” appears when emission-specific features (e.g., myocardial uptake) are being propagated into the attenuation map.

Variations and modifications are contemplated. For example, in one embodiment the following process may be provided:

1. Data is to be collected in at least two energy windows (E0 and E1).E0 is the “peak energy window”, defined as “peak energy+/−dE wherein dE is usually a few percent;E1 is energy window within the scatter energy range wherein the energy range of the “scatter energy window” is below the “peak energy window”, and in various embodiments is wider than the “peak energy window”;Optionally more scatter energy windows are defined such as E2 (or more), that is, the “scatter energy window” is divided to two or more sub-windows, which in some embodiments are non-overlapping, contiguous sub-windows. 2. Defining the system response function—a function that allows estimating the data provided that the object is known. The system response depends on the collimator, the detector, etc.
B. Data Acquisition:Collecting emission data “e(P,E)”, the number of photons collected at detector position “P” at energy window “E”.Here, “E” are the energy windows (E=0, 1, 2, . . . )and “P” is the general designation of the detector position.1. In a multi-pinhole camera P={x,y,p} wherein x, y are the pixel indexes and p is the pinhole index.2. In a rotating SPECT camera, P={x,y,p} wherein x, y are the pixel indexes and p is the projection index (the gantry angle “f”).3. In a camera with a plurality of rotating, collimated heads (such as inFIG. 5), P={x,y,p, f} wherein x, y are the pixel indexes and p is the head index, and f is the angle of head p. Optionally, f is defined by f={fx,fy} if the heads can rotate in 2D.
It should be noted that the dimension of the dataset is x*y*p*e (e=2 for two energy windows: “e0=peak”; and “e1=scatter”, optionally e>2 if the scatter window is subdivided). It also should be noted that the solution is of the same (or lower) dimensionality.
C. Assumption and DefinitionsIt should be noted that a single isotope having a single emission peak is assumed, (e.g. Tc having a peak at 140 keV), but extension to multi-peak isotope or multi-isotope may be provided.1. The Object (patient)The object O(X)={S,D}(X) is defined by both its:“Source Concentration” S(x,y,z) S(X); (S is in “Curie per cc”, x,y,z=X are the voxel index in 3D). Typically, the dimensionality of x=y=z=64 or 128, the source concentration must be a non-negative numberand“body density” D(x,y,z)=D(X); (D in grams per cc). It should be noted that D(x,y,z) may be translated to “Absorption & Compton scattering coefficients” by the energy dependent nuclear parameters cross section parameters. Typically, a linear transformation is applied, but non-linear transformations (that take into account bone chemical composition) are also known.
Tissue Segmentation.To reduce the complexity of computation, the following algorithm may be applied to D(X):1. Air—For D˜0 (or below a threshold)—the tissue is assumed to be “Air”. Specifically if it can be located as being outside the patient boundaries. D is then set to D=02. Lungs—For D between “minimal lung density” and “Maximum lung density”—The tissue a mixture of air and soft tissue, with the required mixture percentage. Specifically if tissue can be located as being inside the patient lungs boundaries. D is left as a variable,3. Soft tissue—For D between “maximal lung density” and “minimal spongy bone density”—The tissue is assumed soft tissue. D is set to the soft tissue average value4. Spongy Bone—For D between “minimal spongy bone density” and “maximal spongy bone density”—The tissue a mixture of soft tissue and bone, with the required mixture percentage. D is left as a variable and the chemical composition is assumed to be the appropriate mixture of hard bone and soft tissue having density D.5. Hard bone—For D between “minimal hard bone density” and “maximal hard bone density”—The tissue a hard bone. D is set to density and hard bone.6. Metal—For D above “maximal hard bone density”—The tissue an implant or foreign object—operator intervention is required.This process is called “segmentation of the patient tissue” (the segments are: air, lungs, soft tissue, spongy bones, and hard bone. However, in some cases, only air, lungs, and soft tissue are considered and bones are ignored and replaced with soft tissue or with “dense soft tissue having same artificially high density).Segmentation reduces the computation as most of the volume comprises air or soft tissue.The simple (linear) transformation from D to absorption and scattering coefficients is:Absorption A(X,E)=a(E)*D(X); wherein “a” is the average absorption coefficient of tissue for energy E (ignoring minor variation due to tissue types)Compton scattering C(X,E,E′)=c(E,E′)*D(X) wherein c(E,E′) is the average is the average Compton scattering coefficient of tissue from energy E to energy E′ (which also define the scattering angle, so in fact we could define C(X,E,t)=c(E,t)*D(X) wherein “t” is the scattering angle)Also note that the total attenuation U(X,E) is:U(X,E)=A(X,E) Sum[{E′}, C(X,E,E′)] wherein Sum[{E′}, C(X,E,E′)] is the summation for all E′ (or angles t) of C(X,E,E′).2. The acquired data setThe emission data set is defined as e(P,E)={e0,e1,e2 . . . }(P) wherein:e0 is the number of photons detected having energy e0, etc. (generally, e0 will be the un-scattered emission energy, while e1,e2 . . . eN are the Compton scattered energies e1,e2 . . . eN<e0.P is generalized pixel number.In a “stationary pinhole camera” P=(p,x,y) wherein p is the pinhole number, x,y, are the pixel indexes associated with the pixel.In a SPECT (rotating) camera, P=(f,x,y) wherein f is the projection angleIn “Spectrum Dynamic” camera, P=(h,f,x,y) wherein h is the head number, f is the head′ angle.

Forward Projection—Estimation of “Peak” Data e0′(P) from S(X) and D(X)Whatever the system is, P indicates the camera configuration. Each P, is associated with a matrix element in the “response function Matrix” M(P,X) which associate (gives the response, or the sensitivity of the detector) a detected photon at generalized pixel P with a radioactive source in location (voxel) X in the body.In a general sense, M(P,X) is the response function, define as the probability of detecting a photon which was randomly emitted from location X at generalized datum pixel P given the camera geometry (including for example collimation, detector sensitivity, etc. This may make M(P,X)=>M(P,X.E) for example given energy dependence of septa penetration and/or detector response), but excluding patient attenuation.With absence of attenuation (D(X)=0), a source S(X) will (statistically speaking—it is estimated to most likely to) produce an acquired data e(P,E)={e0,e1,e2 . . . }(P)={M(P,X)*S(X), 0, 0, . . . } or e0(P)=M(P,X)*S(X). Note that M is very sparse (most of M elements are zero)However, in presence of a real patient, the equation is modified by attenuating each photon by the integrated (summed) total attenuation U(X, e0) (as derived from D(X) along the traveling line from the corresponding X to P.

Forward Scatter Estimation—Estimation of “Scattered” Data e1′(P) from S(X) and D(X) (for Multi Scatter Windows—Also e2′(P), Etc)There are several known methods disclosed in the art to calculate e1′(P) when S′(X) and D′(X) are assumed known.For example, “Monte Carlo simulations” may be used. Generally, these are computational intensive (takes long time). Thus, according to one embodiment, estimation of scattered data by “accelerated Monte Carlo” may comprise the following simplifications:1. Define only one scatter window e1(P) at high energy range of the scatter spectrum.2. As a consequence:A. Only scattering into small angles (“alpha max” which produces less than the max energy loss) needs to be computed, andB. second order scattering can be ignored.C. Attenuation coefficient of the scattered photon may be taken as a single average value U(X,e1′), where e1′ is the average energy of the photons in window e1.The “Monte Carlo” simulation is accelerated if the following steps are ignored:1. For emitted photons:A. ignore propagation direction in which there are no possible valid scattering from them to the detector (taking into account the limited “alpha max” and the collimator acceptance).B. If the path reaches pure air (out of the patient)—terminate the photon.C. Optionally ignore locations where S(X) us below a threshold value.2. Compute only the first scattering process (taking into account D(X)): In this compute only:A. Scattering a less than “alpha max”; andB. Only into valid acceptance directions of the collimator(s)3. Adjust for absorption by (taking into account D′(X)):A. U(X, e0) for the path from X0 to X1 (wherein X0 is the origin of the photon and X1 is the location of the scatter event); andB. U(X, e1′) for the path from X1 to P (wherein X1 is the location of the scatter event, and P is the detection pixel); and

Thus, various embodiments find an accurate estimation S′(X) which is as close as possible to the true source S(X) distribution.

For finding S′(X) an accurate estimation D′(X) which is as close as possible to the true density D(X) is determined. D′(X) is used for attenuation correction of the emission image. D′(X) may be useful for the operator for orienting S′(X) within the patient's body; and for ability to register the image with anatomical images such as CT or MRI.

In some embodiments, a combined reconstruction algorithm may be provided that includes:1. measure data {e0(P), e1(P)}2. start with initial guess {S′(X), D′(X)}3. estimation of {e0′(P), e1′(P)} by forward projection of guess {S′(X), D′(X)}4. compare estimation {e0′(P), e1′(P)} to measured data {e0(P), e1(P)}5. update guess {S′(X), D′(X)} in view of #46. decision if to repeat steps #3 to #5, if not:7. post processing and display the last updated guess {S′(X), D′(X)}. (post processing may comprise filtering and image analysis as known in the art)8. stop

An operational alternative includes: Steps #3, #4, and #5 can be repeated a few times for S′(X) and e0′(P) only (which is more important, and less time consuming), then performing the steps for D′(X) and e1′(P). However, it should be noted that the combined problem is harder in several aspects:1. The object is more complex: O(X)={S,D} (X)—there are two unknown to find per “X”2. The data set e(P) is more complex e(P)={e0, e1}(P)—there are two measured value per “P”3. Estimation of data has two parts:a. Estimation of e0(P) as known in the art; andb. Estimation of e1(P) as disclosed in the accelerated Monte Carlo (above).4. Updating the guess has two parts:a. Updating the source S′(X) as known in the art (keeping for example the following limitations: 1. positivity of S(X), and 2. S(X)=0 outside the patient boundary); andb. Updating the density D′(X) (keeping for example the following limitations: 1. positivity of D′(X), 2. D′(X)=0 outside the patient boundary, 3. segmentation of D′(X) as disclosed above).According to one embodiment, updating the guess D′(X) may comprise the following:1. Subtract the estimated scattered data from measured scattered data to obtain error function:
ERR(P)=e1(P)−e1′(P).2. Reconstruct ERR(P) to produce a suggested change to guessed density deltaD′(X). Reconstruction may be done by methods known in the art. The reconstructed deltaD′(X) may be positive or negative, but it may limited by the following requirements:Positivity of D′(X)—that is deltaD′(X) may not be larger than D′(X)Limited range of D(X)—that is deltaD′(X)+D′(X) must not be larger than “density of hard bone”.deltaD′(X)=0 outside the patient boundary3. Update the guess of the density to be
D′(X)=>D′(X)+deltaD′(X)
Decision to stop (#6) may depend on the number of iterations and/or requiring a minimal match between both e0′ and e0 AND e1′ and e1
Initial Approximation
According to various embodiments:1. Body boundary is calculated by reconstruction of e0(P)+e1(P) (without attenuation correction), and defining the outer perimeter outside which the reconstructed value is less than a threshold value. Optionally a one or two “outward voxels expansion” is performed to ensure that body parts are all included within the boundary.3. Initial guess for D′(X) may be taken as one of:A. Defining all the volume within the body boundary (calculated in #1) as “soft tissue”B. Defining all the volume within the body boundary (calculated in #1) as “soft tissue”, then introducing “average lungs” by “morphing” lungs from “average patient atlas”C. Calculating D′(X) by reconstruction of e1(P) without attenuation correction, but using the segmentation limitation (disclosed above) and body boundary.4. Initial guess for S′(X) is calculated by reconstruction of e0(P) using methods of the art, using optionally:A. Attenuation correction using D′(X) from above; andB. Body boundaries from above.

A technical effect of various embodiments described herein include attenuation compensation using only emission data.