Hypersonic target positioning method based on GNSS ionosphere transient disturbance

Abstract

The invention belongs to the technical field of satellite navigation and ionosphere monitoring, relates to a hypersonic target positioning technology, and particularly provides a hypersonic target positioning method based on GNSS ionosphere transient disturbance, which is used for meeting the reliability requirement and the real-time requirement of hypersonic speed detection. According to the hypersonic target positioning method, firstly, an inclined path TEC (STEC) is calculated through a pseudo-range difference and a carrier phase difference, then the STEC is projected to an ionospheric puncture point (IPP) to obtain a vertical TEC (VTEC), and finally, the detected VTEC change rate is analyzed to obtain track information of a hypersonic target, so that hypersonic target positioning is realized. According to the invention, 24-hour continuous detection of a hypersonic target in an ionized layer can be realized by using a global distribution GNSS receiving station network, the real-time performance is relatively high, the detection precision can be effectively improved, and detection errors caused by the influence of a plasma sheath in a traditional detection mode are also avoided.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of satellite navigation and ionosphere monitoring, and relates to hypersonic target positioning technology, and particularly provides a hypersonic target positioning method based on GNSS ionosphere transient disturbance. BACKGROUND

[0002] The existing space-based infrared early warning technology captures the thermal barrier effect of the thermal signal generated by the hypersonic target due to atmospheric friction through a high-orbit satellite, and the process is divided into two stages: first, the infrared detector on the satellite captures the super-high temperature thermal signal generated by the hypersonic target due to atmospheric friction; then, the target coordinate position is calculated through time difference technology, and the trajectory data is transmitted to the ground base station to trigger the ground-based radar relay tracking. However, the existing technology has two major defects: first, the plasma sheath formed on the surface of the hypersonic target when it passes through the atmosphere will absorb or scatter radar waves, resulting in a high failure probability of ground-based radar relay, and the overall early warning link is prone to breakage; second, the space-based infrared early warning needs to go through multiple links such as "space-based scanning-positioning-ground-based radar relay" cooperation, and the whole process takes a long time. When the ground-based radar receives the plasma sheath interference, the effective early warning window of the system for hypersonic weapons will be further compressed, and the real-time requirement cannot be met. In summary, the detection reliability and real-time performance of the existing technology do not meet the requirements of hypersonic detection, and an innovative technical solution that is resistant to plasma interference and has strong real-time performance is urgently needed. SUMMARY

[0003] The purpose of the present application is to provide a hypersonic target positioning method based on GNSS ionosphere transient disturbance, which meets the reliability and real-time requirements of hypersonic detection. First, the slant path TEC (STEC) is calculated by using the pseudo-range difference and carrier phase difference, then the STEC is projected to the ionospheric piercing point (IPP) to obtain the vertical TEC (VTEC), and finally the trajectory information of the hypersonic target is obtained by analyzing the detected VTEC rate, thereby realizing the positioning of the hypersonic target.

[0004] To achieve the above purpose, the technical solution adopted by the present application is as follows:

[0005] A hypersonic target positioning method based on GNSS ionosphere transient disturbance, characterized in that it comprises the following steps:

[0006] Step 1. Use a GNSS receiver or a fixed receiving base station to detect the GNSS satellite signal and obtain an observation file;

[0007] Step 2. Perform difference calculation based on dual-frequency observation values to obtain the slant path total electron content of the carrier phase observation and the slant path total electron content of the pseudo-range observation, respectively, and simultaneously solve the position information of the ionospheric piercing point according to the coordinate position information of the receiver.

[0008] Step 3. Cycle slip repair and smoothing processing of the slant path total electron content of the carrier phase observation and the pseudorange observation;

[0009] Step 4. Smoothing the slant path total electron content of the pseudorange observation using the slant path total electron content of the carrier phase observation to obtain the absolute value of the ionospheric slant path electron total content;

[0010] Step 5. Projecting the absolute value of the ionospheric slant path electron total content at the piercing point to the ionospheric vertical direction electron total content at the piercing point;

[0011] Step 6. Calculating the change gradient of the ionospheric vertical direction electron total content at the piercing point, and when the change gradient exceeds a preset threshold, determining that the piercing point is an abnormal piercing point, and taking the abnormal piercing point as a hypersonic target trajectory point.

