Method, System, Storage Medium and Device for Determining Volcanic Ionospheric Anomalous Disturbance

Abstract

An embodiment of the present invention discloses a method for determining volcanic ionospheric anomalous disturbances. The method includes: obtaining ground-based GNSS raw observation data; preprocessing the ground-based GNSS raw observation data to obtain an ionospheric TEC data set; removing ionospheric TEC data with a high correlation with the space weather environment according to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC data set in the time domain and the frequency domain, and obtaining an optimal ionospheric TEC data set that is not affected by the space weather environment; improving the LSTM algorithm through the Prophet algorithm, and determining the upper and lower boundaries of volcanic ionospheric TEC anomalous disturbances based on the optimal ionospheric TEC data set; determining whether an anomalous disturbance occurs in the volcanic ionosphere according to the comparison relationship between the ionospheric TEC observation data and the upper and lower boundaries of the volcanic ionospheric TEC anomalous disturbances. The present invention realizes anomalous disturbance identification by data processing and removing space weather interference, and using the Prophet algorithm to improve the LSTM algorithm, and can efficiently monitor volcanic eruptions.

Description

Technical Field

[0001] The present invention relates to the technical field of ionospheric anomaly disturbance detection, and particularly to a method, system, storage medium and device for determining volcanic ionospheric anomaly disturbance. Background Art

[0002] Volcanic eruption is a geological process in which magma, volcanic gases and debris inside the earth are ejected out through openings in the earth's crust. This process not only shapes new geological structures, such as spectacular volcanic cones and vast lava plains, but also has a profound impact on the environment and human society. Volcanic activities can trigger climate change, damage ecosystems, and directly threaten human safety. Therefore, monitoring volcanic activities and accurately predicting volcanic eruptions have become the key to reducing their potential hazards.

[0003] Volcanic eruptions are not limited to surface activities. Their intense energy release can also disturb the earth's space atmosphere environment. Specifically, the impact of volcanic eruptions can trigger mesospheric gravity waves (mesospheric airglow waves), and these waves further trigger co-volcanic ionospheric disturbances CVID, which usually occur within 10 to 45 minutes after a volcanic eruption and appear in a quasi-periodic form. These ionospheric disturbances have a significant impact on the total electron content (TEC) in the ionosphere, providing an important basis for the monitoring and prediction of volcanic activities.

[0004] The TEC in the ionosphere is a key parameter reflecting the state of the ionosphere, and its change is closely related to volcanic activities. The ionospheric disturbances generated during and after a volcanic eruption can, to a certain extent, reflect the characteristics and intensity of the volcanic eruption. More importantly, the abnormal disturbances that appear in the ionosphere before a volcanic eruption often indicate the approach of volcanic activities and are of great significance for the early warning of volcanic eruptions. Summary of the Invention

[0005] Based on this, it is necessary to propose a method for determining volcanic ionospheric anomaly disturbance in view of the above problems.

[0006] A method for determining volcanic ionospheric anomaly disturbance, the method comprising the following steps:

[0007] Obtain ground-based GNSS raw observation data;

[0008] Preprocess the ground-based GNSS raw observation data to obtain an ionospheric TEC data set;

[0009] According to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC data set in the time domain and the frequency domain, remove the ionospheric TEC data with high correlation with the space weather environment to obtain an optimal ionospheric TEC data set that is not affected by the space weather environment;

[0010] Improve the LSTM algorithm through the Prophet algorithm, and determine the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC based on the optimal ionospheric TEC dataset;

[0011] Determine whether there is an abnormal disturbance in the volcanic ionosphere according to the comparison relationship between the ionospheric TEC observation data and the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC.

[0012] In the above solution, the preprocessing of the ground-based GNSS original observation data to obtain the ionospheric TEC dataset specifically includes:

[0013] Based on the GNSS original pseudorange and carrier phase observation equations, determine the slant ionospheric delay Iono, where Ion0 = STEC + DCB, STEC is the slant TEC, and VTEC is the vertical TEC;

[0014] Determine the differential code bias DCB according to the non-differenced non-combined precise point positioning model and the generalized trigonometric series model;

[0015] Determine the slant TEC and vertical TEC of the satellite piercing point of the GNSS station according to the slant ionospheric delay Iono and the differential code bias DCB, and obtain the ionospheric TEC dataset.

[0016] In the above solution, the determination of the slant ionospheric delay Iono based on the GNSS original pseudorange and carrier phase observation equations specifically includes:

[0017] The GNSS original pseudorange and carrier phase observation equations are:

[0018]

[0019]

[0020] Among them, the superscript s and the subscript r respectively represent the satellite and the receiver, and the subscript f represents the L1 or L2 frequency band of the carrier signal; represents the pseudorange observation value, represents the carrier phase observation value; represents the geometric distance between the satellite and the receiver; c represents the speed of light; dt r and dt s respectively represent the clock biases of the receiver and the satellite; Trop represents the tropospheric slant delay, represents the ionospheric slant delay; b r,f and respectively represent the pseudorange hardware delays at the receiver end and the satellite end; λ f represents the carrier phase wavelength, represents the carrier phase integer ambiguity; B r,fand respectively represent the phase hardware delays at the receiver end and the satellite end; ε P and ε L respectively represent the sum of multipath effects, observation noise, and other unmodeled errors in pseudorange and carrier observations;

[0021] Solve the GNSS raw pseudorange and carrier phase observation equations according to the least squares method to determine the slant ionospheric delay Iono.

[0022] In the above solution, the determination of the differential code bias DCB according to the non-differential non-combination precise point positioning model and the generalized trigonometric series model specifically includes:

[0023] Introduce the precise ephemeris and clock offset products released by IGS and use the ionosphere-free combination to estimate the satellite clock offset:

[0024]

[0025] Convert STEC to VTEC, and the conversion formula is as follows:

[0026] VTEC = MF·STEC

[0027]

[0028] where, R represents the radius of the earth; H represents the height of the ionospheric piercing point, α represents a coefficient; z represents the satellite elevation angle;

[0029] Estimate the DCB of the GNSS station using the generalized trigonometric series, and the generalized trigonometric series function model is as follows:

[0030]

