A proton track and dose simulation method based on stopping power and random scattering

Abstract

The application provides a proton track and dose simulation method based on stopping power and random scattering, comprising the following steps: S1. reading proton stopping power, continuous slowing-down approximation range and data of the target material; S2. using a basic radiation parameter table to generate an adaptive scattering step length suitable for the current energy; S3. taking the bypass factor as a criterion to optimize the empirical parameters in the scattering angle correction factor; S4. calculating the standard deviation of the scattering angle to generate a single scattering angle conforming to the Gaussian distribution; S5. calculating the single scattering track and energy loss of the proton; S6. iteratively simulating the single scattering and energy loss calculation until the kinetic energy of the proton is zero or the proton passes through the target material; and S7. outputting the complete track of the proton and the dose distribution curve. The application combines the proton stopping power curve and the scattering theory, takes the bypass factor as a criterion to ensure the correctness of the statistical behavior of particle scattering, and realizes the rapid simulation and accurate calculation of the proton transport track and radiation dose.

Description

Technical Field

[0001] This invention relates to the field of ion radiation and matter interaction, specifically to a method for simulating the radiation track of energetic protons in materials and calculating radiation dose, which is one of the foundations of various proton irradiation applications. Background Technology

[0002] In recent years, proton radiation has become an important method for cancer treatment due to its strong penetrating power and excellent energy deposition characteristics. Meanwhile, with the development of aerospace and space technology, the potential hazards of proton radiation from cosmic rays and solar activity to space equipment and astronauts are becoming increasingly prominent. In addition, proton irradiation also has wide applications in the microelectronics industry, such as semiconductor material modification. Rapid and accurate calculation of proton transport tracks and radiation dose distribution in target materials is a prerequisite for ensuring the accuracy and reliability of important applications such as proton therapy for cancer, space radiation damage protection, and semiconductor material modification.

[0003] Currently, mainstream methods for proton radiation calculations include complex calculations based on Monte Carlo simulations and simplified empirical models based on theoretical models and experimental observations. While the Monte Carlo method offers high computational accuracy, it is computationally complex and time-consuming, and ensuring statistical accuracy often requires enormous computational resources. This limits its feasibility in many real-time applications. Theoretically, the scattering models for protons at different energies differ. Proton acceleration and irradiation experiments are costly, and proton scattering and energy deposition data observed in various energy regions exhibit fragmentation. Therefore, multi-scale simulations of proton radiation often rely on numerous empirical parameters. Simplified empirical models, when dealing with real-world proton multiple scattering problems across energy regions, typically only serve as preliminary estimates. On the other hand, in proton irradiation simulations, accurately characterizing the natural randomness of multiple scattering of protons across different energy regions in the radiation effect is crucial for constructing observational cross-energy region proton tracks and, consequently, accurately calculating three-dimensional energy deposition.

[0004] Therefore, making full use of the accumulated data on the proton stopping ability of various materials, constructing criteria that can simulate the multiple scattering distribution of protons across energy regions, realizing the rapid and accurate simulation of the transport track of protons in the target material, and accurately calculating the radiation dose have become one of the key issues of concern in the proton irradiation-related industries.

[0005] To address the problems of high computational cost and time complexity in existing Monte Carlo methods, and the insufficient accuracy of simplified empirical models in accurately describing the randomness of multiple scattering across energy regions, this invention proposes a proton trajectory and dose simulation method based on stopping power and random scattering. This method combines proton stopping power data with scattering theory, using the detour factor as a global physical criterion, and adaptively adjusts the scattering step size and low-energy scattering corrections to achieve rapid and accurate simulation of proton transport trajectories and radiation dose. This avoids the limitations of complex multi-scale simulations or high-cost experimental observations common in traditional methods.

[0006] In this invention, unless otherwise stated, the following terms have the following meanings:

[0007] "Stopping ability" refers to the mass-stopping ability of a proton in a target material;

[0008] "Continuously slowed approximation range" refers to the proton range calculated based on the continuously slowed approximation model, and its standard notation is [symbol missing]. ;

[0009] The "detour factor" refers to the ratio of the projected range of a proton in the actual medium to the total path length. Summary of the Invention

