阅读说明：本技术 三频信标机数据处理方法及系统 (Data processing method and system for three-frequency beacon machine ) 是由 鲁恒新 于 2021-07-30 设计创作，主要内容包括：本发明提供了一种三频信标机数据处理方法及系统。该方法包括：三频信标机0级数据进行差分相位计算和信号强度计算,生成1级数据产品；生成2级数据产品；获得3级数据产品。本发明提供的三频信标机数据处理方法及系统能够将三频信标机的原始探测数据进行处理,生成实际可用的数据产品。(The invention provides a data processing method and a data processing system of a three-frequency beacon machine. The method comprises the following steps: performing differential phase calculation and signal intensity calculation on the 0-level data of the three-frequency beacon machine to generate a 1-level data product; generating a 2-level data product; and obtaining a 3-level data product. The method and the system for processing the data of the three-frequency beacon machine can process the original detection data of the three-frequency beacon machine to generate an actually usable data product.)

1. A data processing method of a three-frequency beacon is characterized by comprising the following steps:

after the ground three-frequency beacon receiver tracks and locks the satellite-borne three-frequency beacon signal, orthogonal component Q and homodromous component I observation data and signal-to-noise ratio data of the three frequency band signals are output;

differential phase calculation and signal intensity calculation are respectively carried out on the I/Q observation data of the three-frequency beacon receiver, and the differential phase values of VHF, UHF and L and the signal intensity values of the three frequency bands can be calculated;

carrying out phase connection processing on the calculated differential phase value, calculating a group of TEC values of the total electron content of the relative ionized layer every second, calculating an amplitude flicker index of the calculated signal intensity value, calculating a group of flicker indexes every 1 second, carrying out space-time matching, and calculating to obtain a space position corresponding to each group of data at observation time;

by carrying out centralized processing on observation data of a plurality of stations, calculating absolute TEC data of an ionized layer by a multi-station method every second, then carrying out inversion on a single link process by a tomography algorithm to obtain two-dimensional electron density profile data, and extracting the NmF2 anomaly.

2. The method of claim 1 in which the observation data and signal-to-noise ratio data for the quadrature component Q and the common component I of the three band signals are referred to as level 0 data;

the differential phase values of VHF & UHF, UHF & L and the respective signal intensity values of the three frequency bands are referred to as level 1 data;

the spatial position corresponding to each group of data at the observation time is called 2-level data;

two-dimensional electron density profile data obtained by inverting a single link process by using a tomography algorithm and extracted NmF anomalies are called 3-level data.

3. The method of claim 2 in which the level 0 data comprises: the original observation data with time information is generated by the ground station receiver;

the level 1 data includes: carrying out orthogonal solution on the 0-level data of each single station to obtain phase data and amplitude data which are arranged according to time;

the level 2 data includes: carrying out phase connection and scintillation index calculation on the 1-level data of each single station to obtain data with track information;

the level 3 data includes: and on the basis of the 2-level data, inverting the multi-station data to generate an absolute TEC and a two-dimensional profile product of an observation link.

4. The method of claim 2 in which the level 0 data comprises: a signal-to-noise ratio data file;

the level 1 data includes: scientific data, image products and data processing reports;

the level 2 data includes: scientific data, image products and data processing reports;

the level 3 data includes: scientific data and its images and data processing reports.

5. The method for processing the data of the tri-band beacon set of claim 1, wherein the steps of calculating the differential phase and the signal strength of the I/Q observation data of the tri-band beacon receiver to obtain the differential phase values of VHF & UHF and UHF & L and the signal strength values of the three bands comprise: differential phase calculation, and signal strength calculation.

6. The method of claim 5 in which the differential phase calculation comprises:

reading the data of the I path and the Q path of the 0-level VHF frequency band;

calculating a VHF/UHF differential phase;

reading the data of the I path and the Q path of the 0-level L frequency band;

the L/UHF differential phase is calculated.

7. The method of claim 1 in which the signal strength calculation comprises:

reading I path data and Q path data of 0-level VHF, UHF and L frequency bands;

and calculating the signal intensity of the VHF, UHF and L frequency bands.

8. The method of claim 1 in which the phase joining of the calculated differential phase values, calculating a set of TEC values for total electron content of the relative ionosphere per second, calculating a flicker index of amplitude for the calculated signal intensity values, calculating a set of flicker indices every 1 second, performing space-time matching, and calculating the spatial location of each set of data at the observation time comprises: relative TEC calculation, and scintillation index calculation.

9. The method of claim 1 in which the centralized processing of observation data from a plurality of stations is used to compute an ionosphere absolute TEC data per second using a multi-station method, and then a tomographic algorithm is used to invert a single link process to obtain two-dimensional electron density profile data and extract NmF2 anomalies, including:

acquiring relative TEC data of all receiving stations corresponding to the station chain, calculating an integral constant of each receiving station by using a multi-station method, and adding the integral constant to the relative TEC obtained by the receiver to obtain an absolute TEC of each receiving station;

dividing the range of 100 km-500 km above the station chain at intervals according to the height interval of 20km and the latitude interval of 0.5 degrees, and calculating the electron density profile of the divided region through an inversion process.

