Medical image processing device, medical image processing method and program

Abstract

To improve image quality.SOLUTION: A medical image processing apparatus comprises an acquisition unit, a calculation unit, an estimation unit, and a reconstruction unit. The acquisition unit acquires a radiation map 103 and an attenuation map corresponding to initial image reconstruction of a PET scan. The calculation unit calculates, by using a radiative transfer equation (RTE) method, a scatter source map of a subject of the PET scan based on the radiation map and the attenuation map. The estimation unit estimates scatter by using the RTE method and based on the radiation map, the attenuation map, and the scatter source map. The reconstruction unit executes repetitive image reconstruction of the PET scan based on the estimated scatter and data from the PET scan of the subject.SELECTED DRAWING: Figure 4

Description

[Technical field]

[0001] The embodiments disclosed in this specification and the drawings relate to a medical image processing apparatus, a medical image processing method, and a program. [Background technology]

[0002] Positron Emission Tomography (PET) is an imaging method in nuclear medicine based on the use of weakly radioactive marked drugs (i.e., tracers) to image specific features of the body. PET images display the spatial distribution of radioactive drugs over time, allowing the physician or clinician to draw conclusions about, for example, metabolic activity or blood flow.

[0003] In PET imaging, a tracer drug is introduced into the patient being imaged (e.g., via injection, inhalation, or ingestion). After administration, the drug's physical and biomolecular properties cause the drug to concentrate at specific locations within the patient's body. The actual spatial distribution of the drug, the intensity of the areas of drug accumulation, and the speed of the process from administration to the eventual elimination of the drug are all factors that may have clinical importance.

[0004] During this process, the tracer attached to the drug emits a positron, which is the antimatter equivalent of an electron. When the emitted positron collides with an electron, the electron and positron annihilate, resulting in the emission of a pair of gamma rays, each with an energy of 511 keV, which travel approximately 180° apart.

[0005] The spatiotemporal distribution of the tracer is reconstructed via tomographic reconstruction principles, for example, by characterizing each detection event for its energy (i.e., amount of light generated), location, and timing. When two gamma rays are detected within a coincidence time window, they are likely due to the same positron annihilation event and are therefore identified as being a coincidence pair. A line (i.e., Line-Of-Response (LOR)) can be drawn between their locations to determine the likely location of the positron annihilation event. The timing information can also be used to determine the distribution of statistics along the LOR for annihilation based on the time-of-flight information of the two gamma rays. By accumulating a large number of LORs, tomographic reconstruction can be performed to determine a volumetric image of the spatial distribution of radioactivity (e.g., tracer concentration) within the patient.

[0006] Due to health concerns regarding exposure to radiation, clinicians in medical imaging strive to keep radiation doses as low as reasonably achievable. This effort to keep radiation doses as low as reasonably achievable motivates continued improvements in reconstructed image quality while reducing radiation dose and the signal-to-noise ratio of the measured signal.

[0007] For example, scattering is a major degradation factor in PET image reconstruction. The two main methods for scatter correction, namely Monte Carlo simulation and model-based Single Scatter Simulation (SSS), each have drawbacks. The discrete nature of Monte Carlo simulation makes it inherently noisy, but this can be overcome by increasing the number of simulations. However, running a large number of Monte Carlo simulations requires significant time and computational processing, making this approach too slow for commercial implementations in many cases. As an alternative, SSS is relatively fast, making it the preferred approach for commercial applications of PET scatter correction. However, SSS is inherently inaccurate because it only includes single-event scatter estimates and therefore ignores higher-order scatter.

[0008] Therefore, improved methods are desired for performing PET scatter correction and thereby improving the image quality of PET images. Compared to current methods, these improved methods should have either improved accuracy or computational efficiency, or both. [Prior art documents] [Patent documents]

[0009] [Patent Document 1] US Patent Application Publication No. 2018 / 0204356 [Patent Document 2] US Patent Application Publication No. 2018 / 0014806 [Patent Document 3] US Patent Application Publication No. 2015 / 0286785 Summary of the Invention [Problem to be solved by the invention]

[0010] One of the problems to be solved by the embodiments disclosed in this specification and the drawings is to improve image quality. However, the problems to be solved by the embodiments disclosed in this specification and the drawings are not limited to the above problem. Problems corresponding to the effects of each configuration shown in the embodiments described later can also be positioned as other problems. [Means for solving the problem]

[0011] A medical image processing apparatus according to an embodiment includes an acquisition unit, a calculation unit, an estimation unit, and a reconstruction unit. The acquisition unit acquires an emission map and an attenuation map corresponding to an initial image reconstruction of a PET scan. The calculation unit calculates a scatter map of a target of the PET scan based on the emission map and the attenuation map using a Radiative Transfer Equation (RTE) method. The estimation unit estimates scattering based on the emission map, the attenuation map, and the scatter map using the RTE method. A reconstruction unit performs an iterative image reconstruction of the PET scan based on the estimated scattering and data from the PET scan of the target. [Brief description of the drawings]

