A detection index for identifying ionospheric traveling wave disturbances and equatorial ionospheric bubbles

Abstract

The application discloses a detection index capable of identifying ionosphere traveling wave disturbance and equatorial plasma bubble. The index is based on a traditional ionosphere disturbance index ROTI, combines the advantages of a detrended total electron content (DTEC), can effectively detect TIDs and EPBs, and distinguishes the two through a threshold.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of GNSS data processing and ionospheric disturbance detection, and particularly relates to a detection index capable of identifying ionospheric traveling ionospheric disturbance and equatorial plasma bubble. BACKGROUND

[0002] As an important part of the space environment, the ionosphere is a typical non-stationary, multi-scale, and nonlinear complex system. Due to the influence of solar radiation, geomagnetic activity, and other factors, the electron density of the ionosphere often experiences sharp irregular changes, i.e., ionospheric disturbances. Ionospheric disturbance phenomena can significantly affect the radio propagation path through the ionosphere region, and in severe cases, can even cause satellite signal loss, spacecraft malfunction, and even damage. In order to prevent the serious influence of ionospheric disturbances on satellite navigation and positioning and radio communication, monitoring ionospheric disturbance phenomena is of great significance to improve navigation and positioning accuracy and ensure the safety of space activities.

[0003] At middle and low latitudes, traveling ionospheric disturbances (TIDs) and equatorial plasma bubbles (EPBs) are two typical ionospheric disturbance phenomena. TIDs are quasi-periodic fluctuations in ionospheric electron density caused by solar flares, geomagnetic storms, solar eclipses, earthquakes, hurricanes, volcanic eruptions, and other events. EPBs are nighttime ionospheric plasma density depletions related to generalized Rayleigh-Taylor instability triggered by disturbances at the bottom of the ionospheric F layer. The generated plasma bubbles can reach low and even middle latitudes and can spread from a few meters to hundreds of kilometers. Both TIDs and EPBs can cause rapid changes in the amplitude and phase of radio signals passing through the region, thereby affecting the safety and functionality of global satellite navigation systems (GNSS) and possibly causing satellite communication interruption in some extreme cases. Moreover, these two types of events can occur in the same space-time, coexist, couple, or even transform into each other. Therefore, it is crucial to quickly detect and accurately identify TIDs and EPBs.

[0004] At present, the means of detecting ionospheric disturbances such as TIDs and EPBs include airglow imagers, altimeters, high-frequency radars, etc. However, such devices are sparse in distribution, fixed in observation area and expensive in price, which is difficult to meet the needs of large-scale monitoring and early warning. GNSS provides a solution for large-scale and high-precision monitoring of ionospheric disturbances due to its global coverage and all-weather characteristics. However, there are still some new problems and challenges in the detection method of ionospheric disturbances based on GNSS data. The most commonly used GNSS-based ionospheric space environment disturbance detection methods are mainly divided into two categories: two-dimensional detrended total electron content map (DTECMap) and single-station exponential method. The former can accurately represent the characteristics of TIDs, but it is difficult to detect the propagation of TIDs in the absence of sufficient dense reference station data. The latter ROTI index method can detect EPBs, but it is not sensitive to TIDs and cannot effectively identify TIDs. On the other hand, the distinction between ionospheric disturbance types (such as TIDs and EPBs) is usually made by manual visual judgment of images, and it is currently difficult to automatically distinguish and identify the two disturbance types through quantitative indicators. SUMMARY

[0005] In view of the deficiencies in the prior art, the present application provides a detection index IROTI capable of identifying ionospheric traveling wave disturbances and equatorial plasma bubbles, which can effectively solve many difficulties encountered in detecting ionospheric disturbances based on GNSS, and provide a method for detecting and distinguishing TIDs and EPBs for those skilled in the art.

[0006] In order to achieve the above-mentioned purpose, the technical scheme provided by the present application is a detection index capable of identifying ionospheric traveling wave disturbances and equatorial plasma bubbles, comprising the following steps:

[0007] Step 1: obtaining GNSS original observation files and broadcast ephemeris files of an observation area;

[0008] Step 2: detecting cycle slips of GNSS original observation data to obtain epoch numbers and satellite PRN numbers of continuous observation arcs;

[0009] Step 3: constructing a geometry-free phase observation equation based on GNSS original observation data to obtain ionospheric total electron content STEC containing pseudorange and phase hardware delay bias;