[0012] Further, in step 1, the observation file contains: observable satellite number, pseudorange observation value, carrier phase observation value, receiver position information, satellite position information.

[0013] Further, in step 2, the slant path total electron content of the carrier phase observation and the slant path total electron content of the pseudorange observation are respectively:

[0014]

[0015] Wherein, the subscripts 1 and 2 respectively represent the frequency labels of the 1st frequency and the 2nd frequency in the differential calculation; TEC p represents the slant path total electron content of the pseudorange observation, TEC φ represents the slant path total electron content of the carrier phase observation, P represents the pseudorange measurement value, f is the frequency of the incident GNSS signal, and λ represents the carrier wavelength, represents the carrier phase measurement value, and N represents the integer ambiguity in the carrier phase measurement process.

[0016] Further, in step 2, the position information of the ionospheric piercing point is:

[0017] φ IPP = arcsin(sinφ r cosψ IPP + cosφ r sinψ IPP cosA)

[0018]

[0019] Wherein, φ IPP , λ IPP are the geographic latitude and geographic longitude of the ionospheric piercing point, respectively, φ r , λr respectively, latitude and longitude coordinates of the receiver;

[0020] A is the azimuth angle from the station to the satellite, specifically:

[0021]

[0022] wherein x r,O is the distance vector from the receiver to the center of the earth O; x r,s is the distance vector from the receiver to the satellite;

[0023] ψ IPP is the central angle between the station point and the ionospheric piercing point, specifically:

[0024]

[0025] wherein R is the average radius of the earth, E is the elevation angle of the satellite, and H is the height between the piercing point and the ground surface.

[0026] Further, in step 3, the specific process of cycle slip repair is:

[0027] The total electron content of the oblique path is differentiated, and the standard deviation after differentiation is obtained. When the standard deviation exceeds the preset threshold, it is determined that there is a cycle slip.

[0028] The total electron content sequence of the oblique path is segmented according to the cycle slip point. The deviation between the end value of the previous segment and the starting value of the next segment is used as the overall offset, and the previous segment data is translated to make it continuous with the next segment data, specifically represented as:

[0029] TEC_part = bound1-(bound1(end)-bound2(1))

[0030] Wherein TEC_part represents the result after translation of the previous segment data, bound1 represents the previous segment data, bound1(end) represents the end value of the previous segment data, and bound2(1) represents the starting value of the next segment data.

[0031] After the cycle slip repair is completed, the difference between the total electron content of the oblique path under the highest elevation angle before and after the repair is used as the basis to correct the total electron content sequence of the oblique path after the repair, thereby realizing smoothing processing.

[0032] Further, in step 4, the absolute value of the ionospheric oblique path total electron content is:

[0033]

[0034] Wherein, STEC represents absolute value of total electron content of ionosphere oblique path, t represents time of epoch, N represents number of time of epoch in smoothing process, TEC p represents total electron content of oblique path of pseudo-range observation, TEC P,t represents total electron content of oblique path of carrier phase observation. φ,t represents total electron content of oblique path of pseudo-range observation and total electron content of oblique path of carrier phase observation at time of epoch t.

[0035] Further, in step 5, total electron content of ionosphere vertical direction at puncture point is:

[0036]

[0037] Wherein, VTEC represents total electron content of ionosphere vertical direction, STEC represents absolute value of total electron content of ionosphere oblique path, R represents average radius of earth, E represents elevation angle of satellite, H represents height between puncture point and earth surface.

[0038] Further, step 6 is specifically:

[0039] Calculate residual error ΔTEC(t) of total electron content of ionosphere vertical direction at puncture point and total electron content of reference ionosphere model:

[0040] ΔTEC(t) = VTEC(t) - TEC(t) model

[0041] Wherein, TEC(t) model represents total electron content of reference ionosphere model;

[0042] When ΔTEC(t) is greater than preset threshold value, determine as abnormal event caused by hypersonic target passing, and take relative longitude and latitude of puncture point as trajectory coordinates of hypersonic target passing.

