Compton camera image reconstruction method and electronic equipment

Abstract

The invention relates to the field of nuclear radiation detection, and provides a Compton camera image reconstruction method and electronic equipment, and the method comprises the steps: collecting multiple event data generated by a target radioactive source in multiple events; for each event in the multiple events, reconstructing the back-projection cone based on a back-projection algorithm according to the event data of the event to obtain coordinate representation of the spatial position of the back-projection cone; determining an initial coordinate of the target radioactive source according to the coordinate representation of the plurality of back projection cones of the multiple events in the plane where the target radioactive source is located; for each event in the multiple events, determining the weight of the event based on the initial coordinate of the target radioactive source and the event data of the event; and utilizing a model based on a maximum likelihood expectation maximum algorithm to perform iterative reconstruction based on the weight of multiple events to obtain a reconstructed image. By means of the method, the spatial resolution and reconstruction efficiency of the image formed by the Compton camera system can be improved.

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, in the process of image reconstruction using Compton cameras, due to the energy and position uncertainty caused by the measurement error of the detection equipment, the reconstructed image appears blurred to varying degrees. This blurring not only affects the clarity of the image, but also limits the ability of professionals to accurately determine the distribution location of the radiation source 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 event in the multiple events, reconstructing a back-projection cone corresponding to the event based on a back-projection algorithm according to the event data of the event, and obtaining a coordinate representation of the spatial position of the back-projection cone; determining the initial coordinates of the target radioactive source according to 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 event in the multiple events, determining the weight of the event based on the initial coordinates of the target radioactive source and the event data of the event; and using a model based on a maximum likelihood expectation maximization algorithm, iteratively reconstructing based on the weights of the multiple events to obtain a reconstructed image.

[0006] In one embodiment, for each event of 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: 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; 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 angle 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; 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 event in the multiple events, determining the weight of the event according to the difference between the first scattering angle and the second scattering angle includes: determining the standard deviation of multiple second scattering angles corresponding to the multiple events; determining the weight of the event based on the angle difference and the standard deviation, and the formula used includes:

[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 radiation source according to the coordinate representations of the multiple back-projection cones corresponding to the multiple events in the plane where the target radiation source is located includes: constructing an initial density distribution map, the initial density distribution map is located in the plane where the target radiation source is located, and each pixel in the initial density distribution map corresponds to a voxel in the plane; for each event in the multiple events, determining the projection of the back-projection cone corresponding to the event in the initial density distribution map according to the coordinate representation; obtaining a target density distribution map according to the multiple projections of the multiple back-projection cones in the initial density distribution map, each pixel in the target density distribution map corresponds to a density distribution value; and determining the initial coordinates of the target radiation source according to 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, and selecting multiple density distribution values ​​with the top positions from the sorted sequence; and determining multiple initial coordinates of the multiple radiation sources according to the positions of the multiple pixel points corresponding to the selected multiple density distribution values.

[0013] In one embodiment, the method of using a model based on the maximum likelihood expectation maximum algorithm to iteratively reconstruct based on the weights of the multiple events to obtain a reconstructed image includes: determining weight mapping values ​​of the multiple events at each pixel point of 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 maximum 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 iteration expression of the maximum likelihood expectation maximum algorithm based on the weight mapping matrix, the iteration factor matrix and the distribution matrix, and obtaining the reconstructed image by iteratively calculating the iteration 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 its value is 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; Denotes 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 a back projection algorithm according to the event data of the event, including: determining a first unit vector of the axis of the back projection cone corresponding to the event and determining a Compton scattering angle corresponding to the target radiation source in the event according to the event data of the event; determining a 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, and using the position expression of the any point as a coordinate representation of the spatial position of the any back projection cone.

[0018] An embodiment of the present application provides a Compton camera image reconstruction device, the device comprising: a data acquisition module, used 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; a back-projection reconstruction module, used to reconstruct a back-projection cone corresponding to the event based on a back-projection algorithm for each event in the multiple events according to the event data of the event, and obtain a coordinate representation of the spatial position of the back-projection cone; a preliminary determination module, used to determine the initial coordinates of the target radiation source according to the coordinate representations of multiple back-projection cones corresponding to the multiple events on the plane where the target radiation source is located; a weight adjustment module, used to determine the weight of the event for each event in the multiple events based on the initial coordinates of the target radiation source and the event data of the event; and an iterative reconstruction module, used to perform 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.

[0019] An embodiment of the present application provides an electronic device, which includes 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, on which a computer program is stored. 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 the present application can reconstruct the back-projection cone based on a simple back-projection algorithm of the analytical method; preliminarily determine the initial coordinates of the radiation source according to the back-projection cone coordinates; assign different weight values ​​to the detected Compton scattering events with back-projection guidance; further perform iterative calculations based on the weight-corrected back-projection guidance maximum likelihood expectation maximization algorithm to improve the resolution and quality of the image. It can make up for the defects caused by using a single type of reconstruction algorithm, improve the reconstruction speed and noise suppression effect, and significantly improve the spatial resolution of the 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 It 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 imaging principle of a Compton camera provided in an embodiment of the present application.

[0026] Figure 5 It 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 an embodiment of the present application.

[0029] Figure 8 It is a flowchart of an event weight adjustment method provided in one embodiment of the present application.

[0030] Fig. 9 is an example diagram of the original radiation estimation distribution provided by an embodiment of the present application.

[0031] Fig.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] Fig.11 This is an example diagram of an iterative factor matrix image provided by an embodiment of the present application.