[0010] Step 4: screening the STEC obtained in step 3 using the continuous observation arc epoch numbers and satellite PRN numbers obtained in step 2 to obtain STEC sequences in the same continuous arc;

[0011] Step 5: calculating ionospheric disturbance index ROT index using the STEC sequences in the same continuous arc to obtain a time sequence of ROTI index;

[0012] Step 6: Use the moving average method to detrend STEC to obtain the time series of detrended total electron quantity DTEC;

[0013] Step 7: Combine the ROTI index obtained in Step 5 and the DTEC obtained in Step 6 to calculate the IROTI index.

[0014] Furthermore, in step 2, the Turbedit method is used to perform cycle slip detection on the raw observation data obtained in step 1 to obtain the epoch number and satellite PRN number with continuous observation arcs, and the shortest continuous arc is set to epoch μ1 to ensure that the subsequent moving average acquisition of detrended total electron quantity DTEC can proceed normally.

[0015] Furthermore, the GNSS signal in step 3 consists of two parts: pseudorange and carrier phase observations, which can be expressed as:

[0016]

[0017] In the formula, the subscript r and superscript s represent the serial numbers of the GNSS receiver and satellite, respectively, and the subscript i represents the frequency. and These are the pseudorange and phase observations along the line-of-sight (LoS) direction between the satellite and the receiver. dt represents the geometric distance between the satellite and the receiver, c is the speed of light in a vacuum, and dt is the distance between the satellite and the receiver. r and dt s These represent the clock biases of the GNSS receiver and the satellite, respectively. For tropospheric delay along the LoS direction, B is the ionospheric delay at frequency i. r,i and Let b represent the pseudorange hardware delay of the GNSS receiver and the satellite at frequency i, respectively. r,i and λ represents the fractional phase deviation of the GNSS receiver and the satellite at frequency i, respectively. i For wavelength, For phase integer ambiguity, and It is a combination of observation noise and various other errors.

[0018] The geometric phase-free observation equation is obtained by the difference between carrier phase observations at different frequencies, i.e.:

[0019]

[0020] In the formula, These are observations without geometric phase combination. and These are the phase observations along the LoS direction at frequencies f1 and f2, respectively. For the ionospheric delay along a given LoS direction, γ i =(f0 / f i ) 2 λ is a scaling factor, where f0 = 10.23 MHz, and λ1 and λ2 are the wavelengths of frequencies f1 and f2, respectively. and These are the integer ambiguities at frequencies f1 and f2, respectively. This indicates multipath effects, noise errors, and other types of errors.

[0021] The total electron content of the ionosphere (STEC) is defined as the integral of the electron density along the LoS direction from the satellite to the receiver. The relationship between the STEC at frequency f and the ionospheric delay I in units of length can be described as follows:

[0022]

[0023] In the formula, the STEC value is expressed in TECU, where 1 TECU = 10 16 el / m 2 ; Observations without geometric phase combination; γ i =(f0 / f i ) 2 It is a scaling factor, where f0 = 10.23 MHz; This represents multipath effects, noise errors, and other types of errors. The cutoff elevation angle of STEC is set to μ2 to mitigate the effects of multipath effects.

[0024] Furthermore, in step 5, the ROT exponent is the time derivative of the STEC value over two consecutive time periods, i.e.:

[0025]

[0026] In the formula, the ROT value is used as a parameter for measuring phase fluctuation activity and is expressed in TECU / min; These represent the total ionospheric electron content of satellite k at times i-1 and i, respectively, and t i-1 t i is the epoch time, and k is the satellite sequence number.

[0027] ROTI is the standard deviation of ROT over a specific period, and its calculation method is as follows:

[0028]

[0029] In the formula, angle brackets represent the average value over a specific period.

[0030] Furthermore, the formula for calculating the detrended total electron quantity (DTEC) using the moving average method in step 6 is as follows:

[0031]

[0032] In the formula, STEC(i) is the STEC value at time i, N is the sliding window, STEC(j) is the STEC value at time j, and the value of j is in the range of [iN / 2, i+N / 2].

[0033] Furthermore, in step 7, the integral of the DTEC value over a specific period is used as the index IDTEC for determining the occurrence of a disturbance:

[0034]