[0031] In the formula, and λ respectively represent the geographical latitudes of the ionospheric piercing point and the modeling center point; h represents the function related to the local time at the ionospheric piercing point; nmax, mmax, and kmax respectively represent the maximum orders of the polynomial function and the trigonometric series function; Enm, Ck, and Sk respectively represent the model parameters to be estimated of the generalized trigonometric series function.

[0032] In the above solution, the removal of the ionospheric TEC data with high correlation with the space weather environment according to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC data set in the time domain and the frequency domain to obtain the optimal ionospheric TEC data set not affected by the space weather environment specifically includes:

[0033] Obtain the daughter wavelet family by simply stretching and translating the mother wavelet φ(t):

[0034]

[0035] Among them, a represents the translation parameter that controls the position of the wavelet, b represents the scale parameter that controls the width scaling or dilation of the wavelet, and t represents time;

[0036] By giving the ionospheric TEC time series or the solar activity parameter time series, the continuous wavelet transform formula thereof is obtained:

[0037]

[0038] In the formula, φ * (t) represents the complex conjugate function of φ(t), and x represents the ionospheric TEC time series or the solar activity parameter time series;

[0039] Determine the variance of the wavelet energy spectrum measure time series of the ionospheric TEC data set in the time and frequency domains:

[0040] |(WPS) x (a, b)| = |W x (a, b)| 2

[0041]

[0042] Among them, W x (a, b) and W y (a, b) respectively represent the wavelet transforms of x(t) and y(t), represents the complex conjugate of W y (a, b).

[0043] In the above solution, before improving the LSTM algorithm through the Prophet algorithm and determining the upper and lower boundaries of the volcanic ionospheric TEC anomaly perturbation based on the optimal ionospheric TEC data set, the method further includes: obtaining the TEC prediction values of the Prophet algorithm and the LSTM algorithm:

[0044] y(t) = g(t) + s(t) + h(t) + ε t

[0045] In the formula, g(t) represents the trend term; s(t) represents the periodic term; h(t) represents holidays, indicating the influence of holidays that may occur irregularly within one day or multiple days; ε t represents the error term; y(t) represents the TEC prediction value of the Prophet algorithm;

[0046] Construct the relevant parameters of the LSTM neural network:

[0047]

[0048]

[0049]

[0050]

[0051]

[0052] In the formula, C t-1 represents the cell state at time t-1, h t-1 represents the output at the previous moment, x t represents the current input, h t represents the current output, that is, the TEC prediction value of the LSTM algorithm, C t represents the cell state at the current time t.

[0053] In the above solution, improving the LSTM algorithm through the Prophet algorithm and determining the upper and lower boundaries of the volcanic ionospheric TEC anomaly perturbation based on the optimal ionospheric TEC data set specifically include:

[0054] Respectively assign the prediction value of the Prophet algorithm and the prediction value of the LSTM model the initial weights ω1 and ω2 of the ionospheric TEC at time i;

[0055] Through error accuracy analysis with the measured ionospheric TEC, determine the objective function:

[0056]

[0057] Among them, N represents the total number of detection periods; yi represents the measured ionospheric TEC; represents the TEC prediction value of the Prophet algorithm at time i, represents the TEC prediction value of the LSTM algorithm at time i;

[0058] Preset constraint conditions: ω1 + ω2 = 1, ω1 ≥ 0, ω2 ≥ 0;

[0059] Use the least squares method to determine the optimal weights ω1 and ω2, and determine the upper and lower boundaries of the volcanic ionospheric TEC anomaly perturbation.

[0060] In the above solution, the use of the least squares method to determine the optimal weights ω1 and ω2, and determine the upper and lower boundaries of the volcanic ionospheric TEC anomaly perturbation specifically includes:

[0061] Error = TEC actual -TEC predicted

[0062] UB = Mean(Error) + 2×Std(Error)

[0063] LB = Mean(Error) - 2×Std(Error)

[0064] Wherein, TEC actual represents the measured ionospheric TEC; TEC predicted represents the ionospheric TEC predicted by the new algorithm; Mean represents taking the average value; Std represents calculating the standard deviation; UB represents the upper boundary of the ionospheric TEC anomaly; LB represents the lower boundary of the ionospheric TEC anomaly.

[0065] This application also proposes a volcanic ionospheric anomaly disturbance detection system, which includes: a data acquisition unit, a judgment condition presetting unit, and a determination unit;

[0066] The data acquisition unit is used to acquire the raw ground-based GNSS observation data; preprocess the raw ground-based GNSS observation data to obtain an ionospheric TEC data set; remove the ionospheric TEC data with high correlation with the space weather environment according to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC data set in the time domain and the frequency domain, and obtain an optimal ionospheric TEC data set that is not affected by the space weather environment;

[0067] The judgment condition presetting unit is used to improve the LSTM algorithm through the Prophet algorithm, and determine the upper boundary and the lower boundary of the volcanic ionospheric TEC anomaly disturbance based on the optimal ionospheric TEC data set;

[0068] The determination unit is used to determine whether there is an abnormal disturbance in the volcanic ionosphere according to the comparison relationship between the ionospheric TEC data set and the upper boundary and the lower boundary of the volcanic ionospheric TEC anomaly disturbance.

[0069] This application also proposes a readable storage medium, storing a computer program, which when executed by a processor, causes the processor to perform the following steps:

[0070] Acquire the raw ground-based GNSS observation data;

[0071] Preprocess the raw ground-based GNSS observation data to obtain an ionospheric TEC data set;

[0072] Remove the ionospheric TEC data with high correlation with the space weather environment according to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC data set in the time domain and the frequency domain, and obtain an optimal ionospheric TEC data set that is not affected by the space weather environment;

[0073] Improve the LSTM algorithm through the Prophet algorithm, and determine the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC based on the optimal ionospheric TEC dataset;

[0074] Determine whether there is an abnormal disturbance in the volcanic ionosphere according to the comparison relationship between the ionospheric TEC observation data and the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC.

[0075] This application also proposes a computer device, including a memory and a processor. The memory stores a computer program, and the computer program is executed by the processor as follows:

[0076] Obtain the ground-based GNSS raw observation data;

[0077] Preprocess the ground-based GNSS raw observation data to obtain an ionospheric TEC dataset;

[0078] According to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC dataset in the time domain and the frequency domain, remove the ionospheric TEC data with high correlation with the space weather environment, and obtain an optimal ionospheric TEC dataset that is not affected by the space weather environment;