[0012] [Figure 1] FIG. 1 is a diagram showing an example of a PET device according to an embodiment. [Diagram 2] FIG. 2 is a diagram showing an example of a PET device according to an embodiment. [Figure 3A] FIG. 3A is a diagram illustrating gamma ray generation and detection in the absence of scatter events, according to an embodiment. [Figure 3B] FIG. 3B is a diagram illustrating gamma ray generation and detection of a first order scattering event according to an embodiment. [Figure 3C] FIG. 3C is a diagram illustrating gamma ray generation and detection of secondary scattering events according to an embodiment. [Figure 3D] FIG. 3D is a diagram illustrating gamma ray generation and detection of third order scattering events according to an embodiment. [Figure 3E] FIG. 3E is a diagram illustrating gamma ray generation and detection of two-sided first order scatter events according to an embodiment. [Figure 4] FIG. 4 is a flow chart of a method for estimating scatter and reconstructing a PET image, according to an embodiment. [Figure 5A] FIG. 5A is a flow chart of a sub-process of a method for reconstructing a PET image, according to an embodiment. [Figure 5B] FIG. 5B is a flow chart of a sub-process of a method for reconstructing a PET image, according to an embodiment. [Figure 5C] FIG. 5C is a flow chart of a sub-process of a method for reconstructing a PET image, according to an embodiment. [Figure 6A] FIG. 6A is an illustration of a calculated scatterer map, according to an embodiment. [Figure 6B] FIG. 6B is an illustration of estimated scattering flux at detector B, according to an embodiment. [Figure 7] FIG. 7 is an illustration of possible response lines from one detector A that may use an associated scatterer map, according to an embodiment. [Figure 8A] FIG. 8A is an illustration of a calculated isotropic source scattering map, according to an embodiment. [Figure 8B] FIG. 8B is a diagram of estimated scattered flux at detector A and detector B, according to an embodiment. DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS

[0013] Hereinafter, embodiments of a medical image processing apparatus, a medical image processing method, and a program will be described in detail with reference to the drawings.

[0014] The term "plurality," as used herein, is defined as two or more. The term "another," as used herein, is defined as at least a second or more. The terms "including" and / or "having," as used herein, are defined as comprising (i.e., open language). References throughout this document to "one embodiment," "an embodiment," "embodiment," "implementation," "example," or similar terms mean that a particular feature, structure, or characteristic described in connection with an embodiment is included in at least one embodiment of the disclosure. Thus, the appearances of such phrases or in various places throughout this specification are not necessarily all references to the same embodiment. Furthermore, particular features, structures, or characteristics may be combined in any suitable manner in one or more embodiments, without limitation.

[0015] In PET, the measured coincidences include both true coincidences and background signals (e.g., accidental coincidences). To improve the image quality of the reconstructed PET signal, it is desirable to estimate and account for this background signal. For example, the background signal can be accounted for by correcting the data using baseline subtraction based on the estimated background signal, or when a log-likelihood objective function is used to iteratively reconstruct the PET image, the log-likelihood formula can include a background signal term based on the estimated background signal. The background signal includes counts due to accidental events and scatter events. In PET, the background signal is mainly composed of accidental coincidences, also known as accidental and scatter.

[0016] For many annihilation events, only one photon of a pair is detected because the other photon is absorbed or scattered outside the plane of the PET detector ring. Additionally, some photons that reach the scintillation detectors of the PET detector ring go undetected due to the detector's quantum efficiency being less than unity. A detection event in which only one photon of a pair is detected may be called a "single." When two singles from separate annihilations are detected within a coincidence timing window, they are erroneously registered as originating from the same annihilation. This is called an Accidental Coincidence (AC) event, also known as a chance event. Stated differently, an AC event occurs when two unrelated singles are detected within a coincidence timing window.

[0017] Although most scattered photons in the body are not detected at the detector plane, some scattered photons are still detected and registered, resulting in inaccurate LOR. In some embodiments, some of these scattering events, resulting in inaccurate LOR, can be removed by energy discrimination because photons lose a small amount of their energy during the Compton interaction that causes the scattering event. Even so, some scattered photons (scattering) and some accidental coincidences (accidental) are inevitably recorded, and thus the background signal includes accidental and scattering.

[0018] Various corrections can be performed on the PET data, resulting in better image quality. For example, gamma rays are attenuated as they propagate through the patient, detector elements change in their detection efficiency, and accidental and scattered coincidences are recorded along with true coincidence events. Correcting for these effects improves image quality, resulting in images that are clinically useful and quantitative information from PET studies is accurate.

[0019] Gamma rays that strike more or denser matter on their path from the annihilation site to the detector are more likely to be absorbed or scattered (i.e., attenuated) than gamma rays traveling through less dense parts of the body. For this reason, dual energy X-ray computed tomography (CT) can be used to generate three-dimensional maps of material composition and concentration as a function of positrons. With the scattering and attenuation characteristics of various material components known, the reconstructed CT images can be used to determine gamma ray attenuation and scattering cross section as a function of positrons (i.e., attenuation map / image and scattering cross section map / image).

[0020] When images are reconstructed from sinograms without attenuation correction, gamma rays from less dense areas (e.g., the lungs) will be over-represented (e.g., they will appear to emit more gamma rays than they actually do) and more dense tissues will be under-represented. In the absence of attenuation correction, the reconstructed images will be susceptible to artifacts that not only mar the visual appearance of the image but also lead to inaccurate quantification of tracer uptake.