[0033] Fig.12 This is a flowchart of a Compton camera image reconstruction method provided by another embodiment of the present application.

[0034] Fig.13 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 purposes, 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 in the embodiments can be combined with each other without 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 belongs. The terms used herein in the specification of this application are only for the purpose of describing the embodiments in one embodiment 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 way. The following embodiments and features in the embodiments may be combined with each other without conflict.

[0039] In one embodiment, a Compton camera is a device for gamma-ray imaging. The Compton camera uses the Compton scattering principle of gamma rays to reconstruct images and is a non-mechanically collimated detector. When detecting a ray source, this detector does not rely on physical collimation to select a specific direction of radiation, but traces its source by analyzing the energy and angle of scattered photons. Compared with traditional collimated detectors, the Compton camera has a wider field of view, can 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 Figure 1, in a Compton scattering event, an incident photon interacts with a target particle (usually an electron), causing the photon to scatter in different directions and lose some energy in the process. There is a specific relationship between the energy and scattering angle of the scattered photon, and its geometric relationship is the Compton cone. The vertex of the Compton cone represents the location where the scattering event occurred, while the surface of the cone represents all possible directions of the scattered photons.

[0041] In practical applications, the Compton detector can only record the final position of the scattered photon and the position where the photon was scattered, but cannot directly measure the starting position of the incident photon. Therefore, the intersection of the Compton cone is usually calculated through multiple scattering events to infer the position of the source. Calculating and drawing the Compton cone has important applications in 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 event and thus infer the position and properties of the radiation source.

[0042] The image reconstruction algorithms of the relevant Compton camera system are mainly divided into two categories: analytical reconstruction and iterative reconstruction. The analytical reconstruction method uses the method of geometric analysis to reconstruct the back-projection cone of the Compton scattering event, and determines the spatial coordinates of the back-projection cone through numerical calculation. The geometric analysis method does not need to consider the statistical model, and directly performs analytical projection one by one according to the order of event occurrence to complete the image reconstruction. In contrast, the iterative reconstruction method relies on the establishment of a statistical model, and through repeated iterative calculations, the results are gradually approached to the real spatial distribution of the radiation source.

[0043] Among the related geometric analysis methods, the direct backprojection method improves image quality by directly superimposing the backprojection cone. Its main drawback is that it simply superimposes the backprojection cone of each scattering event into the image space. This method does not take into account the statistical model or measurement errors, resulting in artifacts and noise in the image. These problems are particularly evident when the data is incomplete or the scattering angle measurement is inaccurate, making the image blurry. In addition, since the backprojection process is linear, any system noise or background radiation will be directly superimposed on the final image, further reducing the image quality.