[0079] Improve the LSTM algorithm through the Prophet algorithm, and determine the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC based on the optimal ionospheric TEC dataset;

[0080] Determine whether there is an abnormal disturbance in the volcanic ionosphere according to the comparison relationship between the ionospheric TEC observation data and the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC.

[0081] Adopting the embodiment of the present invention has the following beneficial effects: First, obtain the ground-based GNSS raw observation data; preprocess the ground-based GNSS raw observation data to obtain an ionospheric TEC dataset; according to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC dataset in the time domain and the frequency domain, remove the ionospheric TEC data with high correlation with the space weather environment, and obtain an optimal ionospheric TEC dataset that is not affected by the space weather environment; improve the LSTM algorithm through the Prophet algorithm, and determine the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC based on the optimal ionospheric TEC dataset; determine whether there is an abnormal disturbance in the volcanic ionosphere according to the comparison relationship between the ionospheric TEC observation data and the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC. The present invention realizes abnormal disturbance identification by data processing and removing space weather interference, and using the Prophet algorithm combined with the LSTM algorithm, and can efficiently monitor volcanic eruptions. Description of the Drawings

[0082] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the following will briefly introduce the accompanying drawings required in the description of the embodiments or the prior art. Obviously, the accompanying drawings in the following description are only some embodiments of the present invention. For those of ordinary skill in the art, without creative efforts, other accompanying drawings can be obtained based on these drawings.

[0083] Wherein:

[0084] Figure 1 It is a schematic diagram of the process for determining volcanic ionospheric anomaly disturbances in an embodiment;

[0085] Figure 2 It is a schematic diagram of the wavelet coherence spectrum of the TEC time series and space weather parameters at the TONG station in an embodiment;

[0086] Figure 3 It is a schematic diagram of the wavelet coherence spectrum of the TEC time series and space weather parameters at the TONG station in an embodiment;

[0087] Figure 4 It is a schematic diagram of the change of the TEC error sequence at the SAMO station in an embodiment;

[0088] Figure 5 It is a schematic diagram of the TEC sequence variation at the SAMO station in an embodiment. Detailed implementation manners

[0089] The following will clearly and completely describe the technical solutions in the embodiments of the present invention with reference to the accompanying drawings in the embodiments of the present invention. Obviously, the described embodiments are only some embodiments of the present invention, rather than all embodiments; based on the embodiments of the present invention, all other embodiments obtained by those of ordinary skill in the art without creative efforts belong to the scope of protection of the present invention.

[0090] In the following description, a large number of specific details are given to provide a more thorough understanding of the present invention; however, it is obvious to those skilled in the art that the present invention can be implemented without one or more of these details; in other examples, in order to avoid confusion with the present invention, some well-known technical features in the art are not described. It should be understood that the present invention can be implemented in different forms and should not be construed as limited to the embodiments presented here; on the contrary, providing these embodiments will make the disclosure thorough and complete, and will fully convey the scope of the present invention to those skilled in the art.

[0091] The purpose of the terms used herein is only to describe specific embodiments and is not a limitation of the present invention. As used herein, the singular forms "a", "an" and "the" are also intended to include the plural forms unless the context clearly dictates otherwise. It should also be understood that the terms "comprising" and / or "including", when used in this specification, identify the presence of the stated features, integers, steps, operations, elements and / or components, but do not preclude the presence or addition of one or more other features, integers, steps, operations, elements, components and / or groups. As used herein, the term "and / or" includes any and all combinations of the related listed items.

[0092] For ease of understanding, the relevant terms involved in this application are introduced below.

[0093] (1) Ionospheric TEC (Total Electron Content) is a key parameter in ionospheric physics research. It describes the total number of all electrons in a column of unit area in the ionosphere and includes two main forms: VTEC (Vertical Total Electron Content) and STEC (Slant Total Electron Content);

[0094] (2) The Prophet algorithm is a tool for time series prediction. It can efficiently process and predict time series data with complex trends and seasonal variations. This algorithm realizes the rapid iterative optimization of the time series model by combining multiple highly interpretable model components.

[0095] In the prior art, with the rapid development of global navigation satellite systems (GNSS) and space exploration technologies such as occultation, the detection and research of ionospheric anomalies have become a hot topic in the scientific community at home and abroad. However, although these technologies provide rich observational data, in the field of detecting ionospheric abnormal perturbations before volcanic eruptions, there is still a lack of a systematic and efficient method. Currently, most methods focus on the analysis of ionospheric perturbations after volcanic eruptions, while the early warning research before volcanic eruptions is relatively weak.

[0096] Therefore, developing a method that can accurately detect ionospheric abnormal perturbations before volcanic eruptions is of great significance for improving the early warning ability of volcanic eruptions and reducing their potential hazards. The method proposed in this paper is based on this background and aims to effectively detect and give early warnings of ionospheric abnormal perturbations before volcanic eruptions by comprehensively using ground-based GNSS data, space weather data, and advanced data processing and analysis techniques.

[0097] To thoroughly understand the present invention, detailed structures will be presented in the following description to illustrate the technical solutions proposed by the present invention; alternative embodiments of the present invention are described in detail below. However, in addition to these detailed descriptions, the present invention may also have other embodiments.

[0098] As Figure 1 shown, in one embodiment, a method for determining volcanic ionospheric anomaly disturbances is provided. The method for determining volcanic ionospheric anomaly disturbances includes steps S101 to S105, which are described in detail as follows:

[0099] S101. Obtain the original ground-based GNSS observation data;

[0100] Due to its global coverage and high-precision characteristics, ground-based GNSS (Global Navigation Satellite System) is used as a data source to provide a reliable data basis for the monitoring of ionospheric TEC (Total Electron Content).

[0101] S102. Preprocess the original ground-based GNSS observation data to obtain an ionospheric TEC data set;

[0102] By preprocessing the original observation data, noise and outliers can be effectively removed, and the data quality can be improved, thereby obtaining a more accurate ionospheric TEC data set.

[0103] In some embodiments, preprocessing the original ground-based GNSS observation data to obtain an ionospheric TEC data set specifically includes:

[0104] Based on the GNSS original pseudorange and carrier phase observation equations, determine the slant ionospheric delay Iono, where Ion0 = STEC + DCB, STEC is the slant TEC, and VTEC is the vertical TEC;