10. A three-frequency beacon data processing system, comprising:

one or more processors;

a storage device for storing one or more programs,

when executed by said one or more processors, cause said one or more processors to implement a method of three-frequency beacon data processing according to any one of claims 1 to 9.

Technical Field

The invention relates to the technical field of satellite remote sensing, in particular to a data processing method and system of a three-frequency beacon machine.

Background

The Zhangheng-I satellite is the first autonomously developed geophysical field satellite in China, and eight loads are carried on the satellite. The Zhanghenyi satellite is successfully transmitted in 2018, 2.2.8.A global navigation satellite system radio occultation (GRO) receiver carried by the satellite has the occultation observation function of a Global Positioning System (GPS) satellite and a Beidou navigation system (BD) satellite.

Disclosure of Invention

The invention aims to provide a data processing method and a data processing system of a three-frequency beacon machine, which can process original detection data of the three-frequency beacon machine to generate an actually usable data product.

In order to solve the technical problem, the invention provides a data processing method of a three-frequency beacon machine, which comprises the following steps: after the ground three-frequency beacon receiver tracks and locks the satellite-borne three-frequency beacon signal, orthogonal component Q and homodromous component I observation data and signal-to-noise ratio data of the three frequency band signals are output; differential phase calculation and signal intensity calculation are respectively carried out on the I/Q observation data of the three-frequency beacon receiver, and the differential phase values of VHF, UHF and L and the signal intensity values of the three frequency bands can be calculated; carrying out phase connection processing on the calculated differential phase value, calculating a group of TEC values of the total electron content of the relative ionized layer every second, calculating an amplitude flicker index of the calculated signal intensity value, calculating a group of flicker indexes every 1 second, carrying out space-time matching, and calculating to obtain a space position corresponding to each group of data at observation time; by carrying out centralized processing on observation data of a plurality of stations, calculating absolute TEC data of an ionized layer by a multi-station method every second, then carrying out inversion on a single link process by a tomography algorithm to obtain two-dimensional electron density profile data, and extracting the NmF2 anomaly.

In some embodiments, the quadrature component Q and the co-directional component I observation data and signal-to-noise ratio data of the three frequency band signals are referred to as 0-level data; the differential phase values of VHF & UHF, UHF & L and the respective signal intensity values of the three frequency bands are referred to as level 1 data; carrying out phase connection processing on the calculated differential phase value, wherein the calculated relative ionosphere total electron content value (relative TEC) and a scintillation index obtained by carrying out amplitude scintillation index calculation according to the signal intensity value are called as 2-level data; two-dimensional electron density profile data obtained by inverting a single link process by using a tomography algorithm and extracted NmF anomalies are called 3-level data.

In some embodiments, the level 0 data comprises: the original observation data with time information is generated by the ground station receiver; the level 1 data includes: carrying out orthogonal solution on the 0-level data of each single station to obtain phase data and amplitude data which are arranged according to time; the level 2 data includes: carrying out phase connection and scintillation index calculation on the 1-level data of each single station to obtain data with track information; the level 3 data includes: and on the basis of the 2-level data, inverting the multi-station data to generate an absolute TEC and a two-dimensional profile product of an observation link.

In some embodiments, the level 0 data comprises: I/Q observation data and signal-to-noise ratio data files of three frequency bands of VHF, UHF and L; the level 1 data includes: scientific data, image products and data processing reports; the level 2 data includes: scientific data, image products and data processing reports; the level 3 data includes: scientific data and its images and data processing reports.

In some embodiments, the calculating of the differential phase and the signal strength of the I/Q observation data of the triple-band beacon receiver may be performed separately, and the calculating of the differential phase values of VHF & UHF and UHF & L and the signal strength values of the three bands may be performed, including: differential phase calculation, and signal strength calculation.

In some embodiments, the differential phase calculation comprises: reading the data of the I path and the Q path of the 0-level VHF frequency band; calculating a VHF/UHF differential phase; reading the data of the I path and the Q path of the 0-level L frequency band; the L/UHF differential phase is calculated.

In some embodiments, the signal strength calculation comprises: reading I path data and Q path data of 0-level VHF, UHF and L frequency bands; and calculating the signal intensity of the VHF, UHF and L frequency bands.

In some embodiments, performing phase connection processing on the calculated differential phase value, calculating a set of TEC values of total electron content of the relative ionosphere every second, performing amplitude scintillation index calculation on the calculated signal intensity value, calculating a set of scintillation indexes every 1 second, performing space-time matching, and calculating a spatial position corresponding to each set of data at the observation time includes: relative TEC calculation, and scintillation index calculation.

In some embodiments, the centralized processing of observation data of a plurality of stations, the calculation of one ionosphere absolute TEC data per second by using a multi-station method, the inversion of a single link process by using a tomography algorithm to obtain two-dimensional electron density profile data, and the extraction of NmF2 anomalies includes: acquiring relative TEC data of all receiving stations corresponding to the station chain, calculating an integral constant of each receiving station by using a multi-station method, and adding the integral constant to the relative TEC obtained by the receiver to obtain an absolute TEC of each receiving station; dividing the range of 100 km-500 km above the station chain at intervals according to the height interval of 20km and the latitude interval of 0.5 degrees, and calculating the electron density profile of the divided region through an inversion process.

