Proton Imaging Methods and Imaging Equipment

Abstract

This invention provides a proton imaging method and imaging device, applicable to the fields of proton imaging technology and signal processing technology, particularly in the field of particle radiation imaging technology. The method includes: determining a virtual source position based on an initial two-dimensional dose distribution of a proton beam acquired without an imaging object; obtaining multiple sub-beam propagation paths from the virtual source position to the multiple predetermined sub-beam center positions based on the virtual source position and multiple predetermined sub-beam center positions on a predetermined plane; obtaining sub-beam data pairs corresponding to the sub-beam propagation paths and the two-dimensional dose distribution based on the two-dimensional dose distribution and the lateral distribution components determined along the sub-beam propagation paths; obtaining multiple water equivalent path lengths corresponding to at least one projection angle through an iterative algorithm based on the multiple sub-beam data pairs corresponding to at least one projection angle; and performing image reconstruction based on the multiple water equivalent path lengths corresponding to at least one projection angle to obtain a target image of the imaging object.

Description

Technical Field

[0001] This invention relates to the fields of proton imaging technology and signal processing technology, particularly to the field of particle radiation imaging technology, and more specifically to a proton imaging method and imaging device. Background Technology

[0002] After a proton beam passes through an imaging object, the dose distribution it forms in the detector is altered by the object's internal structure and material composition. Based on this change in dose distribution, the equivalent length of the protons within the imaging object can be calculated. Furthermore, by obtaining the equivalent length distribution of the protons within the imaging object at multiple angles and locations, the relative stopping power distribution of the imaging object can be reconstructed and used for range calculation in proton radiotherapy planning, thereby improving the accuracy of proton radiotherapy.

[0003] However, in the existing technology, due to the limited accuracy of the algorithm for calculating the equivalent length of water based on dose distribution changes and the limited dose distribution information obtained by the detector, the final reconstruction result still has a large error and poor imaging quality. Summary of the Invention

[0004] In view of the above problems, the present invention provides a proton imaging method and imaging device.

[0005] According to a first aspect of the present invention, a proton imaging method is provided, comprising: determining a virtual source position based on an initial two-dimensional dose distribution of a proton beam acquired without an imaging object, wherein the initial two-dimensional dose distribution characterizes the dose deposition distribution along the lateral and depth directions formed in a detector after the proton beam passes through free space; obtaining multiple sub-beam propagation paths from the virtual source position to the multiple predetermined sub-beam center positions based on the virtual source position and multiple predetermined sub-beam center positions on a predetermined plane, wherein the predetermined plane is perpendicular to the plane containing the fan-shaped region indicated by the proton beam and passes through the rotational central axis of the imaging object; and, for any sub-beam propagation path and the initial two-dimensional dose distribution... The remaining two-dimensional dose distribution, from the distribution and at least one projection viewpoint, is used to obtain a sub-beam data pair corresponding to the sub-beam propagation path and the two-dimensional dose distribution. The remaining two-dimensional dose distribution represents the dose deposition distribution along the lateral and depth directions formed in the detector after the proton beam passes through the imaging object. Based on the multiple sub-beam data pairs corresponding to at least one projection angle, multiple water equivalent path lengths corresponding to at least one projection angle are obtained through an iterative algorithm. Based on the multiple water equivalent path lengths corresponding to at least one projection angle, image reconstruction is performed to obtain the target image of the imaging object.

[0006] According to an embodiment of the present invention, determining the virtual source location based on the initial two-dimensional dose of the proton beam acquired in the absence of an imaging object includes: constructing a first objective function for the lateral broadening of the proton beam as a function of depth based on the initial two-dimensional dose distribution of the proton beam; and determining the virtual source location based on the first objective function.

[0007] According to an embodiment of the present invention, determining the virtual source location based on a first objective function includes: using the depth obtained by extrapolating the first objective function to a predetermined lateral width as the virtual source location.

[0008] According to an embodiment of the present invention, obtaining multiple sub-beam propagation paths from the virtual source location to the multiple predetermined sub-beam center locations based on the virtual source location and multiple predetermined sub-beam center locations on a predetermined plane includes: determining the central axis of the proton beam based on the initial two-dimensional dose distribution of the proton beam, wherein the central axis of the proton beam is perpendicular to the predetermined plane; determining multiple sub-beam propagation lines based on the virtual source location and the multiple predetermined sub-beam center locations on the predetermined plane, wherein a predetermined interval is maintained between the predetermined sub-beam center locations and a preset point on the predetermined plane, and the preset point is the intersection of the central axis of the proton beam and the predetermined plane; determining the propagation direction of each of the multiple sub-beam propagation lines based on the angle between each of the multiple sub-beam propagation lines and the central axis; and obtaining multiple sub-beam propagation paths from the virtual source location to the multiple predetermined sub-beam center locations based on the propagation directions of each of the multiple sub-beam propagation lines.

[0009] According to an embodiment of the present invention, based on the two-dimensional dose distribution and the lateral distribution component determined along the sub-beam propagation path, multiple sub-beam data pairs corresponding to the sub-beam propagation path and the two-dimensional dose distribution are obtained, including: for any sub-beam propagation path in the sub-beam propagation path, determining the center position of the sub-beam at different depths; determining the lateral broadening of the sub-beam at different depths based on a scattering algorithm; determining the lateral distribution component of the sub-beam along the sub-beam propagation path based on the center position and lateral broadening of the sub-beam at different depths; and obtaining the sub-beam data pairs corresponding to the sub-beam propagation path and the two-dimensional dose distribution based on the two-dimensional dose distribution and the lateral distribution component.

[0010] According to an embodiment of the present invention, obtaining multiple sub-bundle data pairs corresponding to multiple sub-bundle propagation paths and two-dimensional dose distributions based on a two-dimensional dose distribution and a lateral distribution component along the sub-bundle propagation path includes: convolving an initial two-dimensional dose distribution and a lateral distribution component to obtain multiple initial sub-bundle data corresponding to multiple sub-bundle propagation paths and an initial two-dimensional dose distribution; convolving a remaining two-dimensional dose distribution and a lateral distribution component to obtain multiple remaining sub-bundle data corresponding to multiple sub-bundle propagation paths and remaining two-dimensional dose distributions; and determining multiple sub-bundle data pairs corresponding to multiple sub-bundle propagation paths and two-dimensional dose distributions based on the multiple initial sub-bundle data and the multiple remaining sub-bundle data.

[0011] According to an embodiment of the present invention, based on multiple sub-bundle data pairs corresponding to at least one projection angle, an iterative algorithm is used to obtain multiple water equivalent path lengths corresponding to at least one projection angle. This includes: for any projection viewpoint among at least one projection viewpoint, based on the difference between the total remaining sub-bundle data and the total translation sub-bundle data of each of the multiple iteration rounds corresponding to the projection viewpoint, determining multiple target water equivalent path lengths for a target round from multiple candidate water equivalent path lengths for each of the multiple rounds, and using the multiple target water equivalent path lengths as the multiple water equivalent path lengths corresponding to the projection angle; wherein the total remaining sub-bundle data is determined based on the remaining sub-bundle data of each of the multiple sub-bundle data pairs corresponding to the projection viewpoint, and the total translation sub-bundle data of the iteration rounds corresponding to the projection viewpoint is determined based on the translation sub-data of each of the multiple sub-bundle data pairs corresponding to the projection viewpoint, and the translation sub-data is obtained by translating the initial sub-bundle data in the sub-bundle data pairs of the iteration rounds corresponding to the projection viewpoint along the depth direction to the candidate water equivalent path length corresponding to the iteration round.