[0010] In response to the broad needs of applications in this field, the present invention provides a method for simulating proton tracks and doses based on stopping power and random scattering.

[0011] The present invention adopts the following technical solution:

[0012] A method for simulating proton tracks and dose based on stopping power and random scattering includes the following steps:

[0013] S1: Read the proton stopping power, continuous slowing approximate range and deflection factor data of the target material; establish its interpolation function with respect to proton energy; combine the density and radiation length of the target material to construct a basic radiation parameter table covering the incident energy range; and construct a uniform three-dimensional target geometric model.

[0014] S2: Using the basic radiation parameter table constructed in step S1, an adaptive scattering step size suitable for the current energy is generated based on the current energy of the proton, stopping power, and preset energy loss upper limit, providing input parameters for the scattering and energy loss calculation in subsequent steps S4 and S5;

[0015] S3: Using the circumference factor in the basic radiation parameter table as the criterion, calculate the root mean square error between the simulated circumference factor and the theoretical value, and use an optimization algorithm to automatically adjust the values ​​of the empirical parameters A and B in the scattering angle correction factor to minimize the root mean square error.

[0016] S4: Calculate the initial standard deviation of the scattering angle based on the current proton energy, the adaptive scattering step size generated in step S2, and the target radiation length. The initial standard deviation of the scattering angle is adjusted by introducing correction factors that include the empirical parameters A and B determined through optimization in step S3, to obtain the final standard deviation of the scattering angle for the current step. Using the final standard deviation of the scattering angle as the standard deviation of the Gaussian distribution, a single scattering angle conforming to the Gaussian distribution is sampled and generated.

[0017] S5: Update the proton's motion direction based on the scattering angle generated in step S4. Calculate the energy loss within the scattering step based on the stopping power, target density, and current step size in the basic radiation parameter table. Record the energy loss in the spatial grid traversed by the scattering path according to the grid size to complete the recording of energy deposition.

[0018] S6: Repeat steps S2, S4 and S5 in sequence to simulate the continuous scattering and energy loss process of each proton in the target material until the proton energy is exhausted to zero or passes through the target material, thus completing the tracking simulation of a single proton track.

[0019] S7: Repeat the continuous scattering and energy loss process of more than 10,000 protons with the same initial energy in the target material, count the track and energy deposition data of all protons, and output the complete three-dimensional map of proton track, scattering angle-energy density map and dose-depth distribution curve.

[0020] The preferred method, step S1, involves establishing the interpolation function, which typically uses publicly available databases containing proton stopping power, continuously slowed approximate range, and bypass factor data. Cubic spline interpolation is then applied to the original data to construct an interpolation function covering the range of proton incident energy where each parameter varies with energy.

[0021] The optimal adaptive scattering step size for the current energy, as described in step S2, includes: using the current proton energy... And the quality of the target material to prevent the ability As input, the theoretical step size is calculated backward according to the principle that "energy loss does not exceed a preset proportion":

[0022] ,in This is the preset upper limit for the proportion of energy loss. To prevent the quality of the target material from being compromised. The target density is used. The theoretical step size obtained by back-calculation is compared with the preset maximum step size upper limit, and the smaller value is taken as the actual step size, which ensures both computational efficiency in the high-energy range and sufficient accuracy in the low-energy range.

[0023] The optimal method for determining the empirical parameters A and B mentioned in step S3 includes: statistically analyzing the projected ranges of all protons after the simulation. With total path length Calculate the average detour factor ; based on the theoretical value of the circumference factor in the aforementioned basic radiation parameter table. To achieve the objective, construct the root mean square error. Where N is the number of sample data points involved in the calculation; the Nelder-Mead simplex algorithm is used to... Let be the variable, RMSE be the loss function, and the error be minimized iteratively. The condition is met when RMSE remains constant for 20 consecutive iterations or is less than 1 × 10⁻⁶. -4 Stop at a certain time to achieve self-consistent calibration of empirical parameters A and B.

[0024] For the optimal method, step S4, which involves sampling a single scattering angle based on the final standard deviation of the scattering angle, includes: first calculating the initial standard deviation of the scattering angle using the Highland-Lynch-Dahl formula. :