[0105] Determine the differential code bias DCB according to the non-differential non-combined precise point positioning model and the generalized trigonometric series model;

[0106] Determine the slant TEC and vertical TEC of the satellite piercing point of the GNSS station according to the slant ionospheric delay Iono and the differential code bias DCB, and obtain the ionospheric TEC data set.

[0107] Preferably, using GNSS raw pseudorange and carrier phase observation data (RINEX format files) as data sources, the precise ephemeris and clock error products use the post-precise ephemeris and clock errors released by the International GNSS Service, and the Differential Code Bias (DCB) uses the DCB products released by the Chinese Academy of Sciences (CAS). The TEQC data processing software is used to perform data format standardization processing and quality check analysis on the original measurement data. Then, the slant TEC (Slant TEC, STEC) and vertical TEC (Vertical TEC, VTEC) data of the satellite piercing points of the GNSS stations are obtained by solving using the non-differential and non-combination precise point positioning method and the generalized trigonometric series method.

[0108] In some embodiments, based on the GNSS raw pseudorange and carrier phase observation equations, the slant ionospheric delay Iono is determined, specifically including:

[0109] The GNSS raw pseudorange and carrier phase observation equations are:

[0110]

[0111]

[0112] Among them, the superscript s and the subscript r represent the satellite and the receiver respectively, and the subscript f represents the L1 or L2 frequency band of the carrier signal; represents the pseudorange observation value, represents the carrier phase observation value; represents the geometric distance between the satellite and the receiver; c represents the speed of light; dt r and dt s represent the clock errors of the receiver and the satellite respectively; Trop represents the slant tropospheric delay, represents the slant ionospheric delay; b r,f and represent the pseudorange hardware delays at the receiver end and the satellite end respectively; λ f represents the carrier phase wavelength, represents the carrier phase integer ambiguity; B r,f and represent the phase hardware delays at the receiver end and the satellite end respectively; ε P and ε L represent the sum of multipath effects, observation noise and other unmodeled errors in pseudorange and carrier observations respectively;

[0113] According to the least squares method, the GNSS raw pseudorange and carrier phase observation equations are solved to determine the slant ionospheric delay Iono.

[0114] Preferably, the above formulas (1) and (2) are processed to determine the corresponding linearized equations:

[0115]

[0116]

[0117] In the formula, u represents the coefficient matrix of the unknowns after linearization, and x represents the receiver position parameters.

[0118] Specifically, the ionospheric slant delay in the above GNSS raw pseudorange and carrier phase observation equations, that is, the slant ionospheric delay Iono.

[0119] In some embodiments, the differential code bias DCB is determined according to the non-differenced and non-combined precise point positioning model and the generalized trigonometric series model, specifically including:

[0120] Introduce the precise ephemeris and clock offset products released by IGS and use the ionosphere-free combination to estimate the satellite clock offset:

[0121]

[0122] Among them, IGS uses the ionosphere-free combination to estimate the satellite clock offset, and the final clock offset product absorbs the satellite pseudorange hardware delay.

[0123] Substitute Equation (5) into Equations (3) and (4) to obtain the final non-differenced and non-combined precise point positioning model. The estimated slant ionospheric delay (Iono) by this method includes the receiver DCB, that is, Ion0 = STEC + DCB. Therefore, convert STEC to VTEC, and the conversion formula is as follows:

[0124] VTEC = MF·STEC (6)

[0125]

[0126] Among them, R represents the radius of the earth, taking 6371 km; H represents the height of the ionospheric pierce point, with a value of 450 km, α = 0.9782; z represents the satellite elevation angle;

[0127] Use the generalized trigonometric series to estimate the DCB of the GNSS station. The generalized trigonometric series function model is as follows:

[0128]

[0129] In the formula, λ and φ respectively represent the geographical latitudes of the ionospheric piercing point and the modeling center point; h represents the function related to local time at the ionospheric piercing point; nmax, mmax, and kmax respectively represent the maximum orders of the polynomial function and the trigonometric series function; Enm, Ck, and Sk respectively represent the model parameters to be estimated of the generalized trigonometric series function.

[0130] Specifically, h = 2π(t - 14) / 24, where t is the local time.

[0131] Furthermore, based on the parameters related to the space weather environment before volcanic eruptions, including the 10.7-cm radiation flux (F10.7), extreme ultraviolet radiation (EUV), solar wind speed (Vsw) in solar activity parameters, the Bz component in the interplanetary magnetic field IMF, and the Kp index, Dst index, Ap index, AE index, etc. in geomagnetic activity indices, combined with the solved ionospheric TEC, the space weather conditions before volcanic eruptions are accurately analyzed using continuous wavelet transform (CWT) to exclude the interference brought by the space weather environment.

[0132] S103. According to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC dataset in the time domain and the frequency domain, remove the ionospheric TEC data with high correlation with the space weather environment, and obtain the optimal ionospheric TEC dataset that is not affected by the space weather environment;

[0133] Utilizing the characteristics of wavelet analysis in the time domain and the frequency domain, combined with variance analysis, can effectively identify and remove the ionospheric TEC data highly correlated with the space weather environment (such as solar wind, geomagnetic storms, etc.). This step significantly improves the accuracy and reliability of subsequent analysis and avoids misjudgment of external factors on the abnormal disturbances of the volcanic ionosphere.

[0134] In some embodiments, according to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC dataset in the time domain and the frequency domain, remove the ionospheric TEC data with high correlation with the space weather environment, and obtain the optimal ionospheric TEC dataset that is not affected by the space weather environment, which specifically includes:

[0135] The mother wavelet φ(t) is obtained through simple stretching and translation to get the family of daughter wavelets:

[0136]

[0137] where a represents the translation parameter controlling the wavelet position, b represents the scale parameter controlling the width scaling or expansion of the wavelet, and t represents time;

[0138] By giving the ionospheric TEC time series or the solar activity parameter time series, its continuous wavelet transform formula is obtained:

[0139]

[0140] In the formula, φ * (t) represents the complex conjugate function of φ(t), and x represents the ionospheric TEC time series or the solar activity parameter time series;

[0141] Among them, the position of the wavelet in the time domain is given by a, and the position in the frequency domain is given by b. The wavelet transform maps the original sequence into a function of a and b, and at the same time obtains the information of time and frequency.