In addition, the invention also provides a data processing system of the three-frequency beacon machine, which comprises: one or more processors; storage means for storing one or more programs which, when executed by said one or more processors, cause said one or more processors to implement a method of three-beacon data processing according to the foregoing.

According to the data processing method and system of the tri-frequency beacon machine, provided by the invention, the data product with the data precision and format meeting the requirements is generated through a series of operations such as differential phase calculation and signal intensity calculation of original detection data.

Drawings

The foregoing is only an overview of the technical solutions of the present invention, and in order to make the technical solutions of the present invention more clearly understood, the present invention is further described in detail below with reference to the accompanying drawings and the detailed description.

FIG. 1 is a block diagram of a three-frequency beacon system;

FIG. 2 is a triple-frequency beacon observation system;

FIG. 3 is a schematic diagram of a level 0 data naming convention;

FIG. 4 is a schematic diagram of a class 1-4 data naming convention;

FIG. 5 is a schematic illustration of scientific data processing report naming;

FIG. 6 is a schematic illustration of scientific data product image naming;

FIG. 7 is a schematic diagram of a level 1 image product of a tri-frequency beacon;

FIG. 8 is a schematic diagram of a 2-level image product of a tri-frequency beacon;

fig. 9 is a schematic diagram of tri-band beacon 3-level data;

fig. 10 is a diagram of a triple-band beacon level 3 NmF2 dynamic variation;

FIG. 11 is a flow chart of data processing of 0-3 levels of three-frequency beacons;

fig. 12 is a flow chart of triple-band beacon level 0-1 data processing;

fig. 13 is a flow chart of tri-band beacon level 1-level 2 data processing;

FIG. 14 is a relative TEC calculation process flow;

FIG. 15 is a flicker index calculation flow;

FIG. 16 is a level 2-level 3 data processing flow;

FIG. 17 is a schematic of an integration constant calculation;

FIG. 18 is a flow chart of the k2 calculation;

FIG. 19 is a schematic of an electron density inversion grid.

Detailed Description

The preferred embodiments of the present invention will be described in conjunction with the accompanying drawings, and it will be understood that they are described herein for the purpose of illustration and explanation and not limitation.

First, abbreviations to be used hereinafter are introduced. Abbreviations used hereinafter in this application are shown in table 1.

TABLE 1 abbreviations

1 Tri-frequency beacon brief introduction

1.1 System working composition and Block diagram

The tri-frequency beacon observation system can scan the ionized layer in a large range with high time resolution and high horizontal resolution by using a satellite platform, can measure the state of the ionized layer above the F2 layer and the small-scale ionized layer inhomogeneous structure, and can obtain information such as high-precision ionized layer TEC, ionized layer electron density and the like.

The three-frequency beacon observation system consists of a satellite-borne transmitting system and a ground receiving system, and is shown in figure 1. The satellite-borne transmitting system consists of a three-grade beacon transmitter, a three-frequency beacon transmitting antenna and a radio frequency cable between the three-grade beacon transmitter and the three-frequency beacon transmitting antenna; the ground receiving system consists of four parts, namely a receiving antenna, a receiving host, GPS satellite time service and display control.

1.2 principle of operation

The observation principle of the three-frequency beacon system is shown in fig. 2. Through the omnibearing tomography scanning of the ionization layer region to be detected, the scanning projection data is further utilized to reconstruct the spatial distribution of parameters such as the electron density of the ionization layer in the region. The three-frequency beacon observation system transmits coherent signals of three frequencies to the ground through the satellite-borne transmission system, and the three-frequency beacon ground receiving system receives the signals, performs correlation analysis and calculation, and performs inversion on electron density distribution of corresponding positions, thereby realizing tomography of an ionized layer.

A constant-temperature crystal oscillator inside the satellite-borne transmitter provides a high-stability reference signal for the transmitter. After the transmitter is operated, the frequency of the transmitter is converted by the frequency source module to generate coherent radio frequency signals in a VHF frequency band (central frequency: 150.012MHz), a UHF frequency band (central frequency: 400.032MHz) and an L frequency band (central frequency: 1066.752 MHz). The three radio frequency signals are respectively passed through an orthogonal bridge and shunted to generate 6 paths of orthogonal signals with equal amplitude, the orthogonal signals are respectively accessed to an I path input end and a Q path input end of a transmitting antenna after being amplified by power, circularly polarized waves are radiated to the space, and the task of receiving the signals after the signals pass through an ionosphere is completed by a ground three-frequency beacon measuring station. And each receiver on the ground station chain receives signals through an antenna, and generates original observation data after operations such as three-band beacon signal coherent processing and the like.

