Compton camera image reconstruction method and electronic device
By combining the back-projection algorithm with the maximum likelihood expectation maximization algorithm, the ambiguity problem caused by measurement error in Compton camera image reconstruction was solved, achieving clearer radiation source distribution judgment and higher imaging resolution.
Patent Information
- Application Number
- CN202411764278.6
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2024-11-29
- Publication Date
- 2025-09-30
- Estimated Expiration
- 2044-11-29
AI Technical Summary
The existing Compton camera image reconstruction method causes energy and position uncertainty due to measurement errors of the detection equipment, resulting in blurred reconstructed images, affecting image clarity and the accurate judgment of the distribution location of the radiation source.
The back-projection cone is reconstructed based on the back-projection algorithm to preliminarily determine the initial coordinates of the radiation source. Different weight values are assigned according to the back-projection cone coordinates. Iterative reconstruction is performed in combination with the maximum likelihood expectation maximization algorithm to improve image resolution and quality.
The imaging spatial resolution of the Compton camera system is significantly improved, the reconstruction speed and noise suppression effect are enhanced, and the defects of a single reconstruction algorithm are compensated.
Smart Images

Figure CN119722840B_ABST
Abstract
Description
Technical Field
[0001] The present application relates to the field of camera imaging technology, and in particular to a Compton camera image reconstruction method and electronic equipment. Background Art
[0002] As a device designed specifically for gamma-ray imaging, the Compton camera operates based on the physical principle of Compton scattering. It locates the position and energy of the radiation source by capturing and analyzing the energy of scattered and absorbed photons and their position information on the detector.
[0003] However, during image reconstruction using a Compton camera, energy and position uncertainties caused by measurement errors in the detection equipment result in varying degrees of blurring in the reconstructed image. This blurring not only affects image clarity but also limits the ability of professionals to accurately determine the distribution of radioactive sources based on the reconstructed image. Summary of the Invention
[0004] In view of the above, it is necessary to propose a Compton camera image reconstruction method and electronic equipment, which can solve the problem that the distribution of radiation sources cannot be accurately determined based on directly collected data.
[0005] An embodiment of the present application provides a Compton camera image reconstruction method, the method comprising: collecting multiple event data generated by a target radioactive source in multiple events, wherein an event represents Compton scattering of high-energy rays from the target radioactive source; for each of the multiple events, reconstructing a back-projection cone corresponding to the event based on the event data using a back-projection algorithm to obtain a coordinate representation of the spatial position of the back-projection cone; determining the initial coordinates of the target radioactive source based on the coordinate representations of multiple back-projection cones corresponding to the multiple events on the plane where the target radioactive source is located; for each of the multiple events, determining a weight of the event based on the initial coordinates of the target radioactive source and the event data of the event; and performing iterative reconstruction based on the weights of the multiple events using a model based on a maximum likelihood expectation maximization algorithm to obtain a reconstructed image.
[0006] In one embodiment, determining a weight of each event in the plurality of events based on the initial coordinates of the target radiation source and the event data of the event includes: determining an incident direction and an absorption direction of the target radiation source according to the initial coordinates of the target radiation source and the event data of the event; determining a first scattering angle according to an angle between the incident direction and the absorption direction; determining a second scattering angle according to the event data of the event; and determining the weight of the event according to an angular difference between the first scattering angle and the second scattering angle.
[0007] In one embodiment, for each of the multiple events, the event data includes: the first position coordinates of the first action point of the high-energy ray on the first layer detector of the Compton camera in the event, and the second position coordinates of the second action point of the high-energy ray on the second layer detector of the Compton camera; the method also includes: for each of the multiple events, determining the incident direction according to the initial coordinates and the first position coordinates of the event; determining the absorption direction according to the first position coordinates and the second position coordinates of the event; and determining the second scattering angle according to the first deposition energy of the event at the first action point and the second deposition energy at the second action point.
[0008] In one embodiment, for each of the multiple events, determining the weight of the event based on the difference between the first scattering angle and the second scattering angle includes: determining a standard deviation of multiple second scattering angles corresponding to the multiple events; and determining the weight of the event based on the angle difference and the standard deviation, using a formula including:
[0009]
[0010] Among them, w i represents the weight of the i-th event, θ ri represents the first scattering angle, θ i represents the second scattering angle, (θ ri -θ i ) represents the difference, σ θ represents the standard deviation.
[0011] In one embodiment, determining the initial coordinates of the target radioactive source based on the coordinate representations of the multiple back-projection cones corresponding to the multiple events in the plane where the target radioactive source is located includes: constructing an initial density distribution map, the initial density distribution map being located in the plane where the target radioactive source is located, each pixel in the initial density distribution map corresponding to a voxel in the plane; for each event in the multiple events, determining, based on the coordinate representation, a projection of the back-projection cone corresponding to the event in the initial density distribution map; obtaining a target density distribution map based on the multiple projections of the multiple back-projection cones in the initial density distribution map, each pixel in the target density distribution map corresponding to a density distribution value; and determining the initial coordinates of the target radioactive source based on the coordinates of the pixel with the highest density distribution value in the target density distribution map.
[0012] In one embodiment, the method further includes: if the target radiation source includes multiple radiation sources, sorting the density distribution values of the pixel points in the target density distribution map in descending order, selecting multiple density distribution values with the highest ranking from the sorted sequence; and determining multiple initial coordinates of the multiple radiation sources based on the positions of the multiple pixel points corresponding to the selected multiple density distribution values.
[0013] In one embodiment, the iterative reconstruction based on the weights of the multiple events using a model based on the maximum likelihood expectation maximization algorithm to obtain a reconstructed image includes: determining weight mapping values of the multiple events at each pixel point in the image based on the weights of the multiple events and the coordinate representation of the back-projection cone on the plane where the target radiation source is located, and constructing a weight mapping matrix based on the weight mapping values; determining a distribution matrix of the target radiation source in the image based on the number of iterations of the maximum likelihood expectation maximization algorithm; constructing a system matrix based on the weight mapping matrix and the distribution matrix; constructing an iteration factor matrix based on the weight mapping matrix and the system matrix; constructing an iterative expression of the maximum likelihood expectation maximization algorithm based on the weight mapping matrix, the iteration factor matrix, and the distribution matrix, and obtaining the reconstructed image by iteratively calculating the iterative expression.
[0014] In one embodiment, the iterative expression is expressed as:
[0015]
[0016] Where k represents the number of iterations; represents the intensity of the target radiation source at pixel j in the image at the kth iteration; t ij =w i b ij represents the contribution of the i-th event at pixel j; w ij represents the weight of event i at pixel j in the image, w i represents the weight of the i-th event, b ij It represents the sign function determined by the coordinates of the back-projection cone corresponding to the i-th event, and takes the value of 0 or 1; I represents the total number of events, J represents the total number of pixels; T = {t ij} represents the weight mapping matrix; represents the distribution matrix; represents the system matrix; represents the iteration factor matrix.
[0017] In one embodiment, for each event in the multiple events, a back projection cone corresponding to the event is reconstructed based on the event data of the event and a back projection algorithm, including: determining a first unit vector of the axis of the back projection cone corresponding to the event and a Compton scattering angle corresponding to the target radiation source in the event based on the event data; determining a position expression for any point in any back projection cone corresponding to any event based on the first unit vector and the Compton scattering angle, and using the position expression of any point as a coordinate representation of the spatial position of any back projection cone.
[0018] An embodiment of the present application provides a Compton camera image reconstruction device, comprising: a data acquisition module for acquiring multiple event data generated by a target radioactive source in multiple events, wherein an event represents Compton scattering of high-energy rays from the target radioactive source; a back-projection reconstruction module for reconstructing a back-projection cone corresponding to each of the multiple events based on the event data and a back-projection algorithm to obtain a coordinate representation of the spatial position of the back-projection cone; a preliminary determination module for determining the initial coordinates of the target radioactive source based on the coordinate representations of multiple back-projection cones corresponding to the multiple events on the plane where the target radioactive source is located; a weight adjustment module for determining the weight of each of the multiple events based on the initial coordinates of the target radioactive source and the event data; and an iterative reconstruction module for iteratively reconstructing the weights of the multiple events using a model based on a maximum likelihood expectation maximization algorithm to obtain a reconstructed image.
[0019] An embodiment of the present application provides an electronic device, comprising a processor and a memory, wherein the processor is configured to implement the Compton camera image reconstruction method when executing a computer program stored in the memory.
[0020] An embodiment of the present application provides a computer-readable storage medium having a computer program stored thereon. When the computer program is executed by a processor, the Compton camera image reconstruction method is implemented.
[0021] In summary, the Compton camera image reconstruction method described in this application can reconstruct the back-projection cone using a simple analytical back-projection algorithm; preliminarily determine the initial coordinates of the radioactive source based on the back-projection cone coordinates; assign different weights to detected Compton scattering events using back-projection guidance; and further iteratively calculate using a weight-corrected back-projection-guided maximum likelihood expectation maximization algorithm to improve image resolution and quality. This method overcomes the shortcomings of using a single reconstruction algorithm, improving reconstruction speed and noise suppression, and significantly enhancing the spatial resolution of Compton camera system imaging. BRIEF DESCRIPTION OF THE DRAWINGS
[0022] Figure 1 This is an example diagram of the Compton scattering process provided by an embodiment of the present application.
[0023] Figure 2 This is a structural diagram of an electronic device provided in one embodiment of the present application.
[0024] Figure 3 4 is a flowchart of a Compton camera image reconstruction method provided in one embodiment of the present application.
[0025] Figure 4 This is an example diagram of the Compton camera imaging principle provided by an embodiment of the present application.
[0026] Figure 5 This is a flowchart of a back-projection algorithm provided in one embodiment of the present application.
[0027] Figure 6 This is an example diagram of the main variables when reconstructing the back-projection cone provided by an embodiment of the present application.
[0028] Figure 7 This is an example diagram of an image reconstructed using the position of a radiation source on a projection plane obtained by a back-projection algorithm provided in one embodiment of the present application.
[0029] Figure 8 This is a flowchart of an event weight adjustment method provided in one embodiment of the present application.
[0030] Figure 9 This is an example diagram of the original radiation estimation distribution provided by an embodiment of the present application.
[0031] Figure 10 This is an example diagram of the iterative evolution of the radiation source estimation distribution matrix provided in one embodiment of the present application.
[0032] Figure 11 This is an example diagram of an iterative factor matrix image provided in one embodiment of the present application.
[0033] Figure 12 This is a flowchart of a Compton camera image reconstruction method provided by another embodiment of the present application.
[0034] Figure 13 1 is a structural diagram of a Compton camera image reconstruction device provided in one embodiment of the present application. DETAILED DESCRIPTION
[0035] In order to more clearly understand the above-mentioned objectives, features and advantages of the present application, the present application is described in detail below in conjunction with the accompanying drawings and specific embodiments. It should be noted that the embodiments of the present application and the features therein can be combined with each other in the absence of conflict.
[0036] Unless otherwise defined, all technical and scientific terms used herein have the same meaning as those commonly understood by those skilled in the art to which this application pertains. The terms used herein in the specification of this application are for the purpose of describing embodiments in one embodiment only and are not intended to limit this application.
[0037] It should be noted that in this application, "at least one" means one or more, and "more than one" means two or more than two. "And / or" describes the association relationship of associated objects, indicating that three relationships may exist. For example, A and / or B can mean: A exists alone, A and B exist at the same time, and B exists alone, where A and B can be singular or plural. The terms "first", "second", "third", "fourth", etc. (if any) in the specification, claims and drawings of this application are used to distinguish similar objects, rather than to describe a specific order or sequence.
[0038] In the embodiments of the present application, words such as "exemplary" or "for example" are used to indicate examples, illustrations, or descriptions. Any embodiment or design described as "exemplary" or "for example" in the embodiments of the present application should not be interpreted as being more preferred or more advantageous than other embodiments or designs. Specifically, the use of words such as "exemplary" or "for example" is intended to present related concepts in a concrete manner. The following embodiments and features in the embodiments may be combined with each other unless there is a conflict.
[0039] In one embodiment, a Compton camera is a device used for gamma-ray imaging. It utilizes the Compton scattering principle of gamma rays to reconstruct images and is a non-mechanically collimated detector. When detecting radiation sources, this detector does not rely on physical collimation to select a specific direction of radiation, but instead analyzes the energy and angle of scattered photons to trace their origin. Compared to traditional collimated detectors, Compton cameras have a wider field of view, reduce energy loss, and significantly improve detection efficiency.
[0040] In the Compton scattering process, gamma rays collide with free electrons in matter, causing the scattering direction and energy of the gamma rays to change, thus achieving the transfer of energy and momentum. This process mainly occurs on the extranuclear electrons of atoms and does not involve the nucleus. Specifically, refer to Figure 1 As shown in the figure, in a Compton scattering event, an incident photon interacts with a target particle (usually an electron), causing the photon to scatter in different directions, losing some energy in the process. There is a specific relationship between the energy and scattering angle of the scattered photon, which is represented by the Compton cone. The vertex of the Compton cone represents the location of the scattering event, while the surface of the cone represents all possible scattered photon directions.
[0041] In practical applications, Compton detectors can only record the final position of scattered photons and the location where the scattering occurred, but cannot directly measure the initial position of the incident photons. Therefore, the source's location is often inferred by calculating the intersection of the Compton cones from multiple scattering events. Calculating and mapping the Compton cone has important applications in fields such as medical imaging, astrophysics, and nuclear physics. By analyzing the interaction between the Compton cone and the detector, we can better understand the physical properties of the scattering events and infer the location and properties of the radiation source.
[0042] Image reconstruction algorithms for Compton camera systems are primarily categorized into two main types: analytical and iterative reconstruction. Analytical reconstruction methods utilize geometric analysis to reconstruct the back-projection cone of Compton scattering events, determining the spatial coordinates of the back-projection cone through numerical calculations. Geometric analysis methods eliminate the need for statistical models and directly perform analytical projections one by one in the order in which the events occur to complete image reconstruction. In contrast, iterative reconstruction methods rely on establishing a statistical model and, through repeated iterative calculations, gradually approximate the actual spatial distribution of the radiation source.
[0043] Among related geometric analysis methods, direct backprojection improves image quality by directly superimposing backprojection cones. However, its main drawback is that it simply superimposes the backprojection cones of each scattering event into the image space. This approach does not account for statistical models or measurement errors, resulting in artifacts and noise in the image. These problems are particularly pronounced in cases of incomplete data or imprecise scattering angle measurements, resulting in blurred images. Furthermore, because the backprojection process is linear, any system noise or background radiation is directly superimposed on the final image, further degrading image quality.
[0044] Furthermore, the Maximum Likelihood Expectation Maximization (MLEM) algorithm used in related iterative reconstruction methods, while capable of gradually approximating the true distribution of radioactive sources through an iterative process, also has significant drawbacks. As the number of iterations increases, the MLEM algorithm gradually amplifies statistical fluctuations and measurement errors, causing noise in areas with low statistical counts to be gradually amplified, leading to increased noise in the reconstructed image. Furthermore, as iterations proceed, overfitting may occur, meaning the algorithm not only enhances the actual signal but also over-optimizes the noise, ultimately making the noise in the image more prominent. Inaccuracies in the initial estimate can also affect the iterative results. If the initial value deviates significantly from the true value, the error will gradually accumulate during the iteration process, further affecting image quality.
[0045] In complex radiation environments, iterative reconstruction methods are computationally complex and time-consuming. Because multiple iterations are required to approximate the true distribution of radiation sources, each iteration not only requires extensive computational effort but also must account for statistical fluctuations and the accumulation of measurement errors. Furthermore, complex radiation environments may contain multiple radiation sources and background radiation, increasing the complexity and computational burden of the iterative model.
[0046] To address the aforementioned issues, embodiments of the present application provide a Compton camera image reconstruction method that reconstructs a back-projection cone based on a simple analytical back-projection algorithm. Initial coordinates of the radioactive source are preliminarily determined based on the back-projection cone coordinates. Back-projection is used to guide the assignment of weights to detected Compton scattering events. Further iterative calculations are performed using a weight-corrected back-projection-guided maximum likelihood expectation maximization algorithm to improve image resolution and quality. This method overcomes the drawbacks of using a single reconstruction algorithm, improving reconstruction speed and noise suppression, while significantly enhancing the spatial resolution of Compton camera system imaging.
[0047] Figure 2 The electronic device 10 can be a computer, server, mobile phone, tablet computer, notebook computer, or other electronic device. The specific type of electronic device is not limited in the embodiment of the present application.
[0048] like Figure 2 As shown, the electronic device 10 may include a communication module 101, a memory 102, a processor 103, an input / output (I / O) interface 104, and a bus 105. The processor 103 is coupled to the communication module 101, the memory 102, and the I / O interface 104 via the bus 105.
[0049] The communication module 101 may include a wired communication module and / or a wireless communication module. The wired communication module may provide one or more wired communication solutions such as universal serial bus (USB) and controller area network (CAN). The wireless communication module may provide one or more wireless communication solutions such as wireless fidelity (Wi-Fi), Bluetooth (BT), mobile communication network, frequency modulation (FM), near field communication (NFC), infrared technology (IR), etc.
[0050] Memory 102 may include one or more random access memories (RAMs) and one or more non-volatile memories (NVMs). The RAM can be directly read and written by the processor 103 and can be used to store executable programs (e.g., machine instructions) of the operating system or other running programs, as well as user and application data. RAM may include static random-access memory (SRAM), dynamic random access memory (DRAM), synchronous dynamic random access memory (SDRAM), double data rate synchronous dynamic random access memory (DDR SDRAM), etc.
[0051] The non-volatile memory can also store executable programs and user and application data, etc., which can be pre-loaded into the random access memory for direct reading and writing by the processor 103. The non-volatile memory can include disk storage devices and flash memory.
[0052] The memory 102 is configured to store one or more computer programs. The one or more computer programs are configured to be executed by the processor 103. The one or more computer programs include a plurality of instructions. When executed by the processor 103, the plurality of instructions can implement the Compton camera image reconstruction method executed on the electronic device 10.
[0053] In other embodiments, the electronic device 10 further includes an external memory interface for connecting to an external memory to expand the storage capacity of the electronic device 10 .
[0054] The processor 103 may include one or more processing units. For example, the processor 103 may include an application processor (AP), a modem processor, a graphics processing unit (GPU), an image signal processor (ISP), a controller, a video codec, a digital signal processor (DSP), a baseband processor, and / or a neural-network processing unit (NPU). The different processing units may be independent devices or integrated into one or more processors.
[0055] The processor 103 provides computing and control capabilities. For example, the processor 103 is configured to execute a computer program stored in the memory 102 to implement the above-mentioned Compton camera image reconstruction method.
[0056] I / O interface 104 provides a channel for user input and output. For example, I / O interface 104 can be used to connect to various input and output devices, such as a mouse, keyboard, touch screen device, and display screen, allowing users to enter or visualize information. I / O interface 104 can also be used to connect to Compton camera 20 to enable data exchange with Compton camera 20.
[0057] The bus 105 is at least used to provide a channel for mutual communication among the communication module 101 , the memory 102 , the processor 103 , and the I / O interface 104 in the electronic device 10 .
[0058] It should be understood that the structures illustrated in the embodiments of the present application do not constitute a specific limitation on the electronic device 10. In other embodiments of the present application, the electronic device 10 may include more or fewer components than shown, or may combine or separate certain components, or arrange the components differently. The illustrated components may be implemented in hardware, software, or a combination of software and hardware.
[0059] Figure 3 This is a flow chart of a Compton camera image reconstruction method provided by an embodiment of the present application. The Compton camera image reconstruction method is applied to electronic devices, such as Figure 2 The electronic device 10 specifically includes the following steps. According to different requirements, the order of the steps in the flowchart can be changed, and some steps can be omitted.
[0060] Step S301 : collecting multiple event data generated by a target radioactive source in multiple events.
[0061] In one embodiment, the target radiation source may include one or more target radiation sources (or radiation sources), each of which corresponds to multiple events, each of which represents Compton scattering of high-energy radiation from the target radiation source. Specifically, taking a two-layer split Compton camera comprising a scattering layer and an absorption layer as an example, each event includes, but is not limited to: Compton scattering of high-energy radiation (e.g., gamma rays) in the scattering layer, followed by deposition of the Compton-scattered high-energy radiation in the absorption layer after the photoelectric effect.
[0062] In one embodiment, each event in the multiple events corresponds to an event data, and the event data includes two interaction positions and two deposited energies recorded by the scattering layer and the absorption layer.
[0063] The principle of Compton camera image reconstruction based on the event data includes: using a corresponding image reconstruction algorithm, a Compton back-projection cone can be obtained in the imaging space, with the line connecting the two interaction positions as the cone axis and the Compton scattering angle calculated from the energy deposited by the two interactions as the cone vertex angle. The possible location of the radiation source emitting high-energy rays in space is on this Compton back-projection cone. As the number of Compton scattering events continues to increase, the Compton back-projection cones are continuously superimposed, and the true position of the radiation source is continuously strengthened, so that the true distribution of the radiation source in the projection plane in three-dimensional space can be reconstructed.
[0064] refer to Figure 4 As shown in the figure, the incident photons of high-energy rays (such as gamma rays) from the radiation source undergo Compton scattering in the scattering layer, depositing part of the energy and transferring it to the electrons in the medium. The energy measured at this time is the first deposition energy E1 of the recoil electron, and the interaction position is the first action point P1; the high-energy rays that undergo Compton scattering are then completely deposited after the photoelectric effect occurs in the absorption layer. At this time, the energy of the scattered photons measured is the second deposition energy E2, and the interaction position is the second action point P2.
[0065] The radiation source is located on a conical surface with vertex P1, axis P1P2, and vertex angle θ. After multiple events that meet the conditions, many conical surfaces with vertices at different positions and different vertex angles (θ1, θ2, θ3, ...) are obtained. The intersection of these conical surfaces is the location of the radiation source. The Compton scattering angle can be expressed as:
[0066]
[0067] Wherein, θ represents the Compton scattering angle, E1 represents the first deposition energy, E2 represents the second deposition energy, and m e c 2 is the rest energy of the scattered electrons from the radiation source (e.g. 511 keV).
[0068] Step S302 : for each event among the multiple events, reconstruct the back-projection cone corresponding to the event based on the event data and the back-projection algorithm to obtain the coordinate representation of the spatial position of the back-projection cone.
[0069] The back-projection algorithm provided in the embodiment of the present application innovatively determines the coordinate representation of the spatial position of any back-projection cone directly based on a mathematical model. Compared with the related simple back-projection algorithm, the method provided in the embodiment of the present application can reduce the difficulty of reconstructing the back-projection cone and improve the reconstruction speed of determining the back-projection cone.
[0070] In one embodiment, Figure 5 As shown, for each event in the multiple events, reconstructing the back-projection cone corresponding to each event is performed based on the event data of the event and the back-projection algorithm, including the following steps:
[0071] Step S501 : determining, based on the event data, a first unit vector of the axis of the back-projection cone corresponding to the event, and determining the Compton scattering angle corresponding to the target radiation source in the event.
[0072] In one embodiment, the event data includes: a first position coordinate of a first action point of the high-energy ray on a first layer detector of a Compton camera in the event, and a second position coordinate of a second action point of the high-energy ray on a second layer detector of the Compton camera.
[0073] Among them, a three-dimensional rectangular coordinate system oxyz can be established according to the position of the first layer detector or the second layer detector to determine the first position coordinate and the second position coordinate, wherein the unit length in the coordinate system oxyz can be determined according to the size of the imaging pixel point of the projection plane where the radiation source is located, for example, the unit length in oxyz is defined as the length of an imaging pixel point.
[0074] For example Figure 6 As shown, the first center point of the first layer detector or the second center point of the second layer detector can be used as the coordinate origin o; the direction of the line connecting the first center point and the second layer center point is used as the z-axis direction, and the direction from the second center point to the first center point is defined as the positive direction of the z-axis; then the x-axis and y-axis directions are determined according to the right-hand rule; wherein the xoy plane is parallel to the first plane where the first layer detector is located or the second plane where the second layer detector is located.
[0075] In one embodiment, the first unit vector of the axis of the back-projection cone corresponding to the event is determined based on the event data, using the formula comprising:
[0076]
[0077] n x =(X1-X2) / D,
[0078] n y =(Y1-Y2) / D,
[0079] n z =(Z1-Z2) / D,
[0080] Wherein, P1(X1, Y1, Z1) represents the first position coordinate, and P2(X2, Y2, Z2) represents the second position coordinate, Denotes the first unit vector, and D denotes the distance between P1 and P2. Wherein, the length of the first unit vector is 1, and the absolute value of the component in any direction of the first unit vector (for example, |n x |) is in the range [0,1].
[0081] In one embodiment, the event data may also include, but is not limited to: a first deposition energy of the high-energy ray at a first action point of a first layer detector of a Compton camera in the event, and a second deposition energy of the high-energy ray at a second action point of a second layer detector of the Compton camera.
[0082] In one embodiment, for each of the multiple events, the Compton scattering angle corresponding to the target radiation source in the event is determined based on the event data of the event, using a formula including:
[0083]
[0084] Wherein, θ represents the Compton scattering angle, E1 represents the first deposition energy, E represents the second deposition energy, and m e c is the rest energy of the scattered electrons from the radiation source.
[0085] Step S502: Based on the first unit vector and the Compton scattering angle, determine the position expression of any point in any back-projection cone corresponding to any event, and use the position expression of any point as the coordinate representation of the spatial position of any back-projection cone.
[0086] In one embodiment, the determining of the position expression of any point in any back-projection cone corresponding to any event based on the first unit vector and the Compton scattering angle includes: determining the expression of any second unit vector in any back-projection cone with the cone vertex as the starting point based on the length of the component vector in the specified direction in the first unit vector; and determining the position expression of any point in any back-projection cone based on the expression of the any second unit vector.
[0087] In one embodiment, determining an expression of any second unit vector starting from the cone vertex in any back-projected cone based on the length of the component vector in the specified direction in the first unit vector includes:
[0088] When determining the first unit vector The vector n in z When the length of belongs to the numerical range [0,1), the second unit vector is expressed as: in:
[0089]
[0090] Where a=cosθ, b=sinθ, θ represents the Compton scattering angle, represents the azimuth angle between the second unit vector and the first unit vector,
[0091] or,
[0092] When determining the first unit vector The vector n in z When the length is 1, the second unit vector is expressed as: in:
[0093] S x =bc,
[0094] S y =bd,
[0095] S z =a,
[0096] Where a=cosθ, b=sinθ, θ represents the Compton scattering angle, represents the azimuth angle between the second unit vector and the first unit vector,
[0097] In one embodiment, when the first unit vector The vector n in z When the length of is in the range [0,1), 1-n z 2 ≠0, indicating the first unit vector Not parallel to the z-axis. When the first unit vector The vector n in z When the length of is 1, Represents the first unit vector Parallel to the z-axis.
[0098] In one embodiment, the position expression of any point in any back-projection cone is determined based on the expression of any second unit vector, and the formula used includes:
[0099] X=X1+rS x ,
[0100] Y=Y1+rS y ,
[0101] Z=Z1+rS z ,
[0102] Wherein, (X, Y, Z) represents the position expression of any point, (X1, Y1, Z1) represents the coordinates of the vertex of any back-projection cone, r represents the distance between any point and the vertex of any back-projection cone, r∈[0, r1], wherein r1 represents the distance between the vertex of any back-projection cone and the projection plane where the radiation source is located.
[0103] In one embodiment, when r changes from 0 to r1 at a rate of Δr, and Likewise As the rate of change of r changes from 0 to 2π, the coordinates of any point on the back-projection cone can be calculated, thereby obtaining the spatial position of the back-projection cone. Here, r1 is the length of the outer line of the Compton cone, which can be determined by the distance from the projection surface to the first interaction point (for example, P1(X1, Y1, Z1)). As the value of r changes from 0 to r1, each value of r corresponds to a ring on the two-dimensional cross-section of the back-projection cone at that distance. By calculating the ring at each distance, the coordinates of the entire back-projection cone in space can be obtained.
[0104] Among them, Δr and The relationship can be calculated by the following formula:
[0105]
[0106] Through the above embodiments, the problems of slow calculation speed and poor noise suppression effect of the analytical reconstruction method of the relevant Compton camera system can be solved. Based on the simple and fast Compton back-projection cone analytical reconstruction method designed by ourselves, the Compton back-projection cone can be reconstructed quickly and accurately, which significantly improves the reconstruction speed and effectively suppresses noise.
[0107] Step S303 : determining the initial coordinates of the target radiation source according to the coordinate representations of the multiple back-projection cones corresponding to the multiple events on the plane where the target radiation source is located.
[0108] In one embodiment, taking a single point target radiation source as an example, a simple back-projection algorithm yields a series of superimposed Compton back-projection cones in the imaging space. However, based on the Compton scattering principle, only one point on each Compton back-projection cone represents the true location of the radiation source; the remaining portion is unwanted noise. Therefore, the preliminary coordinates of the radiation source can be obtained by superimposing the statistical results of the back-projection data.
[0109] For example Figure 7 The figure below shows an example of a reconstructed image of the position of a radioactive source on a projection plane, obtained using the back-projection algorithm provided in an embodiment of the present application. The horizontal axis represents the length of the projection plane where the radioactive source is located, and the vertical axis represents the width of the projection plane where the radioactive source is located. Each circular ring represents a cross-section of the calculated back-projection cone, and the blank space in the center surrounded by all the circular rings can be used to indicate the approximate position of the radioactive source on the projection plane.
[0110] In one embodiment, determining the initial coordinates of the target radioactive source based on the coordinate representations of multiple back-projection cones corresponding to the multiple events in the plane where the target radioactive source is located includes: constructing an initial density distribution map, the initial density distribution map being located in the plane where the target radioactive source is located, each pixel in the initial density distribution map corresponding to a voxel in the plane; determining, based on the coordinate representations, the projection of the back-projection cone corresponding to each event in the initial density distribution map; obtaining a target density distribution map based on multiple projections of the multiple back-projection cones in the initial density distribution map, each pixel in the target density distribution map corresponding to a density distribution value; and determining the initial coordinates of the target radioactive source based on the coordinates of the pixel with the highest density distribution value in the target density distribution map.
[0111] In one example, taking the target radiation source as a single point radiation source, the initial density distribution map can be regarded as a density distribution map to be filled in which all data are empty values (for example, 0, indicating white). The plane where the initial density distribution map is located is located in the plane where the target radiation source is located, and each pixel corresponds to a voxel in the plane. Based on the coordinate representation of the back-projection cone in the plane where the target radiation source is located, multiple projections of the back-projection cone in the initial density distribution map can be determined, wherein each projection corresponds to a back-projection cone, for example, the projection can be a Compton projection ellipse. Based on the coordinates of each voxel projected in the plane, the pixel points corresponding to each voxel in the initial density distribution map can be filled, for example, with 1 (indicating black). In this way, the target density distribution map (for example) can be obtained based on multiple projections. Figure 7As shown). The density distribution value of each pixel in the target density distribution map can be traversed, and the coordinates of the pixel with the highest density distribution value in the traversal result are selected as the initial coordinates of the target radiation source. For example, if pixel A is the intersection of k (e.g., 5) projections, then pixel A has been filled k times 1, and the density distribution value of pixel A is k; if the density distribution value k in the traversal result is the maximum value, the coordinates of pixel A in the target density distribution map can be used as the initial coordinates of the target radiation source.
[0112] In one embodiment, during Compton imaging, simply counting the highest value in the density distribution cannot accurately determine the location of multiple point-like radioactive sources. Therefore, a more sophisticated method is required to locate multiple radioactive sources and ensure that each event can be properly attributed to the correct source.
[0113] In one embodiment, if the target radiation source includes multiple radiation sources, the density distribution values of the pixel points in the target density distribution map are sorted in descending order, and multiple density distribution values with the highest ranking are selected from the sorted sequence; and multiple initial coordinates of the multiple radiation sources are determined based on the positions of the multiple pixel points corresponding to the selected multiple density distribution values.
[0114] In one example, if the target radiation source includes multiple radiation sources, the method for determining the initial coordinates of a single point radiation source can be referred to. Multiple largest density distribution values can be selected from the density distribution values of the pixels in the target density distribution map, in descending order. Each selected density distribution value corresponds to a radiation source, and the initial coordinates of the corresponding radiation source are determined based on the position of the pixel corresponding to each selected density distribution value.
[0115] Through the above embodiment, multiple Compton cones can be drawn based on the detection of multiple events, and the initial position of the radiation source can be determined simply and quickly based on the positions of multiple intersections of Compton ellipses formed by multiple detections and projections.
[0116] Step S304 : For each event in the multiple events, determine a weight of the event based on the initial coordinates of the target radiation source and the event data of the event.
[0117] In one embodiment, a simple back-projection method based on an analytical approach only considers the distribution probability of the position of its results in mathematical statistics. It does not take into account that during the reconstruction process, due to the physical limitations of the detector itself, the data generated by the detector detection will inevitably have a certain degree of uncertainty, and the Compton cone is a single event. The probability distribution of the radiation source position is not an accurate representation of the actual distribution of the radiation source. Therefore, the analytical method reflects the statistical distribution of the radiation source, and the quality of the reconstructed image needs to be improved.
[0118] In one example, when reconstructing a Compton image, the uncertainty generated by detector detection can be divided into two components: energy uncertainty and position uncertainty. Energy uncertainty manifests itself as an error in the electron volt values measured by the two detectors, leading to errors in the calculated theoretical Compton scattering angle. The better the energy resolution of the detector crystals, the lower the energy uncertainty. Position uncertainty manifests itself as an error in the position coordinates of the Compton scattering interaction detected by the two detectors, using the centroid coordinates of the crystal array. Consequently, there is an error between the centroid coordinates and the actual position of the Compton interaction. Smaller individual crystals in the array result in higher position resolution.
[0119] Therefore, although the simple back-projection algorithm proposed by the present invention is simple and fast, it will produce many artifacts, and the position resolution of the image needs to be improved due to energy and position uncertainty. To obtain better image quality, it can be combined with iterative reconstruction method for optimization.
[0120] The iterative reconstruction method provided in the embodiments of this application innovatively compares the back-projected cone obtained by preliminary reconstruction using the analytical reconstruction method. This method accounts for the error caused by uncertainty propagation, assigning different weights to detected Compton scattering events and adjusting the weights accordingly. This allows for iterative reconstruction in subsequent processes, incorporating event weights into the reconstruction, to significantly improve image resolution and quality.
[0121] In one embodiment, Figure 8 As shown, for each event in the multiple events, determining the weight of the event based on the initial coordinates of the target radiation source and the event data of the event includes the following steps:
[0122] Step S801 : determining the incident direction and the absorption direction of the target radiation source according to the initial coordinates of the target radiation source and the event data of the event.
[0123] In one embodiment, an axis can be determined by calculating the line connecting the first interaction position and the second interaction position. A Compton cone is drawn with the axis as the axis center. Each direction on the cone surface may be the incident direction of the radiation source, for example Figure 1 shown.
[0124] In one embodiment, for each of the multiple events, the event data includes: a first position coordinate of a first impact point of the high-energy ray on a first layer of detectors in a Compton camera during the event, and a second position coordinate of a second impact point of the high-energy ray on a second layer of detectors in the Compton camera. The method further includes: determining the incident direction based on the initial coordinates and the first position coordinates of the event; and determining the absorption direction based on the first position coordinates and the second position coordinates of the event.
[0125] Through the above embodiment, the scattering and absorption position coordinates of each event detected by the detector can be used to connect the radiation source position and the scattering position to obtain the detected incident direction, and connect the scattering position and the absorption position to obtain the absorption direction.
[0126] Step S802: determining a first scattering angle according to the angle between the incident direction and the absorption direction.
[0127] In one embodiment, the angle θ between the incident direction and the absorption direction can be ri as the first scattering angle of event i.
[0128] Step S803: determining a second scattering angle according to the event data of the event.
[0129] In one embodiment, the second scattering angle is determined based on a first deposition energy of the event at the first action point and a second deposition energy at the second action point, using a formula including:
[0130]
[0131] Among them, θ i represents the Compton scattering angle of event i, E1 represents the first deposition energy of event i, E2 represents the second deposition energy of event i, m e c 2 is the rest energy of the scattered electrons from the radiation source of event i (e.g. 511 keV).
[0132] Step S804: determining a weight of the event according to an angle difference between the first scattering angle and the second scattering angle.
[0133] In one embodiment, determining the weight of the event based on the difference between the first scattering angle and the second scattering angle includes: determining a standard deviation of multiple second scattering angles corresponding to the multiple events; and determining the weight of the event based on the angle difference and the standard deviation, using a formula including:
[0134]
[0135] Among them, w i represents the weight of the i-th event (or time i), θ ri represents the first scattering angle, θ i represents the second scattering angle, (θ ri -θ i ) represents the difference, σ θ represents the standard deviation.
[0136] In one example, σ θ It reflects the influence of energy uncertainty on the weight. The two scattering angles θ ri and θ i The difference between the two scattering angles reflects the impact of positional uncertainty on the weight. A smaller difference between the two scattering angles indicates that the Compton cone is closer to the radiation source, and the Compton scattering event corresponding to this Compton cone can be assigned a larger weight. A larger difference between the two scattering angles indicates that the Compton cone is farther away from the radiation source, and the Compton scattering event corresponding to this Compton cone can be assigned a smaller weight.
[0137] Through the above embodiment, the difference between the incident and absorption directions can be calculated to determine the actually measured Compton scattering angle, and the Compton scattering angle calculated based on the energy information is compared with the actually measured scattering angle. Different weights are assigned to different events based on the difference between the two, thereby preliminarily screening out scattering events with higher credibility.
[0138] In one example, if the target radiation source is a multi-point source, the coordinates of the highest density distribution value can be determined using a preliminary analytical reconstruction method. The source coordinates corresponding to this location are identified, and the Compton scattering events passing through the pixel range corresponding to this coordinate are determined. This coordinate is provisionally used as the source location of this Compton event, and the angular differences of the events are calculated and weighted. The locations with the second highest density distribution value and other higher values are processed in sequence. For each coordinate, the angular difference of the Compton scattering events passing through this location is calculated and assigned a corresponding weight.
[0139] In one example, during the weight assignment process, if the direct back-projection path of a Compton scattering event passes through multiple radiation source locations simultaneously, multiple weight calculations are performed on that event, and the result with the highest weight is retained as the final radiation source location for that event. Through this embodiment, each event can be appropriately weighted, significantly improving the accuracy of multi-point radiation source reconstruction.
[0140] Through the above embodiment, by considering energy uncertainty and position uncertainty and assigning different weights to Compton scattering events, the position of the radiation source can be reconstructed more accurately based on the event weights in the subsequent iterative reconstruction process, effectively reducing artifacts, and significantly improving the resolution and quality of the reconstructed image, thereby enhancing the accuracy and reliability of image reconstruction.
[0141] Step S305 : using a model based on the maximum likelihood expectation maximization algorithm, iteratively reconstruct based on the weights of the multiple events to obtain a reconstructed image.
[0142] The weight-corrected maximum likelihood expectation maximization (MLEM) method provided in the embodiment of the present application innovatively divides the maximum likelihood expectation maximization into four parts: a radiation estimation distribution matrix (hereinafter referred to as the distribution matrix), a system response matrix (hereinafter referred to as the system matrix), an iteration factor matrix, and a weight adjustment mapping relationship matrix (hereinafter referred to as the weight mapping matrix). Based on the independence of events, events are allocated to different process accelerations, and iterative calculations are completed with modular calculations, which can effectively improve computing efficiency and suppress noise through weight correction.
[0143] In one embodiment, the iterative reconstruction based on the weights of the multiple events using a model based on the maximum likelihood expectation maximization algorithm to obtain a reconstructed image includes: determining weight mapping values of the multiple events at each pixel point in the image based on the weights of the multiple events and the coordinate representation of the back-projection cone on the plane where the target radiation source is located, and constructing a weight mapping matrix based on the weight mapping values; determining a distribution matrix of the target radiation source in the image based on the number of iterations of the maximum likelihood expectation maximization algorithm; constructing a system matrix based on the weight mapping matrix and the distribution matrix; constructing an iteration factor matrix based on the weight mapping matrix and the system matrix; constructing an iterative expression of the maximum likelihood expectation maximization algorithm based on the weight mapping matrix, the iteration factor matrix, and the distribution matrix, and obtaining the reconstructed image by iteratively calculating the iterative expression.
[0144] In one embodiment, in order to better introduce the matrices of the above four parts of maximum likelihood expectation maximization, the following introduces the maximum likelihood expectation maximization iterative reconstruction method based on weight correction, which includes:
[0145] For each Compton scattering event assigned a different weight, a likelihood function is constructed based on the Maximum Likelihood Expectation (MLE) algorithm in the MLEM algorithm to evaluate the probability of predicting the radiation contribution intensity P under a given radiation source distribution parameter λ. This likelihood function can be expressed as:
[0146]
[0147] The parameter λj represents the intensity or activity of the radiation source at the image pixel j, and the observation data includes the total number of events detected by the detector I and the total number of pixels in the image pixel array J.
[0148] Since the events corresponding to the emission of the radioactive source satisfy the Poisson distribution, the likelihood function can be decomposed into the product of multiple independent pixel observations. i Represents the prediction contribution matrix composed of all pixel positions in event i, where the matrix P i The elements in are represented as p ij , p ij represents the predicted contribution value of each event i on pixel j, and pi is the total contribution value of event i to all pixels, that is, all elements p in the matrix Pi ij The harmony.
[0149] The goal of the MLE algorithm is to determine the optimal update direction of the parameter λj by taking the derivative of the likelihood function and setting it to zero, thereby achieving the best estimation of the radiation source distribution and image reconstruction. The derivative of the likelihood function and setting it to zero can be expressed as the following formula (2):
[0150]
[0151] in, Formula (3); t ij =w i b ij represents the contribution of event i at pixel j, w i Indicates the weight of event i. b ij It represents the sign function determined by the coordinates of the back-projection cone corresponding to the i-th event, and takes the value 0 or 1. If the back-projection cone corresponding to event i passes through pixel j, then b ij The value is 1. If the back-projection cone corresponding to event i does not pass through pixel j, then b ij The value is 0.
[0152] The Expectation Maximization (EM) algorithm in the MLEM algorithm updates the radiometric distribution variables in the maximum likelihood estimation through an iterative optimization process. The EM algorithm iterates in two steps: the expectation E-step and the maximization M-step.
[0153] In the E step, given the image contribution estimate of the current radiation source distribution parameter λj, the expectation under the current conditions is calculated. The M step updates the parameter estimate λ by maximizing the expectation function. When the M step is completed through the maximum likelihood (MLE), the expectation maximization (EM) algorithm is promoted to the maximum likelihood expectation maximization (MLEM) algorithm, and X ijFor each event i at pixel j, X i It is the total contribution value of event i to all pixels.
[0154] The total contribution X of event i is known i Equal to p i Under the condition, the contribution of event i to pixel j is X ij Equal to p ij The conditional expected value of the conditional probability is X ij In X i =p i Under the condition, the value is p ij Therefore, the conditional expectation is expressed mathematically as the estimated value X ij is the weighted average of all possible values under the given condition Xi=Pi, which can be simplified as
[0155]
[0156] Because the conditional expectation is given the radiation source distribution λ of the current k-th iteration k , to predict the possible contribution value p ij And j0 belongs to J, which is any pixel position j, so the iterative expression is:
[0157]
[0158] The conditional expected prediction value is used to replace the given contribution value condition of the next iteration to maximize the expectation. The direction of maximization needs to be derived by the maximum likelihood algorithm. Substituting formula (5) into formula (3) to obtain the image distribution iteration relationship, the relationship between the k+1th iteration λ and the kth iteration λ is obtained, that is, the iterative expression of the maximum likelihood expectation maximum algorithm:
[0159]
[0160] In one embodiment, the reconstruction iterative implementation method realizes image reconstruction, and only needs to complete the calculation of formula (6). Formula (6) can be decomposed into four parts: radiation estimation distribution matrix, system response matrix, iteration factor matrix and weight adjustment mapping relationship matrix.
[0161] In one embodiment, the weight adjustment mapping relationship matrix is the contribution t of event i in j after weight adjustment. ij , is the event weight w i and back projection b ij Composition, where t ij =w i b ij , Because b ijIt is calculated by analytical method. For example, each event i is a Compton projection ellipse composed of 360 pixel positions. The event analytical value of each position is 1, so the contribution value of event i at pixel j is w ij With w i The values are the same.
[0162] In one embodiment, the response position of each event (e.g., the position of each pixel of the back-projected conic ellipse) requires a large number of computational steps, which poses significant challenges to the algorithm's memory requirements, algorithm convergence speed, and parallel processing. Since many positions in the system's response space have essentially no response, resulting in very sparse data in this vast response space, this algorithm uses a list-based approach with each event independent. Each event is mapped individually to pi = 1 × w i .
[0163] Using lists to store system response space data significantly reduces memory consumption compared to traditional large matrix storage. In lists, the data for each event is independent, eliminating the need to pre-allocate large amounts of memory for the entire system matrix. Only the data needed for actual calculations can be stored.
[0164] The calculation of the weight in the embodiment of the present application is based on the b obtained by analytical calculation. ij Assign weight value w to the position i , rather than traversing every position in the imaging space, it can achieve more accurate calculations while reducing computational complexity.
[0165] In one embodiment, the radiation estimate distribution matrix is represented as It expresses the image reconstruction and estimation of the distribution of radioactive sources. The initial value of k = 0 is the analytical method for weight correction (for example Figure 9 As the iteration proceeds, the result of the radiation source estimation distribution matrix gradually evolves into the reconstructed image, such as Figure 10 shown.
[0166] In one embodiment, the system response matrix is a matrix of event lists, with the matrix size being the number of events, and expressing the sum of the contributions of the events to the entire image space. The system response matrix contains the sum of the Compton projection position weight of each event multiplied by the radiation source estimate corresponding to each position, and its expression is:
[0167] In one embodiment, reference Figure 11As shown in the figure, taking the pixel coordinate (0,3) of the radioactive source in the image as an example, an iteration factor matrix consisting of 1500 Compton scattering events is plotted, reflecting the sensitivity image of each position in the current iteration image. The iteration factor matrix is the matrix of the reconstructed image space. The iteration factor matrix size is the pixel size of the reconstructed image. It is the ratio of the event weight calculated for each event distributed at the current position to the system response of the event, and the expression is: As the iterations progress, the reconstruction algorithm continuously adjusts the distribution of radioactive sources in the image based on the current reconstructed image and the iteration factor matrix. In each iteration, the reconstructed image is updated to better match the detected event data and the system response model. At the same time, the iteration factor matrix is also updated with each iteration to reflect the sensitivity of the updated radioactive source distribution image to the detected event. Figure 11 As shown, the iterative factor matrix is significantly more sensitive at the source location than at locations further away. As iterations proceed, the sensitivity difference between the source and its surroundings gradually decreases, allowing the final reconstructed image to better reflect the actual source distribution and avoid overfitting and reconstruction bias. The iterative factor matrix is adjusted based on the current reconstructed image at each iteration, ensuring that the update at each pixel position more accurately reflects the detected event data.
[0168] In one embodiment, the final expression of the weight-corrected back-projection-guided maximum likelihood expectation maximization algorithm image reconstruction is as follows:
[0169]
[0170] Where k represents the number of iterations; represents the intensity of the target radiation source at pixel j in the image at the kth iteration; t ij =w i b ij represents the contribution of the i-th event at pixel j; w ij represents the weight of event i at pixel j in the image, represents the weight of the i-th event, b ij It represents the sign function determined by the coordinates of the back-projection cone corresponding to the i-th event, and takes the value of 0 or 1; I represents the total number of events, J represents the total number of pixels; T = {t ij} represents the weight mapping matrix; represents the distribution matrix; represents the system matrix; represents the iteration factor matrix.
[0171] In one embodiment, the specific iterative process of the weight-corrected back-projection-guided maximum likelihood expectation maximization algorithm includes: first, setting the initial radiation source distribution parameter λ(0) by the direct back-projection value after weight assignment, and then in each iteration, obtaining the system response matrix by calculating the sum of the contributions of each event, obtaining the iteration factor matrix based on the system response matrix and the current radiation estimation distribution matrix, and calculating the product of the iteration factor matrix and the current radiation estimation distribution matrix and then adjusting the mapping relationship matrix ratio. The calculation of the next radiation estimation distribution is completed. As the iteration proceeds, the radiation source estimation distribution reconstruction diagram is shown as follows: Figure 10 , the true position of the radiation source is gradually enhanced based on the sensitivity. When the sensitivity of position j is high, the estimate of this position is increased in the next iteration. Then, the high reconstructed estimate of the j-pixel position in this iteration will suppress the further adjustment of the sensitivity of this position, gradually reaching convergence. As the iteration factor matrix converges, the radiation source reconstruction also approaches convergence. Because the calculation of each event is independent, these events can be calculated in parallel on different processors at the same time. Each independent event or data set is assigned to a process in the process pool. The weight distribution value of the specific Compton event, as well as the system response and sensitivity, are calculated and the results are stored back in the shared structure. After all processes complete the calculation tasks, the main process will collect the calculation results of all subprocesses and merge them into the final data set. The parallel computing strategy can effectively reduce the total calculation time and improve the speed and efficiency of image reconstruction by optimizing resource utilization.
[0172] For example Figure 12 FIG. 1 is a flowchart of a Compton camera image reconstruction method according to another embodiment of the present application. The Compton camera image reconstruction method according to the embodiment of the present application includes:
[0173] Based on a self-designed analytical method, the Compton back-projection cone is reconstructed quickly and accurately. Taking a single point radioactive source as an example, the preliminary coordinates of the radioactive source are obtained by superposition statistics of the back-projection data.
[0174] Using the coordinates of the scattering and absorption positions of each event detected by the detector, the source position and the scattering position are connected to obtain the detected incident direction. The scattering position and the absorption position are then connected to obtain the absorption direction. The detected Compton scattering angle is determined by calculating the difference between the incident and absorption directions. The Compton scattering angle calculated using energy information is compared with the actual measured scattering angle, and different weights are assigned based on the difference, thereby preliminarily screening out scattering events with high credibility.
[0175] These weights are used as inputs to the Maximum Likelihood Expectation Maximization (MLEM) algorithm for iterative optimization. The MLEM algorithm iterates through two main steps: in the expectation step, the expected value of each scattering event under the current parameters is calculated, and the contribution of each event to the image is adjusted based on the weights. In the maximization step, based on the results of the expectation calculation, the radiation source distribution parameters are updated to maximize the likelihood function. The iterative process breaks down the Maximum Likelihood Expectation Maximization algorithm into four parts: the radiation estimation distribution matrix, the system response matrix, the iteration factor matrix, and the weight adjustment mapping relationship matrix. Events are assigned to different processes for acceleration based on their independence, and the iterative calculations are completed using modular computing, effectively improving computational efficiency and suppressing noise through weight correction.
[0176] Through this iterative process, the algorithm gradually approaches the true distribution of the radioactive source, ultimately achieving a more accurate reconstructed image. Weight adjustment plays a crucial role in this process. By assigning greater weight to high-confidence events, it effectively suppresses the impact of noise and prevents it from interfering with the reconstructed image. Furthermore, weight adjustment accelerates the convergence of the MLEM algorithm and reduces unnecessary computation. Combining the analytical preliminary reconstruction with the iterative optimization of the MLEM algorithm not only increases reconstruction speed but also significantly improves the spatial resolution and quality of the image, overcoming the shortcomings of a single method and achieving a more comprehensive and accurate Compton imaging reconstruction.
[0177] The method provided in the embodiments of this application addresses the slow computational speed and poor noise suppression of current analytical reconstruction methods for Compton camera systems. By developing a simple and fast analytical reconstruction method using Compton back-projection cones, the method significantly improves reconstruction speed and effectively suppresses noise. To address the issue of insufficient spatial resolution, the method combines the advantages of both analytical and iterative algorithms, using back-projection guidance to assign different weights to detected Compton scattering events before iterative reconstruction. This method first uses the developed analytical reconstruction method to perform a preliminary reconstruction of the back-projection cone, then compares the actual detected scattering angles with the angles calculated based on the Compton scattering principle, adjusts the weights, and significantly improves image resolution and quality through iterative reconstruction. This combined advantage of both methods not only improves reconstruction speed and noise suppression, but also significantly enhances the spatial resolution of Compton camera system imaging.
[0178] Figure 13 1 is a structural diagram of a Compton camera image reconstruction device provided in one embodiment of the present application.
[0179] In some embodiments, the Compton camera image reconstruction device 40 may include a plurality of functional modules composed of computer program segments. The computer program of each program segment in the Compton camera image reconstruction device 40 may be stored in a memory of an electronic device and executed by at least one processor to perform (see Figure 3Description) Function for image reconstruction from a Compton camera.
[0180] In this embodiment, the Compton camera image reconstruction device 40 can be divided into multiple functional modules based on the functions they perform. These functional modules may include: a data acquisition module 401, a backprojection reconstruction module 402, a preliminary determination module 403, a weight adjustment module 404, and an iterative reconstruction module 405. As used herein, a module refers to a series of computer program segments that can be executed by at least one processor and perform a fixed function, and are stored in a memory. The functional implementation of each module in this embodiment can be found in the definition of the Compton camera image reconstruction method above, and will not be repeated here.
[0181] The data acquisition module 401 is configured to acquire multiple event data generated by a target radiation source in multiple events, wherein an event represents Compton scattering of high-energy rays from the target radiation source.
[0182] The back projection reconstruction module 402 is used to reconstruct the back projection cone corresponding to each event of the multiple events based on the event data of the event and the back projection algorithm to obtain the coordinate representation of the spatial position of the back projection cone.
[0183] The preliminary determination module 403 is configured to determine the initial coordinates of the target radiation source according to the coordinate representations of the multiple back-projection cones corresponding to the multiple events on the plane where the target radiation source is located.
[0184] The weight adjustment module 404 is configured to determine a weight of each event among the multiple events based on the initial coordinates of the target radiation source and the event data of the event.
[0185] The iterative reconstruction module 405 is configured to utilize a model based on a maximum likelihood expectation maximization algorithm to perform iterative reconstruction based on the weights of the multiple events to obtain a reconstructed image.
[0186] An embodiment of the present application further provides a computer-readable storage medium, on which a computer program is stored. The computer program includes program instructions. The method implemented when the program instructions are executed can refer to the methods in the above-mentioned embodiments of the present application.
[0187] The computer-readable storage medium may be an internal memory of the electronic device described in the above embodiment, such as a hard disk or memory of the electronic device. The computer-readable storage medium may also be an external storage device of the electronic device, such as a plug-in hard disk, a smart memory card (SMC), a secure digital (SD) card, a flash memory card, etc. equipped on the electronic device.
[0188] In some embodiments, the computer-readable storage medium may include a program storage area and a data storage area, wherein the program storage area may store an operating system, applications required for at least one function, etc.; the data storage area may store data created according to the use of the electronic device, etc.
[0189] In the above embodiments, the description of each embodiment has its own focus. For parts that are not described or recorded in detail in a certain embodiment, reference can be made to the relevant description of other embodiments.
[0190] Those skilled in the art will appreciate that the units and algorithm steps of each example described in conjunction with the embodiments disclosed herein can be implemented in electronic hardware, or a combination of computer software and electronic hardware. Whether these functions are performed in hardware or software depends on the specific application and design constraints of the technical solution. Professional and technical personnel can use different methods to implement the described functions for each specific application, but such implementation should not be considered beyond the scope of this application.
[0191] In the embodiments provided in this application, it should be understood that the disclosed devices / terminal equipment and methods can be implemented in other ways. For example, the device / terminal equipment embodiments described above are merely illustrative. For example, the division of the modules or units is merely a logical function division. In actual implementation, there may be other division methods, such as multiple units or components can be combined or integrated into another system, or some features can be ignored or not executed. Another point is that the mutual coupling or direct coupling or communication connection shown or discussed can be through some interfaces, indirect coupling or communication connection of devices or units, which can be electrical, mechanical or other forms.
[0192] The units described as separate components may or may not be physically separate, and the components shown as units may or may not be physical units, that is, they may be located in one place or distributed across multiple network units. Some or all of these units may be selected to achieve the purpose of this embodiment according to actual needs.
[0193] The above-described embodiments are only used to illustrate the technical solutions of the present application, rather than to limit them. Although the present application has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the aforementioned embodiments, or make equivalent replacements for some of the technical features therein. These modifications or replacements do not deviate the essence of the corresponding technical solutions from the spirit and scope of the technical solutions of the various embodiments of the present application, and should all be included in the scope of protection of the present application.
Claims
1. A Compton camera image reconstruction method, characterized in that: The method comprises: collecting multiple event data generated by a target radioactive source in multiple events, wherein the events represent Compton scattering of high-energy rays from the target radioactive source; For each event in the plurality of events, reconstructing a back-projection cone corresponding to the event based on the event data of the event and a back-projection algorithm to obtain a coordinate representation of a spatial position of the back-projection cone; determining the initial coordinates of the target radioactive source based on the coordinate representations of the multiple back-projection cones corresponding to the multiple events on the plane where the target radioactive source is located; For each of the multiple events, a weight of the event is determined based on the initial coordinates of the target radiation source and the event data of the event, including: determining an incident direction and an absorption direction of the target radiation source based on the initial coordinates of the target radiation source and the event data of the event; determining a first scattering angle based on an angle between the incident direction and the absorption direction; determining a second scattering angle based on the event data of the event; and determining a weight of the event based on an angular difference between the first scattering angle and the second scattering angle, including: determining a standard deviation of multiple second scattering angles corresponding to the multiple events; and determining the weight of the event based on the angular difference and the standard deviation, using a formula including: ,in, represents the weight of the i-th event, represents the first scattering angle, represents the second scattering angle, Indicates the difference, represents the standard deviation; A model based on a maximum likelihood expectation maximization algorithm is used to iteratively reconstruct the weights of the multiple events to obtain a reconstructed image, including: determining weight mapping values of the multiple events at each pixel point in the image based on the weights of the multiple events and the coordinate representation of the back-projection cone on the plane where the target radiation source is located, and constructing a weight mapping matrix based on the weight mapping values; determining a distribution matrix of the target radiation source in the image based on the number of iterations of the maximum likelihood expectation maximization algorithm; constructing a system matrix based on the weight mapping matrix and the distribution matrix; constructing an iteration factor matrix based on the weight mapping matrix and the system matrix; constructing an iterative expression of the maximum likelihood expectation maximization algorithm based on the weight mapping matrix, the iteration factor matrix and the distribution matrix, and obtaining the reconstructed image by iteratively calculating the iterative expression.
2. The Compton camera image reconstruction method according to claim 1, wherein: For each event in the plurality of events, the event data includes: a first position coordinate of a first point of action of the high-energy ray on a first layer of detectors of the Compton camera in the event, and a second position coordinate of a second point of action of the high-energy ray on a second layer of detectors of the Compton camera; The method further comprises: For each event in the plurality of events, determining the incident direction according to the initial coordinates and the first position coordinates of the event; determining the absorption direction according to the first position coordinates and the second position coordinates of the event; The second scattering angle is determined according to a first deposition energy of the event at the first action point and a second deposition energy of the event at the second action point.
3. The Compton camera image reconstruction method according to claim 1, wherein: Determining the initial coordinates of the target radiation source according to the coordinate representations of the multiple back-projection cones corresponding to the multiple events on the plane where the target radiation source is located includes: Constructing an initial density distribution map, wherein the initial density distribution map is located in the plane where the target radiation source is located, and each pixel point in the initial density distribution map corresponds to a voxel in the plane; For each event in the plurality of events, determining, according to the coordinate representation, a projection of a back-projection cone corresponding to the event in the initial density distribution map; Obtaining a target density distribution map according to multiple projections of the multiple back-projection cones in the initial density distribution map, wherein each pixel point in the target density distribution map corresponds to a density distribution value; The initial coordinates of the target radiation source are determined according to the coordinates of the pixel point with the highest density distribution value in the target density distribution map.
4. The Compton camera image reconstruction method according to claim 3, wherein: The method further comprises: If the target radiation source includes multiple radiation sources, the density distribution values of the pixels in the target density distribution map are sorted in descending order, and the multiple density distribution values with the highest ranking are selected from the sorted sequence; The multiple initial coordinates of the multiple radiation sources are determined according to the positions of the multiple pixel points corresponding to the selected multiple density distribution values.
5. The Compton camera image reconstruction method according to claim 1, wherein: The iterative expression is expressed as: , Where k represents the number of iterations; represents the intensity of the target radiation source at pixel j in the image at the kth iteration; represents the contribution of the i-th event at pixel j; Representing an event At the pixel point of the image The weight of represents the weight of the i-th event, It represents the sign function determined by the coordinates of the back-projection cone corresponding to the i-th event, and takes the value of 0 or 1; , I represents the total number of events, J represents the total number of pixels; T= represents the weight mapping matrix; represents the distribution matrix; represents the system matrix; represents the iteration factor matrix.
6. The Compton camera image reconstruction method according to claim 1, wherein: For each event of the multiple events, reconstructing a back-projection cone corresponding to the event based on the event data of the event and a back-projection algorithm includes: determining, based on event data of the event, a first unit vector of an axis of a back-projection cone corresponding to the event, and determining a Compton scattering angle corresponding to the target radiation source in the event; Based on the first unit vector and the Compton scattering angle, a position expression of any point in any back-projection cone corresponding to any event is determined, and the position expression of the any point is used as the coordinate representation of the spatial position of the any back-projection cone.
7. An electronic device, characterized in that: The electronic device includes a processor and a memory, and the processor is configured to implement the Compton camera image reconstruction method according to any one of claims 1 to 6 when executing a computer program stored in the memory.