[0142] Determine the variance of the wavelet energy spectrum measure time series of the ionospheric TEC data set in the time and frequency domains:

[0143] |(WPS) x (a,b)| = |W x (a,b)| 2 (11)

[0144]

[0145] Among them, W x (a,b) and W y (a,b) respectively represent the wavelet transforms of x(t) and y(t), represents the complex conjugate of W y (a,b).

[0146] Preferably, in the application, to quantify the relationship between TEC and related parameters, the cross-wavelet transform of two time series x(t) and y(t) is defined as shown in Equation (4), where W x (a,b) and W y (a,b) are respectively the wavelet transforms of x(t) and y(t), is the complex conjugate of W y (a,b):

[0147] Such as Figure 2 、 Figure 3 shows an example of using wavelet transform to exclude the interference of the solar-terrestrial space environment on the detection of ionospheric anomalies before the Tonga volcano eruption in 2022.

[0148] By applying cross-wavelet transform and wavelet coherence spectrum, the correlation between TEC over the TONG station and six solar-terrestrial space environment parameters from January 1st to 16th was analyzed, as shown in Figure 2 and Figure 3 ; Figure 2The color bar on the right represents the cross-wavelet power spectral density, and the arrow direction represents the phase relationship between the two: a rightward arrow indicates that the two sequences are in phase, a leftward arrow indicates an anti-phase, a vertically downward arrow indicates that the former sequence changes 1 / 4 of a period ahead of the latter sequence, and a vertically upward arrow indicates that the latter sequence changes 1 / 4 of a period ahead of the former sequence; Figure 3 The wavelet coherence spectrum of TEC and other parameters is given, and the color bar on the right is the wavelet coherence value, which characterizes the strength of coherence.

[0149] Such as Figure 2 As shown in the cross-wavelet transform of the TEC time series and space weather parameters at the TONG station from January 1 to 16, 2022. Among them, the closed area of the thick black solid line passed the standard red noise test at the 95% confidence level, indicating the significance of the period; the conical area below the thin black solid line is the Cone of Influence (COI) area, which is the area where the edge effect of the wavelet transform data has a greater impact. The thick red solid line represents the volcanic eruption time, as Figure 3 As shown in the schematic diagram of the wavelet coherence spectrum of the TEC time series and space weather parameters at the TONG station from January 1 to 16, 2022.

[0150] In summary, by comparing the resonance period, phase relationship, and coherence between TEC and other parameters, the spatial environment interference on the 6th, 8th, 9th, 14th, and 15th can be effectively excluded.

[0151] S104. Improve the LSTM algorithm through the Prophet algorithm, and determine the upper and lower boundaries of the abnormal perturbation of volcanic ionospheric TEC based on the optimal ionospheric TEC dataset;

[0152] The Prophet algorithm is used to improve the LSTM (Long Short-Term Memory) algorithm, which combines the accuracy of the Prophet in time series prediction and the advantage of the LSTM in dealing with long-term dependence relationships. This improvement enables the model to more accurately capture the non-linear changes and long-term trends in the ionospheric TEC data, thereby more precisely determining the upper and lower boundaries of the abnormal perturbation of volcanic ionospheric TEC.

[0153] In some embodiments, before improving the LSTM algorithm through the Prophet algorithm and determining the upper and lower boundaries of the abnormal perturbation of volcanic ionospheric TEC based on the optimal ionospheric TEC dataset, the method further includes:

[0154] Obtain the TEC prediction values of the Prophet algorithm and the LSTM algorithm:

[0155] y(t) = g(t) + s(t) + h(t) + ε t (14)

[0156] where \(g(t)\) represents the trend term; \(s(t)\) represents the periodic term; \(h(t)\) represents holidays, indicating the impact of holidays that may occur irregularly within one or more days; \(\epsilon\) t represents the error term; \(y(t)\) represents the TEC prediction value of the Prophet algorithm;

[0157] Relevant parameters for constructing the LSTM neural network:

[0158]

[0159]

[0160]

[0161]

[0162]

[0163] where \(C\) t-1 represents the cell state at time \(t - 1\), \(h\) t-1 represents the output at the previous moment, \(x\) t represents the current input, \(h\) t represents the current output, i.e., the TEC prediction value of the LSTM algorithm, \(C\) t represents the cell state at the current time \(t\).

[0164] The Prophet algorithm is a tool for time series prediction, which realizes the rapid iterative optimization of the time series model. For trend change points, the algorithm will automatically detect. The LSTM method has a recursive network structure, which consists of several LSTM modules. The LSTM module includes a memory unit, an input gate, an output gate, and a forget gate, thus having the ability of long-term memory. Incorporating the Prophet algorithm into the LSTM network makes the newly constructed method have higher sensitivity and better applicability when detecting ionospheric anomalies.

[0165] Above, \(g(t)\) is used to simulate the non-periodic changes in the time series values; \(s(t)\) represents the periodic term, such as weekly changes and seasonal changes; \(h(t)\) represents holidays, indicating the impact of holidays that may occur irregularly within one or more days; \(\epsilon\) t represents the error term, assuming it follows a normal distribution.

[0166] When applying the Prophet algorithm, the principal component analysis method is used to determine the applicability of the three models of the trend term g(t), the periodic term s(t), and the holiday term h(t) in the Prophet algorithm for detecting ionospheric TEC anomalies before volcanic eruptions. At the same time, four types of hyperparameters (degree of trend change, seasonal flexibility, holiday flexibility, and seasonal multiplication) are adjusted. The initial default values of the four types of hyperparameters are set to 0.05, 10, 10, and additive. The adjustment is made to minimize the RMSE, thereby determining the optimal detection combination model.

[0167] Preferably, the memory units in the LSTM network tend to retain the useful information in the sequence model, and the most critical forget gate will cause the information selected previously to be forgotten. Its input is composed of the state value c at time t-1 t-1 , the input value x at time t t and the output value h at the previous moment. t-1 It consists of.

[0168] The established LSTM model sequence includes a bidirectional LSTM layer (units = 64), a Dropout layer (rate = 0.2), an LSTM layer (units = 64), a Dropout layer (rate = 0.2), and a Dense layer (units = 1). The optimizer is selected as Adam, the loss function is mse. When training the model, batch = 64, epochs = 12, validation_split = 0.3 are selected, and the data of the previous 27 days are used for rolling prediction of the TEC of the next day to obtain the LSTM model for best detecting ionospheric TEC anomalies.