The three-frequency beacon receiver system mainly comprises a receiving antenna, a receiving host, a GPS satellite time service part and a display control part. The receiving antenna mainly completes the receiving and preliminary preprocessing of the satellite beacon signals. The receiving antenna uses orthogonal dipoles which are fed orthogonally, and receives signals in a circular polarization mode. The antenna processing unit mainly completes the functions of filtering, combining and power amplifying of the beacon signals. The receiver mainly completes the core processing process of the beacon signals, and the core processing process mainly comprises the acquisition of the satellite beacon signals, the high-precision tracking of the satellite beacon signals, the coherent processing of the three-band beacon signals and the like. The receiver mainly comprises a radio frequency front end and an intermediate frequency signal processing part, wherein a radio frequency front end unit completes coherent down-conversion of beacon signals, and an intermediate frequency signal processing unit completes realization of a core algorithm. The GPS satellite time service part completes accurate time service to the receivers and can realize the time synchronization function among a plurality of receivers.

1.3 data output

The observation data of the three-frequency beacon machine is output by the ground station receiver. The receiver calculates the time of the satellite flying over the ephemeris file, and starts the machine in advance to prepare for receiving the satellite signal. When the satellite flies over the ground station, a decimal I/Q observation data file and a decimal signal-to-noise ratio data file can be generated by the receiver in each orbit.

2 standard data products

2.1 hierarchical definition

According to classification rules of satellite-to-ground observation data product classification (GB/T32453-2015) and classification and definition of electromagnetic monitoring satellite data products, the observation data of the triple-frequency beacon comprises 0-level, 1-level, 2-level and 3-level data products, and the definition of the standard data products of the triple-frequency beacon is shown in Table 2.

TABLE 2 TBB data product hierarchy definition and data product List

2.2 data product introduction at level

2.2.21-level data product

The level 1 data is single station phase data and amplitude data obtained by calculating I/Q data of three frequency points of 0 level data VHF, UHF and L. The level 1 data product comprises: scientific data, image products, and data processing reports.

(1) Scientific data

The 1-level scientific data file of the three-frequency beacon machine comprises the following physical quantity data:

differential phase data of VHF & UHF bands

Differential phase data of UHF & UHF bands

Differential phase data of L & UHF bands

Signal amplitude data of VHF band

Signal amplitude data of UHF band

Signal amplitude data of L band

The format of the data file of the three-frequency beacon 1 level is shown in the table 3 and the table 4.

Table 3 three-frequency beacon 1-level file attribute specification table

Serial number	Attribute name	Attribute content	Remarks for note
				1	SOFTVERSION	Program version number
2	PAYLOADID	Instrument code
				3	ORBITNUM	Track number
4	ORBITTYPE	Track classification	1: lifting rail, 0: falling rail

Table 4 three-frequency beacon 1 level data description table

(2) Level 1 image product

A level 1 image product of the three-frequency beacon receiver records the variation trend of the differential phase and the signal amplitude of VHF, UHF and L frequency bands of each track along with time, and the variation trend comprises 6 groups of data curves:

VHF & UHF differential phase-t

UHF & UHF differential phase-t

L & UHF differential phase-t

VHF signal amplitude-t

UHF Signal amplitude-t

L signal amplitude-t.

A schematic diagram of a level 1 image product of a tri-band beacon receiver is shown in FIG. 7.

(3) Level 1 data processing report

The three-frequency beacon machine 1 level data processing report comprises the following components:

(1) processing the version number: v0.1

(2) The starting time is yyymmdd HH and MM is SS

(3) Inputting data: 0 level data file name

(4) Assistance data

Power amplification temperature: normal/abnormal; marking an abnormal time period;

voltages of 12V, 9.5V and +/-5V: normal/abnormal; marking an abnormal time period;

(5) process of treatment

Is normal

② lack of numbers

③ data file corruption

(6) End time: yyymmdd HH MM: SS.ZZZ

(7) Outputting a product: level 1 data product filename

2.2.32-level data product

The level 2 data is relative TEC data obtained by performing mathematical processing such as phase connection on the phase data of the signals in the VHF, UHF, and L bands of the level 1 data, and ionospheric scintillation index data calculated using the amplitude data of the signals in the VHF, UHF, and L bands. The level 2 data product comprises: scientific data, image products, and data processing reports.

(1) Scientific data

The three-frequency beacon machine 2-level scientific data file comprises the following contents:

relative TEC values in VHF and UHF bands

Relative TEC values of L & UHF bands

Flicker index value of VHF band

Flicker index value in UHF frequency band

Flicker index value of L frequency band

The data format of the three-frequency beacon 2 level is shown in table 5 and table 6.

TABLE 5 three-frequency Beacon 2-level data File Attribute Specification

Serial number	Attribute name	Attribute content	Remarks for note
				1	SOFTVERSION	Program version number
2	PAYLOADID	Instrument code
				3	ORBITNUM	Track number
4	ORBITTYPE	Track classification	1: lifting rail, 0: falling rail

TABLE 6 THREE-FREE SIGNAL 2-STAGE DATA TABLE FORMAT DEFINATION

(2) Level 2 image product

The three-frequency beacon machine level 2 image product comprises the following 5 groups of curves:

relative TEC value-t of VHF & UHF frequency band

L & UHF band relative TEC value-t

In VHF band scintillation index-t

Flicker index-t in UHF band

Scintillation index-t in L band

A schematic diagram of a 2-level image product of a three-frequency beacon is shown in fig. 8.