[0044] In addition, the Maximum Likelihood Expectation Maximization (MLEM) algorithm used in the related iterative reconstruction method, although it can gradually approach the true distribution of the radiation source through the iterative process, also has significant defects. As the number of iterations increases, the MLEM algorithm will gradually amplify statistical fluctuations and measurement errors, so that the noise in the low statistical count area will be gradually amplified, resulting in increased noise in the reconstructed image. In addition, as the iterations proceed, overfitting may occur, that is, the algorithm not only strengthens the actual signal, but also over-optimizes the noise, ultimately making the noise in the image more significant. The inaccuracy of the initial estimate will also affect the iteration results. If the initial value deviates far from the true value, the error will gradually accumulate during the iteration process, further affecting the image quality.

[0045] In complex radiation environments, iterative reconstruction methods are relatively complex and time-consuming. Since multiple iterations are required to approximate the true distribution of radiation sources, each iteration not only has to deal with a large number of calculations, but also has to consider the accumulation of statistical fluctuations and measurement errors. In addition, there may be multiple radiation sources and background radiation in a complex radiation environment, which increases the complexity and computational burden of the iterative model.

[0046] To solve the above problems, the embodiment of the present application provides a Compton camera image reconstruction method, which can reconstruct the back-projection cone based on a simple back-projection algorithm of the analytical method; preliminarily determine the initial coordinates of the radiation source according to the back-projection cone coordinates; assign different weight values ​​to the detected Compton scattering events with back-projection guidance; further perform iterative calculations based on the weight-corrected back-projection guidance maximum likelihood expectation maximization algorithm to improve the resolution and quality of the image. It can make up for the defects caused by using a single type of reconstruction algorithm, improve the reconstruction speed and noise suppression effect, and significantly improve the spatial resolution of the Compton camera system imaging.

[0047] Figure 2 The electronic device 10 can be a computer, a server, a mobile phone, a tablet computer, a notebook computer, or other electronic devices. The specific type of the electronic device is not limited in any way 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 through 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), controller area network bus (CAN), etc. 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 technology (NFC), infrared technology (IR), etc.

[0050] The memory 102 may include one or more random access memories (RAM) and one or more non-volatile memories (NVM). The random access memory can be directly read and written by the processor 103, and can be used to store executable programs (such as machine instructions) of the operating system or other running programs, and can also be used to store user and application data. The random access memory 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 may also store executable programs and user and application data, etc., and may be loaded into the random access memory in advance for direct reading and writing by the processor 103. The non-volatile memory may include a disk storage device and a flash memory.

[0052] The memory 102 is used 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 multiple instructions. When the multiple instructions are executed by the processor 103, the Compton camera image reconstruction method executed on the electronic device 10 can be implemented.