[0012] According to an embodiment of the present invention, the method further includes: for any one of the multiple rounds, based on a second objective function, determining the difference between the total remaining sub-bundle data and the total translation sub-bundle data of the round corresponding to the projection viewpoint, according to the total remaining sub-bundle data, the total translation sub-bundle data, and the candidate water equivalent path length of the iteration round corresponding to the projection viewpoint; wherein, the second objective function includes a data fitting term, a prior constraint term, and a sparsity term, the data fitting term is used to determine the difference between the total remaining sub-bundle data and the total translation sub-bundle data of the iteration round corresponding to the projection viewpoint, the prior constraint term is used to constrain the multiple candidate water equivalent path lengths of the multiple iteration rounds corresponding to the projection viewpoint according to the multiple initial water equivalent path lengths corresponding to the projection viewpoint, and the sparsity term is used to perform sparse constraints on the candidate water equivalent path lengths.

[0013] According to an embodiment of the present invention, image reconstruction is performed based on multiple water equivalent path lengths corresponding to at least one projection angle to obtain a target image of the imaging object, including: processing and back-projecting the multiple water equivalent path lengths corresponding to at least one projection angle based on a filtered back-projection algorithm to obtain a target image of the imaging object.

[0014] A second aspect of the present invention provides an imaging apparatus, comprising: a detector configured to acquire an initial two-dimensional dose distribution of a proton beam in the absence of an imaging object, and to acquire residual two-dimensional dose distributions at at least one projection viewpoint, wherein the initial two-dimensional dose distribution characterizes the dose deposition distribution along the lateral and depth directions formed in the detector after the proton beam passes through free space, and the residual two-dimensional dose distribution characterizes the dose deposition distribution along the lateral and depth directions formed in the detector after the proton beam passes through an imaging object; and a processor configured to: obtain a distribution from the virtual source position to a plurality of predetermined sub-beam center positions on a predetermined plane of the imaging object. Multiple sub-beam propagation paths at multiple predetermined sub-beam center positions; for any propagation path in the propagation path, and any two-dimensional dose distribution in the initial two-dimensional dose distribution and the remaining two-dimensional dose distribution under at least one projection viewpoint, sub-beam data pairs corresponding to the sub-beam propagation path and the two-dimensional dose distribution are obtained based on the two-dimensional dose distribution and the lateral distribution component determined along the sub-beam propagation path; multiple water equivalent path lengths corresponding to at least one projection angle are obtained based on the multiple sub-beam data pairs corresponding to at least one projection angle; and image reconstruction is performed on the imaging object based on the multiple water equivalent path lengths corresponding to at least one projection angle to obtain the target image of the imaging object.

[0015] A third aspect of the present invention provides an electronic device, comprising: one or more processors; and a memory for storing one or more computer programs, wherein the one or more processors execute the one or more computer programs to implement the steps of the method described above.

[0016] A fourth aspect of the present invention also provides a computer-readable storage medium having a computer program or instructions stored thereon, wherein the computer program or instructions, when executed by a processor, implement the steps of the above-described method.

[0017] A fifth aspect of the present invention also provides a computer program product, including a computer program or instructions, which, when executed by a processor, implement the steps of the above-described method.

[0018] According to embodiments of the present invention, by acquiring the initial two-dimensional dose distribution of the proton beam, the virtual source position can be accurately determined based on the initial two-dimensional dose distribution. This virtual source position is then used to divide the initial two-dimensional dose distribution and the remaining two-dimensional dose distribution, resulting in multiple sub-bundle data pairs composed of initial and remaining sub-bundle data. This achieves information separation of the proton beam along different propagation paths. Furthermore, calculating the water equivalent path length corresponding to each sub-bundle based on multiple sub-bundle data pairs can significantly reduce the multiple Coulomb scattering coupling caused by lateral broadening and superposition between different sub-bundles, thereby reducing the impact of scattering contamination on resolution. Based on this, when reconstructing the target image, higher precision can be achieved in distinguishing different materials within the imaged object, effectively improving the spatial resolution and accuracy of the final image. Attached Figure Description

[0019] The above-described features, other objects, and advantages of the present invention will become clearer from the following description of embodiments of the invention with reference to the accompanying drawings, in which:

[0020] Figure 1 An application scenario diagram of the proton imaging method and device according to an embodiment of the present invention is shown.

[0021] Figure 2 A flowchart of a proton imaging method according to an embodiment of the present invention is shown.

[0022] Figure 3 A schematic diagram illustrating the determination of the propagation directions of multiple sub-beams according to an embodiment of the present invention is shown.

[0023] Figure 4 A schematic diagram of proton imaging according to an embodiment of the present invention is shown.

[0024] Figure 5 A block diagram of an electronic device suitable for implementing a proton imaging method according to an embodiment of the present invention is shown. Detailed Implementation

[0025] Hereinafter, embodiments of the present invention will be described with reference to the accompanying drawings. However, it should be understood that these descriptions are exemplary only and are not intended to limit the scope of the invention. In the following detailed description, numerous specific details are set forth to provide a thorough understanding of the embodiments of the invention for ease of explanation. However, it will be apparent that one or more embodiments may be practiced without these specific details. Furthermore, descriptions of well-known structures and techniques are omitted in the following description to avoid unnecessarily obscuring the concept of the invention.

[0026] The terminology used herein is for the purpose of describing particular embodiments only and is not intended to limit the invention. The terms “comprising,” “including,” etc., as used herein indicate the presence of the stated features, steps, operations, and / or components, but do not exclude the presence or addition of one or more other features, steps, operations, or components.

[0027] All terms used herein (including technical and scientific terms) have the meanings commonly understood by those skilled in the art, unless otherwise defined. It should be noted that the terms used herein are to be interpreted in a manner consistent with the context of this specification, and not in an idealized or overly rigid way.

[0028] When using expressions such as "at least one of A, B and C", they should generally be interpreted in accordance with the meaning that is commonly understood by those skilled in the art (e.g., "a system having at least one of A, B and C" should include, but is not limited to, a system having A alone, a system having B alone, a system having C alone, a system having A and B, a system having A and C, a system having B and C, and / or a system having A, B and C, etc.).

[0029] When a proton beam propagates through a medium, the energy deposition exhibits typical Bragg peak characteristics. Specifically, the energy release reaches its maximum at the end of the propagation path and then rapidly decays.

[0030] Applying proton beams to the medical field, such as cancer treatment, allows most of the proton beam's energy to be concentrated on the tumor area of ​​the patient, reducing the risk of damage to surrounding healthy tissues.

[0031] In clinical proton therapy, the accuracy of proton beam range prediction directly affects the precision of proton therapy. Typically, patients first undergo X-ray computed tomography (XCT) to obtain CT images of their tissues. These images are then converted into a relative stopping power (RSP) distribution for protons, and the proton beam range is calculated based on the RSP. However, in these techniques, the conversion of CT values ​​from CT images to RSP introduces an error of approximately 3.5% to 5%, thus reducing the accuracy of proton beam range prediction.