(3) Level 2 data processing report

The three-frequency beacon machine 2-level data processing report file comprises the following components:

(1) processing software version number: v0.1

(2) Starting time: yyymmdd HH MM: SS.ZZZ

(3) Inputting data: level 1 data file name

(4) Assistance data

(5) Process of treatment

The treatment process is normal

② abnormal phenomena and possible causes

(6) End time: yyymmdd HH MM: SS.ZZZ

(7) Outputting data: 2-level data file name

2.2.43-level data product

The 3-level data is obtained by performing inversion processing on the basis of the 2-level data to generate an absolute TEC and a two-dimensional profile product of an observation link. According to the electron density two-dimensional distribution result, the abnormality of NmE, NmF1 and NmF2 is extracted. The 3-level data product comprises: scientific data and its images and data processing reports.

Scientific data

The 3-stage data product formats are detailed in tables 7 and 8.

TABLE 7 triple Beacon 3-level data File Attribute Specification

Serial number	Attribute name	Attribute content	Remarks for note
				1	PAYLOADID	Instrument code
2	ORBITNUM	Track number
				3	ORBITTYPE	Track classification	Lifting rail and lowering rail
4	SOFTVERSION	Program version number	VR0.1

Table 8 three-beacon 3-level data table format description

(2) Image product

Schematic diagrams of 3-level data products are shown in fig. 9 and 10.

(2) Level 3 data processing report

The type of the 3-level data processing report file is a TXT file. The components are as follows:

(1) processing software version number: v0.1

(2) Starting time: yyymmdd HH MM: SS.ZZZ

(3) Inputting data: (a) a level 2 data file name; (b) level 3 data of the first 5 revisiting tracks, listing their filenames

(4) Assistance data

Firstly, seismic events: seismic records of more than 6 levels in a 100km range within 1 day before the track time

Magnetic emotion index: kp, Dst, etc. of the track corresponding to time

③ IRI model

(5) Process of treatment

The treatment process is normal

② abnormal phenomena and possible causes

(6) End time: yyymmdd HH MM: SS.ZZZ

(7) Outputting data: 3-level data file name

3 data processing flow

After the ground tri-band beacon receiver tracks and locks the satellite-borne tri-band beacon signal, orthogonal component Q and homodromous component I observation data and signal-to-noise ratio data (Level 0 data) of the three frequency band signals are output. And respectively carrying out differential phase calculation and signal intensity calculation on the I/Q observation data of the three-frequency beacon receiver, and calculating to obtain VHF & UHF and UHF & L differential phase values and respective signal intensity values (Level 1 data) of three frequency bands. Carrying out phase connection processing on the calculated differential phase value, and calculating a group of TEC values of the total electron content of the relative ionized layer every second; and calculating the amplitude flicker indexes of the calculated signal intensity values, and calculating a group of flicker indexes every 1 second. Performing space-time matching, and calculating to obtain a space position (Level 2 data) corresponding to each group of data at the observation time; by carrying out centralized processing on observation data of a plurality of stations, calculating absolute TEC data of an ionized layer by a multi-station method every second, then carrying out inversion on a single link process by a tomography algorithm to obtain two-dimensional electron density profile data, and extracting the NmF2 anomaly (Level 3). The flow chart of the data processing at each level of the triple-frequency beacon system is shown in fig. 11.

4 data processing method

4.10-level data Generation of level 1 data

The processing flow of generating the level 1 data from the triple-frequency beacon level 0 data is shown in fig. 12.

4.1.1 differential phase calculation

Step 1: reading 0-level VHF frequency band I path and Q path data

The 0-level data of the three-frequency beacon is mainly two groups of data of an I path (syntropy) and a Q path (orthogonality) which are obtained by pairwise combination of VHF and UHF and combination of L and UHF frequency bands and phase difference processing. UHF is the reference frequency for differential processing when the receiver locks the signal, so the obtained effective differential data only has two combinations of VHF/UHF and L/UHF. 50 sets of data per second.

Step 2: calculating VHF/UHF differential phase

And performing arc tangent processing according to a formula (1) and calculating to obtain a differential phase corresponding to the VHF/UHF frequency band.

Φ(t)＝arctan(Q/I) (1)

And 3, step 3: reading 0-level L frequency band I path and Q path data

And reading the data of the L frequency band I path and the Q path, wherein 50 groups of data are read every second.

And 4, step 4: calculating L/UHF differential phase

And performing arc tangent calculation according to a formula (1) to obtain differential phase data corresponding to the L/UHF frequency band.

4.1.2 Signal Strength calculation

Step 1: reading I-path and Q-path data of 0-level VHF, UHF and L frequency bands

And respectively reading the I path data and the Q path data of the 0-level VHF, UHF and L frequency bands of the three-frequency beacon.

Step 2: calculating the signal strength of VHF, UHF and L frequency bands

The signal strength calculation is performed according to equation (2). And calculating the power values to obtain three groups of signal intensity value data corresponding to the VHF, UHF and L frequency bands respectively, wherein the unit is dBm.

