Parameter limiting Monte Carlo method for determining solar proton event energy spectrum

Abstract

The invention provides a parameter-limited Monte Carlo method for determining a solar proton event energy spectrum, which comprises the following steps of: S1, acquiring a proton time-varying flux of each energy channel of a space proton detector, and judging the beginning and the end of a solar proton event; s2, calculating the total flux of protons in each energy channel after the background in the solar activity calm period is removed during the occurrence period of the solar proton event; s3, performing logarithmic transformation on parameters of the proton double-power-law energy spectrum function according to the total flux obtained in the S2 to determine a sampling range, sampling the parameters in the determined sampling range by using a uniform and random Monte Carlo method, and fitting the energy spectrum by using an integral method; and S4, obtaining a group of optimal parameters of the proton continuous energy spectrum, and drawing an energy spectrum curve in combination with an energy spectrum function. According to the Monte Carlo method provided by the invention, the calculation efficiency is improved in the fitting process, and the fitting precision is improved. By using the integral method, the calibration requirement on the satellite observation energy channel and energy section is avoided, and the problem of error accumulation caused by calibration is solved.

Description

Technical Field

[0001] The invention relates to the field of solar-terrestrial space environment monitoring and data analysis, and in particular to a parameter-limited Monte Carlo method for determining the energy spectrum of solar proton events. Background Art

[0002] Solar proton events are one of the main manifestations of solar activity. They are essentially the release of energy stored in the magnetic field of the solar atmosphere. During the release process, a large number of high-energy particles are ejected. The intensity of high-energy particles observed in the geosynchronous orbit increases by several orders of magnitude. Since most of the high-energy particles ejected from the sun are protons, accounting for more than 90% of the total number of high-energy particles, they are called solar proton events. Solar proton events are a violent solar eruption activity phenomenon that causes a sharp increase in proton flux in a short period of time, thus causing serious radiation hazards to satellites and astronauts in orbit.

[0003] Existing energy spectrum fitting methods usually rely on calibrating the effective energy points of the satellite detection energy channel. The selection of these effective energy points is relatively fixed, but in actual applications, as the spectrum shape changes, the real effective energy points often deviate. In order to overcome this limitation, the present invention proposes an energy spectrum fitting scheme based on the integral method. This method can make the fitted energy spectrum adaptively adjusted to adapt to the actual situation of the detection data, thereby eliminating the need for strict calibration of the energy bands in the satellite observation energy channel. Through this adaptive processing, the fitting process can better reflect the changes in the spectrum shape in the actual observation data, can effectively avoid the error accumulation problem introduced by the calibration process, and significantly improve the accuracy of the fitting. Summary of the invention

[0004] The object of the present invention is to provide a parameter-constrained Monte Carlo method for determining the energy spectrum of solar proton events to solve the problems existing in the above background. The present invention adopts the following technical solutions: A parameter-constrained Monte Carlo method for determining the energy spectrum of a solar proton event, characterized in that it comprises the following steps: S1: Obtain the time-varying proton flux of each energy channel of the space proton detector to determine the beginning and end of the solar proton event; S2: Calculate the total flux of protons in each energy channel during a solar proton event after removing the background of the quiet period of solar activity; S3: Logarithmically transform the parameters of the proton double power-law energy spectrum function according to the total flux obtained in S2 to determine the sampling range, use the uniform random Monte Carlo method to sample the parameters in the determined sampling range, and use the integral method to fit the energy spectrum; S4: Obtain a set of optimal parameters for the proton continuous energy spectrum and draw the energy spectrum curve in combination with the energy spectrum function; Preferably, in step S1, a solar proton event is defined as: an observed energy continuously greater than 10 The total flux of protons is not less than 10 The proton flux enhancement process lasts for no less than 15 minutes. The start of a solar proton event is defined as a proton flux greater than or equal to 10 The first of three consecutive data points, the end of the event is when the flux is greater than or equal to 10 for the last time .