[0043] Based on the above technical scheme, the present application has the advantages that a hypersonic target positioning method based on GNSS ionosphere transient disturbance is provided, and has the following advantages:

[0044] 1) The present application utilizes global distribution GNSS receiving station network, and can realize 24-hour continuous detection on hypersonic target in ionosphere;

[0045] 2) The present application can quickly and effectively identify ionosphere disturbance excited by hypersonic target by real-time inversion of TEC and correlation of synchronous puncture point anomaly, so as to lock position of hypersonic target, and has strong real-time performance;

[0046] 3) This invention utilizes the influence of the plasma sheath formed by the hypersonic target in the ionosphere for position detection, avoiding detection errors caused by the influence of the plasma sheath in traditional detection methods.

[0047] 4) This invention integrates observations from multiple systems such as GPS and BeiDou in the GNSS system, which can effectively improve detection accuracy. Attached Figure Description

[0048] Figure 1 This is a flowchart illustrating the hypersonic target localization method based on GNSS ionospheric transient disturbances in this invention.

[0049] Figure 2 This is a schematic diagram of the ionospheric puncture point in this invention. Detailed Implementation

[0050] To make the objectives, technical solutions, and beneficial effects of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments.

[0051] This embodiment provides a hypersonic target localization method based on GNSS ionospheric transient disturbances, the process of which is as follows: Figure 1 As shown, the specific steps include:

[0052] Step 1. Use a GNSS receiver or fixed receiving base station to detect GNSS satellite signals and obtain a standardized RINEX format observation file containing the detection results, including: observable satellite number, pseudorange observation value, carrier phase observation value, receiver position information, satellite position information, etc.

[0053] Step 2. Perform differential calculations based on dual-frequency observations to remove noise interference such as relativistic effect errors and multipath errors, and obtain the total electron content (STEC) of the slant path calculated by carrier phase and pseudorange respectively; and calculate the latitude and longitude information of the GNSS satellite signal at the ionospheric puncture point based on the coordinate position information of the receiver, that is, simultaneously calculate the position information of the ionospheric puncture point (IPP);

[0054] The equations for calculating pseudorange and carrier phase observations for dual-frequency receivers (such as GPS L1 / L2 and BDS B1 / B2) are as follows:

[0055] P k =ρ+c(δt) s -δt r )+T+I k,p +ε pk

[0056]

[0057] Wherein: subscript k represents the mark of carrier signal of different frequency; represents carrier phase measurement value; P k represents pseudo-range measurement value; λ k represents carrier wavelength; ρ represents geometric distance between receiver and satellite; c represents light speed; N k represents integer ambiguity in carrier phase measurement process, which is a constant; δt s and δt r respectively represent clock error of satellite and receiver; I k,p and respectively represent ionospheric delay error of pseudo-range and ionospheric delay error of carrier phase corresponding to GNSS carrier signal of different frequency; T represents tropospheric delay error; ε pk and respectively represent noise error in pseudo-range measurement value and carrier phase measurement value;

[0058] In the above formula, the difference of carrier signal frequency only affects ionospheric delay error I k , and other errors can be eliminated by double-frequency difference, epoch smoothing and the like; through difference calculation, the combination observation equation of double-frequency GNSS pseudo-range and carrier phase group is respectively obtained as follows:

[0059] P1-P2=(I1-I2)+ε p1 -ε p2

[0060]

[0061] At this time, the difference equation only contains ionospheric delay error, integer ambiguity and noise error, wherein, the integer ambiguity can be removed by subsequent carrier phase smoothing pseudo-range, and the noise error can be ignored;

[0062] Therefore, the ionospheric delay error is solved, and specifically:

[0063] Using Appleton Hartree formula, the refractive index in ionosphere can be approximately obtained as follows:

[0064]

[0065] Wherein, N e is total electron content on propagation path, and f is frequency of incident GNSS signal;

[0066] Because ionosphere is a dispersive medium, the propagation speed of each frequency component exists difference when GNSS signal propagates in the ionosphere, which leads to difference of propagation speed of electromagnetic wave packet energy (group velocity v g) less than light speed c, the ionospheric delay of pseudorange of GNSS signal is positive; while the propagation speed (phase velocity v p ) greater than light speed c, the ionospheric delay of carrier phase is negative; thus, the phase velocity v p and group velocity v g of GNSS signal with carrier frequency f are respectively:

[0067]

[0068] From the above formula, on the fixed signal propagation path, the ionospheric delay of pseudorange and phase observation of GNSS signal is affected by the frequency f and the total electron content N e on the propagation path; then the ionospheric delay correction value of pseudorange is:

[0069]

[0070] Where, TEC p represents the total ionospheric electron content measured by pseudorange;

[0071] The ionospheric delay correction value of carrier phase is:

[0072]

[0073] Where, TEC φ represents the total ionospheric electron content measured by carrier phase;

[0074] Finally, the obtained ionospheric delay correction value is brought into the combined observation equation of dual-frequency GNSS pseudorange and carrier phase group respectively:

[0075]

[0076] The total electron content of slant path obtained by pseudorange and carrier phase observation value respectively is:

[0077]

[0078] Ionospheric piercing point (IPP) position information solving:

[0079] The longitude and latitude of piercing point are solved by the following formula:

[0080] φ IPP = arcsin (sinφ r cosψ IPP + cosφ r sinψ IPP cosA)

[0081]

[0082] Where, φ IPP , λ IPP These are the geographical latitude and longitude of the ionospheric puncture point, respectively, φ r and λ r Let A be the latitude and longitude coordinates of the known receiving station, and let ψ be the azimuth angle from the station to the satellite. IPP The geocentric angle between the measuring station and the ionospheric puncture point (IPP);

[0083] Ionospheric puncture point such as Figure 2 As shown, from geometric relationships, we can obtain:

[0084]

[0085] Where, x r,O x is the distance vector from the receiver to the Earth's center O; r,s This is the distance vector from the receiver to the satellite;

[0086] The azimuth angle from the station to the satellite can be obtained:

[0087]

[0088] The geocentric angle ψ between the measuring station and IPP IPP :

[0089]

[0090] Where Z is the angle formed by the distance from the receiver to the puncture point and the distance from the Earth's center to the puncture point;

[0091] From the Law of Sines, we can obtain:

[0092]

[0093] Where R is the average radius of the Earth, E is the elevation angle of the satellite, and H is the height between the puncture point and the Earth's surface;

[0094] Therefore, we can conclude that:

[0095]

[0096] Step 3. Perform cycle slip repair and smoothing on the total electron content of the slant path obtained from carrier phase and pseudorange due to signal interruption and other reasons to reduce the error caused by satellite signal interruption;

[0097] The slant path electron total content calculated by carrier phase is differentiated, and the standard deviation after differentiation is calculated to express the "average distance" of the deviation from the average value. When the standard deviation exceeds the pre-set threshold value, it can be determined that there is a cycle slip. Then the STEC sequence is segmented according to the cycle slip point, and the deviation of the end value (bound1 (end)) of the previous segment and the start value (bound2 (1)) of the next segment is taken as the overall offset, and the previous segment data is translated to make it continuous with the next segment data:

[0098] TEC_part = bound1 - (bound1 (end) - bound2 (1))

[0099] Wherein, TEC_part represents the result after the previous segment data is translated;

[0100] After the splicing corrected segmented data is eliminated, the difference between the original STEC before and after the repair is taken as the basis to translate the repaired STEC sequence as a whole, so as to ensure that the repair time is consistent with the original value;

[0101] Step 4. The STEC obtained by the smooth pseudo-range observation of the carrier phase observation with high observation accuracy is used to remove the integer ambiguity in the carrier phase observation value, and the ionospheric absolute electron total content is obtained;

[0102] The pseudo-range observation accuracy is not high, and the obtained TEC p The accuracy is about 1 TECU; the observation accuracy of the carrier phase is much higher than that of the pseudo-range observation, and the obtained TEC φ The accuracy can reach 0.01 TECU, but the calculation value contains integer ambiguity, and only the relative TEC is obtained; therefore, in order to obtain high-precision TEC measurement value, the advantages of the two observation methods should be complementary. The TEC φ The TEC p obtained by the smooth pseudo-range observation, the TEC obtained by the pseudo-range and carrier phase at the time of epoch t is subtracted, and the mean value difference is obtained by N epochs, and then the initial TEC calculated by the differential pseudo-range is added, so as to obtain the ionospheric absolute TEC value:

[0103]

[0104] Step 5. The ionospheric slant path electron total content STEC at the piercing point is converted into the ionospheric vertical direction electron total content VTEC at the piercing point by using the projection function, so as to more accurately reflect the TEC change in the horizontal direction of the ionosphere;

[0105] The ionospheric vertical direction TEC at the piercing point is obtained from the geometric relationship shown in Figure 2

[0106] ​VTEC = STEC x cosZ

[0107] That is:

[0108]

[0109] Step 6. The change gradient of the ionospheric total electron content is obtained using the calculated ionospheric vertical path total electron content VTEC at the piercing point, and when the threshold is exceeded, the abnormal disturbance caused by the hypersonic target is determined, and the abnormal piercing point is the trajectory point of the hypersonic target;

[0110] The calculated ionospheric vertical TEC at the piercing point is one-to-one corresponding to the longitude and latitude coordinates at the piercing point, and the residual of the model TEC calculated using the International Reference Ionosphere (IRI) model or other regional ionospheric maps is obtained:

[0111] ΔTEC(t) = VTEC(t) - TEC(t) model

[0112] Wherein, TEC(t) model represents the total electron content of the reference ionospheric model;

[0113] When ΔTEC(t) is greater than the pre-set threshold σ, it is determined that the abnormal event caused by the hypersonic target passing through, and the longitude and latitude of the piercing point corresponding to the abnormal point are the trajectory coordinates of the hypersonic target passing through.

[0114] The beneficial effects of the present application will be described in detail in combination with simulation tests.

[0115] The observation files of the three observation stations of South Korea OSN3, SEJN and SUWN on February 7, 2016 were used to detect the trajectory information of the Hwansong-4 launched by North Korea on February 7, 2016; the average radius of the earth was set to 6371.009 x 10 3 m, the ionospheric height was 350 x 10 3 m, the speed of light was 299792458 m / s, and the threshold for cycle slip repair was 0.5.

[0116] The TEC at the ionospheric piercing point was inversed using the GPS satellite signals received on the same day, wherein the STEC obtained by the pseudo-range was calculated using the C1 and P2 band signals, and the STEC obtained by the carrier phase was calculated using the L1 and L2 band signals.

[0117] After the cycle slip repair and smoothing processing, the obtained carrier phase smoothed pseudo-range is used, and finally the projection function is used to obtain the final VTEC; the residual obtained by subtracting the model TEC calculated by using the international reference ionosphere (IRI) model from the set threshold 0.5TECu / s is compared, and the abnormal point positions of (9.15°N, 89.85°E) and (24.54°N, 120°E) are obtained, which are judged as the position information of the satellite passing through the ionosphere, which is consistent with the true trajectory information.

[0118] The above is only a specific embodiment of the present application, any feature disclosed in the specification can be replaced by other equivalent or similar purpose alternative features unless specifically described; all features disclosed, or steps in all methods or processes, except for mutually exclusive features and / or steps, can be combined in any way.

Claims

1. A method for hypersonic target positioning based on GNSS ionospheric transient disturbances, characterized in that, The method comprises the following steps: Step 1: GNSS satellite signals are detected by using a GNSS receiver or a fixed receiving base station to obtain an observation file; Step 2: difference calculation is performed based on double-frequency observation values to obtain slant-path total electron content of carrier phase observation and slant-path total electron content of pseudo-range observation respectively, and position information of an ionospheric piercing point is obtained by synchronous solution according to coordinate position information of the receiver; Step 3: cycle slip repair and smoothing processing are performed on the slant-path total electron content of carrier phase observation and pseudo-range observation; Step 4: the slant-path total electron content of carrier phase observation is used to smooth the slant-path total electron content of pseudo-range observation to obtain absolute value of ionospheric slant-path electron total content; Step 5: the absolute value of ionospheric slant-path electron total content at the piercing point is projected into ionospheric vertical direction electron total content at the piercing point; Step 6: a change gradient of the ionospheric vertical direction electron total content at the piercing point is calculated, and when the change gradient exceeds a preset threshold, the piercing point is determined as an abnormal piercing point, and the abnormal piercing point is taken as a hypersonic target trajectory point.