[0021] To apply attenuation correction, an attenuation map is generated, by which the attenuation through the patient for all LORs can be determined. For example, in a stand-alone PET device, this can be done in a transmission scan, in which an external positron source rotates around the patient and the attenuation of the transmitted gamma rays is determined. In a combined PET / CT scanner, the acquired CT images can be used for PET attenuation correction. That is, the acquisition unit described below acquires the attenuation map. The attenuation map is based on, for example, a CT scan of the object.

[0022] Furthermore, in certain embodiments, joint estimation methods can be used to simultaneously obtain both the attenuation map and the activity distribution for the PET data.

[0023] Like the attenuation map, the map of scattering cross section as a function of positrons can also be used to improve the image quality of the PET image by taking background signal into account in the PET reconstruction process. For example, the PET image can be reconstructed using a Maximum-Likelihood Expectation Maximization (ML-EM) method by iteratively updating the reconstructed image f according to the following equation (1): That is, the reconstruction unit described below uses the Maximum-Likelihood Expectation Maximization method to perform iterative image reconstruction of the PET scan.

[0024]

number

[0025] Here, g i is the measured count at the i-th LOR, H is the system matrix, and f j k is the estimated activity in the jth voxel at the kth iteration, and the average background signal is s i The method described herein uses the RTE to calculate s _i To provide a better (e.g., faster and / or less noisy) estimate of

[0026] As introduced above, scatter is a major degradation factor in PET image reconstruction. Two related methods for scatter correction, namely Monte Carlo simulation and model-based SSS, each have drawbacks. The discrete nature of Monte Carlo simulation makes it inherently noisy, and running a large number of Monte Carlo simulations requires significant time and computational processing making it infeasible for general industrial applications. Although SSS is relatively fast, making commercial applications of PET scatter correction the preferred approach, it is inherently inaccurate because it only involves a single-scatter estimate, thereby ignoring higher order scatter.

[0027] The method described herein uses the RTE approach to estimate scattering (i.e., the scattering term s_ i ), which has advantages over Monte Carlo simulation and SSS methods. Compared to the SSS approach, the method described herein provides more accurate scattering estimation, including both the total scattering for one detector in each LOR and the total scattering for both two detectors in each LOR. Compared to the MC approach, the method described herein achieves faster scattering estimation and lower noise.

[0028] Referring now to the figures, in which like reference numerals indicate identical or corresponding parts throughout the several views, Figures 1 and 2 show a PET device 800 including a number of Gamma-Ray Detectors (GRDs) (e.g., GRD1, GRD2 through GRDN), each configured as a rectangular detector module. According to one embodiment, the detector ring includes 40 GRDs. In another embodiment, there are 48 GRDs, with a greater number of GRDs being used to create a larger aperture size for the PET device 800. The GRDs include scintillator crystal arrays for converting gamma rays (e.g., at optical, infrared, and ultraviolet wavelengths) into scintillation photons that are detected by a photodetector.

[0029] Each GRD may include a two-dimensional array of individual detector crystals that absorb gamma radiation and emit scintillation photons. The scintillation photons may be detected by a two-dimensional array of photomultiplier tubes (PMTs), which are also disposed in the GRD. A waveguide may be disposed between the array of detector crystals and the PMTs. Furthermore, each GRD may include multiple PMTs of various sizes, each PMT disposed to receive scintillation photons from multiple detector crystals. Each PMT may generate an analog signal indicating when a scintillation event occurs, and the energy of the gamma ray that generates the detection event. Furthermore, photons emitted from one detector crystal may be detected by multiple PMTs, and based on the analog signals generated by each PMT, the detector crystal corresponding to the detection event may be determined using, for example, Anger logic and crystal decoding. It should be noted that the Anger operation is not necessarily required if there is a one-to-one correspondence between the crystals and the photon detectors.

[0030] FIG. 2 shows a schematic diagram of a PET system having a GRD arranged to detect gamma rays emitted from an object OBJ. The GRD can measure the timing, location, and energy corresponding to each gamma ray detection. In one embodiment, the gamma ray detectors are arranged in a ring as shown in FIG. 1 and FIG. 2. The detector crystals can be scintillator crystals, which have individual scintillator elements arranged in a two-dimensional array, and the scintillator elements can be any known scintillation material. The PMTs can be arranged such that light from each scintillator element is detected by multiple PMTs to enable Anger calculation and crystal decoding of scintillation events.

[0031] 2 shows an example of the arrangement of a PET device 800, where an object OBJ to be imaged is placed on a table 816, and GRD modules GRD1 to GRDN are arranged around the object OBJ and the table 816. The GRDs may be fixedly connected to a circular component 820 that is fixedly connected to a gantry 840. The gantry 840 houses many parts of the PET imaging device. The gantry 840 of the PET imaging device also includes an open aperture through which the object OBJ and the table 816 can pass, and gamma rays emitted in opposite directions from the object OBJ due to an annihilation event can be detected by the GRDs, and the timing and energy information can be used to determine coincidence counts for gamma ray pairs.