Wherein, 10lg (I)²+Q²) The processing of (1) is to convert the relative voltage value produced by the receiver processing into a power value, and 231(dB) is the channel gain amount (including antenna gain, rf gain and algorithm gain).

4.21-level data Generation 2-level data

The processing flow of generating the 2-level data from the three-frequency beacon 1-level data is shown in fig. 13.

4.2.1 relative TEC calculation

TEC refers to the path integral of the ionospheric electron density along the satellite-receiver link. The differential doppler phase recorded by the coherent beacon receiver is limited to 0-360 deg., so that the data must be concatenated to become the corresponding relative TEC. The phase connection method is to take a proper threshold value from 0-360 degrees of data phase to judge the phase and connect the data, and the specific implementation steps are shown in fig. 14.

The specific treatment method comprises the following steps:

step 1: reading differential phase data

And respectively reading the differential phase data of the VHF & UHF and UHF & L frequency bands.

Step 2: taking an initial value and setting a phase threshold

Let t be 1, k (1) be 1, and D be 300.

And 3, step 3: differential calculation

For the phase data having n points, the difference calculation is performed on the phase data as shown in equation (3).

Δ(t)＝Φ(t+1)-Φ(t)，t＝1,2,....,n-1 (3)

And 4, step 4: phase flip number calculation

And (3) comparing the delta value of each data point with the phase threshold value D according to the relation of the formula (4) to obtain the turnover number k of 2 pi taken by each phase data during connection.

And 5, step 5: phase connection processing

The differential phase sequence of each data point is subjected to phase value increase and decrease, i.e. phase concatenation, in turn according to equation (5).

Φ'_c(t)＝Φ(t)+2kπ (5)

And 6, step 6: calculating a minimum phase value

After the connection is completed, the minimum phase value Z in the phase sequence is calculated as in equation (6).