[0025] .

[0026] Then introduce a low-energy empirical correction factor The final standard deviation of the scattering angle is obtained. ;by The variance is derived from a Gaussian distribution. Polar angle obtained by sampling azimuth exist Uniform sampling rotates the proton motion direction to a new direction, achieving physical consistency of the single scattering angle and low-energy enhancement effect.

[0027] The preferred method, the spatially resolved energy deposition recording method described in step S5, includes: discretizing each proton path into depth slices according to a set thickness, and calculating the spatial overlap length between the path and each slice. The energy loss in this step is distributed proportionally to the length. The corresponding slice is dynamically expanded and the energy is accumulated to finally output the dose-depth distribution curve, which can be directly used for proton therapy depth dose curve comparison.

[0028] The optimal multi-dimensional result output and visualization described in step S7 includes: a 3D track plot: recording the 3D coordinate string of each proton trajectory, supporting Matplotlib 3D line plot output; a scattering angle-energy density plot: displaying the relationship between the single-step scattering angle and the current energy in a log-log 2D histogram, with color mapping of event density; a dose-depth distribution curve: normalizing the slice energy deposition to 100% of the maximum dose; and a detour factor-energy curve: outputting a comparison chart of the optimized simulated detour factor and detour factor data in a public database to verify the reliability of the model.

[0029] By adopting the above technical solution, the beneficial effects of the present invention are:

[0030] This invention provides a method for rapidly constructing target material physical parameters based on a public database and interpolation encapsulation. Compared with traditional manual table lookup or experimental measurement, it can complete the continuous invocation of energy-stopping ability, energy-continuously slowed approximate range, and energy-bypass factor in milliseconds, significantly shortening the simulation preprocessing cycle and eliminating interpolation errors caused by parameter discretization.

[0031] This invention proposes a coupled strategy of "adaptive step size + low-energy scattering correction": on the one hand, by adjusting the adaptive step size in step S2, the low-energy step size is automatically reduced by setting an upper limit for energy loss; on the other hand, by correcting the low-energy scattering in step S3, a correction factor is introduced. The Highland-Lynch-Dahl formula is enhanced at low energies, with parameters A and B automatically optimized through deflection factor verification. These two steps work synergistically to ensure that the proton scattering angle distribution within the proton incident energy range simultaneously satisfies both microscopic cross-sectional statistics and macroscopic deflection factor verification, thus resolving the problem of systematically small scattering angles in the low-energy region encountered by traditional methods.

[0032] This invention embeds the detour factor as a global physical criterion into the parameter inversion process. Using the theoretical value of the detour factor in the basic radiation parameter table as the target, it employs the Nelder-Mead algorithm to automatically optimize parameters A and B, achieving a closed loop of "simulation-verification-correction." Compared to manual parameter tuning, the optimized model maintains computational efficiency while controlling the detour factor error to less than 1 × 10⁻⁶. -4 This ensures that the position of the Bragg peak and the dose drop after the peak are highly consistent with the experimental data. Attached Figure Description

[0033] Figure 1 A flowchart of a proton trajectory and dose simulation method based on stopping power and random scattering provided by the present invention;

[0034] Figure 2 This is a density diagram showing the relationship between the single-step scattering angle and the current energy provided in an embodiment of the present invention;

[0035] Figure 3 A flowchart of the optimization iterative empirical parameter algorithm provided in the embodiments of the present invention;

[0036] Figure 4 This is an example of a parameter inversion RMSE iteration diagram based on detour factor as a criterion provided in an embodiment of the present invention.

[0037] Figure 5 This is a comparison chart of the simulated detour factor and the theoretical value provided in the embodiments of the present invention;

[0038] Figure 6 This is a comparison diagram of the simulated dose-depth distribution and the literature distribution curves provided in this embodiment of the invention on the same coordinate system;

[0039] Figure 7 This is a three-dimensional image of the complete proton track provided in an embodiment of the present invention. Detailed Implementation

[0040] To make the above-mentioned objects, features, and advantages of the present invention more apparent and understandable, specific embodiments of the present invention are described in detail below. These specific embodiments are for illustrative purposes only and are not intended to limit the scope of the invention.