[0032] Also shown in FIG. 2 are circuits and hardware for acquiring, storing, processing, and distributing gamma ray detection data. The circuits and hardware include a processor 870, a network controller 874, a memory 878, and a data acquisition system (DAS) 876. The PET imaging device also includes a data channel that transmits detection measurements from the GRD to the DAS 876, the processor 870, the memory 878, and the network controller 874. The data acquisition system 876 can control the acquisition, digitization, and routing of detection data from the detector. In one embodiment, the DAS 876 controls the movement of the table 816. The processor 870 performs functions including reconstructing an image from the detection data according to the method 100 (shown in FIG. 4 and described below), pre-reconstruction processing of the detection data, and post-reconstruction processing of the image data, as discussed herein.

[0033] According to an embodiment, the processor 870 of the PET device 800 of Figures 1 and 2 may be configured to perform the method 100 as described herein. The processor 870 may include a CPU that may be implemented as a separate logic gate, such as an Application Specific Integrated Circuit (ASIC), Field Programmable Gate Array (FPGA), or other Complex Programmable Logic Device (CPLD). The FPGA or CPLD implementation may be coded in VHDL, Verilog, or any other hardware description language, and the code may be stored directly in electronic memory within the FPGA or CPLD or as a separate electronic memory.

[0034] The processor 870 has each of the functions of the acquisition function, the calculation function, the estimation function, and the reconstruction function. In the embodiment, each processing function performed by the acquisition function, the calculation function, the estimation function, and the reconstruction function is stored in the memory 878 in the form of a program executable by a computer. The processor 870 is a processor that realizes the function corresponding to each program by reading and executing the program from the memory 878. In other words, the processor in the state in which each program is read has each function shown in the processor 870 in FIG. 2. Note that in FIG. 2, the processing functions performed by the acquisition function, the calculation function, the estimation function, and the reconstruction function are described as being realized by a single processor 870, but it is also possible to configure the processor 870 by combining multiple independent processors, and to realize the function by each processor executing a program. In other words, each of the above-mentioned functions may be configured as a program, and one processor 870 may execute each program. As another example, a specific function may be implemented in a dedicated independent program execution circuit. Note that in FIG. 2, the acquisition function, the calculation function, the estimation function, and the reconstruction function are examples of an acquisition unit, a calculation unit, an estimation unit, and a reconstruction unit, respectively.

[0035] Additionally, memory 878 may be a hard disk drive, CD-ROM drive, DVD drive, Flash drive, RAM, ROM, or any other electronic storage known in the art. Memory 878 may be non-volatile, such as ROM, EPROM, EEPROM, or FLASH memory. Memory 878 may also be volatile, such as static or dynamic RAM, and a processor, such as a microcontroller or microprocessor, may be provided to manage the electronic memory as well as the interaction between the FPGA or CPLD and the memory.

[0036] Alternatively, the CPU in the processor 870 may execute a computer program including a set of computer readable instructions to perform the method 100 described herein, the program being stored in any of the non-volatile electronic memories and / or hard disk drives, CDs, DVDs, Flash drives, or any other known storage media as described above. Furthermore, the computer readable instructions may be provided as a utility application, background daemon, or operating system component, or combination thereof, executing in conjunction with a processor, such as an Intel® Zenon® processor or an AMD® Opteron® processor, and an operating system, such as Microsoft® Vista®, Unix®, Solaris®, Linux®, Apple®, macOS®, and other operating systems known to those skilled in the art. Furthermore, the CPU may be implemented as multiple processors that function cooperatively to execute instructions in parallel.

[0037] In one embodiment, the reconstructed image may be displayed on a display, which may be an LCD display, a CRT display, a plasma display, an OLED, an LED, or any other display known in the art.

[0038] A network controller 874, such as an Intel Ethernet Pro Network Interface Card from Intel Corporation of America, can interface between various components of the PET imaging device. Additionally, the network controller 874 can also interface with an external network. As can be appreciated, the external network can be a public network, such as the Internet, or a private network, such as a LAN or WAN network, or any combination thereof, and can also include PSTN or ISDN sub-networks. The external network can also be wired, such as an Ethernet network, or wireless, such as a cellular network, including EDGE, 3G, and 4G wireless cellular systems. The wireless network can also be WiFi, Bluetooth, or any other known wireless form of communication.

[0039] 3A through 3E provide a visual representation of gamma ray production, scattering, and detection from zero order events through higher order events (eg, two or more scatterings).

[0040] 3A shows a diagram of gamma ray detection in the absence of scattering. Positron annihilation 510 generates two gamma rays emitted approximately 180° apart along a first trajectory 512(A) and a second trajectory 512(B). The two gamma rays are detected by a first detector 540(A) and a second detector 540(B), respectively, so that the LOR can be measured / obtained between the two detectors.

[0041] In contrast, Figure 3B shows an example when a scattering event occurs along the path of the second trajectory 512 (B). Such a scattering event may be an example of one-sided scattering herein. In this case, the gamma ray is scattered at an angle θ at a point 520, and the scattered gamma ray is detected at another second detector 540 (B), such that the resulting LOR does not pass through the point of positron annihilation 510.