Z＝min(Φ'_c(t)) (6)

And 7, step 7: minimum phase return to zero calculation

Φ_c(t)＝Φ'_c(t)-Z (7)

Wherein phi'_c(t) is the phase sequence connected, and the minimum phase value Z is subtracted in turn according to the formula (7), and finally the differential Doppler phase sequence phi which is connected is obtained_c(t)。

And 8, step 8: calculating a correlation coefficient C_D

Calculating the correlation coefficient C according to the formula (8)_D. Wherein, c is 299792458m/s, f_r＝16.668MHz，m_iRepresenting the frequency multiplication, if the combination of VHF/UHF is adopted, m is adopted₁＝9，m₂24; if the UHF/L combination is taken, then m is taken₁＝24，m₂＝64。

Step 9: calculating the relative TEC value

And calculating to obtain the relative TEC value of the ionized layer according to the formula (4-9). TEC values are reported in TECU (1 TECU: 1016el. m-2).

TEC_r＝C_D×Φ_c(t) (9)

4.2.2 flicker index calculation

The process of the three-frequency beacon amplitude flicker index S4 is shown in fig. 15.

Step 1: reading signal strength values

And respectively reading the data of the signal intensity values of the VHF frequency band, the UHF frequency band and the L frequency band.

Step 2: data packet

Grouping was performed in 1 group every 50 numbers (1 second length).

And 3, step 3: calculating the mean value of the squares of each group of data

The average of the squares of the individual values of each set of data was calculated according to equation (10).

And 4, step 4: calculating the square of the mean value of each group of data

The square of the mean value of each set of data was calculated according to equation (11).

And 5, step 5: calculating a flicker index value

A flicker index value is calculated every 1 second according to equation (12). Generally 0.1. ltoreq. S4<0.3 is considered ionospheric weak scintillation, 0.3. ltoreq. S4<0.6 indicates ionospheric medium scintillation, and S4>0.6 indicates ionospheric strong scintillation.

4.32-level data Generation of 3-level data

The process flow of generating 3-level data from three-frequency beacon 2-level data is shown in fig. 16.

(1) And acquiring relative TEC data of all receiving stations corresponding to the station chain, calculating an integral constant of each receiving station by using a multi-station method, and adding the integral constant to the relative TEC obtained by the receiver to obtain the absolute TEC of each receiving station.

The specific calculation steps of the absolute TEC are as follows:

the absolute TEC is first calculated using the dual frequency beacon.

As shown in FIG. 17, let the skew TEC from receiver A to satellite S be I_s1Projected vertical TEC is I_v1Then the two have the following relationship:

wherein D is₁Referred to as the transfer function of full oblique to full vertical electron content, χ, on path S1_APerpendicular to the path S1 and the receiving station AAnd (4) an included angle. Oblique total electron content I during satellite transit_s1The relationship with the differential doppler phase is expressed as:

wherein I_A(t_k) For receiver A at t_k(k 1, …, n) observing the given differential relative TEC; phi₁(0) Is the unknown integration constant of the receiving station a to be solved.Its value is equal to 1.3X 10,3, where c is the speed of light, f₀＝16.668MHz，n₁＝3,n₂＝8；Φ_d1(t_k) The receiver A observes a corresponding differential Doppler phase at each moment; formula (14) can be represented as

Wherein D is_1kAnd I_v1kAnd respectively the transfer function of the ray corresponding to the k-th observation point A and the vertical TEC. Similarly, for receiving station B there are:

the formulas (15), (16) can be changed to

Approximately, it is assumed that for two receiving stations at different locations on the ground, the vertical TECs of the two receiving stations are equal. Thus for the two-station method I_v1k≈I_v2kUsing the least squares error bound I_v1kAnd I_v2kThe error between them, the objective function of the two-station method, can be expressed as follows:

by finding phi from E₁(0) And phi₂(0) The partial derivative value of the partial derivative is made to be zero, and the partial derivative value is matched with a determinant to be solved. Simultaneous equations The unknown integral constant phi of the two stations can be obtained₁(0) And phi₂(0)。

Substituting equation (14) may then derive the dual-frequency absolute TEC for each station from the observed differential relative TEC as shown below:

generalizing to m receiving stations, the initial phases of multiple stations can also be solved according to the same concept, i.e. the so-called multi-station method. The objective function of the multi-station method is formula (4-20).

The unknown phase integral constant phi of each station can be obtained by solving the minimum value of E_i(0) (i ═ 1, ·, m), the dual-frequency absolute TEC of each station is obtained.

Second, the dual-frequency auxiliary three-frequency solution is used for solving the absolute TEC of the three-frequency beacon

Three-frequency beacon has one more carrier compared with double-frequency beaconThe carriers of the two frequencies can be differentiated as follows:

among the three-frequency beacons are:

f₁＝9f₀＝150.012MHz

f₂＝24f₀＝400.032MHz

f₃＝64f₀＝1066.752MHz

therefore, three groups of carriers can obtain three differential phase expressions (23) by using the differential Doppler technology:

wherein f0 is 16.668MHz, n1 is 9, n2 is 24, n3 is 64, Δ P₁₂、ΔP₁₃、ΔP₂₃The absolute value of the differential phase comprises the differential Doppler phase and a phase integration constant given by a receiver.

Then:

from equation (24):

ΔP₁₂(n₁ ²+n₂ ²)＝ΔP₁₃n₂ ² (25)

the above formula can be rewritten as:

(K₁₂+Δφ₁₂)(n₁ ²+n₂ ²)＝(k₁₃+Δφ₁₃)n₂ ²or

K₁₂(n₁ ²+n₂ ²)-K₁₃n₂ ²＝Δφ₁₃n₂ ²-Δφ₁₂(n₁ ²+n₂ ²) (26)

Wherein, K₁₂And K₁₃Is Δ P₁₂And Δ P₁₃Integer part, Δ φ₁₂And delta phi₁₃Is the fractional part thereof, given by the receiver data. Based on prime number principle, there are:

x₁₂(n₁ ²+n₂ ²)-x₁₃n₂ ²＝1 (27)

wherein n is₁And n₂The coprime, simultaneous equations (26), (27) can solve K₁₂And K₁₃The expression is as follows:

K₁₂＝[Δφ₁₃n₂ ²-Δφ₁₂(n₁ ²+n₂ ²)]x₁₂+kn₂ ²

K₁₃＝[Δφ₁₃n₂ ²-Δφ₁₂(n₁ ²+n₂ ²)]x₁₃+k(n₁ ²+n₂ ²) (28)

where K is an integer, then from K₁₂And K₁₃Can find out Delta P₁₂And Δ P₁₃Then converting Δ P₁₃Or Δ P₁₂Carrying in formula (23) and finishing to obtain:

in triple-band beacons n₁＝3，n₂＝8，f₀16.668 MHz. Brought into the above formula to obtain

TEC＝8.3165×10¹⁶×[(Δφ₁₃x₁₂-Δφ₁₂x₁₃)_mod1+k₂] (30)

Wherein 8.3165X 1016 (unit: el/m2) is called fuzzy coefficient, k, of TEC measured by tri-band beacon₂Is a positive integer which is undetermined as the phase integration constant of the dual-frequency beacon.

From equations (19), (30), we can simplify the solution for dual and tri-band TEC as:

TEC＝k₁×1.2995×10¹⁵+φ₁ (31)

TEC＝k₂×8.3165×10¹⁶+φ₂ (32)

wherein k is₁And k₂Respectively representing dual-frequency and tri-frequency phase integral constants, the values of which are integers phi₁And phi₂Representing the fractional parts of the dual and tri-band differential phases, respectively. Equation (31) can calculate the phase integral constant k including the error by the multi-station method₁And then combining the measured differential phase fraction phi of the ground receiver₁The absolute TEC containing the error can be calculated. Because the TECs of the same inversion region at the same time are the same, substituting the absolute TEC calculated by the formula (31) into the formula (32) to obtain k₂The calculation formula of (2):

since the absolute TEC calculated by equation (31) has an error, k calculated by equation (33)₂There is also typically an error that we need to correct accurately with the tri-band phase measurement. Research shows that the error of the dual-frequency measurement TEC is in the order of 1-3 TECU (1TECU is 1.0 × 1016el/m2), and the triple-frequency phase ambiguity coefficient is 8.3TECU (namely, the three-frequency phase ambiguity coefficient is grouped by 8.3 TECU), and is combined with k₂Is a basic premise of an integer, for k₂The flow chart of the algorithm for performing the exact solution is shown in fig. 18.

(2) Dividing the range of 100 km-500 km above a station chain at intervals according to the height interval of 20km and the latitude interval of 0.5 degrees (fig. 19), wherein in the ionospheric tomography process by satellite signal observation, the obtained total ionospheric electron content TEC is the integral of the ionospheric electron density along the path from a satellite to a receiver, and can be represented as:

TEC＝∫_pN_eds (34)

where Ne represents the electron density distribution and p is the path between the receiver and the satellite.

The electron density inversion calculation process is as follows:

firstly, ionospheric CT is carried out by using a function basis model method of spherical harmonic analysis and empirical orthogonal function combination. The spherical harmonics and empirical orthogonal functions are selected for the representation of electron density:

wherein θ andrespectively, latitude and longitude, and h, altitude.Is a normalized legendre function, and m, n, k are the order of the spherical harmonics and Empirical Orthogonal Functions (EOF), where the EOF is obtained by empirical mode or measured electron density data by singular value decomposition.

After an inversion area is determined according to the satellite orbit and the receiver position, the electron density continuously distributed in the area to be inverted is discretized into N different grids, and the electron density in each grid is assumed to be constant. Then, effective rays in the inversion region are determined, wherein the effective rays refer to rays of which the intersection points with upper and lower height boundaries (spherical shells concentric with the earth surface) of the inversion demarcated region fall in the inversion region. The projection length of the path traversed by the ray in each grid can be calculated as follows: firstly, the spherical coordinates of the satellite orbit and the ground receiving stationIs converted intoIs rectangular coordinates (X, Y, Z), and the conversion formula is:

wherein R ═ R_e+h，R_eRepresenting the radius of the earth. And then solving the coordinates of each ray passing through the intersection point on each grid in a rectangular coordinate system by using a linear equation. The equation of a straight line for a ray between a ground receiving station and a satellite can be expressed as:

wherein (X)_sat,Y_sat,Z_sat) And (X)_sta,Y_sta,Z_sta) Respectively, representing the position coordinates of the satellite and the tri-band beacon receiver (station). (X, Y, Z) represents any point on the line, K represents a proportionality constant, and K is a determinable constant for each fixed point (X, Y, Z). We can derive from the above formula:

if the intersection point is located on a plane at a distance R from the earth's center point, R is known, so there are:

X²+Y²+Z²＝R² (39)

substituting the above formula into the formula, and obtaining a quadratic equation AK²Solving the + BK + C as 0 to obtain K, replacing the K in the formula to obtain the position coordinate of the intersection point, and then calculating

The projection length Δ l between the grids is obtained. The set of projection lengths of all rays among grids in the inversion region can be represented as a projection matrix of ionosphere electron density inversion. After the projection matrix is determined, the problem of ionospheric electron density inversion can be simplified by conversion to:

y＝Hx (41)

where x is the weight coefficient of the basis function and the absolute TEC can be represented by a column of y H is an operator consisting of an empirical orthogonal function and the spherical harmonic and geometric path projection length. Since the satellite beacon multi-station method has a large error in estimating the absolute TEC, we correct the above matrix, and we adopt a method of mutually differentiating the corresponding receiver TECs to influence the multi-station method error, which is expressed as:

b＝y-y₀＝(H-H₀)x＝Ax (42)

where b is the difference of the TEC data and y0 is a reference measurement for a given receiver where forward operator a includes H and reference operator H0. Solving matrix equation by using a method of Truncated Singular Value Decomposition (TSVD) regularization, and decomposing the matrix A into a matrix A according to the singular value decomposition theory of the matrix

In the formula: u ═ U₁,u₂,...,u_n)，u_iIs a left singular value vector; v ═ V (V)₁,v₂,...,v_n)，v_iIs the right singular value vector; sigma ═ σ [ ([ sigma ] ])₁,σ₂,...,σ_n)，σ_iBeing singular values, σ₁≥σ₂≥…≥σ_nAnd k is more than or equal to 0 and is a truncation coefficient. After the basis function weight x is obtained by solving the formula (44), the initial distribution value of the electron density can be obtained by re-bringing the formula (41).

Secondly, taking the inversion result of the function base model as an iteration initial value, and performing second iteration reconstruction on the formula (42) by adopting an ART algorithm of the pixel base model to obtain an ionospheric disturbance structure with a smaller scale, wherein the calculation method comprises the following steps:

where q denotes the iteration round, the specific order may be generated sequentially or by random number method, i ═ mod (q, L) +1, L is the number of rows in matrix y, t_iThe i-th row element, R, of the matrix y_iDenotes the ith element of R, and λ is the relaxation factor. The ART algorithm is fast in convergence, generally, about 10-20 iterations are taken, and a proper iteration termination threshold value can be set in the specific execution process. And obtaining the final ionized layer electron density distribution through the inversion of the two steps.

The above description is only for the purpose of illustrating the preferred embodiments of the present invention and is not to be construed as limiting the present invention in any way, and it will be apparent to those skilled in the art that the above description of the present invention can be applied to various modifications, equivalent variations or modifications without departing from the spirit and scope of the present invention.