[0041] Please see Figure 1 This is a flowchart illustrating a proton trajectory and dose simulation method based on stopping power and random scattering provided by the present invention, comprising the following steps:

[0042] S1: Read the proton stopping power, continuous slowing approximate range and deflection factor data of the target material; establish its interpolation function with respect to proton energy; combine the density and radiation length of the target material to construct a basic radiation parameter table covering the incident energy range; and construct a uniform three-dimensional target geometric model.

[0043] S2: Using the basic radiation parameter table constructed in step S1, calculate an adaptive scattering step size suitable for the current energy based on the current energy of the proton, stopping power, and preset energy loss upper limit, providing input parameters for the scattering and energy loss calculations in subsequent steps S4 and S5;

[0044] S3: Using the circumference factor in the basic radiation parameter table as the criterion, calculate the root mean square error between the simulated circumference factor and the theoretical value, and use an optimization algorithm to automatically adjust the values ​​of the empirical parameters A and B in the scattering angle correction factor to minimize the root mean square error.

[0045] S4: Calculate the initial standard deviation of the scattering angle based on the current proton energy, the adaptive scattering step size generated in step S2, and the target radiation length. The initial standard deviation of the scattering angle is adjusted by introducing correction factors that include the empirical parameters A and B determined through optimization in step S3, to obtain the final standard deviation of the scattering angle for the current step. Using the final standard deviation of the scattering angle as the standard deviation of the Gaussian distribution, a single scattering angle conforming to the Gaussian distribution is sampled and generated.

[0046] S5: Update the proton's motion direction based on the scattering angle generated in step S4. Calculate the energy loss within the scattering step based on the stopping power, target density, and current step size in the basic radiation parameter table. Record the energy loss in the spatial grid traversed by the scattering path according to the grid size to complete the recording of energy deposition.

[0047] S6: Repeat steps S2, S4 and S5 in sequence to simulate the continuous scattering and energy loss process of a single proton in the target material until the proton energy is exhausted to zero or passes through the target material, thus completing the simulation of tracking the track of a single proton.

[0048] S7: Repeat the continuous scattering and energy loss process of more than 10,000 protons with the same initial energy in the target material, count the track and energy deposition data of all protons, and output the complete three-dimensional map of proton track, scattering angle-energy density map and dose-depth distribution curve.

[0049] This embodiment uses liquid water (H2O) as an example, calculating the incident energy as 160 MeV and 10 5 The continuous slowing and multiple scattering processes of a proton in water are described, and the Bragg peak and dose distribution are given.

[0050] Step S1 involves establishing the interpolation function, which typically uses publicly available databases for proton stopping power, continuous slowing approximate range, and deflection factor data. A cubic spline interpolation method is used on the original data, combined with the target density and radiation length, to construct an interpolation function covering the range of proton incident energy where each parameter varies with energy.

[0051] Specifically, this embodiment obtains the data file "pstar_water.txt" of proton incident liquid water from the NIST PSTAR public database, which includes the proton mass stopping power S(E) in the range of 0.01–200 MeV and the continuous slowing approximate range. (E) and detour factor The original data is then subjected to cubic spline interpolation to obtain callable functions get_stopping_power(E), get_csda_range(E), and get_detour_factor(E), which are then combined with the density of liquid water. and radiation length The program encapsulates the data into a unified interface "waterDB" for direct use in subsequent adaptive step size, scattering angle sampling, and detour factor verification. If the incident energy E0 is between two original data points, the program automatically uses the interpolation result, eliminating the need for additional experimental measurements and enabling simulation to start with "zero parameters".

[0052] Step S1 involves constructing a uniform three-dimensional geometric model with a thickness along the z-axis in a Cartesian coordinate system. The lateral dimension is larger than the maximum lateral diffusion range of protons, giving it a uniform density, and the slice thickness is set for subsequent energy deposition recording.

[0053] Specifically, this embodiment sets up The cubic liquid water target material imparts uniform density. and radiation length and set The slice thickness is used for subsequent energy deposition recording.