2. The method of claim 1, wherein, In step 1, the observation file comprises: observable satellite number, pseudo-range observation value, carrier phase observation value, receiver position information and satellite position information.

3. The method of claim 1, wherein, In step 2, the slant-path total electron content of carrier phase observation and the slant-path total electron content of pseudo-range observation are respectively: where the subscripts 1 and 2 represent the first and second frequency labels in the differential calculation; TEC p represents the slant-path total electron content (TEC) for the pseudorange observation φ represents the slant-path total electron content (TEC) for the carrier phase observation, P represents the pseudorange measurement, f is the frequency of the incoming GNSS signal, and λ represents the carrier wavelength, represents the carrier phase measurement, and N represents the integer ambiguity in the carrier phase measurement process.

4. The method of claim 1, wherein, In step 2, the position information of the ionospheric piercing point is: φ IPP = arcsin(sinφ r cosψ IPP + cosφ er sinψ IPP cosA) where φ IPP , λ IPP are the geographic latitude and longitude of the ionospheric pierce point, respectively, and φ r , λ r are the latitude and longitude coordinates of the receiver, respectively. A is an azimuth angle from the station to the satellite, and specifically is: where x r,O is the receiver-to-geocenter distance vector; x r,s is the receiver-to-satellite distance vector; Ψ IPP The central angle between the station point and the ionospheric piercing point, specifically: Wherein, R is the average radius of the earth, E is the elevation angle of the satellite, and H is the height between the piercing point and the ground.

5. The method of claim 1, wherein, In step 3, the specific process of cycle slip repair is: The slant-path electron total content is differentiated to obtain a standard deviation after differentiation, and when the standard deviation exceeds a preset threshold, it is determined that there is cycle slip; The slant-path electron total content sequence is segmented according to the cycle slip point, and the deviation of the end value of the previous segment and the start value of the next segment is taken as the overall offset, and the previous segment data is translated to make it continuous with the next segment data, and specifically represented as: TEC_part = bound1-(bound1(end)-bound2(1)) Wherein, TEC_part represents the result after translation of the previous segment data, bound1 represents the previous segment data, bound1(end) represents the end value of the previous segment data, and bound2(1) represents the start value of the next segment data; After the cycle slip repair is completed, the difference between the slant-path electron total content under the highest elevation angle observation before and after the repair is taken as the benchmark to perform overall translation correction on the slant-path electron total content sequence after the repair, so as to realize smoothing processing.

6. The method of claim 1, wherein, In step 4, the absolute value of ionospheric slant-path electron total content is: Wherein, STEC represents the absolute value of the ionosphere slant path total electron content, t represents the time of the epoch, N represents the number of epoch time in the smoothing process, TEC p represents the slant path total electron content of the pseudo-range observation P,t represents the slant path total electron content of the pseudo-range observation φ,t represents the slant path total electron content of the pseudo-range observation and the slant path total electron content of the carrier phase observation at the time of the epoch t.

7. The method of claim 1, wherein, In step 5, the ionospheric vertical direction electron total content at the piercing point is: Wherein, VTEC represents the ionospheric vertical direction electron total content, STEC represents the absolute value of ionospheric slant-path electron total content, R represents the average radius of the earth, E represents the elevation angle of the satellite, and H represents the height between the piercing point and the ground.

8. The method of claim 1, wherein, Step 6 is specifically: The residual error ΔTEC(t) of the ionospheric vertical direction electron total content at the piercing point and the electron total content of a reference ionospheric model is calculated: ATEC(t) = VTEC(t) - TEC(t) model TEC(t) = TEC(t - 1) + K(t) + D(t) + E(t) + N(t) (1) model represents the total electron content of the reference ionospheric When the ΔTEC(t) is greater than the preset threshold value, the abnormal event caused by the hypersonic target passing is determined, and the longitude and latitude of the relative piercing point are taken as the track coordinates of the hypersonic target passing.

Images