[0169] The Prophet algorithm has good applicability for detecting significant ionospheric anomalies, but it cannot accurately detect relatively weak ionospheric anomalies; while the LSTM method solves this defect, but there is a deviation in the stability of detecting anomalies. Therefore, the LSTM algorithm is improved through the Prophet algorithm, and the upper and lower boundaries of the volcanic ionospheric TEC anomaly perturbation are determined through the optimal ionospheric TEC dataset.

[0170] In some embodiments, the LSTM algorithm is improved through the Prophet algorithm, and the upper and lower boundaries of the volcanic ionospheric TEC anomaly perturbation are determined based on the optimal ionospheric TEC dataset, specifically including:

[0171] The predicted values of the Prophet algorithm are respectively given and the predicted values of the LSTM model The initial weights ω1 and ω2 of the ionospheric TEC at time i.

[0172] Through error accuracy analysis with the measured ionospheric TEC, the objective function is determined:

[0173]

[0174] Among them, N represents the total number of detection periods; yi represents the measured ionospheric TEC; represents the TEC prediction value of the Prophet algorithm at time i, represents the TEC prediction value of the LSTM algorithm at time i.

[0175] Preset constraint conditions: ω1 + ω2 = 1, ω1 ≥ 0, ω2 ≥ 0;

[0176] Use the least squares method to determine the optimal weights ω1 and ω2, and determine the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC.

[0177] Specifically, y(t) in formula (14) is the prediction value of the Prophet algorithm h in formula (15) t is the prediction value of the LSTM model

[0178] In some embodiments, use the least squares method to determine the optimal weights ω1 and ω2, and determine the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC, specifically including:

[0179] Error = TEC actual -TEC predicted (17)

[0180] UB = Mean(Error) + 2 × Std(Error) (18)

[0181] LB = Mean(Error) - 2 × Std(Error) (19)

[0182] In the formula, TEC actual represents the measured ionospheric TEC; TEC predicted represents the ionospheric TEC predicted by the new algorithm; Mean represents taking the average value; Std represents calculating the standard deviation; UB represents the upper boundary of the abnormal ionospheric TEC; LB represents the lower boundary of the abnormal ionospheric TEC.

[0183] S105. Determine whether there is an abnormal disturbance in the volcanic ionosphere according to the comparison relationship between the ionospheric TEC observation data and the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC.

[0184] By comparing the ionospheric TEC observation data with the determined upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC, it is possible to clearly identify whether volcanic activities have caused abnormal disturbances in the ionosphere. This method is not only scientific and efficient but also can provide important references for volcanic monitoring and early warning.

[0185] In some embodiments, the values of ω1 and ω2 are set to 0.6 and 0.4 respectively, Figure 4 and Figure 5 shows the ionospheric anomalies before the Tonga volcano eruption in 2022 detected by this method. As Figure 4 shown: the ionospheric TEC changes over SAMO station from January 1st to 16th, 2022 (the red line represents the true TEC value, the blue line represents the algorithm prediction value, and the thick red solid line represents the volcanic eruption time). The thick red solid line marks the specific time of the volcanic eruption. By comparing the red and blue lines, it can be clearly seen that there are significant fluctuations in the ionospheric TEC before and after the volcanic eruption. Before the volcanic eruption, the algorithm prediction value (blue line) and the true value (red line) show a relatively close trend. However, at the time of the volcanic eruption and for some time after that, the prediction value shows a large deviation, reflecting the strong disturbance of the volcanic eruption on the ionosphere.

[0186] As Figure 5 shown: the ionospheric TEC anomalies over SAMO station from January 1st to 16th, 2022 (the black line and the blue line represent the upper and lower boundaries of the TEC anomaly respectively, and the thick red solid line represents the volcanic eruption time); the thick red solid line also marks the time of the volcanic eruption. It can be seen from the figure that before and after the volcanic eruption, the TEC value quickly crosses the anomaly boundary, indicating that volcanic activity has a significant abnormal disturbance on the ionosphere. By comparing the deviation degree of the TEC value from the anomaly boundary, the specific impact degree and duration of the volcanic eruption on the ionosphere can be quantitatively evaluated.

[0187] In summary, before the volcanic eruption, the method proposed in this application can better predict the change trend of the ionospheric TEC. By setting the upper and lower boundaries of the TEC anomaly, the ionospheric abnormal disturbance caused by the volcanic eruption is effectively identified, proving that this scheme has high sensitivity and accuracy in determining the volcanic ionospheric abnormal disturbance. By real-time monitoring the change of the ionospheric TEC, signs of volcanic activity can be detected in advance, providing valuable warning time for relevant departments.

[0188] This application also proposes a volcanic ionospheric abnormal disturbance detection system, which includes: a data acquisition unit, a judgment condition preset unit, and a determination unit;

[0189] The data acquisition unit is used to obtain the raw ground-based GNSS observation data; preprocess the raw ground-based GNSS observation data to obtain the ionospheric TEC data set; remove the ionospheric TEC data with high correlation with the space weather environment according to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC data set in the time domain and the frequency domain, and obtain the optimal ionospheric TEC data set that is not affected by the space weather environment;

[0190] A judgment condition preset unit is used to improve the LSTM algorithm through the Prophet algorithm and determine the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC based on the optimal ionospheric TEC dataset.

[0191] A determination unit is used to determine whether an abnormal disturbance occurs in the volcanic ionosphere according to the comparison relationship between the ionospheric TEC dataset and the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC.

[0192] This application also proposes a readable storage medium storing a computer program, which, when executed by a processor, causes the processor to perform the following steps:

[0193] Obtain the original ground-based GNSS observation data;

[0194] Preprocess the original ground-based GNSS observation data to obtain an ionospheric TEC dataset;

[0195] According to the variances of the time series of the wavelet energy spectrum measure of the ionospheric TEC dataset in the time domain and the frequency domain, remove the ionospheric TEC data with a high correlation with the space weather environment to obtain an optimal ionospheric TEC dataset that is not affected by the space weather environment;

[0196] Improve the LSTM algorithm through the Prophet algorithm and determine the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC based on the optimal ionospheric TEC dataset;