[0035] In the formula, |DTEC(i)| is the absolute value of DTEC at time i, Δt is the time interval of discrete DTEC, n is the number of DTECs in a specific period, and the IDTEC value is expressed in TECU·min.

[0036] The IROTI index, combining IDTEC and ROTI, is proposed. This index considers changes in amplitude and phase to detect ionospheric disturbances and thus distinguish between TIDs and EPBs. The specific calculation formula can be expressed as:

[0037] IROTI = IDTEC·ROTI(8)

[0038] In the formula, the IROTI value is expressed in TECU. 2 Expressed in units.

[0039] When IROTI < μ3TECU 2 When IROTI > μ4TECU, the ionospheric disturbance type can be determined to be TIDs. 2 When this occurs, the disturbance type can be considered as EPBs.

[0040] Compared with the prior art, the present invention has the following advantages:

[0041] Based on the traditional ionospheric disturbance index ROTI, the IROTI index is proposed by combining the advantages of detrended total electron quantity DTEC. It can effectively solve the problem that the current traditional ionospheric disturbance detection indices TIDs and EPBs are difficult to detect at the same time, and the two can be identified and distinguished by a threshold. Attached Figure Description

[0042] Figure 1 This is the overall process framework for calculating the IROTI index in this invention.

[0043] Figure 2This is a DTEC time series diagram of artificially simulated TIDs and EPBs signals according to an embodiment of the present invention, wherein, Figure 2 (a) is the analog signal of TIDs. Figure 2 (b) is the analog signal of EPBs.

[0044] Figure 3 The values ​​of ROTI, IROTI, and DTEC calculated by the HKWS station on June 22, 2015, according to this embodiment of the invention, are as follows: Figure 3 (a) ROTI calculated for the HKWS station. Figure 3 (b) IROTI calculated for the HKWS station. Figure 3 (c) DTEC calculated for HKWS station.

[0045] Figure 4 The values ​​of ROTI, IROTI, and DTEC calculated by the METU station on July 27, 2012, according to this embodiment of the invention, are as follows: Figure 4 (a) ROTI calculated for METU station. Figure 4 (b) IROTI calculated for the TIDs period Figure 4 (c) IROTI calculated for the EPBs period Figure 4 (d) DTEC calculated for the TIDs period and Figure 4 (e) is the DTEC calculated for the EPBs period.

[0046] Figure 5 The above are statistical results of ROTI and IROTI during the TIDs period and the ionospheric calm period in embodiments of the present invention, wherein... Figure 5 (a) Statistical results of ROTI during TIDs and ionospheric calm periods. Figure 5 (b) Statistical results of IROTI during TIDs and ionospheric calm periods.

[0047] Figure 6 The statistical results of ROTI and IROTI for TIDs and EPBs events in embodiments of the present invention are presented. Detailed Implementation

[0048] This invention provides a detection index that can identify ionospheric traveling wave disturbances and equatorial plasma bubbles. Based on the traditional ionospheric disturbance index ROTI, this index combines the advantages of detrended total electron quantity (DTEC) and can effectively detect TIDs and EPBs, and distinguish between them through thresholds.

[0049] The technical solution of the present invention will be further described below with reference to the accompanying drawings and embodiments.

[0050] like Figure 1As shown, the present invention provides a detection index capable of identifying ionospheric traveling wave disturbances and equatorial plasma bubbles, comprising the following steps:

[0051] Step 1: Obtain the raw GNSS observation files and broadcast ephemeris files for the observation area.

[0052] Step 2: Perform cycle slip detection on the raw GNSS observation data to obtain the epoch number and satellite PRN number with continuous observation arcs.

[0053] The Turbedit method is used to perform cycle slip detection on the raw observation data obtained in step 1 to obtain the epoch number and satellite PRN number with continuous observation arcs. The shortest continuous arc is set to 300 epochs to ensure that the subsequent moving average acquisition of detrended total electron quantity DTEC can be carried out normally.

[0054] Step 3: Construct a geometric phase-free observation equation based on the raw GNSS observation data to obtain the total electron content (STEC) of the ionosphere, which includes pseudorange and phase hardware delay bias.

[0055] GNSS signals consist of two parts: pseudorange and carrier phase observations, which can be represented as:

[0056]