[0042] Figures 3C and 3D show the case of higher order scattering. Figure 3C shows second order scattering and Figure 3D shows third order scattering. Scattering point 520 is the last point to which the gamma ray scatters before being detected by the second detector 540(B). Scattering point 530 in Figure 3C and scattering points 530(1) and 5320(2) in Figure 3D are scattering points that occur before the last scattering point.

[0043] 3E shows the case of two-sided scattering, where the gamma ray generated by the positron annihilation is scattered at least once at either point 520(A) or point 520(B), respectively. The scattered gamma ray can then be detected at detector 540(A) and detector 540(B) in the case of first order scattering, and an LOR is formed therebetween. It can be seen that for two-sided scattering events, higher order scattering events are possible.

[0044] The PET device 800 described above with respect to Figures 1 and 2 can be configured in view of Figures 3A through 3E to carry out the methods described above and below in connection with Figures 4, 5A, 5B, and 5C.

[0045] 4 is a flow chart of a method 100 for estimating scatter using an RTE-based approach and reconstructing a PET image using the estimated scatter. It will be appreciated that the method 100 described in connection with FIG. 4 for a given LOR may be applied for a given number of detectors and all LORs of a given volume of interest. In this manner, iteratively reconstructed images may be iteratively generated based on the estimated scatter for one LOR of the volume of interest, a given region of the volume of interest, or the entire volume of the volume of interest.

[0046] For a given LOR (eg, the LOR between detector A and detector B), the scatter may be calculated as follows:

[0047]

number

[0048] Here, equation (3) below describes a one-sided scattering event for detector B:

[0049]

number

[0050] Also, the following equation (4) represents a one-sided scattering event for detector A:

[0051]

number

[0052] For each of the above equations, referring to FIG. 6A, A and B are detector A and detector B, respectively, O is the vanishing point, and ε A (ε A ' ) is the detector efficiency (i.e., the scattering energy) of detector A at 511 keV, and q O (V O ) is the radiation concentration (i.e., radiation volume) at O, and σ AO is the detector cross-sectional area perpendicular to the photon beam between O and detector A, and R AO is the distance from detector A to O, and μ is the linear extinction coefficient at 511 keV. A is the scattering flux calculated from the RTE for detector A, and R B is the scattering flux calculated from the RTE for detector B, and I ABrepresents a two-sided scattering event for detector A and detector B from O, where each of the two gamma rays generated at O ​​and received at the detectors is scattered at least once.

[0053] With this understanding, to determine the scatter for a given LOR between detector A and detector B, the scatter flux at detectors A and B must be determined by the RTE. In this way, the scatter event I A , I B , and I AB may be determined.

[0054] Accurate scattering solutions, including simulation of first order and multi-order scattered fluxes with accurate physical models, can be obtained using the RTE as shown in equation (5) below.

[0055]

number

[0056] where ψ(r(vector),E,Ω(hat)) is the relative intensity of the photon flux at point r(vector), energy E, and direction Ω(hat), E'(Ω') and E(Ω(hat)) are the incoming and outgoing energies (angles) of the flux, q(r(vector),E,Ω(hat)) is the emission map, n(hat) is the normal direction of the interface, f(r(vector),E,E',Ω(hat),Ω') is the scattering cross section (for PET, only Compton scattering is considered), and μ(r(vector),E) is the total attenuation coefficient of various tissues from the μ map. In an embodiment, the μ map can be obtained from a radiographic image, such as a computed tomography image, or a similar probing scan of the object. That is, the acquisition unit acquires an attenuation map corresponding to the initial image reconstruction of the PET scan.

[0057] Here, the radiation boundary condition given by the following equation (6) can be applied:

[0058]

number

[0059] Also, in one embodiment of the RTE-based approach, first order scattering can be simulated as follows:

[0060]

number

[0061] Here, the subscript 1 in ψ1(r,E,Ω) denotes first order scattering at point r, E is the energy, Ω is the unit vector in the direction of propagation for the photon flux, and the subscript O in ψ0(r',E',Ω') denotes zeroth order scattering (i.e., the photon beam in the absence of scattering). Furthermore, the vanishing point O denotes the starting point for the photon path, and f(r',E,E',Ω,Ω') is the scattering cross section including both Compton and Rayleigh scattering for photon PET. Finally, the variable μ(r,E) denotes the total attenuation coefficient for a photon at point r and energy E. This integral equation may be solved by discretizing the coordinates r (vector symbol), E, ​​Ω (hat symbol), r' (vector symbol), E', Ω', and r'' (vector symbol) and then solving numerically.

[0062] In addition to the first order scattering ψ(r,E,Ω), the multiple scattering ψ(r,E,Ω) can include higher order scattering. For example, the multiple scattering can be calculated iteratively as in Equation (8):

[0063]

number

[0064] Here, the first term on the right hand side is first order scattering and the second term on the right hand side (i.e., the integral) represents higher order scattering.

[0065] And at position r D The scattered photon field at the detector, located at (vector symbol), can be expressed as follows:

[0066]

number

[0067] In certain embodiments, the scattering cross section term f(r′, E, E′, Ω, Ω′) in the integrand and the photon flux term ψ s (r′(vector),E′,Ω′(hat symbol)) and ψ0(r′(vector),E′,Ω′(hat symbol)) may be expanded and expressed using the lowest order spherical harmonic terms in the series expansion to simplify the calculation. For example, the first scattered flux may be calculated by the discretized integral formula given by the following equation (10):