[0005] Preferably, in step S2, within 30 days before the solar proton event starting point determined in step S1, the proton flux of each observed energy channel is sliding averaged in units of 6 hours, and the minimum value of the 6-hour sliding average is selected as the background proton flux of each energy channel. The function form for sliding average is as follows: ,in, is the current time point, is the time point corresponding to the observation value, It is at the moment The 6-hour moving average of It is a window The number of data points in It is in the window The data point values ​​within represents the time window from hour t-6 to time t.

[0006] Preferably, the background proton flux is removed from the time-varying flux of the detector to obtain a proton flux increment, and time integration is performed to obtain the total proton flux of each energy channel during the solar proton event.

[0007] Preferably, in the step S3, the double power-law energy spectrum function of the solar proton event is fitted using four free parameters to describe the behavior characteristics of the solar proton spectrum line in a wide energy range, and the energy spectrum function is in the form of: ,in, is the total flux of protons, is the amplitude of the continuous energy spectrum, is the proton energy, It is the energy turning point between the low energy segment and the high energy segment. is the spectral index of the low energy band, is the spectral index of the high energy band.

[0008] Preferably, based on the total proton flux obtained in step S2, the parameters of the double power-law energy spectrum function are logarithmically transformed to determine the parameters of the double power-law energy spectrum function: , and sampling range.

[0009] Preferably, in step S3, the fitting effect is evaluated using the integration method and the logarithmic root mean square error formula: .

[0010] The detection flux is obtained by multiplying the channel width by the flux value of the channel, and the fitting flux is calculated by integrating the energy within the channel range, corresponding to represents the logarithmic root mean square error, represents the number of GOES satellite observation channels used for fitting, Indicates the actual detection results. represents the fitting result, is the detector channel width, and Corresponding to The upper and lower limits of is the flux increment of the corresponding energy channel.

[0011] Preferably, the fitting accuracy is controlled by setting the random number seed and the sampling times, and a Monte Carlo method of multi-parameter uniform random sampling is used to uniformly and randomly extract a set of energy spectrum parameter values ​​in the parameter sampling range given in step S3, and calculate the corresponding .

[0012] Preferably, in the step S4, after completing all sampling times, select A set of energy spectrum parameters corresponding to the minimum value, and the energy spectrum function is drawn in combination with the energy spectrum function.

[0013] In summary, due to the adoption of the above technical solution, the beneficial effects of the present invention are: 1. The present invention effectively achieves the accurate fitting of the continuous energy spectrum of solar proton events by introducing the Monte Carlo method. Specifically, the time-varying flux data of protons in each energy channel is obtained through a space proton detector. The original data is firstly subjected to background removal and preprocessing to ensure the accuracy and reliability of the input data. Subsequently, the continuous energy spectrum objective function of the solar proton event is fitted in combination with the Monte Carlo method. Compared with traditional methods, the present invention significantly improves the computational efficiency and the accuracy of the fitting results during the fitting process, providing a more efficient technical path for the study of the energy spectrum of solar proton events.

[0014] 2. In view of the limitation that the existing energy spectrum fitting method relies on the calibration of the effective energy points of the satellite detector energy channel, the present invention proposes a new energy spectrum fitting method based on the integration method, which can realize adaptive energy spectrum fitting. In the traditional method, the selection of effective energy points is relatively fixed, but in actual measurement, these energy points will be offset due to the change of the energy spectrum shape, resulting in the error accumulation problem in the calibration result. The present invention adopts the integration method to avoid the calibration step, directly process the detection data in the whole range, so that the fitting results can reflect the changing characteristics of the energy spectrum shape in real time, thereby accurately obtaining a more realistic energy spectrum distribution and avoiding the potential system errors caused by the energy point calibration requirements.