[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 processor (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), etc. 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 used to execute a computer program stored in the memory 102 to implement the above-mentioned Compton camera image reconstruction method.

[0056] The I / O interface 104 is used to provide a channel for user input or output. For example, the I / O interface 104 can be used to connect various input and output devices, such as a mouse, a keyboard, a touch device, a display screen, etc., so that the user can enter information or visualize information. The I / O interface 104 can also be used to connect to the Compton camera 20 to achieve data interaction with the 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 is to be understood that the structure illustrated in the embodiment of the present application does 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 in the figure, or combine some components, or split some components, or arrange the components differently. The components shown in the figure may be implemented in hardware, software, or a combination of software and hardware.

[0059] Figure 3 : 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 in the flowchart 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 radiation source in multiple events.

[0061] In one embodiment, the target radiation source may include one or more target radiation sources (or referred to as radiation sources), each of which corresponds to multiple events, and the event indicates that high-energy rays from the target radiation source undergo Compton scattering. Specifically, taking a double-layer separated Compton camera including a scattering layer and an absorption layer as an example, each event includes, but is not limited to: high-energy rays (such as gamma rays) undergo Compton scattering in the scattering layer, and the high-energy rays after Compton scattering are all deposited after the photoelectric effect occurs in the absorption layer.

[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 reconstructing the Compton camera image based on the event data includes: using the 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 based on 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 real position of the radiation source is continuously strengthened, so that the real distribution of the radiation source in the projection plane in the three-dimensional space can be reconstructed.

[0064] refer to Figure 4 As shown, the incident photons of high-energy rays (such as gamma rays) from the radiation source are Compton scattered 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 deposited energy E1 of the recoil electrons, and the interaction position is the first action point P1; the high-energy rays that are subsequently Compton scattered are completely deposited after the photoelectric effect occurs in the absorption layer. At this time, the energy of the scattered photons measured is the second deposited energy E2, and the interaction position is the second action point P2.

[0065] The location of the radiation source is on a conical surface with P1 as the vertex, the line connecting P1P2 as the axis, and the vertex angle θ. After many events that meet the conditions, many conical surfaces with vertices at different positions and different vertex angles (θ1, θ2, θ3, ...) can be obtained. The intersection of these conical surfaces is the location of the radiation source. Among them, 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, m e c 2 is the rest energy of the scattered electrons of 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, based on the event data of the event, the back-projection cone corresponding to each event is reconstructed based on the back-projection algorithm, including the following steps:

[0071] Step S501, determining, according to 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 of detectors or the second layer of detectors to determine the first position coordinates and the second position coordinates, wherein the unit length in the coordinate system oxyz can be determined according to the size of the imaging pixel point on 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, and the formula used includes:

[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 coordinates, P2(X2, Y2, Z2) represents the second position coordinates, represents the first unit vector, and D represents the distance between P1 and P2. The length of the first unit vector is 1, and the absolute value of the component of any direction of the first unit vector (for example |n x |) is in the range of [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 event in the multiple events, the Compton scattering angle corresponding to the target radiation source in the event is determined according to the event data of the event, and the formula used includes:

[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 of 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 method of determining 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 a 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 second unit vector.

[0087] In one embodiment, the expression of determining any second unit vector starting from the cone vertex in any back-projected cone according to 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 value range [0,1), 1-n z 2 ≠0, indicating the first unit vector is not parallel to the z-axis. 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 according to 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 When the rate of change of changes from 0 to 2π, the coordinates of any point on the back-projection cone can be calculated to obtain the spatial position of the back-projection cone. Among them, 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)). In the process of r changing 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 by using a simple and fast Compton back-projection cone analytical reconstruction method designed by itself, so as to quickly and accurately reconstruct the Compton back-projection cone, significantly improve the reconstruction speed and effectively suppress 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 the target radiation source as a single point radiation source as an example, after using a simple back-projection algorithm, what is obtained is a series of superimposed Compton back-projection cones in the imaging space. According to the Compton scattering principle, only one point on each Compton back-projection cone is the real position of the radiation source, and the other parts are unnecessary noise interference. Therefore, the preliminary coordinates of the radiation source can be obtained by the superposition statistics of the back-projection data.

[0109] For example Figure 7 As shown, it is an example of an image reconstructed from the position of the radiation source in the projection plane obtained by the back-projection algorithm provided in an embodiment of the present application. The horizontal axis represents the length of the projection plane where the radiation source is located, and the vertical axis represents the width of the projection plane where the radiation source is located. Each ring represents a cross-sectional ring of the calculated back-projection cone, and the central blank space surrounded by all the rings can be used to indicate the approximate position of the radiation source in the projection plane.

[0110] In one embodiment, 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 in the plane where the target radiation source is located includes: constructing an initial density distribution map, 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; determining the projection of the back-projection cone corresponding to each event in the initial density distribution map according to the coordinate representation; obtaining a target density distribution map according to the multiple projections of the multiple back-projection cones in the initial density distribution map, and each pixel point in the target density distribution map corresponds to a density distribution value; and determining the initial coordinates of the target radiation source according to the coordinates of the pixel point 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 as an example, 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), and 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 point 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, and the projection can be, for example, 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), so that the target density distribution map (for example) can be obtained based on multiple projections. Figure 7As shown). The density distribution value of each pixel point in the target density distribution map can be traversed, and the coordinates of the pixel point 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 point A is the intersection of k (for example, 5) projections, then pixel point A has been filled k times 1 cumulatively, and the density distribution value of pixel point A is k; if the density distribution value k in the traversal result is the maximum value, the coordinates of pixel point 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, the method of simply counting the highest value of the density distribution value cannot accurately determine the position of multiple point-shaped radiation sources. Therefore, a more sophisticated method is needed to find the positions of multiple radiation sources to ensure that each event can be reasonably attributed to the correct radiation 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 first positions are selected from the sorted sequence; and 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.