[0068]

number

[0069] Here, Y lm *(Ω) is the complex conjugate of the spherical harmonic function of order l and m, and ψ(r,E,l,m) is the intensity of the first scattered flux in the spherical harmonic domain. The spherical harmonic function can be given by the following equation (11):

[0070]

number

[0071] Here, Y lm(Ω) is the spherical harmonic function of angle l and degree m, and P _l m are the associated Legendre polynomials, N is a normalization constant, and θ and φ represent the colatitude and longitude, respectively. The number of spherical harmonics used to approximate the first scattered flux may depend on the material composition and scattering cross section at point r′ (vector symbol) and the type of scattering (e.g., Compton and Rayleigh scattering).

[0072] Furthermore, the multiple scattering fluxes can be calculated using the discretized integral spherical harmonics formula given by the following equation (12):

[0073]

number

[0074] Here, ψ s (k+1) (r (vector symbol), E, ​​l, m) and ψ s k where (r',E',l,) are the flux intensities for multiple scattering including all scattering events up to and including the k+1th order, and f(r',E,E',l) is the lth coefficient when the scattering cross section is expanded using the Legendre polynomials.

[0075]

number

[0076] By defining the above equation (13), the iteration equation defined above is ψ s (k+1) = Aψ s k +ψ1.

[0077] Therefore, scattering can be accurately simulated by including both the first scattered flux and the multiple scattered flux to represent an accurate physical model.

[0078] 4, as described above, scattering may be estimated in subprocess 110 of method 100. In other words, for a given LOR 115, R A (i.e., the total one-sided scattering flux calculated from the RTE for detector A), R B (total one-sided scattering flux calculated from RTE for detector B), and I AB (i.e., the total two-sided scatter events for detector A and detector B from O) can be determined.

[0079] In other words, in an embodiment, the calculator uses the RTE technique to calculate a scatter map (scatter flux) of the PET scanned object based on the emission map and the attenuation map, and the estimator uses the RTE technique to estimate scatter based on the emission map, the attenuation map and the calculated scatter map.

[0080] Thus, in sub-process 110 of method 100, the first order and higher order scatter fluxes for a given LOR 115 may be estimated by the RTE. That is, the estimated scatter includes contributions from first order and higher order scatter flux terms. The estimation of scatter in sub-process 110 of method 100 is performed using the radiation map obtained in step 103 of method 100. The radiation map obtained in step 103 of method 100 may be an initial reconstruction of the PET scan, or a "coarse" reconstruction of the PET scan. That is, the acquisition unit may acquire a radiation map corresponding to an initial image reconstruction of the PET scan. The radiation map forms the basis from which the scatter may be estimated in sub-process 110 of method 100. Furthermore, for example, R BIt will be appreciated that a scatter estimation may be performed for each LOR of associated detector A to determine R. For example, the multiple LORs may include detector B, detector B', and detector B''. A similar process is performed to determine R. A In other words, the estimator estimates the scatter for each LOR of the detector. The scatter estimation, as performed in sub-process 110 of method 100, is described in further detail in relation to FIG.

[0081] In step 120 of method 100, the estimated scatter is utilized for iterative or "fine" image reconstruction, taking into account the raw PET scan data. During "fine" image reconstruction, the raw PET scan data is obtained in step 107 of method 100 and may be modified by the ML-EM algorithm for PET. That is, the reconstructor performs iterative image reconstruction of the PET scan based on the estimated scatter and data from the PET scan of the object.

[0082] Essentially, the iterative update equation for the ML-EM algorithm for PET is:

[0083]

number

[0084] Here, g i is the measured count at the i-th LOR, and f j k is the estimated activity in the jth voxel at the kth iteration.

[0085] In addition, the following equation (15) holds true.

[0086]

number

[0087] The scatter correction in the reconstruction of the iterations is given by S ABOr, the estimated scattering s considering the above formula i , the update equation for the iteration with scatter correction is:

[0088]

number

[0089] Here, g i is the measured count at the i-th LOR, H is the system matrix, and f j k is the estimated activity in the jth voxel at the kth iteration, S AB s i which includes counts due to random and scattered events.

[0090] The above-described calculations of the scatter estimation performed in sub-process 110 of method 100 for a given LOR 115, and the iterative (i.e., "fine") image reconstruction in step 120 of method 100, may be iteratively performed for each given LOR 115 of a particular detector. With appropriate modifications, the above-described method may also be iteratively performed for each LOR 115 of another of the multiple detectors of the detector ring. In this manner, the entire scan volume of the patient may be evaluated and the reconstructed image therefrom may be iteratively updated to provide a more accurate image free of scatter.

[0091] The above-mentioned iterative image reconstruction of step 120 of method 100 may be performed, and the results may be used to iteratively generate a reconstructed image in step 125 of method 100 .

[0092] As mentioned above, scattering is, for a given LOR115, R A (i.e., the total one-sided scattering flux calculated from the RTE for detector A), R B (total one-sided scattering flux calculated from RTE for detector B), and I AB(i.e., the total two-sided scattering events for detector A and detector B from O). That is, scattering is estimated based on annihilation events where only one of the two generated gamma rays is scattered, and based on annihilation events where both of the generated gamma ray pairs are scattered.