[0032] To address the aforementioned issues, proton imaging technology has been proposed, which can directly measure the RSP distribution of the imaged object, reduce the conversion error from CT value to RSP, and thus improve the accuracy of range calculation.

[0033] The aforementioned proton imaging technique utilizes a detector to acquire the one-dimensional dose distribution of a proton beam after it passes through an imaging object. This one-dimensional dose distribution is then used to calculate the proton range. By comparing the range difference between the proton beam with and without the imaging object, the water-equivalent length of the imaging object is obtained, and the RSP distribution of the imaging object is reconstructed. However, this proton imaging technique only acquires one-dimensional dose distribution information in the depth direction. The lateral broadening of the proton beam is typically on the order of millimeters to centimeters, causing protons to pass through multiple regions simultaneously when incident on the imaging object, resulting in overlapping range information from different paths. Furthermore, multiple Coulomb scattering of protons within the object further amplifies the lateral broadening of the proton beam, leading to a decrease in the spatial resolution of the reconstructed image and making it difficult to accurately reconstruct internal structural details.

[0034] In view of this, this invention proposes a proton imaging method. Addressing the problem in related technologies where one-dimensional dose distribution imaging only acquires information in the depth direction, leading to poor imaging resolution, this invention collects two-dimensional dose deposition distributions of the proton beam along both the lateral and depth directions. This yields the initial two-dimensional dose distribution of the proton beam before passing through the imaging object and the remaining two-dimensional dose distribution after passing through the object, thereby introducing two-dimensional spatial information to extract lateral scattering characteristics. Furthermore, to accurately extract the water equivalent length from the two-dimensional dose distribution, this invention divides the initial and remaining two-dimensional dose distributions based on the virtual source location, obtaining multiple sub-beam data pairs that include initial and remaining sub-beam data. This achieves separation of the proton beam propagation information at different locations. Based on these multiple sub-beam data pairs, multiple water equivalent path lengths are obtained, reducing resolution loss caused by proton beam broadening and multiple Coulomb scattering. This improves the spatial resolution of proton imaging when reconstructing images based on multiple water equivalent path lengths.

[0035] Figure 1 An application scenario diagram of the proton imaging method and device according to an embodiment of the present invention is shown.

[0036] like Figure 1 As shown, the application scenarios according to this embodiment may include proton sources, collimators, imaging objects, detectors, and processors.

[0037] The proton source is used to emit a proton beam. After being collimated by a collimator, the proton beam passes through the imaging object. As the proton beam passes through the imaging object, the dose of the proton beam is absorbed by the imaging object, resulting in a difference between the initial two-dimensional dose distribution of the proton beam before passing through the imaging object and the remaining two-dimensional dose distribution after passing through the imaging object.

[0038] The detector acquires the initial two-dimensional dose distribution and the remaining two-dimensional dose distribution of the proton beam, and transmits the acquired initial two-dimensional dose distribution and the remaining two-dimensional dose distribution to the processor.

[0039] The processor processes the received initial two-dimensional dose distribution and the remaining two-dimensional dose distribution to complete the image reconstruction of the imaging object and obtain the target image.

[0040] Figure 2 A flowchart of an imaging method according to an embodiment of the present invention is shown.

[0041] like Figure 2 As shown, the proton imaging method of this embodiment includes operations S210 to S250.

[0042] In operation S210, the virtual source location is determined based on the initial two-dimensional dose distribution of the proton beam acquired in the absence of an imaging object.

[0043] The initial two-dimensional dose distribution characterizes the dose deposition distribution along the lateral and depth directions formed in the detector after the proton beam passes through free space.

[0044] In operation S220, based on the virtual source position and the multiple predetermined sub-bundle center positions of the predetermined plane, multiple sub-bundle propagation paths from the virtual source position to the multiple predetermined sub-bundle center positions are obtained.

[0045] The predetermined plane is perpendicular to the plane containing the sector region indicated by the proton beam, and the predetermined plane passes through the rotation center axis of the imaging object.

[0046] In operation S230, for any sub-beam propagation path in the sub-beam propagation path, and any two-dimensional dose distribution in the initial two-dimensional dose distribution and the remaining two-dimensional dose distribution under at least one projection viewpoint, a sub-beam data pair corresponding to the sub-beam propagation path and the two-dimensional dose distribution is obtained based on the two-dimensional dose distribution and the lateral distribution component determined along the sub-beam propagation path.

[0047] The remaining two-dimensional dose distribution characterizes the dose deposition distribution along the lateral and depth directions formed in the detector after the proton beam passes through the imaging object.

[0048] In operation S240, based on multiple sub-bundle data pairs corresponding to at least one projection angle, multiple water equivalent path lengths corresponding to at least one projection angle are obtained through an iterative algorithm.

[0049] In operation S250, image reconstruction is performed based on the multiple water equivalent path lengths corresponding to at least one projection angle to obtain the target image of the imaging object.

[0050] The imaging object can be any shape. By collecting the remaining two-dimensional dose distribution formed in the detector after the proton beam passes through the imaging object and absorbs part of its energy, the imaging object can be imaged.

[0051] The proton beam used to pass through the imaging object is emitted by a proton source. The proton source is fixed, and by rotating the imaging object, the two-dimensional dose distribution of the proton beam is collected from different projection angles, thereby realizing the acquisition of multi-angle imaging data.

[0052] The initial two-dimensional dose distribution can be the dose deposition distribution along the lateral and depth directions formed in the detector before the proton beam passes through the imaging object. The depth direction can be the direction of the proton beam's incident axis, used to describe the dose variation with penetration distance; the lateral direction is perpendicular to the depth direction, used to describe the dose distribution of the proton beam in the lateral direction; the dose distribution can be the energy distribution of the proton beam in the detector. The remaining two-dimensional dose distribution can be the remaining two-dimensional dose distribution formed in the detector after the proton beam passes through the imaging object and some energy is absorbed by the imaging object.

[0053] The initial two-dimensional dose distribution and the remaining two-dimensional dose distribution can be acquired by a detector placed at the rear of the imaging object. Specifically, when the proton beam is incident on the detector's scintillator, the scintillator deposits the proton beam's energy and generates scintillation light. The detector's camera then begins exposure, acquiring the two-dimensional visible light distribution image generated by the scintillator and transmitting it to the processor. The processor performs optical distortion correction, filtering, and noise reduction on the two-dimensional visible light distribution image to restore the two-dimensional dose distribution of the proton beam.

[0054] The virtual source position is not the actual proton source. After the proton source emits the proton beam, it is adjusted by the collimator and then diverges. The backward extensions of the proton beam at different positions in the divergent state intersect at a certain point in space, which is the virtual source position of the proton beam.

[0055] The predetermined plane can be determined based on the rotational central axis of the imaging object and the plane containing the sector region indicated by the proton beam. The plane passing through the rotational central axis of the imaging object and perpendicular to the proton beam is used as the predetermined plane. During propagation, the proton beam forms a geometric sector on the designated plane; this geometric sector is the sector region indicated by the proton beam, and the designated plane is the plane containing the sector region indicated by the proton beam.

[0056] The rotation center axis of the imaging object can be the central axis used to rotate the imaging object when detecting data of the imaging object at different projection angles.

[0057] Sub-beam data pairs can be dose deposition distributions along a proton beam in a certain sub-beam propagation direction.