[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, and multiple largest density distribution values ​​can be selected from the density distribution values ​​of the pixel points 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 point of each selected density distribution value.

[0115] Through the above embodiments, multiple Compton cones can be drawn based on the detection of multiple events, and the initial position of the radiation source can be simply and quickly determined 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 the 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, the simple back-projection method based on the analytical method only considers the distribution probability of the position of its result in mathematical statistics, and does not take into account that due to the physical limitations of the detector itself during the reconstruction process, 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 of the distribution of the radiation source position is not accurately the true distribution of the radiation source. Therefore, the analytical method reflects the distribution of radiation sources in a statistical sense, and the quality of the estimated reconstructed image needs to be improved.

[0118] In one example, when reconstructing the Compton image, the uncertainty generated by the detector detection can be divided into two parts: energy uncertainty and position uncertainty. Among them, the energy uncertainty is manifested as: there is an error in the electron volt value measured by the two detectors, so there is an error in the calculation of the theoretical Compton scattering angle. The better the energy resolution of the detector crystal, the lower the energy uncertainty. The position uncertainty is manifested as: when the two detectors detect the coordinates of the position of the Compton scattering interaction, the center of gravity coordinate position of the crystal array is used. Therefore, there is an error between the center of gravity coordinate position and the actual position of the Compton effect. The smaller the size of a single crystal in the crystal array, the higher the position resolution.

[0119] Therefore, although the simple back-projection algorithm proposed by the present invention is simple and has a fast imaging speed, it will have more artifacts, and the position resolution of the image imaging needs to be improved due to energy and position uncertainty. In order to obtain better image quality, it can be optimized in combination with the iterative reconstruction method.

[0120] The iterative reconstruction method provided in the embodiment of the present application innovatively uses the back-projection cone obtained by preliminary reconstruction based on the analytical reconstruction method, compares the actually detected scattering angle with the angle calculated based on the Compton scattering principle, expresses the error caused by the uncertainty transmission, assigns different weight values ​​to the detected Compton scattering events and adjusts the weights. In this way, iterative reconstruction can be performed in combination with the event weights in the subsequent process to significantly improve the resolution and quality of the image.

[0121] In one embodiment, Figure 8 As shown, for each event in the multiple events, based on the initial coordinates of the target radiation source and the event data of the event, determining the weight 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 event in the multiple events, 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. The method further includes: determining the incident direction according to the initial coordinates and the first position coordinate of the event; and determining the absorption direction according to the first position coordinate and the second position coordinate of the event.

[0125] Through the above embodiments, 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 an 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 according to a first deposition energy of the event at the first action point and a second deposition energy at the second action point, and the formula used includes:

[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, and 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 according to the difference between the first scattering angle and the second scattering angle includes: determining the standard deviation of multiple second scattering angles corresponding to the multiple events; determining the weight of the event based on the angle difference and the standard deviation, and the formula used includes:

[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 reflects the influence of position uncertainty on the weight. The smaller the difference between the two scattering angles, the less the Compton cone deviates from the radiation source, and a larger weight can be assigned to the Compton scattering event corresponding to this Compton cone. The larger the difference between the two scattering angles, the greater the distance that the Compton cone deviates from the radiation source, and a smaller weight can be assigned to the Compton scattering event corresponding to this Compton cone.

[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 by 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 radiation source, the position coordinates of the highest value of the density distribution value can be determined based on the preliminary analytical reconstruction method, the radiation source coordinates corresponding to the position can be identified, and the Compton scattering events passing through the pixel range corresponding to the coordinates can be determined. The coordinates are temporarily used as the radiation source position of the Compton event, and the angle difference of the event is calculated and weighted. The second highest value and other higher value positions of the density distribution value are processed in turn, and the angle difference of the Compton scattering events passing through the position is calculated for each coordinate, and the corresponding weight is assigned.

[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 positions at the same time, multiple weight calculations are performed on the event, and the result with the highest weight is finally retained as the final radiation source position of the event. Through the above embodiment, each event can be reasonably assigned a weight, thereby significantly improving the accuracy of multi-point radiation source reconstruction.

[0140] Through the above embodiments, 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, artifacts can be effectively reduced, and the resolution and quality of the reconstructed image can be significantly improved, thereby improving the accuracy and reliability of image reconstruction.

[0141] Step S305 , using a model based on the maximum likelihood expectation maximization algorithm, iteratively reconstructs 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). Events are allocated to different process accelerations based on the independence of events, 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 method of using a model based on the maximum likelihood expectation maximum algorithm to iteratively reconstruct based on the weights of the multiple events to obtain a reconstructed image includes: determining weight mapping values ​​of the multiple events at each pixel point of 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 maximum 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 iteration expression of the maximum likelihood expectation maximum algorithm based on the weight mapping matrix, the iteration factor matrix and the distribution matrix, and obtaining the reconstructed image by iteratively calculating the iteration 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 with different weights, a likelihood function is constructed based on the Maximum Likelihood Expectation (MLE) algorithm in the MLEM algorithm to evaluate the possibility of predicting the radiation contribution intensity P under a given radiation source distribution parameter λ. The likelihood function can be expressed as:

[0146]

[0147] The parameter λj represents the intensity or activity of the radiation source at the position of 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 radiation 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 sum.

[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. Taking 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 represents the weight of event i. ij 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. 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 completes the update of the radiometric distribution variables in the maximum likelihood estimation through an iterative optimization process. The EM algorithm iterates through two steps: the expectation E step and the maximization M step.

[0153] In the E step, given the image contribution estimate of the current 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 by 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 X ij In X i =p i Under the condition, the value is p ij The conditional expectation is thus mathematically expressed 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 calculation to maximize the expectation. The maximum likelihood algorithm is used to derive the maximum direction. 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 iteration expression of the maximum likelihood expectation maximum algorithm:

[0159]

[0160] In one embodiment, the reconstruction iteration 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 ijis 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 calculation steps, which brings significant challenges to the algorithm memory requirements, algorithm convergence speed, and parallel processing. Since many positions in the system response space have basically no response, resulting in very sparse data in this huge response space, this algorithm uses list processing and each event is independent. Each event is mapped separately as pi = 1×w i .

[0163] Using a list to store system response space data can significantly reduce memory consumption compared to traditional large matrix storage. In the list, the data of each event is independent, and there is no need to pre-allocate a large amount of memory space for the entire system matrix. Only the data required for actual calculations can be saved.

[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 radiation sources. When the initial value of k = 0, the analytical method for weight correction is constructed (for example Fig. 9 As the iterations proceed, the result of the radiation source estimation distribution matrix gradually evolves into the reconstructed image, such as Fig.10 shown.

[0166] In one embodiment, the system response matrix is ​​a matrix of event lists, the matrix size is the number of events, and it expresses the total contribution of the event to the entire image space. The system response matrix contains the sum of the weight value of the Compton projection position of each event multiplied by the estimated value of the radiation source corresponding to each position, and its expression is:

[0167] In one embodiment, reference Fig.11As shown in the figure, taking the position of the radiation source at the pixel coordinate (0,3) of the image as an example, an iteration factor matrix composed of 1500 Compton scattering events is drawn, reflecting the sensitivity image of each position of the current iteration number 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 at each position. The expression is As the iterations proceed, the reconstruction algorithm will continuously adjust 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 will be updated to better match the detected event data and the system response model. At the same time, the iteration factor matrix will also be updated with the iterations to reflect the sensitivity of the updated radioactive source distribution image to the detected event after each iteration. Fig.11 As shown in the figure, the sensitivity of the iterative factor matrix at the radiation source position is significantly higher than that at the position far from the radiation source. As the iteration proceeds, the sensitivity difference between the radiation source and its surroundings of the iterative factor matrix will gradually decrease, so that the final reconstructed image better reflects the actual distribution of the radiation source and will not overfit to produce reconstruction bias. The iterative factor matrix will be adjusted according to the current reconstructed image in each iteration to ensure that the update of each pixel position can more accurately reflect the detected event data.

[0168] In one embodiment, the final expression of image reconstruction based on weight-corrected back-projection guided maximum likelihood expectation maximization algorithm 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 of 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 its value is 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; Denotes 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, the direct back-projection value after weight allocation is set as the initial radiation source distribution parameter λ(0), and then in each iteration, the system response matrix is ​​obtained by calculating the sum of contributions of each event, and the iteration factor matrix is ​​obtained based on the system response matrix and the current radiation estimation distribution matrix, and the product of the iteration factor matrix and the current radiation estimation distribution matrix is ​​calculated and then the ratio of the weight adjustment mapping relationship matrix is ​​calculated. 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: Fig.10 , the true position of the radiation source will be gradually strengthened according to the sensitivity. When the sensitivity of position j is high, the estimation of this position will be increased in the next iteration. Then the reconstruction estimate of the j pixel position in this iteration is high, which will inhibit the continued adjustment of the sensitivity of this position and gradually reach convergence. With the convergence of the iteration factor matrix, the reconstruction of the radiation source will also tend to converge. Since the calculation of each event is independent, these events can be calculated in parallel on different processors at the same time. Assign each independent event or data set to a process in the process pool. Calculate the weight distribution value of a specific Compton event as well as the system response and sensitivity, and store the results back to 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 computing time and improve the speed and efficiency of image reconstruction by optimizing resource utilization.

[0172] For example Fig.12 FIG. 2 is a flowchart of a Compton camera image reconstruction method provided by another embodiment of the present application. The Compton camera image reconstruction method provided by the embodiment of the present application includes:

[0173] Based on the 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 scattering and absorption position coordinates of each event detected by the detector, the radiation source position and the scattering position can be connected to obtain the detected incident direction, and then the scattering position and the absorption position can be 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, so as to preliminarily screen out scattering events with higher credibility.

[0175] These weights are used as inputs, and the Maximum Likelihood Expectation Maximization (MLEM) algorithm is introduced 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 in combination with 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 splits the maximum likelihood expectation maximization into four parts: the radiation estimation distribution matrix, the system response matrix, the iteration factor matrix, and the weight adjustment mapping relationship matrix, and distributes the events to different processes for acceleration based on the independence of the events, completing the iterative calculation with modular calculations, effectively improving the computational efficiency and suppressing noise through weight correction.

[0176] Through this iterative process, the algorithm gradually approaches the real distribution of the radiation source, and finally obtains a more accurate reconstructed image. Weight adjustment plays an important role in this process. By giving higher weights to high-confidence events, the influence of noise is effectively suppressed, and the interference of noise on the reconstructed image is avoided. At the same time, weight adjustment also speeds up the convergence speed of the MLEM algorithm and reduces unnecessary calculations. Combining the preliminary reconstruction of the analytical method with the iterative optimization of the MLEM algorithm not only improves the reconstruction speed, but also significantly improves the spatial resolution and quality of the image, making up for the shortcomings of a single method and achieving a more comprehensive and accurate Compton imaging reconstruction.

[0177] The method provided in the embodiment of the present application develops a simple and fast Compton back-projection cone analytical reconstruction method to address the problems of slow calculation speed and poor noise suppression effect of the current analytical reconstruction method of the Compton camera system, which significantly improves the reconstruction speed and effectively suppresses noise. In order to address the problem that the spatial resolution of the image is not high enough, the advantages of the analytical and iterative algorithms are combined to propose that different weight values ​​are assigned to the detected Compton scattering events under the guidance of back-projection and then iterated. The back-projection cone is initially reconstructed using the developed analytical reconstruction method, and the actual detected scattering angle is compared with the angle calculated based on the Compton scattering principle, and the weight is adjusted. The resolution and quality of the image are significantly improved through iterative reconstruction. In this way, the advantages of the two methods are combined, which not only improves the reconstruction speed and noise suppression effect, but also significantly improves the spatial resolution of the Compton camera system imaging.

[0178] Fig.13 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 programs of the various program segments 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 according to the functions it performs. The functional modules may include: a data acquisition module 401, a back-projection reconstruction module 402, a preliminary determination module 403, a weight adjustment module 404, and an iterative reconstruction module 405. The module referred to in this application refers to a series of computer program segments that can be executed by at least one processor and can complete fixed functions, which are stored in a memory. In this embodiment, the functional implementation of each module in the Compton camera image reconstruction device 40 can refer to the above definition of the Compton camera image reconstruction method, and will not be repeated here.

[0181] The data acquisition module 401 is used to acquire multiple event data generated by the 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 used to determine the initial coordinates of the target radiation source according to the coordinate representation of 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 used to use a model based on the 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 (Smart Media Card, SMC), a secure digital (Secure Digital, SD) card, a flash card (Flash 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, an application 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 emphasis. For parts that are not described or recorded in detail in a certain embodiment, reference can be made to the relevant descriptions of other embodiments.

[0190] Those of ordinary skill 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 to be beyond the scope of this application.

[0191] In the embodiments provided in the present 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 only schematic. For example, the division of the modules or units is only a logical function division. There may be other division methods in actual implementation, 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 separated, and the components shown as units may or may not be physical units, that is, they may be located in one place or distributed on multiple network units. Some or all of the units may be selected according to actual needs to achieve the purpose of the solution of this embodiment.

[0193] The embodiments described above 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, a person skilled in the art should understand that the technical solutions described in the aforementioned embodiments may still be modified, or some of the technical features may be replaced by equivalents. Such modifications or replacements do not deviate the essence of the corresponding technical solutions from the spirit and scope of the technical solutions of the embodiments of the present application, and should all be included in the protection scope 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 radiation source in multiple events, wherein the event represents Compton scattering of high-energy rays from the target radiation source; For each event of the multiple events, reconstructing a back-projection cone corresponding to the event based on the event data 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 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; For each event of the plurality of events, determining a weight of the event based on the initial coordinates of the target radiation source and the event data of the event; A model based on a maximum likelihood expectation maximization algorithm is used to iteratively reconstruct based on the weights of the multiple events to obtain a reconstructed image.

2. The Compton camera image reconstruction method according to claim 1, characterized in that: For each event of the plurality of events, based on the initial coordinates of the target radiation source and the event data of the event, determining a weight of the event, including: 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; Determining a first scattering angle according to an angle between the incident direction and the absorption direction; determining a second scattering angle based on the event data of the event; A weight of the event is determined according to an angular difference between the first scattering angle and the second scattering angle.

3. The Compton camera image reconstruction method according to claim 2, characterized in that: For each event in the multiple events, the event data includes: a first position coordinate of a first action point 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 action point of the high-energy ray on a second layer of detectors of the Compton camera; The method further comprises: For each event of 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.

4. The Compton camera image reconstruction method according to claim 3, characterized in that: For each event of the multiple events, determining the weight of the event according to a difference between the first scattering angle and the second scattering angle includes: determining a standard deviation of a plurality of second scattering angles corresponding to the plurality of events; Based on the angle difference and the standard deviation, the weight of the event is determined using the formula: 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.

5. The Compton camera image reconstruction method according to claim 1, characterized in that: 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.

6. The Compton camera image reconstruction method according to claim 5, characterized in that: 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 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.

7. The Compton camera image reconstruction method according to claim 1, characterized in that: The step of using a model based on a maximum likelihood expectation maximization algorithm to iteratively reconstruct based on the weights of the multiple events to obtain a reconstructed image comprises: Determine the weight mapping values ​​of the multiple events at each pixel point of 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 construct a weight mapping matrix according to 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; An iterative expression of the maximum likelihood expectation maximization algorithm is constructed according to the weight mapping matrix, the iterative factor matrix and the distribution matrix, and the reconstructed image is obtained by iteratively calculating the iterative expression.

8. The Compton camera image reconstruction method according to claim 7, characterized in that: 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; 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 its value is 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.

9. The Compton camera image reconstruction method according to claim 1, characterized in that: For each event of the multiple events, reconstructing the back-projection cone corresponding to the event based on the event data of the event and based on the back-projection algorithm includes: Determine, according to the event data of the event, a first unit vector of the axis of the back-projection cone corresponding to the event, and determine the 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 any point is used as a coordinate representation of the spatial position of any back-projection cone.

10. An electronic device, characterized in that: The electronic device comprises a processor and a memory, and the processor is used to implement the Compton camera image reconstruction method as described in any one of claims 1 to 9 when executing a computer program stored in the memory.

Images