[0093] 5A, 5B, and 5C are flow charts illustrating a sub-process 110 of a method 100 for estimating scatter using an RTE-based approach. Figures 5A through 5C consider one annihilation event in which two gamma rays, gamma A and gamma B, are generated in opposite directions.

[0094] In an embodiment, FIG. 5A estimates the scattered flux detected at detector B when only gamma B scatters, relative to the LOR between detector A and detector B. In an embodiment, FIG. 5B estimates the scattered flux detected at detector A when only gamma A scatters, relative to the LOR between detector A and detector B. In an embodiment, FIG. 5C estimates the scattered flux detected at both detector A and detector B when both gamma A and gamma B scatter, relative to the LOR between detector A and detector B. In this way, I A , I B , and I AB may be determined.

[0095] With reference to FIG. 5A, a scattering cross section map may be calculated using the RTE method described above for a given LOR when only one gamma ray of an annihilation event is scattered. The scattering cross section map described as a scattering cross section map in the context of the RTE method is hereinafter referred to as a scattering source map. This term is a more visual representation of what is shown in the figures. As introduced above, one-sided scattering can be understood by looking at FIG. 5A and FIG. 5B, where various scattering sources affect the one gamma ray generated at the annihilation point. The one gamma ray of the scattering annihilation event may be detected, for example, at detector B. Thus, as described above, R B may be calculated for detector A by the RTE.

[0096] In step 105A of sub-process 110 of FIG. 5A, for example, a scatter map for detector A may be calculated using the RTE method, performed using the attenuation map obtained in step 102 of sub-process 110 and the emission map obtained in step 103 of sub-process 110, and implemented in estimating the scatter detected at detector B. The scatter map, an example of which is shown in FIG. 6A, is calculated iteratively based on a discretization of the PET scan volume or a cross-section thereof, where each discretized region of the PET scan volume may be considered a scatter source. In an embodiment, the scatter map may be determined in part based on a CT scan. The scatter map (i.e., the scatter cross-section map) may be calculated iteratively with the following equations (17) and (18):

[0097]

number

[0098]

number