[0058] Given a determined virtual source location, the two-dimensional dose distribution of the proton beam can be spatially partitioned based on this location, resulting in multiple sub-beam data pairs. Before passing through the imaging object, the proton beam forms an initial two-dimensional dose distribution in the detector; after passing through the imaging object, the dose distribution formed in the detector is updated to the remaining two-dimensional dose distribution. Because the imaging object has different energy deposition effects on proton beams at different spatial locations, the corresponding sub-beams also change. Therefore, each sub-beam data pair includes initial and remaining sub-beam data along the same propagation path, reflecting the changes in the proton beam's characteristics before and after penetrating the imaging object along that path.

[0059] The propagation path of each sub-bundle data pair is determined by the virtual source location and a specified location on a predetermined plane, which can be a location on the predetermined plane determined according to a predetermined interval.

[0060] The predetermined sub-bundle center position can be determined based on a predetermined interval. Multiple predetermined sub-bundle center positions can be obtained using the predetermined interval, and then the sub-bundle propagation path of each sub-bundle can be determined based on the line connecting the virtual source position and the predetermined sub-bundle center position.

[0061] The lateral distribution component can be a one-dimensional Gaussian distribution of the sub-bundle at different depth positions associated with the predetermined sub-bundle center position and lateral broadening.

[0062] Given multiple sub-beam propagation paths, the initial sub-beam data in the sub-beam data pair can be determined based on the initial two-dimensional dose distribution and the lateral distribution components along the sub-beam propagation paths.

[0063] Given multiple sub-beam propagation paths, the remaining sub-beam data in the sub-beam data pair can be accurately determined based on the remaining two-dimensional dose distribution and the lateral distribution components along the sub-beam propagation paths.

[0064] The water equivalent path length can be a difference parameter between the remaining two-dimensional dose distribution and the initial two-dimensional dose distribution. The water equivalent path length for each of multiple sub-bundle data pairs can be a difference parameter between the initial sub-bundle data and the remaining sub-bundle data for each sub-bundle data pair. The water equivalent path length is used for image reconstruction to determine the relative stopping power distribution of the imaged object, i.e., the target image.

[0065] According to embodiments of the present invention, by acquiring the initial two-dimensional dose distribution of the proton beam, the virtual source position can be accurately determined based on the initial two-dimensional dose distribution. This virtual source position is then used to divide the initial two-dimensional dose distribution and the remaining two-dimensional dose distribution, resulting in multiple sub-bundle data pairs composed of initial and remaining sub-bundle data. This achieves information separation of the proton beam along different propagation paths. Furthermore, calculating the water equivalent path length corresponding to each sub-bundle based on multiple sub-bundle data pairs can significantly reduce the multiple Coulomb scattering coupling caused by lateral broadening and superposition between different sub-bundles, thereby reducing the impact of scattering contamination on resolution. Based on this, when reconstructing the target image, higher precision can be achieved in distinguishing different materials within the imaged object, effectively improving the spatial resolution and accuracy of the final image.

[0066] According to an embodiment of the present invention, determining the virtual source location based on the initial two-dimensional dose of the proton beam acquired in the absence of an imaging object includes: constructing a first objective function for the lateral broadening of the proton beam as a function of depth based on the initial two-dimensional dose distribution of the proton beam; and determining the virtual source location based on the first objective function.

[0067] Based on the initial two-dimensional dose distribution of the proton beam obtained without imaging objects, the half-width at half-maximum (FWHM) (i.e., lateral broadening) of the proton beam at different depth positions is extracted from the initial two-dimensional dose distribution. A first objective function is constructed to describe the lateral broadening of the proton beam as a function of depth, which is used to fit the beam diameter broadening law of the proton beam. The first objective function can be expressed as y=ax+b, where x represents the depth and y represents the FWHM of the proton beam at different depths.

[0068] By fitting the first objective function using the least squares algorithm, parameters a and b are obtained to determine the law of change of the dimension with depth. Then, given the first objective function, the virtual source position is calculated.

[0069] According to an embodiment of the present invention, determining the virtual source location based on a first objective function includes: using the depth obtained by extrapolating the first objective function to a predetermined lateral width as the virtual source location.

[0070] The predetermined dimension can be a position where the lateral width is 0, that is, the half-width of the proton beam is 0.

[0071] Based on the first objective function y=ax+b, given the parameters a and b, the first objective function is extrapolated to a predetermined lateral width, i.e., y = 0, and the corresponding depth x is solved. This x is the virtual source position of the proton beam.

[0072] According to an embodiment of the present invention, by extrapolating the first objective function to a lateral broadening of a predetermined size, the depth position corresponding to the position where the half-width at half-maximum (HWHM) is 0 can be determined based on the beam diameter broadening law of the proton beam, so as to accurately determine the virtual source position.

[0073] According to an embodiment of the present invention, obtaining multiple sub-beam propagation paths from the virtual source location to the multiple predetermined sub-beam center locations based on the virtual source location and multiple predetermined sub-beam center locations on a predetermined plane includes: determining the central axis of the proton beam based on the initial two-dimensional dose distribution of the proton beam, wherein the central axis of the proton beam is perpendicular to the predetermined plane; determining multiple sub-beam propagation lines based on the virtual source location and the multiple predetermined sub-beam center locations on the predetermined plane, wherein a predetermined interval is maintained between the predetermined sub-beam center locations and a preset point on the predetermined plane, and the preset point is the intersection of the central axis of the proton beam and the predetermined plane; determining the propagation direction of each of the multiple sub-beam propagation lines based on the angle between each of the multiple sub-beam propagation lines and the central axis; and obtaining multiple sub-beam propagation paths from the virtual source location to the multiple predetermined sub-beam center locations based on the propagation directions of each of the multiple sub-beam propagation lines.

[0074] The predetermined plane is perpendicular to the plane containing the sector region indicated by the proton beam and passes through the rotational axis of the imaging object.

[0075] Among them, the center positions of multiple predetermined sub-beams are located on the central vertical line of the predetermined plane, which is the intersection line between the plane where the proton beam indicates the fan-shaped region and the predetermined plane.

[0076] The distance between the predetermined sub-beam center position and the preset point on the predetermined plane is a predetermined interval. For example, the distance between the first predetermined sub-beam center position and the preset point is 0.5, the distance between the second predetermined sub-beam center position and the imaging isocenter position is 1.0, and so on, with the interval between every two predetermined sub-beam center positions being consistent.

[0077] Based on the angle between the propagation line and the baseline, the slopes of each propagation line are determined, thereby determining the propagation direction of each propagation line.

[0078] Based on the propagation direction of each of the multiple propagation lines, multiple sub-bundle propagation paths are obtained from the virtual source position to the center positions of multiple predetermined sub-bundles.

[0079] Regarding how to determine the direction of transmission, the following will be combined with... Figure 3 Please provide an explanation.

[0080] Figure 3 A schematic diagram illustrating the determination of the propagation direction according to an embodiment of the present invention is shown.

[0081] like Figure 3As shown, a predetermined plane is determined based on the plane of the sector region indicated by the proton beam and the rotation center axis of the imaging object. Multiple predetermined sub-beam center positions are determined on the predetermined plane based on a predetermined interval. The proton beam center axis is used as a reference line, and the line connecting each predetermined sub-beam center position and the virtual source position is used as a propagation line. The propagation direction of each predetermined sub-beam center position is determined based on the angle between the propagation line and the reference line.