[0054] The adaptive scattering step size suitable for the current energy mentioned in step S2 includes: for each proton, when the current energy is E, according to... The theoretical step size is calculated and compared with the preset upper limit. The smaller value between the two is taken as the actual step size to ensure that the high-energy range moves with large and fast steps and the low-energy range moves with small and precise steps.

[0055] Specifically, in this embodiment, the proton incident energy is The preset energy loss ratio is 5%. To ensure the calculation accuracy in this embodiment, the upper limit of the step size is preferably set to 0.01 cm. The theoretical step size is calculated to be... ,Pick That is, 0.01cm is used as the actual step size for the first simulation step.

[0056] Please see Figure 2 This is a density diagram showing the relationship between the single-step scattering angle of protons in the liquid water target and the current energy in this embodiment.

[0057] The method for determining the empirical parameters A and B in step S3 includes: statistically analyzing the projected range of all protons after the simulation. With total path length Calculate the average detour factor Based on the theoretical value of the circumference factor in the basic radiation parameter table. To achieve the objective, construct the root mean square error. Where N is the number of sample data points involved in the calculation; the Nelder-Mead simplex algorithm is used to... Let be the variable, RMSE be the loss function, and the error be minimized iteratively. The condition is met when RMSE remains constant for 20 consecutive iterations or is less than 1 × 10⁻⁶. -4Stop at a certain time to achieve self-consistent calibration of empirical parameters A and B.

[0058] Please see Figure 3 This is a flowchart of the optimization iterative empirical parameter algorithm used in this embodiment. Specifically, the initial values ​​of empirical parameters A and B are both set to 1. The Nelder-Mead algorithm is used to iteratively call a short simulation of 500 protons. The simulation stops after an average of 20 iterations when the RMSE value remains unchanged for 20 consecutive iterations, resulting in the optimal values ​​of A=1.7013 and B=47.7987. This set of parameters is fixed and used in subsequent formal dose simulations to achieve the optimization of the algorithm. arrive The error between the full-energy-range scattering angle distribution and the theoretical value of the PSTAR bypass factor is <0.1%.

[0059] Please see Figure 4 , which is the RMSE after parameter inversion iteration in this embodiment using the detour factor as the criterion.

[0060] Please see Figure 5 This is a comparison curve of the detour factor simulated using the optimized parameter values ​​in this embodiment and the theoretical value.

[0061] Step S4, which involves sampling a single scattering angle based on the final standard deviation of the scattering angle, includes: first calculating the initial standard deviation of the scattering angle using the Highland-Lynch-Dahl formula. :

[0062]

[0063] Then introduce a low-energy empirical correction factor The final standard deviation of the scattering angle is obtained. ;by The variance is derived from a Gaussian distribution. Polar angle obtained by sampling azimuth exist Uniform sampling rotates the direction of proton motion to a new direction.

[0064] Specifically, in this embodiment, during the simulation of proton motion, the total energy of the proton includes kinetic energy and rest energy. Rest energy refers to the energy corresponding to the resting mass of the proton, with a value of 938.272 MeV. The relative velocity of a proton is calculated as the ratio of momentum to total energy. The momentum of a proton is expressed as the square root of its total energy minus its rest energy. The current step size, The target radiation length is given. The total energy of the first step of the simulation is the initial incident energy plus the static energy: 160 MeV + 938.272 MeV = 1098.272 MeV, where the proton momentum in the first step is... The relative velocity of the protons in the first step Step size of the first step , It is 36.08 Unit processing =36.08 =36.08cm, then the initial standard deviation of the scattering angle is :

[0065]

[0066] Final scattering angle standard deviation From Gaussian distribution Polar angle obtained by sampling azimuth exist Uniform sampling rotates the direction of proton motion to a new direction.

[0067] The method for calculating energy loss in step S5 includes: calling the basic radiation parameter table to obtain the mass stopping power S(E) corresponding to the current energy E, and combining it with the target density. and the adaptive step size determined in step S2 Through formula Calculate the approximate energy loss due to continuous slowing down within this step size. Then, update the current proton energy to... This provides the correct energy input for the next step of scattering and energy loss simulation. Energy loss It will be used for subsequent spatial grid energy deposition records.