[0099] For the first calculation or step of the iterative process, s (r', E', Ω') may be set to 0, and subsequent iterations of the calculation may be determined according to the desired accuracy.

[0100] Using a scattering source map calculated for detector A in step 105A of sub-process 110 and information measured by detector B, the attenuation map obtained in step 102 of sub-process 110 and the emission map obtained in step 103 of sub-process 110, the scattering flux may be calculated in step 112 of sub-process 110 according to the following equation (19):

[0101]

number

[0102] The scattered flux at detector B calculated above includes estimates of first and higher order scattered flux, as shown in FIG. 6B.

[0103] According to an embodiment, the above-mentioned scatter source map for detector A may be used with multiple LORs resulting from vanishing point O in which detector A is included. For example, scatter flux estimates at detector B', detector B'', and detector B''' may be calculated with the above-mentioned scatter source map where LORs may be obtained between A and B', A and B'', and A and B'''. In an example, as shown in FIG. 7, there may be at least two LORs based on the range of detector A, the LORs including detector B and detector B'.

[0104] With respect to FIG. 5B, a scatter map similar to that above in FIG. 5A may be calculated for a given LOR when only one gamma ray of an annihilation event is scattered. Such one-sided scattering can be understood from FIG. 6A and FIG. 6B, where various scatter sources affect the one gamma ray generated at the annihilation point. The one gamma ray of the scattered annihilation event may be detected, for example, by detector A. Thus, as described above, R A may be calculated for detector B by the RTE.

[0105] In step 105B of sub-process 110 of FIG. 5B, for example, a scatter map for detector B may be calculated using the RTE and implemented in estimating the scatter detected by detector A based on information from the attenuation map obtained in step 102 of sub-process 110 and the emission map obtained in step 103 of sub-process 110. The scatter map may be iteratively calculated based on a discretization of the PET scan volume or a cross-section thereof, where each discretized region of the PET scan volume may be considered a scatter source. In an embodiment, the scatter map may be determined in part based on a CT scan. The scatter map (i.e., the scatter cross-section) may be iteratively calculated according to the following equations (20) and (21):

[0106]

number

[0107]

number

[0108] For the first calculation or step of the iterative process, s (r', E', Ω') may be set to 0, and the next iteration of the calculation may be determined according to the desired accuracy.

[0109] Once the scattering source map for detector B has been calculated in step 105B of sub-process 110, based on information measured by detector A and obtained in step 102 of sub-process 110 and the emission map obtained in step 103 of sub-process 110, the scattering flux may be calculated in step 113 of sub-process 110 according to the following equation (22):

[0110]

number

[0111] The scattered flux at detector A calculated above includes estimates of first and higher order scattered flux.

[0112] According to an embodiment, the above-mentioned scatter source map for detector B may be used with multiple LORs resulting from vanishing point O in which detector B is included. For example, scatter flux estimates at detector A', detector A'', and detector A''' may be calculated with the above-mentioned scatter source map where LORs may be obtained between B and A', B and A'', and B and A'''. In an example, there may be at least two LORs based on the range of detector B, the LORs including detector A and detector A'.

[0113] In FIG. 5A and FIG. 5B, I A and I B After calculating each of the one-sided scatter events detected by detector A and detector B along the LOR, which represent I,, we must consider the two-sided scatter events scattered by both gamma A and gamma B. For this purpose, we use I, which represents the scatter events of detector A and detector B, respectively. AB may be calculated. Initially, a scatter map may be calculated in step 105AB of sub-process 110 and may be informed by the attenuation map obtained in step 102 of sub-process 110 and the emission map obtained in step 103 of sub-process 110. The scatter map (i.e., the scatter cross section map) calculated in step 105AB of sub-process 110 may be an isotropic source scatter cross section map calculated using the following equation (23): i.e., the scatter map of the subject of the PET scan is an isotropic source scatter map.

[0114]

number

[0115] An example of an isotropic source scattering cross section map is shown in Figure 8A. The scattered flux detected at each of detectors A and B, shown in Figure 8B, in the isotropic source scattering cross section map calculated in step 105AB of sub-process 110, may be estimated in step 114 of sub-process 110 and informed by the attenuation map obtained in step 102 of sub-process 110 and the emission map obtained in step 103 of sub-process 110. The scattering may be estimated by using the following equations (24) and (25), respectively:

[0116]

number

[0117]

number

[0118] The scattered fluxes calculated above for each of detector A and detector B include estimates of first and higher order scattered fluxes.

[0119] According to an embodiment, the isotropic source scattering cross section map described above may be used with multiple LORs resulting from the vanishing point O.

[0120] Obviously, many modifications and variations are possible in light of the above teachings, and it is therefore to be understood that, within the scope of the appended claims, the invention may be practiced otherwise than as specifically described herein.

[0121] Thus, the foregoing discussion discloses and describes merely exemplary embodiments of the present invention. As will be understood by those skilled in the art, the present invention may be embodied in other specific forms without departing from its spirit or essential characteristics. Accordingly, the disclosure of the present invention is intended to be illustrative, but not limiting, of the scope of the present invention as well as the remaining claims. This disclosure, including any readily identifiable variations of the teachings herein, in part, defines the scope of the preceding claim terms such that the subject matter of the invention is not available to the public.

[0122] According to at least one of the embodiments described above, image quality can be improved.

[0123] Although some embodiments have been described, these embodiments are presented as examples and are not intended to limit the scope of the invention. These embodiments can be implemented in various other forms, and various omissions, substitutions, modifications, and combinations of embodiments can be made without departing from the spirit of the invention. These embodiments and their modifications are included in the scope of the invention and its equivalents as described in the claims, as well as in the scope and spirit of the invention. [Explanation of symbols]

[0124] 870 Processor 874 Network Controller 876 Data Acquisition System 878 Memory

Claims

1. an acquisition unit for acquiring an emission map and an attenuation map corresponding to an initial image reconstruction of a PET scan; a calculation unit for calculating a scatter map of the object of the PET scan based on the emission map and the attenuation map using a Radiative Transfer Equation (RTE) method; an estimator for estimating scattering based on the emission map, the attenuation map, and the scattering source map using the RTE method; a reconstruction unit that performs an iterative image reconstruction of the PET scan based on the estimated scattering and data from the PET scan of the object; and A medical image processing device comprising:

2. The medical imaging system of claim 1 , wherein the estimated scatter includes contributions from first order scatter flux terms and higher order scatter flux terms.

3. The medical imaging system of claim 1 , wherein the attenuation map is based on a dual-energy X-ray CT scan of the object.

4. The medical imaging system of claim 1 , wherein the estimated scatter is estimated based on annihilation events in which only one of two generated gamma rays is scattered.

5. The medical imaging system of claim 1 , wherein the estimated scatter is estimated based on annihilation events in which both of a generated gamma ray pair are scattered.

6. The medical imaging system of claim 5 , wherein the scatter map of the object of the PET scan is an isotropic scatter map.

7. The medical image processing apparatus of claim 1 , wherein the reconstructor is configured to perform the iterative image reconstruction of the PET scan using a maximum likelihood expectation maximization method.

8. The medical image processing apparatus according to claim 1 , wherein the estimation unit estimates the scattering for each LOR (Line Of Response) of a detector.

9. A medical image processing method performed by a medical image processing device, comprising: obtaining emission and attenuation maps corresponding to an initial image reconstruction of the PET scan; calculating a scatter map of the object of the PET scan based on the emission map and the attenuation map using a Radiative Transfer Equation (RTE) method; estimating scattering based on the emission map, the attenuation map, and the scattering source map using the RTE method; performing an iterative image reconstruction of the PET scan based on the estimated scattering and data from the PET scan of the object. Medical image processing method.

10. obtaining emission and attenuation maps corresponding to an initial image reconstruction of the PET scan; calculating a scatter map of the object of the PET scan based on the emission map and the attenuation map using a Radiative Transfer Equation (RTE) method; estimating scattering based on the emission map, the attenuation map, and the scattering source map using the RTE method; A program that causes a computer to perform a process for performing an iterative image reconstruction of the PET scan based on the estimated scattering and data from the PET scan of the object.

Images