[0082] According to an embodiment of the present invention, by determining the propagation direction of each predetermined sub-beam center position based on the angle between the propagation line and the reference line, the depth dose component corresponding to each propagation direction can be accurately determined in the two-dimensional dose distribution of the proton beam based on the propagation direction.

[0083] According to an embodiment of the present invention, based on the two-dimensional dose distribution and the lateral distribution component determined along the sub-beam propagation path, multiple sub-beam data pairs corresponding to the sub-beam propagation path and the two-dimensional dose distribution are obtained, including: for any sub-beam propagation path in the sub-beam propagation path, determining the center position of the sub-beam at different depths; determining the lateral broadening of the sub-beam at different depths based on a scattering algorithm; determining the lateral distribution component of the sub-beam along the sub-beam propagation path based on the center position and lateral broadening of the sub-beam at different depths; and obtaining the sub-beam data pairs corresponding to the sub-beam propagation path and the two-dimensional dose distribution based on the two-dimensional dose distribution and the lateral distribution component.

[0084] Sub-bundles at different depths have different full width at half maximum (FWHM) broadenings.

[0085] For each sub-beam, the full width at half maximum (FWHM) at different depths is calculated based on a scattering algorithm to determine the lateral broadening of the sub-beam at different depths.

[0086] The lateral distribution components are determined by the one-dimensional Gaussian distribution based on the sub-bundle center positions and lateral broadening at different depths.

[0087] According to an embodiment of the present invention, by combining the propagation direction of the sub-bundle with its lateral broadening characteristics during propagation, the lateral distribution of the sub-bundle at different depths is accurately calculated, so that the generated sub-bundle data is closer to the propagation and scattering characteristics of the actual proton beam, thereby improving the accuracy of the obtained sub-bundle data pairs.

[0088] According to an embodiment of the present invention, obtaining multiple sub-bundle data pairs corresponding to multiple sub-bundle propagation paths and two-dimensional dose distributions based on a two-dimensional dose distribution and a lateral distribution component along the sub-bundle propagation path includes: convolving an initial two-dimensional dose distribution and a lateral distribution component to obtain multiple initial sub-bundle data corresponding to multiple sub-bundle propagation paths and an initial two-dimensional dose distribution; convolving a remaining two-dimensional dose distribution and a lateral distribution component to obtain multiple remaining sub-bundle data corresponding to multiple sub-bundle propagation paths and remaining two-dimensional dose distributions; and determining multiple sub-bundle data pairs corresponding to multiple sub-bundle propagation paths and two-dimensional dose distributions based on the multiple initial sub-bundle data and the multiple remaining sub-bundle data.

[0089] The initial sub-bundle data can be determined based on the convolution results of the initial two-dimensional dose distribution and the lateral distribution components.

[0090] The remaining sub-bundle data can be determined based on the convolution results of the remaining two-dimensional dose distribution and the lateral distribution components.

[0091] The initial and remaining sub-bundle data for the same sub-bundle propagation path are identified as sub-bundle data pairs corresponding to the sub-bundle propagation path. Then, based on multiple initial and remaining sub-bundle data sets, sub-bundle data pairs for each of the multiple sub-bundle propagation paths are determined.

[0092] According to an embodiment of the present invention, based on multiple sub-bundle data pairs corresponding to at least one projection angle, an iterative algorithm is used to obtain multiple water equivalent path lengths corresponding to at least one projection angle. This includes: for any projection viewpoint among at least one projection viewpoint, based on the difference between the total remaining sub-bundle data and the total translation sub-bundle data of each of the multiple iteration rounds corresponding to the projection viewpoint, determining multiple target water equivalent path lengths for a target round from multiple candidate water equivalent path lengths for each of the multiple rounds, and using the multiple target water equivalent path lengths as the multiple water equivalent path lengths corresponding to the projection angle; wherein, the total remaining sub-bundle data is determined based on the remaining sub-bundle data of each of the multiple sub-bundle data pairs corresponding to the projection viewpoint, the total translation sub-bundle data of the iteration rounds corresponding to the projection viewpoint is determined based on the translation sub-data of each of the multiple sub-bundle data pairs corresponding to the projection viewpoint, and the translation sub-data of the sub-bundle data pairs corresponding to the projection viewpoint is obtained by translating the initial sub-bundle data in the sub-bundle data pairs corresponding to the projection viewpoint along the depth direction to the candidate water equivalent path length corresponding to the iteration round.

[0093] Among them, multiple initial water equivalent path lengths are used as initial estimates for the iterative solution process. Their number corresponds one-to-one with multiple sub-bundle data pairs and can be obtained based on the characteristics of the remaining two-dimensional dose distribution, empirical models, or preset default values.

[0094] Determining the equivalent water path lengths corresponding to the projection viewpoint requires multiple iterations. For each iteration, the candidate equivalent water path lengths for each sub-bundle data pair are adjusted so that shifting the initial sub-bundle data along the depth direction by the corresponding candidate equivalent water path length approximates the corresponding remaining sub-bundle data as closely as possible. By updating the candidate equivalent water path lengths through multiple iterations, the difference between the shifted initial sub-bundle data and the remaining sub-bundle data gradually decreases. When the difference meets the convergence condition or reaches the preset number of iterations, the final target equivalent water path length is determined.

[0095] In this process, due to the mutual influence between sub-bundles caused by multiple Coulomb scattering, the data of any sub-bundle not only contains information about its own propagation path, but is also affected by the coupling interference of the propagation characteristics of neighboring sub-bundles. Therefore, in each iteration, it is necessary to jointly calculate the total remaining sub-bundle data and the total translated sub-bundle data based on the candidate water equivalent path length of each sub-bundle data pair, and use the overall difference between the two as the evaluation index to ensure that the update of the candidate water equivalent path length can reflect the coupling relationship between all sub-bundles, thereby obtaining a more stable solution that is closer to the real physical process.

[0096] Summing multiple remaining sub-bundle data can be done by summing the doses of multiple remaining sub-bundle data to obtain the total remaining sub-bundle data.

[0097] The dose of the initial sub-bundle data is translated along the depth direction based on the candidate water equivalent path length diameter to obtain multiple translated initial sub-bundle data. Furthermore, the doses of the multiple translated initial sub-bundle data are summed to obtain the total translated sub-bundle data.

[0098] According to an embodiment of the present invention, after completing multiple iterations, the target water equivalent path length can be determined from the candidate water equivalent path lengths obtained in each iteration. Specifically, one of the following two methods can be used.

[0099] One approach is to terminate early based on a threshold: when the candidate water equivalent path length obtained in a certain round makes the difference between the total translation sub-bundle and the total remaining sub-bundle satisfy a preset difference threshold, the iteration can be terminated, and the candidate water equivalent path length of that round is determined as the target water equivalent path length.

[0100] Another approach is to select the optimal path based on a limited number of iterations: after the iteration reaches a preset number of iterations, the candidate water equivalent path lengths obtained in each iteration can be compared, and the candidate water equivalent path length that minimizes the difference between the total translation sub-bundle and the total remaining sub-bundle can be selected as the target water equivalent path length.