[0197] Determine whether an abnormal disturbance occurs in the volcanic ionosphere according to the comparison relationship between the ionospheric TEC observation data and the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC.

[0198] This application also proposes a computer device, including a memory and a processor, where the memory stores a computer program, and the computer program is executed by the processor to perform the following steps:

[0199] Obtain the original ground-based GNSS observation data;

[0200] Preprocess the original ground-based GNSS observation data to obtain an ionospheric TEC dataset;

[0201] According to the variances of the time series of the wavelet energy spectrum measure of the ionospheric TEC dataset in the time domain and the frequency domain, remove the ionospheric TEC data with a high correlation with the space weather environment to obtain an optimal ionospheric TEC dataset that is not affected by the space weather environment;

[0202] Improve the LSTM algorithm through the Prophet algorithm and determine the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC based on the optimal ionospheric TEC dataset;

[0203] According to the comparison relationship between the ionospheric TEC observation data and the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC, it is determined whether the volcanic ionosphere has an abnormal disturbance.

[0204] Those of ordinary skill in the art can understand that all or part of the processes in the methods of the above embodiments can be completed by instructing relevant hardware through a computer program. The program can be stored in a non-volatile computer-readable storage medium. When the program is executed, it can include the processes of the embodiments of the above methods. Among them, any reference to a memory, storage, database, or other medium used in the various embodiments provided in the present application can include non-volatile and / or volatile memories. Non-volatile memories can include read-only memory (ROM), programmable ROM (PROM), electrically programmable ROM (EPROM), electrically erasable programmable ROM (EEPROM), or flash memory. Volatile memories can include random access memory (RAM) or external cache memory. By way of illustration and not limitation, RAM is available in many forms, such as static RAM (SRAM), dynamic RAM (DRAM), synchronous DRAM (SDRAM), double data rate SDRAM (DDR SDRAM), enhanced SDRAM (ESDRAM), synchronous link (Synchlink) DRAM (SLDRAM), Rambus direct RAM (RDRAM), direct memory bus dynamic RAM (DRDRAM), and Rambus dynamic RAM (RDRAM), etc.

[0205] The technical features of the above embodiments can be combined arbitrarily. For the sake of brevity of description, not all possible combinations of the technical features in the above embodiments are described. However, as long as there is no contradiction in the combination of these technical features, it should be considered to be within the scope described in this specification.

[0206] The above-described embodiments only represent several implementation manners of the present application. Their descriptions are relatively specific and detailed, but they should not be construed as limiting the patent scope of the present application. It should be noted that for those of ordinary skill in the art, without departing from the concept of the present application, several modifications and improvements can still be made, and these all belong to the protection scope of the present application. The above-disclosed is only the preferred embodiment of the present invention, and of course, it cannot be used to limit the scope of the rights of the present invention. Therefore, the equivalent changes made according to the claims of the present invention still fall within the scope covered by the present invention.

Claims

1. A method for determining volcanic ionospheric anomaly disturbances, characterized in that, The method includes: Obtaining the original ground-based GNSS observation data; Preprocessing the original ground-based GNSS observation data to obtain an ionospheric TEC data set; Removing the ionospheric TEC data with high correlation with the space weather environment according to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC data set in the time domain and the frequency domain, and obtaining an optimal ionospheric TEC data set that is not affected by the space weather environment; Improving the LSTM algorithm through the Prophet algorithm, and determining the upper and lower boundaries of the volcanic ionospheric TEC anomaly perturbation based on the optimal ionospheric TEC data set, specifically including: Respectively assign the predicted values of the Prophet algorithm and the predicted values of the LSTM algorithm The initial weights ω1 and ω2 of the ionospheric TEC at time i; Determining the objective function through error accuracy analysis with the measured ionospheric TEC; Among them, N represents the total number of detection time periods; yi represents the measured ionospheric TEC; represents the TEC prediction value of the Prophet algorithm at the i-th moment; represents the TEC prediction value of the LSTM algorithm at the i-th moment; Presetting the constraint conditions: ω1 + ω2 = 1, ω1 ≥ 0, ω2 ≥ 0; Determine the optimal weights ω1 and ω2 using the least squares method, and determine the upper and lower boundaries of the volcanic ionospheric TEC anomaly perturbation, specifically including: Error = TEC actual -TEC predicted UB = Mean(Error) + 2 × Std(Error) LB = Mean(Error) - 2 × Std(Error) where TEC actual represents the measured ionospheric TEC; TEC predicted represents the ionospheric TEC predicted by the new algorithm; Mean represents taking the average value; Std represents calculating the standard deviation; UB represents the upper boundary of the ionospheric TEC anomaly; LB represents the lower boundary of the ionospheric TEC anomaly; Determining whether there is an abnormal perturbation in the volcanic ionosphere according to the comparison relationship between the ionospheric TEC observation data and the upper and lower boundaries of the volcanic ionospheric TEC anomaly perturbation.

2. The method for determining volcanic ionospheric anomalous disturbances according to claim 1, wherein The preprocessing of the original ground-based GNSS observation data to obtain an ionospheric TEC data set specifically includes: Based on the GNSS original pseudorange and carrier phase observation equations, determining the slant ionospheric delay Iono, where Ion0 = STEC + DCB, STEC is the slant TEC, and VTEC is the vertical TEC; Determining the differential code bias DCB according to the non-differenced non-combined precise point positioning model and the generalized trigonometric series model; Determining the slant TEC and vertical TEC of the satellite piercing point of the GNSS station according to the slant ionospheric delay Iono and the differential code bias DCB, and obtaining an ionospheric TEC data set.

3. The method for determining volcanic ionospheric abnormal disturbances according to claim 2, wherein The determining of the slant ionospheric delay Iono based on the GNSS original pseudorange and carrier phase observation equations specifically includes: The GNSS original pseudorange and carrier phase observation equations are: Wherein, the superscript s and subscript r respectively represent the satellite and the receiver, and the subscript f represents the L1 or L2 frequency band of the carrier signal; represents the pseudorange observation value, represents the carrier phase observation value; represents the geometric distance between the satellite and the receiver; c represents the speed of light; dt r and dt s respectively represent the clock biases of the receiver and the satellite; Trop represents the tropospheric slant delay, represents the ionospheric slant delay; b r,f and respectively represent the pseudorange hardware delays at the receiver end and the satellite end; λ f represents the carrier phase wavelength, represents the carrier phase integer ambiguity; B r,f and respectively represent the phase hardware delays at the receiver end and the satellite end; ε P and ε L respectively represent the sum of multipath effects, observation noise, and other unmodeled errors in pseudorange and carrier observations; Solving the GNSS original pseudorange and carrier phase observation equations according to the least squares method to determine the slant ionospheric delay Iono.