[0015] 3. In order to get closer to the real energy spectrum distribution of solar proton events, the present invention selects the Monte Carlo method for energy spectrum fitting, which further optimizes the calculation efficiency and accuracy while ensuring the retention of real energy spectrum information. Compared with traditional methods, the Monte Carlo fitting method shows stronger flexibility and adaptability, and can be flexibly adjusted according to the changes in actual observation data. It not only effectively improves the fitting accuracy, but also can be applied to the fitting requirements of other objective functions, providing technical support for a wider range of scientific research. BRIEF DESCRIPTION OF THE DRAWINGS

[0016] Figure 1 It is a flow chart of the present invention; Figure 2 It is a time-varying flux diagram of protons in each energy channel of the space proton detector of the present invention; Figure 3 This is a schematic diagram of the time span of solar proton events of the present invention; Figure 4 The solar proton energy spectrum determined by the present invention. DETAILED DESCRIPTION

[0017] In order to make the above-mentioned objects, features and advantages of the present invention more obvious and easy to understand, specific embodiments of the present invention are described in detail below.

[0018] The embodiment of the present invention provides a Monte Carlo method for determining the proton continuum energy spectrum of a solar proton event, such as Figure 1 As shown, the following steps are included: S1: Obtain the time-varying proton flux of each energy channel of the space proton detector to determine the beginning and end of the solar proton event S2: Calculate the total flux of protons in each energy channel during a solar proton event after removing the background of the quiet period of solar activity; S3: Logarithmically transform the parameters of the proton double power-law energy spectrum function based on the total flux obtained in S2 to determine the sampling range, use the uniform random Monte Carlo method to sample the parameters in the determined sampling range, and use the integral method to fit the energy spectrum; S4: Obtain a set of optimal parameters for the proton continuous energy spectrum and draw the energy spectrum curve in combination with the energy spectrum function; A Monte Carlo method for determining the proton continuum energy spectrum of a solar proton event is introduced below in conjunction with a specific embodiment: GOES15 satellite data is used to determine the continuum energy spectrum of the solar proton event on September 10, 2017 as an example; In step S1, the time-varying proton flux of each energy channel of the space proton detector is obtained to determine the beginning and end of the solar proton event.

[0019] Specifically, the solar proton event list is provided by the National Oceanic and Atmospheric Administration (NOAA) of the United States. The time-varying flux data of 9 different energy channels of the space proton detector corresponding to the date of the solar proton event, corresponding to P1-P10 (excluding P7), with a resolution of every 5 minutes are obtained from the GOES15 satellite data website, such as Figure 2 shown.