[0101] According to an embodiment of the present invention, the method further includes: for any one of the multiple rounds, based on a second objective function, determining the difference between the total remaining sub-bundle data and the total translation sub-bundle data of the round corresponding to the projection viewpoint, according to the total remaining sub-bundle data, the total translation sub-bundle data, and the candidate water equivalent path length of the iteration round corresponding to the projection viewpoint; wherein, the second objective function includes a data fitting term, a prior constraint term, and a sparsity term, the data fitting term is used to determine the difference between the total remaining sub-bundle data and the total translation sub-bundle data of the iteration round corresponding to the projection viewpoint, the prior constraint term is used to constrain the multiple candidate water equivalent path lengths of the multiple iteration rounds corresponding to the projection viewpoint according to the multiple initial water equivalent path lengths corresponding to the projection viewpoint, and the sparsity term is used to perform sparse constraints on the candidate water equivalent path lengths.

[0102] The second objective function can be shown in formula (1).

[0103] (1)

[0104] Where z is the depth, For data fitting terms, These are prior constraints. For sparse constraint terms, For multiple target water equivalent path lengths, For multiple candidate water equivalent path lengths, Let be the equivalent path length of the i-th sub-bundle data pair among multiple candidate water equivalent path lengths. Let be the initial water equivalent path length corresponding to the i-th sub-bundle data pair among multiple initial water equivalent path lengths.

[0105] For the data fitting term, the data fitting term is used to determine the total remaining sub-bundle data. Total translational sub-beam data The differences between them, That is, the total remaining sub-bundle data For n remaining sub-bundles of data The summation is obtained as follows: This represents the data of the i-th remaining sub-bundle. The total translation sub-bundle data represents the data from n translation sub-bundles. The summation is obtained as follows: This represents the data of the i-th initial sub-bundle.

[0106] For the prior constraint term, β is the prior weight coefficient, and the prior constraint term is used to determine the equivalent path length of the i-th candidate water. Equivalent path length of initial water The difference between them ensures that the obtained candidate water equivalent path length does not deviate excessively from the initial water equivalent path length. The value can be selected according to actual needs, and no restrictions are imposed here.

[0107] For sparse terms, sparse terms Used to suppress noise while preserving the edges and details of the candidate water equivalent path length, λ is the prior weight coefficient, and W is a sparse transform operator, such as wavelet transform or discrete cosine transform.

[0108] According to an embodiment of the present invention, the difference between the total translation sub-bundle data and the total residual sub-bundle data of the projection viewpoint corresponding to the candidate water equivalent path length is determined based on a second objective function. The data fitting term in the second objective function is used to compare the total translation sub-bundle data with the total residual sub-bundle data, thereby accurately selecting the candidate water equivalent path length that best conforms to the actual propagation law; the prior constraint term is used to limit the candidate water equivalent path length of each sub-bundle data pair from deviating excessively from its corresponding initial water equivalent path length; the sparsity term is used to suppress the interference of noise on the update of the candidate water equivalent length during the iteration process, thereby improving the stability and robustness of the iteration process.

[0109] According to an embodiment of the present invention, image reconstruction is performed based on multiple water equivalent path lengths corresponding to at least one projection angle to obtain a target image of the imaging object. This includes: processing and back-projecting the multiple water equivalent path lengths corresponding to at least one projection angle based on a filtered back-projection algorithm to obtain a target image of the imaging object.

[0110] The filtered back projection algorithm can process multiple water equivalent path lengths, assign weights to multiple water equivalent path lengths corresponding to different projection angles, suppress high-frequency noise, and obtain accurate target images of the imaging object.

[0111] A second aspect of the present invention provides an imaging device, such as... Figure 1 As shown, it includes a detector and a processor.

[0112] The detector is configured to: acquire an initial two-dimensional dose distribution of a proton beam in the absence of an imaging object, and acquire residual two-dimensional dose distributions at at least one projection viewpoint, wherein the initial two-dimensional dose distribution characterizes the dose deposition distribution in the detector along the lateral and depth directions after the proton beam passes through free space, and the residual two-dimensional dose distribution characterizes the dose deposition distribution in the detector along the lateral and depth directions after the proton beam passes through an imaging object.

[0113] The processor is configured to: obtain multiple sub-beam propagation paths from the virtual source location to the multiple predetermined sub-beam center locations based on the virtual source location and multiple predetermined sub-beam center locations on a predetermined plane of the imaging object; for any propagation path in the propagation path, and any two-dimensional dose distribution in the initial two-dimensional dose distribution and the remaining two-dimensional dose distribution under at least one projection viewpoint, obtain sub-beam data pairs corresponding to the sub-beam propagation path and the two-dimensional dose distribution based on the two-dimensional dose distribution and the lateral distribution component determined along the sub-beam propagation path; obtain multiple water equivalent path lengths corresponding to at least one projection angle based on the multiple sub-beam data pairs corresponding to at least one projection angle; and reconstruct the image of the imaging object based on the multiple water equivalent path lengths corresponding to at least one projection angle to obtain the target image of the imaging object.

[0114] In addition, the imaging device also includes a proton source for generating a proton beam, and simultaneously generates a signal that is transmitted to the detector while the proton beam is being generated.

[0115] The imaging device also includes a collimator, which can be made of metals such as brass or nickel, with a 1mm×50mm through-hole in the center to collimate the proton beam width into a flat fan beam.

[0116] The detector comprises a shielding enclosure, a scintillator, and a camera, used to detect the two-dimensional dose distribution information of the proton beam. The detector begins operation after receiving a signal from the proton source. The shielding enclosure consists of a support frame and light-shielding material with a light-blocking rate of over 99%, used to isolate external light sources and prevent them from interfering with the detection process. When the fan-shaped proton beam is incident on the scintillator, energy is deposited in the scintillator, generating scintillation light, and the camera captures the two-dimensional visible light distribution image generated by the scintillator.

[0117] The processor performs optical distortion correction, filtering, and noise reduction on the two-dimensional visible light distribution acquired by the camera, restoring the visible light distribution to a two-dimensional dose distribution of the proton beam. Simultaneously, the two-dimensional proton dose distribution is divided into several sub-beams to obtain multiple sub-beam data pairs, and the equivalent water path length corresponding to at least one projection angle is determined, thereby enabling image reconstruction of the imaged object.

[0118] Figure 4 An imaging schematic diagram according to an embodiment of the present invention is shown.

[0119] like Figure 4 As shown, image reconstruction was performed on imaging object 1 and imaging object 2 respectively. Due to the different materials of imaging object 1 and imaging object 2, their scattering and energy loss effects on the proton beam differ, resulting in differences in their respective residual two-dimensional dose distributions. The water equivalent path lengths extracted based on these different residual two-dimensional dose distributions are different, thus leading to different reconstructed relative stopping power images. From... Figure 4As shown in the transverse and longitudinal cross-sectional views, the RSP images obtained based on the imaging method of this invention can accurately reproduce the geometry of the imaged object, with only slight blurring at the interfaces between different materials. For imaged object 1, the relative error of its acrylic material is within 1%, and the absolute error in the air region is less than 0.02. For imaged object 2, the relative errors of its acrylic, polytetrafluoroethylene, polyethylene, and polyoxymethylene resin are −0.31%, 0.50%, −1.21%, and 0.61%, respectively, and the absolute errors in the air region are all less than 0.02. These results demonstrate that the imaging method of this invention can achieve high-precision image reconstruction.

[0120] Figure 5 A block diagram of an electronic device suitable for implementing an imaging method according to an embodiment of the present invention is shown.