[0068] Specifically, in this embodiment, the mass stopping power corresponding to an incident energy of 160 MeV is... The energy loss in the first step of the simulation is calculated as follows: Then, the current energy of the proton is updated to... .

[0069] Step S4, the spatially resolved energy deposition record, includes: calculating the geometric overlap of each proton path slice at a thickness of 1 mm, and calculating the overlap length. Linear amortized energy loss For the corresponding slice, dynamically expand the slice array and accumulate energy.

[0070] Specifically, in this embodiment After statistical analysis of 160MeV protons, the depth-relative dose distribution map was normalized to 100%. Depth-dose data were read from the "160MeV.csv" file, with depth values ​​as the x-axis and relative dose values ​​as the y-axis. The dose values ​​were also normalized to 100% of their maximum values. The depth-relative dose distribution map in this embodiment was plotted on the same coordinate system and compared with the Bragg peak position deviation of <0.2 cm. Obtaining the document 160MeV.csv: The relative dose and incident depth relationship curve of 160MeV protons in water simulated by the TRIM calculation module of the SRIM program was obtained from the publicly available document "SRIM simulation study on the distribution characteristics of Bragg peak in proton beam therapy" ([1] Tian Lixia, Fang Bingbing, Zhu Guanghao, et al. SRIM simulation study on the distribution characteristics of Bragg peak in proton beam therapy [J]. China Radiation Health, 2021, 30(01)). The relationship curve was digitized by the professional online data extraction tool WebPlotDigitizer, and the data points on the relative dose and incident depth relationship curve of 160MeV protons in water were extracted. The extracted data points were organized into a structured table file with the first column being depth and the second column being relative dose. The file was named "document 160MeV.csv".

[0071] Please see Figure 6 This embodiment is based on the literature. The dose-depth distribution curve was calculated using 160 MeV protons, and the slice energy deposition was normalized to 100% of the maximum dose. The curves were then compared with those in the literature on the same coordinate system.

[0072] The multi-dimensional output results described in step S6 include: dose-depth distribution curve: the slice energy deposition is normalized to 100% of the maximum dose and compared with the literature curve on the same coordinate system; scattering angle-energy density map: a log-log two-dimensional histogram showing the relationship between the single-step scattering angle and the current energy, with color mapping of event density; detour factor-energy curve: outputting a comparison chart of the optimized simulated detour factor and the PSTAR theoretical value to verify the reliability of the model; complete track 3D map: recording the 3D coordinate string of each proton trajectory, supporting Matplotlib 3D visualization, for teaching and model diagnosis.

[0073] Please see Figure 7 To illustrate the three-dimensional track pattern of protons on the liquid water target in this embodiment, 100 proton tracks were drawn to avoid overcrowding.

Claims

1. A method for simulating proton tracks and dose based on stopping power and random scattering, characterized in that, Includes the following steps: S1: Read the proton stopping power, continuous slowing approximate range, and deflection factor data of the target material; Establish its interpolation function with respect to proton energy, combine the density and radiation length of the target material, construct a basic radiation parameter table covering the incident energy range, and construct a uniform three-dimensional target geometric model; S2: Using the basic radiation parameter table constructed in step S1, an adaptive scattering step size suitable for the current energy is generated based on the current energy E of the proton, the stopping power S(E), and the preset energy loss upper limit. S3: Using the scattering factor in the basic radiation parameter table as a criterion, calculate the root mean square error between the simulated scattering factor and the theoretical value, and adjust the scattering angle correction factor using an optimization algorithm. The values ​​of empirical parameters A and B in the equation are chosen to minimize the root mean square error. S4: Calculate the initial standard deviation of the scattering angle based on the current proton energy, the adaptive scattering step size generated in step S2, and the target radiation length. The initial standard deviation of the scattering angle is adjusted by introducing a scattering angle correction factor that includes the empirical parameters A and B determined through optimization in step S3, thus obtaining the final standard deviation of the scattering angle for the current scattering step size. Using the final standard deviation of the scattering angle as the standard deviation of the Gaussian distribution, a single scattering angle conforming to the Gaussian distribution is sampled and generated. S5: Update the proton's motion direction based on the scattering angle generated in step S4. Calculate the energy loss within the scattering step based on the stopping power, target density, and current scattering step size in the basic radiation parameter table. Record the energy loss in the spatial grid traversed by the scattering path according to the grid size to complete the recording of energy deposition. S6: Repeat steps S2, S4 and S5 in sequence to simulate the continuous scattering and energy loss process of a single proton in the target material until the proton energy is exhausted to zero or passes through the target material, thus completing the simulation of tracking the track of a single proton. S7: Repeat the continuous scattering and energy loss process of more than 10,000 protons with the same initial energy in the target material, count the track and energy deposition data of all protons, and output the three-dimensional map of the complete proton track, the scattering angle-energy density map and the dose-depth distribution curve; The method for generating the adaptive scattering step size in S2 is as follows: set an upper limit on the proportion of single-step energy loss to the current energy. Based on the stopping power S(E) and target density And the current energy E, through the formula Calculate the theoretical step size and compare it with the preset upper limit of the step size. Take the smaller of the two values ​​as the adaptive scattering step size.