[0057] In the formula, the subscript r and superscript s represent the serial numbers of the GNSS receiver and satellite, respectively, and the subscript i represents the frequency. and These are the pseudorange and phase observations along the line-of-sight (LoS) direction between the satellite and the receiver. dt represents the geometric distance between the satellite and the receiver, c is the speed of light in a vacuum, and dt is the distance between the satellite and the receiver. r and dt s These represent the clock biases of the GNSS receiver and the satellite, respectively. For tropospheric delay along the LoS direction, B is the ionospheric delay at frequency i. r,i and Let b represent the pseudorange hardware delay of the GNSS receiver and the satellite at frequency i, respectively. r,i and λ represents the fractional phase deviation of the GNSS receiver and the satellite at frequency i, respectively. i For wavelength, For phase integer ambiguity, and It is a combination of observation noise and various other errors (such as the multipath effect of pseudorange and carrier phase observations).

[0058] Because pseudorange observations are more susceptible to multipath effects and observation noise, the accuracy of ionospheric TEC (Transient Electron Diagnostic and Transient Measurement) decreases. Therefore, the geometric-phase-free combination method can extract ionospheric observations along the LoS (LoS) direction, achieving higher accuracy in ionospheric perturbation analysis than the commonly used phase-smoothing pseudorange method. The geometric-phase-free observation equation is obtained by the difference between carrier phase observations at different frequencies, i.e.:

[0059]

[0060] In the formula, These are observations without geometric phase combination. and These are the phase observations along the LoS direction at frequencies f1 and f2, respectively. For the ionospheric delay along a given LoS direction, γ i =(f0 / f i ) 2 λ is a scaling factor, where f0 = 10.23 MHz, and λ1 and λ2 are the wavelengths of frequencies f1 and f2, respectively. and These are the integer ambiguities at frequencies f1 and f2, respectively. This indicates multipath effects, noise errors, and other types of errors.

[0061] The total electron content of the ionosphere (STEC) is defined as the integral of the electron density along the LoS direction from the satellite to the receiver. The relationship between the STEC at frequency f and the ionospheric delay I in units of length can be described as follows:

[0062]

[0063] In the formula, the STEC value is expressed in TECU, where 1 TECU = 10 16 el / m 2 ; Observations without geometric phase combination; γ i =(f0 / f i ) 2 It is a scaling factor, where f0 = 10.23 MHz; This indicates multipath effects, noise errors, and other types of errors. The cutoff elevation angle of STEC is set to 30° to mitigate the effects of multipath effects.

[0064] Without geometric phase combination, the integer ambiguity of ionospheric delay and phase cannot be separated. However, the integer ambiguity of phase observations remains unchanged in continuous observation arcs. The difference in integer ambiguity between signals of different frequencies is a constant, which can be eliminated in the subsequent calculation of disturbance detection index. This is the reason for calculating cycle slip detection in step 2.

[0065] Step 4: Use the epoch number of the continuous observation arc segment obtained in Step 2 and the satellite PRN number to filter the STEC obtained in Step 3 to obtain the STEC sequence in the same continuous arc segment.

[0066] Step 5: Calculate the ionospheric disturbance index ROT index using STEC sequences in the same continuous arc segment, and then obtain the time series of the ROTI index.

[0067] The ROT exponent is the time derivative of the STEC value over two consecutive time periods, i.e.:

[0068]

[0069] In the formula, the ROT value is used as a parameter for measuring phase fluctuation activity and is expressed in TECU / min; These represent the total ionospheric electron content of satellite k at times i-1 and i, respectively, and t i-1 t i is the epoch time, and k is the satellite sequence number.

[0070] ROTI is the standard deviation of ROT over a specific period, and its calculation method is as follows:

[0071]

[0072] In this embodiment, the specific period is set to 5 minutes, and the angle brackets represent the average value within the 5-minute time interval.

[0073] Step 6: Use the moving average method to detrend STEC to obtain the time series of detrended total electron quantity (DTEC).

[0074] The formula for calculating DTEC using the moving average method is as follows:

[0075]

[0076] In the formula, STEC(i) is the STEC value at time i; N is the sliding window, which is set to 30 minutes in this embodiment; STEC(j) is the STEC value at time j, and the value of j is in the range of [iN / 2, i+N / 2].

[0077] Step 7: Combine the ROTI index obtained in Step 5 and the DTEC obtained in Step 6 to calculate the IROTI index.