4. The method for determining volcanic ionospheric anomaly disturbance according to claim 3, wherein The determining of the differential code bias DCB according to the non-differenced non-combined precise point positioning model and the generalized trigonometric series model specifically includes: Introducing the precise ephemeris and clock error products released by IGS and using the ionosphere-free combination to estimate the satellite clock error: Converting STEC to VTEC, and the conversion formula is as follows: VTEC = MF · STEC where, R represents the radius of the earth; H represents the height of the ionospheric piercing point, α represents the coefficient; z represents the satellite elevation angle; Estimating the DCB of the GNSS station using the generalized trigonometric series, and the generalized trigonometric series function model is as follows: Wherein, represents the vertical total electron content at the height h of the ionospheric piercing point where the geographical latitude is ; and represent the geographical latitudes of the ionospheric piercing point and the modeling center point respectively; h represents the height of the ionospheric piercing point; n, m, and k are all natural numbers, representing the orders of the polynomial function and the trigonometric series function respectively; n max , m max and k max represent the maximum orders of the polynomial function and the trigonometric series function respectively; E nm is the model parameter to be estimated for the polynomial function part, and C k and S k represent the model parameters to be estimated for the generalized trigonometric series function respectively.

5. The method for determining volcanic ionospheric abnormal disturbances according to claim 4, characterized in that The removing of the ionospheric TEC data with high correlation with the space weather environment according to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC data set in the time domain and the frequency domain, and obtaining an optimal ionospheric TEC data set that is not affected by the space weather environment specifically includes: Obtaining the sub-wavelet family by simply stretching and translating the mother wavelet φ(t): Among them, a represents the translation parameter for controlling the wavelet position, b represents the scale parameter for controlling the width scaling or dilation of the wavelet, and t represents time; By given ionospheric TEC time series or solar activity parameter time series, obtain its continuous wavelet transform formula: where φ * (t) represents the complex conjugate function of φ(t), and x represents the ionospheric TEC time series or the solar activity parameter time series; Determine the variance of the wavelet energy spectrum measure time series of the ionospheric TEC data set in the time and frequency domains: |(WPS) x (a,b)|=|W x (a,b)| 2 Among them, W x (a, b) and W y (a, b) respectively represent the wavelet transforms of x(t) and y(t), denotes W y (a, b)'s complex conjugate.

6. The method for determining volcanic ionospheric abnormal disturbances according to claim 5, wherein Before determining the upper and lower boundaries of the volcanic ionospheric TEC abnormal perturbation based on the optimal ionospheric TEC data set by improving the LSTM algorithm through the Prophet algorithm, this method further includes: obtaining the TEC prediction values of the Prophet algorithm and the LSTM algorithm: y(t) = g(t) + s(t) + h(t) + ε t where \(g(t)\) represents the trend term; \(s(t)\) represents the periodic term; \(h(t)\) represents holidays, indicating the impact of holidays that may occur irregularly within one or more days; \(\varepsilon\) t represents the error term; \(y(t)\) represents the TEC prediction value of the Prophet algorithm; Construct relevant parameters of the LSTM neural network: where C t-1 represents the cell state at time t-1, h t-1 represents the output at the previous time, x t represents the current input, h t represents the current output, i.e., the TEC prediction value of the LSTM algorithm, C t represents the cell state at the current time t.

7. A volcanic ionospheric anomaly detection system, characterized in that, The system includes: a data acquisition unit, a judgment condition preset unit, and a determination unit; The data acquisition unit is used to obtain the ground-based GNSS raw observation data; preprocess the ground-based GNSS raw observation data to obtain an ionospheric TEC data set; according to the variances of the wavelet energy spectrum measure time series of the ionospheric TEC data set in the time domain and the frequency domain, remove the ionospheric TEC data with high correlation with the space weather environment, and obtain an optimal ionospheric TEC data set that is not affected by the space weather environment; The judgment condition presetting unit is used to improve the LSTM algorithm through the Prophet algorithm, and determine the upper and lower boundaries of the abnormal disturbance of the volcanic ionospheric TEC based on the optimal ionospheric TEC data set, specifically including: respectively assigning the predicted values of the Prophet algorithm and the predicted values of the LSTM algorithm The initial weights ω1 and ω2 of the ionospheric TEC at time i; through error accuracy analysis with the measured ionospheric TEC, determine the objective function: Among them, N represents the total number of detection time periods; yi represents the measured ionospheric TEC; represents the TEC prediction value of the Prophet algorithm at the i-th moment, represents the TEC prediction value of the LSTM algorithm at the i-th moment; Preset the constraint conditions: ω1 + ω2 = 1, ω1 ≥ 0, ω2 ≥ 0; Determine the optimal weights ω1 and ω2 using the least squares method, and determine the upper and lower boundaries of the volcanic ionospheric TEC anomaly perturbation, specifically including: Error = TEC actual -TEC predicted UB = Mean(Error) + 2 × Std(Error) LB = Mean(Error) - 2 × Std(Error) Wherein, TEC actual represents the measured ionospheric TEC; TEC predicted represents the ionospheric TEC predicted by the new algorithm; Mean represents taking the average value; Std represents calculating the standard deviation; UB represents the upper boundary of the ionospheric TEC anomaly; LB represents the lower boundary of the ionospheric TEC anomaly; The determination unit is used to determine whether there is an abnormal perturbation in the volcanic ionosphere according to the comparison relationship between the ionospheric TEC data set and the upper and lower boundaries of the volcanic ionospheric TEC abnormal perturbation.

8. A readable storage medium stores a computer program, and when the computer program is executed by a processor, the processor is caused to execute the steps of the method according to any one of claims 1 to 6.

9. A computer device includes a memory and a processor, the memory stores a computer program, and when the computer program is executed by the processor, the processor is caused to execute the steps of the method according to any one of claims 1 to 6.

Images