[0121] like Figure 5 As shown, an electronic device 500 according to an embodiment of the present invention includes a processor 501, which can perform various appropriate actions and processes according to a program stored in a ROM 502 (read-only memory) or a program loaded from a storage portion 508 into a RAM 503 (random access memory). The processor 501 may include, for example, a general-purpose microprocessor (e.g., a CPU), an instruction set processor and / or an associated chipset and / or a special-purpose microprocessor (e.g., an application-specific integrated circuit (ASIC)), etc. The processor 501 may also include onboard memory for caching purposes. The processor 501 may include a single processing unit or multiple processing units for performing different actions of the method flow according to an embodiment of the present invention.

[0122] RAM 503 stores various programs and data required for the operation of electronic device 500. Processor 501, ROM 502, and RAM 503 are interconnected via bus 504. Processor 501 executes various operations of the method flow according to embodiments of the present invention by executing programs in ROM 502 and / or RAM 503. It should be noted that the programs may also be stored in one or more memories other than ROM 502 and RAM 503. Processor 501 may also execute various operations of the method flow according to embodiments of the present invention by executing programs stored in said one or more memories.

[0123] According to an embodiment of the present invention, the electronic device 500 may further include an input / output (I / O) interface 505, which is also connected to a bus 504. The electronic device 500 may also include one or more of the following components connected to the input / output (I / O) interface 505: an input section 506 including a keyboard, mouse, etc.; an output section 507 including a cathode ray tube (CRT), liquid crystal display (LCD), etc., and a speaker, etc.; a storage section 508 including a hard disk, etc.; and a communication section 509 including a network interface card such as a LAN card, modem, etc. The communication section 509 performs communication processing via a network such as the Internet. A drive 510 is also connected to the input / output (I / O) interface 505 as needed. A removable medium 511, such as a disk, optical disk, magneto-optical disk, semiconductor memory, etc., is installed on the drive 510 as needed so that computer programs read from it can be installed into the storage section 508 as needed.

[0124] The present invention also provides a computer-readable storage medium, which may be included in the device / apparatus / system described in the above embodiments; or it may exist independently and not assembled into the device / apparatus / system. The computer-readable storage medium carries one or more programs, which, when executed, implement the method according to the embodiments of the present invention.

[0125] According to embodiments of the present invention, the computer-readable storage medium may be a non-volatile computer-readable storage medium, such as including, but not limited to: portable computer disks, hard disks, random access memory (RAM), read-only memory (ROM), erasable programmable read-only memory (EPROM or flash memory), portable compact disk read-only memory (CD-ROM), optical storage devices, magnetic storage devices, or any suitable combination thereof. In the present invention, the computer-readable storage medium may be any tangible medium containing or storing a program that can be used by or in conjunction with an instruction execution system, apparatus, or device. For example, according to embodiments of the present invention, the computer-readable storage medium may include ROM 502 and / or RAM 503 and / or one or more memories other than ROM 502 and RAM 503 described above.

[0126] Embodiments of the present invention also include a computer program product comprising a computer program containing program code for performing the methods shown in the flowchart. When the computer program product is run on a computer system, the program code is used to enable the computer system to implement the imaging method provided by the embodiments of the present invention.

[0127] When the computer program is executed by the processor 501, it performs the functions defined in the system / apparatus of this invention. According to embodiments of the invention, the systems, apparatuses, modules, units, etc., described above can be implemented by computer program modules.

[0128] In one embodiment, the computer program may rely on a tangible storage medium such as an optical storage device or a magnetic storage device. In another embodiment, the computer program may also be transmitted and distributed in the form of signals over a network medium, and may be downloaded and installed via the communication section 509, and / or installed from a removable medium 511. The program code contained in the computer program can be transmitted using any suitable network medium, including but not limited to: wireless, wired, etc., or any suitable combination thereof.

[0129] In such an embodiment, the computer program can be downloaded and installed from a network via communication section 509, and / or installed from removable medium 511. When the computer program is executed by processor 501, it performs the functions defined in the system of this embodiment of the invention. According to embodiments of the invention, the systems, devices, apparatuses, modules, units, etc., described above can be implemented by computer program modules.

[0130] According to embodiments of the present invention, program code for executing the computer programs provided in the embodiments of the present invention can be written in any combination of one or more programming languages. Specifically, these computational programs can be implemented using high-level procedural and / or object-oriented programming languages, and / or assembly / machine languages. Programming languages ​​include, but are not limited to, languages ​​such as Java, C++, Python, "C", or similar programming languages. The program code can be executed entirely on the user's computing device, partially on the user's device, partially on a remote computing device, or entirely on a remote computing device or server. In cases involving remote computing devices, the remote computing device can be connected to the user's computing device via any type of network, including a local area network (LAN) or a wide area network (WAN), or it can be connected to an external computing device (e.g., via the Internet using an Internet service provider).

[0131] The flowcharts and block diagrams in the accompanying drawings illustrate the architecture, functionality, and operation of possible implementations of systems, methods, and computer program products according to various embodiments of the present invention. In this regard, each block in a flowchart or block diagram may represent a module, segment, or portion of code containing one or more executable instructions for implementing a specified logical function. It should also be noted that in some alternative implementations, the functions indicated in the blocks may occur in a different order than those indicated in the drawings. For example, two consecutively indicated blocks may actually be executed substantially in parallel, and they may sometimes be executed in reverse order, depending on the functions involved. It should also be noted that each block in a block diagram or flowchart, and combinations of blocks in a block diagram or flowchart, may be implemented using a dedicated hardware-based system that performs the specified function or operation, or using a combination of dedicated hardware and computer instructions.

[0132] Those skilled in the art will understand that the features described in the various embodiments of the present invention can be combined and / or combined in various ways, even if such combinations or combinations are not explicitly described in the present invention. In particular, the features described in the various embodiments of the present invention can be combined and / or combined in various ways without departing from the spirit and teachings of the present invention. All such combinations and / or combinations fall within the scope of the present invention.

[0133] The embodiments of the present invention have been described above. However, these embodiments are merely illustrative and not intended to limit the scope of the invention. Although various embodiments have been described above, this does not mean that the measures in the various embodiments cannot be used advantageously in combination. Various substitutions and modifications can be made by those skilled in the art without departing from the scope of the invention, and all such substitutions and modifications should fall within the scope of the invention.

Claims

1. A method of proton imaging, characterized by, The method comprises the following steps: determining a virtual source position according to an initial two-dimensional dose distribution of a proton beam acquired without an imaging object, wherein the initial two-dimensional dose distribution represents a dose deposition distribution along a transverse direction and a depth direction formed in a detector after the proton beam passes through free space; obtaining a plurality of sub-beam propagation paths from the virtual source position to a plurality of predetermined sub-beam center positions according to the virtual source position and the plurality of predetermined sub-beam center positions on a predetermined plane, wherein the predetermined plane is perpendicular to a plane where a fan-shaped region indicated by the proton beam is located, and the predetermined plane passes through a rotation center axis of the imaging object; for any one of the sub-beam propagation paths and any one of the initial two-dimensional dose distribution and a remaining two-dimensional dose distribution under at least one projection viewing angle, obtaining a sub-beam data pair corresponding to the sub-beam propagation path and the two-dimensional dose distribution according to the two-dimensional dose distribution and a transverse distribution component determined along the sub-beam propagation path, wherein the remaining two-dimensional dose distribution represents a dose deposition distribution along a transverse direction and a depth direction formed in the detector after the proton beam passes through the imaging object; obtaining a plurality of water equivalent path lengths corresponding to the at least one projection angle respectively by an iterative algorithm according to a plurality of sub-beam data pairs corresponding to the at least one projection angle respectively; performing image reconstruction according to the plurality of water equivalent path lengths corresponding to the at least one projection angle respectively to obtain a target image of the imaging object.