[0078] Figure 2 A set of DTEC time series plots using artificially simulated TIDs and EPBs signals is presented. Figure 2 (a) is the analog signal of TIDs, which mainly consists of the wave-like characteristic signal (W1) during the calm period and the TIDs (D1) signal during the disturbance period. Figure 2(b) is the analog signal of EPBs, mainly composed of the wave-like characteristic signal (W2) during the calm period and the EPBs signal (D2) during the disturbance period. Figure 2 (a) It can be seen that for most of the 5-minute time window (D1) during which TIDs occur, the time derivative of DTEC is relatively close. This indicates that the ROT during the TIDs period is basically the same, resulting in a low ROTI level, which is why the traditional ROTI method cannot effectively identify TIDs. Since the DTEC value during the occurrence of both TIDs and EPBs is always greater than the DTEC value during the calm period, this invention calculates the integral of the DTEC value over a specific period as the exponent IDTEC (where DTEC is the time derivative of ... Figure 2 (the gray area in (a)):

[0079]

[0080] In the formula, DTEC(i) is the DTEC value at time i, Δt is the time interval of discrete DTEC, n is the number of DTECs within 5 minutes, and the DTEC value is usually expressed in TECU·min.

[0081] contrast Figure 2 (a) and Figure 2 (b), TIDs and EPBs Figure 2 The IDTEC values ​​in the middle gray area are all at high levels, making it impossible to distinguish between these two types of disturbances. However, from... Figure 3 (b) It can be seen that the DTEC time derivative values ​​differ significantly within the 5-minute time window when EPBs occur, indicating that ROTI is at a high level when EPBs occur. Therefore, combining ROTI with the detection of TIDs and EPBs can further distinguish between the two.

[0082] This invention combines IDTEC and ROTI to propose the IROTI index, which considers changes in amplitude and phase to detect ionospheric disturbances and thus distinguish between TIDs and EPBs. The specific calculation formula can be expressed as follows:

[0083] IROTI = IDTEC·ROTI(8)

[0084] In the formula, the IROTI value is expressed in TECU. 2 Expressed in units.

[0085] When IROTI < 4.91TECU 2 When IROTI > 7.28 TECU, the ionospheric disturbance type can be determined to be TIDs. 2 When this occurs, the disturbance type can be considered as EPBs.

[0086] Figure 3This presents the ROTI, IROTI, and DTEC results calculated by the HKWS station on June 22, 2015. On this day, a large-scale TIDs event caused by a strong geomagnetic storm was observed in the Hong Kong area. Figure 3 As shown in (c), between 22:00UT and 22:30UT, the DTEC value exhibits a wavy structure (TIDs) with varying amplitudes. Generally, when the ROTI value is less than 0.25 TECU / min, the ionosphere can be considered to be in a calm state. However, from... Figure 3 (a) It can be seen that the ROTI value at this station did not change significantly during the TIDs period, and was always less than 0.2 TECU / min, indicating that ROTI cannot be used to detect TIDs. From Figure 4 (b) It can be seen that, compared with ROTI, the IROTI value increases significantly during the occurrence of TIDs and can be clearly distinguished from the ionospheric quiescent period, which indicates that IROTI can be used to detect TIDs.

[0087] Figure 4 This displays the ROTI, IROTI, and DTEC values ​​for TIDs and EPBs observed over METU on July 27, 2012. From... Figure 4 (d) It can be seen that the DTEC value of the G12 satellite exhibited a wave-like characteristic (TIDs) lasting for nearly 1 hour around 11:40. (Comparison) Figure 4 (a) and Figure 4 (b) During the period when TIDs occurred, the IROTI value increased from 0.06 TECU. 2 Increased to 0.65TECU 2 The overall trend of ROTI is less than 0.25 TECU / min, indicating that ROTI cannot detect TIDs. Figure 4 Of the three time points indicated by the arrows on the left, through Figure 4 (d) It can be determined that no significant disturbance occurred at this time, but ROTI suddenly rose sharply to 0.36 TECU / min, while IROTI remained relatively stable, staying below 0.2 TECU. This sudden increase in ROTI value is likely due to local ionospheric instability caused by intense solar activity at midday (LT = 12.4 h). This indicates that IROTI has good anti-interference capabilities in TIDs detection and can effectively avoid the influence of midday ionospheric activity. Figure 4 (e) It can be seen that the DTEC values ​​of all satellites fluctuated drastically to varying degrees around 19:25, and within the corresponding time period Figure 5 The ROTI value in (a) also reached 0.89 TECU / min, and based on the ROTI value, it was determined to be an EPBs. When an EPBs event occurs, Figure 5In (c), the value of IROTI also increased to different degrees at the same time as ROTI, which indicates that IROTI can also be used to detect EPBs.

[0088] Figure 5 The statistical results of ROTI and IROTI during TID events and ionospheric quiescence periods are presented. Solid and hollow dots represent the maximum values ​​of the correlation indices during TID events and quiescence periods, respectively. Figure 5 (a) It can be seen that the maximum value of ROTI during the ionospheric calm period is 0.06–0.31 TECU / min, and the maximum value of ROTI during the disturbed period is 0.08–0.63 TECU / min, which largely overlaps with the "calm" period. Figure 5 (The gray area in (a)). Because most of the ROTI of normal ionospheric changes can also reach this region, TIDs events in the overlapping region cannot be successfully detected using the ROTI index alone. Figure 5 In (a), the vast majority of TIDs events are in overlapping regions, indicating that using ROTI to monitor TIDs has certain limitations. However, as... Figure 6 As shown in (b), the maximum IROTI value during the calm period is between 0.06 and 0.39 TECU. 2 Between these values, the maximum IROTI for detected TIDs events ranged from 0.22 to 7.28 TECU. 2 The overlap between the two is very small. Figure 6 (b) The gray part shows that the IROTI index proposed in this invention has a greater improvement in the identification effect of TIDs compared with the ROTI index.

[0089] ​ The statistical results of ROTI and IROTI for TIDs and EPBs events within the selected time period are presented, where light and dark dots represent the maximum values ​​of the correlation indices with EPBs and TIDs events, respectively. ​ It can be seen that the maximum value of IROTI when EPBs events occur ranges from 4.91 to 72.78 TECU. 2 Within the range, the overall level is higher than that of TIDs events. The gray area represents the overlap between TIDs and EPBs events. TIDs and EPBs events in this overlapping area are difficult to classify using the IROTI index and belong to a fuzzy range. Ionospheric disturbances in this range may involve the simultaneous occurrence of TIDs and EPBs. The upper and lower bounds of this fuzzy range are IROTI = 4.91 TECU. 2 And IROTI = 7.28 TECU 2 This can be considered a threshold range for distinguishing between TIDs and EPBs. When IROTI < 4.91 TECU 2When IROTI > 7.28 TECU, the ionospheric disturbance type can be determined to be TIDs. 2 When this occurs, the disturbance type can be considered as EPBs.

[0090] The specific embodiments described herein are merely illustrative of the spirit of the invention. Those skilled in the art to which this invention pertains may make various modifications or additions to the described specific embodiments or use similar methods to replace them, without departing from the spirit of the invention or exceeding the scope defined by the appended claims.

Claims

1. A detection index for identifying ionospheric traveling wave disturbances and equatorial ionospheric bubbles, characterized in that, The method comprises the following steps: Step 1, obtaining GNSS raw observation files and broadcast ephemeris files of an observation area; Step 2, detecting cycle slips of GNSS raw observation data to obtain ephemeris numbers and satellite PRN numbers of continuous observation arcs; Step 3, constructing a geometry-free phase observation equation based on GNSS raw observation data to obtain ionospheric total electron content (STEC) containing pseudorange and phase hardware delay bias; Step 4, screening the STEC obtained in step 3 by using the ephemeris numbers and satellite PRN numbers of continuous observation arcs obtained in step 2 to obtain STEC sequences in the same continuous arc; Step 5, calculating an ionospheric disturbance index (ROT) index by using the STEC sequences in the same continuous arc to obtain a time sequence of the ROTI index; Step 6, performing detrending processing on the STEC by using a sliding average method to obtain a time sequence of detrended total electron content (DTEC); Step 7, calculating an IROTI index by combining the ROTI index obtained in step 5 and the DTEC obtained in step 6.

2. A detection index for identifying ionospheric traveling wave disturbances and equatorial ionospheric bubbles as claimed in claim 1, characterized in that: The original observation data obtained in step 1 is subjected to cycle slip detection by using the Turbedit method in step 2 to obtain the ephemeris number and satellite PRN number with continuous observation arc segments, and the shortest continuous arc segment is set as The ephemeris is used to ensure that the subsequent sliding average DTEC can be normally obtained.

3. The detection index for identifying ionospheric traveling wave disturbances and equatorial plasma bubbles as described in claim 1, characterized in that: The GNSS signal in step 3 is composed of two parts of pseudorange and carrier phase observation value, and can be expressed as: where the subscripts r and superscripts s denote the serial numbers of the GNSS receiver and the satellite, respectively, the subscript i denotes the frequency, and are the pseudorange and phase observations along the line-of-sight (LoS) direction between the satellite and the receiver, denotes the geometric distance between the satellite and the receiver, c is the speed of light in vacuum, and denote the clock biases of the GNSS receiver and the satellite, respectively, is the tropospheric delay along the LoS direction, is the ionospheric delay at the frequency i , and denote the pseudorange hardware delays of the GNSS receiver and the satellite at the frequency i , and denote the phase fractional biases of the GNSS receiver and the satellite at the frequency i , is the wavelength, is the phase integer ambiguity, and is the combined term of the observation noise and other various errors.

4. The detection index for identifying ionospheric traveling wave disturbances and equatorial plasma bubbles as described in claim 3, characterized in that: The geometry-free phase observation equation in step 3 is obtained by the difference value of carrier phase observation values of different frequencies, that is: wherein is the geometry-free phase combination observation, and are the f 1, f 2 phase observations at the given LoS direction at frequency is the ionospheric delay along the given LoS direction, is a scale factor, where Mhz, denotes the frequency, i takes the value 1 or 2, and are the f 1, f 2 wavelengths at the frequency and are the f 1, f 2 integer ambiguities at the frequency denotes multipath effects, noise errors and other types of errors.

5. The detection index for identifying ionospheric traveling wave disturbances and equatorial plasma bubbles as described in claim 4, characterized in that: The ionospheric total electron content (STEC) is defined in step 3 as the integral of the electron density along the LoS direction from the satellite to the receiver, the frequency being f STEC and ionospheric delay I The relationship in units of length can be described as: wherein the STEC value is expressed in TECU, 1 TECU = 10 16 el / m 2 ; is the geometry-free phase combination observation; is a scale factor, wherein Mhz; denotes multipath effects, noise errors and other types of errors; the cut-off height angle of the STEC is set to in order to weaken the influence of multipath effects.

6. A detection index for identifying ionospheric traveling wave disturbances and equatorial ionospheric bubbles as claimed in claim 1, wherein: ###0002### The ROT index in step 5 is the time derivative of the STEC values in two continuous time periods, that is: In the formula, the ROT value is a parameter for measuring phase fluctuation activity, expressed in TECU / min; , are the total electron content of the ionosphere at the time of the satellite k pass i -1, i moment, , is the epoch time, k is the satellite serial number; The ROTI is the standard deviation of ROT in a specific period, and the specific calculation method is as follows: In the formula, the angle brackets represent the average value in a specific period.

7. The detection index for identifying ionospheric traveling wave disturbances and equatorial plasma bubbles as described in claim 1, characterized in that: The formula for calculating the detrended total electron content (DTEC) in step 6 by using the sliding average method is as follows: In the formula, is the STEC value at the time of the first i is the STEC value at the time of the first is the STEC value at the time of the first j is the STEC value at the time of the first j is the STEC value at the time of the first is the STEC value at the time of the first 8. The detection index for identifying ionospheric traveling wave disturbances and equatorial plasma bubbles as described in claim 1, characterized in that: In step 7, the integral of the DTEC value in a specific period is calculated as an index IDTEC for judging the occurrence of disturbance: wherein is i the absolute value of the DTEC at the moment, is the time interval of the discrete DTEC, n is the number of DTEC in a certain period, the IDTEC value is expressed in TECU·min.

9. A detection index for identifying ionospheric traveling wave disturbances and equatorial ionospheric bubbles as claimed in claim 7, characterized in that: In step 7, the IROTI index is proposed by combining IDTEC and ROTI, which considers the changes of amplitude and phase to detect ionospheric disturbances and then distinguish TIDs and EPBs, and the specific calculation formula can be expressed as: In the formula, the IROTI value is expressed in TECU 2 units; When IROTI < 0 TECU 2 , it can be judged that the ionospheric disturbance type is TIDs, and when IROTI > 0 TECU 2 , it can be considered that the disturbance type is EPBs, , is a constant, and the specific value is set according to actual needs.

Images