2. The method for simulating proton tracks and doses based on stopping power and random scattering according to claim 1, characterized in that, The steps for determining parameters A and B in S3 through iterative optimization are as follows: Step 1: Set the correction factor Initial values ​​for parameters A and B; Step 2: Based on the current parameters A and B, perform a track simulation experiment on the proton set used for parameter optimization: calculate the scattering angle and update the proton direction using correction factors that include the current parameters A and B; calculate the energy loss and iterate the simulation until the proton energy is depleted to zero; Step 3: After the simulation test is completed, the projected range and total path length of all experimental protons are statistically analyzed, and the average detour factor simulation value is calculated. Step 4: Compare the calculated average simulated value of the detour factor with the theoretical value of the detour factor in the basic radiation parameter table. If the error does not meet the preset convergence condition, the values ​​of parameters A and B are automatically adjusted using the optimization algorithm, and steps 2 to 4 are repeated for iterative optimization. Step 5: When the error meets the convergence condition, output the parameters A and B determined at this time as the final optimized empirical parameters for the target material.

3. The method for simulating proton tracks and doses based on stopping power and random scattering according to claim 2, characterized in that, S3 is the offline parameter optimization stage before simulation. It is executed once after steps S1 and S2 are completed, and the empirical parameters A and B are obtained and then not changed. The loop that is repeatedly executed in step S6 is S2, S4 to S5, using the fixed parameters A and B optimized in S3.

4. The method for simulating proton tracks and doses based on stopping power and random scattering according to claim 1, characterized in that, The initial standard deviation of the scattering angle in S4 Calculated using the Highland-Lynch-Dahl formula: ， in The relative velocity of a proton is calculated as the ratio of momentum to total energy. The momentum of a proton is expressed as the square root of its total energy minus its rest energy. The current step size, The radiation length of the target material; The final standard deviation of the scattering angle in S4: .

5. The method for simulating proton tracks and doses based on stopping power and random scattering according to claim 1, characterized in that, The method described in S1 for constructing an energy-dependent basic radiation parameter table covering the incident energy range includes: using cubic spline interpolation to construct interpolation functions for stopping power, continuously slowing approximate range, and circumduction factor with respect to proton energy; and constructing the interpolation functions, target density, and radiation length into an energy-dependent basic radiation parameter table for use in subsequent simulation steps.

6. The method for simulating proton tracks and doses based on stopping power and random scattering according to claim 1, characterized in that, The energy loss calculation method described in S5 includes: based on the stopping power S(E) and target density in the basic radiation parameter table. and current step size Calculate the energy loss within this scattering step. And update the current energy of the proton to .

7. The method for simulating proton tracks and doses based on stopping power and random scattering according to claim 1, characterized in that, The method described in S5 for recording energy loss in the spatial grid traversed by the scattering path according to the grid size includes: dividing the spatial path of each proton step into several depth slices according to a pre-set thickness, calculating the spatial overlap length between the path and each slice, and linearly allocating the energy loss of that step to the corresponding slice according to the overlap length ratio, thereby achieving high-resolution energy deposition recording in the depth direction.

Images