2. The method of claim 1, wherein, The method for determining the virtual source position according to the initial two-dimensional dose of the proton beam acquired without the imaging object comprises the following steps: constructing a first target function of transverse broadening of the proton beam varying with depth according to the initial two-dimensional dose distribution of the proton beam; determining the virtual source position according to the first target function.

3. The method of claim 2, wherein, The method for determining the virtual source position according to the first target function comprises the following steps: taking a depth obtained when the first target function is extrapolated to the transverse broadening being a predetermined size as the virtual source position.

4. The method of claim 1, wherein, The method for obtaining the plurality of sub-beam propagation paths from the virtual source position to the plurality of predetermined sub-beam center positions according to the virtual source position and the plurality of predetermined sub-beam center positions on the predetermined plane comprises the following steps: determining a center axis of the proton beam according to the initial two-dimensional dose distribution of the proton beam, wherein the center axis of the proton beam is perpendicular to the predetermined plane; determining a plurality of sub-beam propagation lines according to the virtual source position and the plurality of predetermined sub-beam center positions on the predetermined plane, wherein the predetermined sub-beam center positions maintain a predetermined interval with preset points on the predetermined plane, and the preset points are intersection points of the center axis of the proton beam and the predetermined plane; determining a propagation direction of each of the plurality of sub-beam propagation lines according to an included angle between each of the plurality of sub-beam propagation lines and the center axis; and obtaining the plurality of sub-beam propagation paths from the virtual source position to the plurality of predetermined sub-beam center positions according to the propagation direction of each of the plurality of sub-beam propagation lines.

5. The method of claim 1, wherein, The process of obtaining multiple sub-beam data pairs corresponding to the sub-beam propagation path and the two-dimensional dose distribution based on the two-dimensional dose distribution and the lateral distribution component determined along the sub-beam propagation path includes: For any of the sub-beam propagation paths, Determine the center positions of the sub-bundles at different depths corresponding to the sub-bundle propagation path; Based on the scattering algorithm, the lateral broadening of the sub-beam at different depths is determined; Based on the center position and lateral widening of the sub-bundle at different depths, the lateral distribution component of the sub-bundle along the propagation path of the sub-bundle is determined; Based on the two-dimensional dose distribution and the lateral distribution component, sub-bundle data pairs corresponding to the sub-bundle propagation path and the two-dimensional dose distribution are obtained.

6. The method of claim 5, wherein, The process of obtaining multiple sub-beam data pairs corresponding to the multiple sub-beam propagation paths and the two-dimensional dose distribution based on the two-dimensional dose distribution and the lateral distribution component along the sub-beam propagation path includes: Convolve the initial two-dimensional dose distribution and the lateral distribution component to obtain multiple initial sub-beam data corresponding to the multiple sub-beam propagation paths and the initial two-dimensional dose distribution; Convolve the remaining two-dimensional dose distribution and the lateral distribution component to obtain multiple remaining sub-bundle data corresponding to the multiple sub-bundle propagation paths and the remaining two-dimensional dose distribution; Based on the multiple initial sub-bundle data and the multiple remaining sub-bundle data, multiple sub-bundle data pairs corresponding to the multiple sub-bundle propagation paths and the two-dimensional dose distribution are determined.

7. The method of any one of claims 1-6, wherein, The step of obtaining multiple water equivalent path lengths corresponding to each of the at least one projection angle through an iterative algorithm based on multiple sub-bundle data pairs corresponding to each of the at least one projection angle includes: For any projection viewpoint among the at least one projection viewpoint, based on the difference between the total remaining sub-bundle data and the total translation sub-bundle data of each of the multiple iteration rounds corresponding to the projection viewpoint, the multiple target water equivalent path lengths of the target round are determined from the multiple candidate water equivalent path lengths of each of the multiple iteration rounds, and the multiple target water equivalent path lengths are taken as the multiple water equivalent path lengths corresponding to the projection angle. The total remaining sub-bundle data is determined based on the remaining sub-bundle data of multiple sub-bundle data corresponding to the projection viewpoint. The total translation sub-bundle data of the iteration round corresponding to the projection viewpoint is determined based on the translation sub-data of multiple sub-bundle data of the iteration round corresponding to the projection viewpoint. The translation sub-data is obtained by translating the initial sub-bundle data in the sub-bundle data pair of the iteration round corresponding to the projection viewpoint along the depth direction by the equivalent path length of the candidate water corresponding to the iteration round.

8. The method of claim 7, wherein, The method further includes: For any of the multiple iteration rounds, the difference between the total remaining sub-bundle data and the total translation sub-bundle data of the iteration round corresponding to the projection viewpoint is determined based on the second objective function, according to the total remaining sub-bundle data, the total translation sub-bundle data, and the candidate water equivalent path length of the iteration round corresponding to the projection viewpoint. The second target function includes a data fitting term, a prior constraint term and a sparse term, the data fitting term is used to determine the difference between the total residual beamlet data and the total translation beamlet data of the iteration round corresponding to the projection view angle, the prior constraint term is used to constrain the candidate water equivalent path lengths of multiple iteration rounds corresponding to the projection view angle according to the multiple initial water equivalent path lengths corresponding to the projection view angle, and the sparse term is used to sparsely constrain the candidate water equivalent path lengths.

9. The method of claim 1, wherein, The image reconstruction according to the multiple water equivalent path lengths corresponding to the at least one projection angle respectively includes: Based on a filtered back-projection algorithm, the multiple water equivalent path lengths corresponding to the at least one projection angle respectively are processed to obtain the target image of the imaging object.

10. An image forming apparatus characterized by comprising: It includes: The detector is configured to: acquire an initial two-dimensional dose distribution of a proton beam in the absence of an imaging object, and acquire a respective residual two-dimensional dose distribution under at least one projection view angle, wherein the initial two-dimensional dose distribution represents the dose deposition distribution along the transverse direction and the depth direction formed in the detector after the proton beam passes through the free space, and the residual two-dimensional dose distribution represents the dose deposition distribution along the transverse direction and the depth direction formed in the detector after the proton beam passes through the imaging object; and The processor is configured to: According to the virtual source position and the multiple predetermined beamlet center positions of the predetermined plane of the imaging object, multiple beamlet propagation paths from the virtual source position to the multiple predetermined beamlet center positions are obtained; For any one of the propagation paths and any one of the initial two-dimensional dose distribution and the residual two-dimensional dose distribution under at least one projection view angle, a beamlet data pair corresponding to the beamlet propagation path and the two-dimensional dose distribution is obtained through an iterative algorithm according to the two-dimensional dose distribution and the transverse distribution component determined along the beamlet propagation path; According to the multiple beamlet data pairs corresponding to the at least one projection angle respectively, the multiple water equivalent path lengths corresponding to the at least one projection angle respectively are obtained; and According to the multiple water equivalent path lengths corresponding to the at least one projection angle respectively, the imaging object is image reconstructed to obtain a target image of the imaging object.

Images