[0020] Specifically, the energy of the solar proton event corresponding to the date is greater than 10 The time-varying flux of protons every 5 minutes, with energies greater than 10 ( Figure 3 The proton flux is defined as the time span of the event. The start of a solar proton event is defined as the time span when the proton flux is greater than or equal to 10 The first of three consecutive data points, the end of the event is when the flux is greater than or equal to 10 for the last time , and the duration of the solar proton event is determined accordingly, such as Figure 3 As shown, in this embodiment, the solar proton event starts at 16:45 on September 10, 2017, and ends at 17:20 on September 14, 2017.

[0021] In step S2, the total flux of protons in each energy channel during a solar proton event is calculated after removing the background of the quiet period of solar activity.

[0022] Specifically, within 30 days before the solar proton event starting point determined in step S2, the proton flux of each observed energy channel is averaged in units of 6 hours, and the minimum value of the 6-hour sliding average is selected as the background proton flux of each energy channel. The function form for sliding average is: ,in, is the current time point, is the time point corresponding to the observation value, is the 6-hour sliding mean at time, It is a window The number of data points in It is in the window The data point values ​​within Indicates from hour Before to time time window.

[0023] Specifically, in this embodiment, the 30-day period is from August 10, 2017 to September 9, 2017, and the proton flux of the P1-P10 (excluding P7) energy channels is calculated by a 6-hour sliding average, and the minimum value of all calculation results is selected to obtain the background fluxes, respectively, P1: 2.547089, P2: 0.081328, P3: 0.043777, P4: 0.017357, P5: 0.014189, P6: 0.002710, P8: 0.002230, P9: 0.000678, P10: 0.000449, all in units of , the background proton flux is removed from the time-varying flux of the detector to obtain the proton flux increment, and the increment of the total proton flux of each energy channel during the solar proton event is obtained by integrating according to the duration of the solar proton event. The specific increments are P1: 374061660, P2: 10483843, P3: 4402597, P4: 2245267, P5: 437784, P6: 33170, P8: 1258, P9: 467, P10: 161, all in units of .

[0024] In step S3, the parameters of the proton double power-law energy spectrum function are logarithmically transformed to determine the sampling range.

[0025] Specifically, the double power-law energy spectrum function is in the form of: , four free parameters are used to describe the behavior of solar proton lines over a wide energy range, among which, represents the total flux of protons, represents the proton energy, represents the amplitude of the continuous energy spectrum, Indicates the energy turning point between the low energy segment and the high energy segment. represents the spectral index of the low energy band, represents the spectral index of the high energy band, and this function is used for the fitting objective function in the embodiment.

[0026] Specifically, based on the total proton flux obtained in step S3, the parameters of the energy spectrum double power law function are determined: , and The sampling range, ,in is the maximum flux observed in all channels, corresponding to the P1 channel, ,in is the increment of each energy channel (corresponding to P1-P10, excluding P7), and the energy channel widths corresponding to P1-P10 (excluding P7) are (unit ( )): P1: 0.74–4.2, P2: 4.2–8.7, P3: 8.7–14.5, P4: 15–40, P5: 38–82, P6: 84–200, P7: 110–900, P8: 330–420, P9: 420–510, P10: 510–700. Since the energy channel range of P7 overlaps with that of P8, P9, and P10, this energy channel is discarded in the present invention, and the energy midpoint of the energy channel is determined according to the energy channel range. , thus we get the corresponding ,according to Get , for The maximum value in is determined according to the fitted energy range , the above logarithmic separation is used to reduce the parameter sampling space, thereby accelerating the fitting speed and improving the fitting accuracy. In this embodiment, the total flux of P1 is 3.47833069e8 Sure The parameter space of is (8,9), and , , At the same time, from the double power law energy spectrum function form, we can know that in the actual fitting process, we must always ensure .

[0027] In step S3, the Monte Carlo method with uniform random sampling is used to obtain the optimal parameters of the proton continuum energy spectrum.

[0028] Specifically, the logarithmic root mean square error formula is used to evaluate the fitting effect during the fitting process using the integral method: ,in, represents the logarithmic root mean square error, Indicates the number of GOES satellite observation channels used for fitting. In this embodiment, the number of channels is , Indicates the actual detection results. represents the fitting result, is the detector energy channel width. The fitting curve is integrated corresponding to the observed energy band to evaluate the error with the actual value. and Corresponding to The upper and lower limits of and can be obtained from the previous step. is the flux increment of the corresponding energy channel, the increments are P1:3.47833069e8, P2:8.82130877e6, P3:2.94803329e6, P4:2.17665950e6, P5:3.76976159e5, P6:3.22259432e4, P8:1.26484759e3, P9:4.32669353e2, P10:1.56672084e2, the units are .

[0029] In step S4, an energy spectrum curve is drawn using a set of optimal parameters for obtaining a proton continuous energy spectrum and combining the energy spectrum function.

[0030] Specifically, the seed of random number generation is set to the current time to ensure the randomness of the Monte Carlo fitting. According to the preprocessing, the sampling times are set to 100,000 times to ensure the reliability of the sampling. The Monte Carlo method of multi-parameter uniform random sampling is used to uniformly and randomly extract a set of energy spectrum parameter values ​​in the parameter sampling space given in step S3, and calculate the corresponding , repeat sampling and record the , until the set number of samplings is completed, and the smallest of all samples is output Corresponding energy spectrum parameters, such as Figure 4 As shown, in this embodiment, the minimum The value is 0.181643, and the corresponding parameter values ​​are =8.75345, =1.777613, =3.723261, =2.141541.

Claims

1. A parameter-constrained Monte Carlo method for determining the energy spectrum of solar proton events, characterized in that: The following steps are involved: S1: Obtain the time-varying proton flux of each energy channel of the space proton detector to determine the beginning and end of the solar proton event; S2: Calculate the total flux of protons in each energy channel during a solar proton event after removing the background of the quiet period of solar activity; S3: Logarithmically transform the parameters of the proton double power-law energy spectrum function according to the total flux obtained in S2 to determine the sampling range, use the uniform random Monte Carlo method to sample the parameters in the determined sampling range, and use the integral method to fit the energy spectrum; S4: Obtain a set of optimal parameters for the proton continuous energy spectrum and draw the energy spectrum curve in combination with the energy spectrum function; The method of using the Monte Carlo method to sample parameters in S3 is to control the fitting accuracy by setting the random number seed and the number of samplings, and use the Monte Carlo method of multi-parameter uniform random sampling to uniformly and randomly extract a set of energy spectrum parameter values ​​in a given parameter sampling range, and use the integral method and the logarithmic root mean square error formula to evaluate the fitting effect. The error formula used to evaluate the fitting effect is as follows: ,in, represents the logarithmic root mean square error, represents the number of GOES satellite observation channels used for fitting, Indicates the actual detection results. Indicates the fitting result. The fitting curve is integrated corresponding to the observed energy segment to evaluate the error with the actual value. is the detector channel width, and The corresponding upper and lower limits are, is the flux increment of the corresponding energy channel.

2. A parameter-constrained Monte Carlo method for determining the energy spectrum of solar proton events according to claim 1, characterized in that: The method described in S3 for evaluating the fitting effect using the integral method and the logarithmic root mean square error formula is: Compare the detection results using the logarithmic root mean square And the fitting results , to quantify the difference between them, the calculation steps include: first take the detection result And the fitting results The logarithmic ratio of these ratios is then averaged and squared. The actual detection result is obtained by measuring the detector channel width. The flux increment of the corresponding energy channel The fitting result is obtained by multiplying the flux within the energy channel.

3. A parameter-constrained Monte Carlo method for determining the energy spectrum of solar proton events according to claim 1, characterized in that: The method described in S1 for determining the start and end of a solar proton event is to obtain the time-varying flux of protons in each energy channel of the GOES satellite proton detector, and use the time-varying flux of protons with energy greater than 10 MeV as the basis for determining the start and end time of the event.

4. A parameter-constrained Monte Carlo method for determining the energy spectrum of solar proton events according to claim 1, characterized in that: The method for calculating the total proton flux of each energy channel after removing the background of the quiet period of solar activity during the solar proton event described in S2 is: in the 30-day period before the solar proton event, the proton flux of each observed energy channel is averaged in units of 6 hours, and the minimum value of the 6-hour sliding average is selected as the background proton flux of each energy channel. The function of the sliding average is as follows: ,in, is the current time point, is the time point corresponding to the observation value, is the 6-hour moving average at time t, It is a window The number of data points in It is in the window The data point values ​​within Represents the time window from t-6 hours ago to time t.

5. A parameter-constrained Monte Carlo method for determining the energy spectrum of solar proton events according to claim 1, characterized in that: The method for determining the sampling range described in S3 is: using a double power-law energy spectrum function to fit solar proton events, using four free parameters to describe the behavior characteristics of solar proton spectral lines in a wide energy range, and performing logarithmic transformation on the parameters of the double power-law energy spectrum function to determine the sampling range of each parameter. The double power-law energy spectrum function is in the following form: ,in, is the total flux of protons, is the amplitude of the continuous energy spectrum, is the proton energy, It is the energy turning point between the low energy segment and the high energy segment. is the spectral index of the low energy band, is the spectral index of the high energy band.

6. A parameter-constrained Monte Carlo method for determining the energy spectrum of solar proton events according to claim 1, characterized in that: The optimal parameters described in S4 are selected after completing all sampling times. A set of energy spectrum parameters corresponding to the minimum value.

Images