Zhang Heng I satellite induction type magnetometer data processing method and system
阅读说明：本技术 张衡一号卫星感应式磁力仪数据处理方法及系统 (Zhang Heng I satellite induction type magnetometer data processing method and system ) 是由 王桥 于 2021-07-30 设计创作，主要内容包括：本发明提供了一种张衡一号卫星感应式磁力仪数据处理方法及系统。该方法包括：感应式磁力仪0级数据进行数据格式转换及参数校正,生成1级数据产品；生成2级数据产品；获得3级数据产品；生成4级数据产品。本发明提供的张衡一号卫星感应式磁力仪数据处理方法及系统能够将卫星上感应式磁力仪的原始探测数据进行处理,生成实际可用的数据产品。(The invention provides a data processing method and a data processing system for a Zhang Heng I satellite induction type magnetometer. The method comprises the following steps: performing data format conversion and parameter correction on the 0-level data of the induction type magnetometer to generate a 1-level data product; generating a 2-level data product; obtaining a 3-level data product; and 4-level data products are generated. The data processing method and the data processing system for the Zhang Heng I satellite induction type magnetometer can process the original detection data of the induction type magnetometer on the satellite to generate an actually usable data product.)
1. A data processing method for a Zhang Heng I satellite induction type magnetometer is characterized by comprising the following steps:
performing data format conversion and parameter correction on the 0-level data of the induction type magnetometer to generate a 1-level data product;
generating three-component waveform and power spectrum physical quantity data with geographic and geomagnetic coordinate information by utilizing coordinate transformation according to information such as the position of a satellite point, the satellite attitude and the like, and labeling auxiliary information such as earthquake, magnetic index and the like to generate a 2-level data product;
selecting characteristic frequency band waveform or power spectrum data in 2-level data products as input, resampling the variable magnetic field data to generate a time sequence product of revisiting track observation data of a global range and a Chinese area, and labeling earthquake and magnetic index information to obtain 3-level data products;
selecting characteristic frequency band waveforms or power spectrum data in 2-level data products as input, and selecting a global range and a Chinese region overhead variable magnetic field waveform and specific frequency band power spectrum data in a sliding manner to obtain a spatial distribution background field and a dynamic variation amplitude to generate 4-level data products;
wherein, 1 level data product includes: the waveform and power spectrum physical quantity data of three components of the magnetometer recorded by each track in ULF, ELF and VLF frequency bands change along with time;
the level 2 data product comprises: the variation trend of three-component waveform and power spectrum physical quantity data of the inductive magnetometer recorded by each track in ULF, ELF and VLF frequency bands along with time and space;
the 3-level data product comprises: on the basis of 2-level data, resampling the variable magnetic field data to generate a time series product of revisit orbit observation data of a global scope and a Chinese area, and labeling earthquake and magnetic index information;
the 4-level data product comprises: file names of images of the global magnetic field spatial distribution and subgraphs contained therein.
2. The data processing method of Zhang Heng I satellite induction type magnetometer of claim 1, wherein the 1-level, 2-level, 3-level and 4-level data products comprise: scientific data, image products, and data processing reports.
3. The data processing method of Zhang Heng I satellite induction type magnetometer of claim 1, wherein the induction type magnetometer level 0 data is subjected to data format conversion and parameter correction to generate a level 1 data product, and the method comprises the following steps:
the method comprises the steps of carrying out format conversion on binary 0-level data, judging whether Fourier transform is needed according to the type of scientific data, carrying out orthogonality correction on spectral lines one by one in a frequency domain range, carrying out correction and physical quantity conversion by using a transfer function of each frequency band of a magnetic sensor, if the waveform data is the waveform data, converting the waveform data into a time domain range through inverse Fourier transform, converting frequency spectrum data into power spectrum data, and finally generating corrected three-component waveform and power spectrum quantity data, corresponding processing reports and fast-view image products.
4. The data processing method of Zhang Heng I satellite induction type magnetometer of claim 3, wherein in the calibration process, the in-orbit calibration data and calibration signal calibration parameters obtained from ground test are extracted for error analysis to judge the working state of the induction type magnetometer.
5. The Zhang Heng I satellite induction type magnetometer data processing method according to claim 1, wherein three-component waveform and power spectrum physical quantity data with geographic and geomagnetic coordinate information are generated by coordinate transformation according to information such as the position of a sub-satellite point, the satellite attitude and the like, auxiliary information such as earthquake and magnetic index is labeled, and a 2-level data product is generated, and the method comprises the following steps:
converting the three-component magnetic field data into a satellite coordinate system by utilizing the relationship between the sensor coordinate system and the satellite coordinate system;
and converting the satellite coordinate system SBC into the geographic coordinate system GEO according to the conversion matrix.
6. The Zhang Heng I satellite induction type magnetometer data processing method according to claim 1, wherein the method comprises the steps of selecting characteristic frequency band waveform or power spectrum data in 2-level data products as input, resampling the data of a changing magnetic field to generate a time series product of revisit orbit observation data of a global scope and a Chinese area, labeling earthquake and magnetic index information, and obtaining 3-level data products, wherein the method comprises the following steps:
according to the CSES satellite orbit revisiting period, by using the current orbit number, the corresponding revisiting orbit data of the previous 30 days can be obtained;
aiming at 2A-level waveform and power spectrum data of a half-orbit, grids are divided according to latitude intervals of 0.05 degrees and 0.1 degrees;
calculating a median value in each grid through time sequence analysis;
and according to the obtained median, drawing and generating a data processing report.
7. The data processing method of Zhang Heng I satellite induction type magnetometer of claim 6, wherein the step of drawing and generating a data processing report according to the obtained median value comprises:
distinguishing three components and three frequency bands from revisit orbit data of the current 5 days and sequences obtained by calculating the revisit orbit data of the previous 30 days, namely resampling time sequences, drawing a picture, and meanwhile, annotating the magnetic emotion index of the current orbit time in an icon;
drawing seismic event sequences of which the magnitude is greater than 5, 6 and more than 7 occurring in the first 30 days of a longitude area range of plus or minus 5 units of each track;
drawing the global position of the current orbit and the seismic events with the magnitude greater than 6 occurring within 7 days before the annotation;
and generating a data processing report.
8. The Zhang Heng I satellite induction type magnetometer data processing method according to claim 1, wherein the method selects a characteristic frequency band waveform or power spectrum data in 2-level data products as input, and slidably selects a global area and a Chinese area overhead variable magnetic field waveform and a specific frequency band power spectrum data to obtain a spatial distribution background field and a dynamic variation amplitude to generate 4-level data products, and comprises the following steps:
dividing the globe into a plurality of grids at spatial intervals of 2.5 ° (latitude) x 5 ° (longitude);
performing a spatial distribution analysis;
and adding the seismic event and the magnetic index to generate a 4-level data processing report.
9. The Zhang Heng I satellite induction type magnetometer data processing method of claim 8, wherein performing spatial distribution analysis comprises:
according to 30-day time interval and 5 ° (longitude) multiplied by 2.5 ° (latitude) space interval, sliding to select 2-grade data of magnetic field over global and Chinese regions, and calculating median B of all tracks in each space intervalmQuartering site ofAnd) And a fractional bit difference IQR and an abnormality discrimination value Bc(ii) a And by applying the median value B in each intervalmCarrying out interpolation to obtain global and national spatial distribution background fields;
calculating the median and the quartile difference of the varying magnetic field in each grid according to the varying magnetic field data of the current 5 days;
and calculating the difference value of the median value in each interval of each day and the background median value of the previous 30 days, normalizing by using the maximum value in each grid, and interpolating the normalized values of the intervals to obtain the global and national daily dynamic change spatial distribution map.
10. A Zhangheng I satellite induction type magnetometer data processing system is characterized by comprising:
one or more processors;
a storage device for storing one or more programs,
when executed by the one or more processors, cause the one or more processors to implement the satellite induced magnetometer data processing method 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 Zhang Heng I satellite induction type magnetometer data processing method and system.
Background
The Zhang Heng I satellite is the first autonomous geophysical field satellite developed in China, and the induction type magnetometer is one of eight loads carried on the satellite.
Disclosure of Invention
The invention aims to provide a data processing method and a data processing system for an induction magnetometer of Zhang Heng I satellite, which can process original detection data of the induction magnetometer on the satellite to generate a practically usable data product.
In order to solve the technical problem, the invention provides a data processing method of a Zhang Heng I satellite induction type magnetometer, which comprises the following steps: performing data format conversion and parameter correction on the 0-level data of the induction type magnetometer to generate a 1-level data product; generating three-component waveform and power spectrum physical quantity data with geographic and geomagnetic coordinate information by utilizing coordinate transformation according to information such as the position of a satellite point, the satellite attitude and the like, and labeling auxiliary information such as earthquake, magnetic index and the like to generate a 2-level data product; selecting characteristic frequency band waveform or power spectrum data in 2-level data products as input, resampling the variable magnetic field data to generate a time sequence product of revisiting track observation data of a global range and a Chinese area, and labeling earthquake and magnetic index information to obtain 3-level data products; selecting characteristic frequency band waveforms or power spectrum data in 2-level data products as input, and selecting a global range and a Chinese region overhead variable magnetic field waveform and specific frequency band power spectrum data in a sliding manner to obtain a spatial distribution background field and a dynamic variation amplitude to generate 4-level data products; wherein, 1 level data product includes: the waveform and power spectrum physical quantity data of three components of the magnetometer recorded by each track in ULF, ELF and VLF frequency bands change along with time; the level 2 data product comprises: the variation trend of three-component waveform and power spectrum physical quantity data of the inductive magnetometer recorded by each track in ULF, ELF and VLF frequency bands along with time and space; the 3-level data product comprises: on the basis of 2-level data, resampling the variable magnetic field data to generate a time series product of revisit orbit observation data of a global scope and a Chinese area, and labeling earthquake and magnetic index information; the 4-level data product comprises: file names of images of the global magnetic field spatial distribution and subgraphs contained therein.
In some embodiments, level 1, level 2, level 3, and level 4 data products each include: scientific data, image products, and data processing reports.
In some embodiments, the data format conversion and parameter correction are performed on the induction type magnetometer level 0 data to generate a level 1 data product, and the method comprises the following steps: the method comprises the steps of carrying out format conversion on binary 0-level data, judging whether Fourier transform is needed according to the type of scientific data, carrying out orthogonality correction on spectral lines one by one in a frequency domain range, carrying out correction and physical quantity conversion by using a transfer function of each frequency band of a magnetic sensor, if the waveform data is the waveform data, converting the waveform data into a time domain range through inverse Fourier transform, converting frequency spectrum data into power spectrum data, and finally generating corrected three-component waveform and power spectrum quantity data, corresponding processing reports and fast-view image products.
In some embodiments, during the calibration process, the on-track calibration data and calibration signal calibration parameters obtained by a ground test are extracted to perform error analysis to judge the working state of the induction type magnetometer.
In some embodiments, three-component waveform and power spectrum physical quantity data with geographic and geomagnetic coordinate information are generated by coordinate transformation according to information such as the position of a sub-satellite point, the satellite attitude and the like, and auxiliary information such as earthquake, magnetic index and the like is labeled to generate a 2-level data product, which comprises: converting the three-component magnetic field data into a satellite coordinate system by utilizing the relationship between the sensor coordinate system and the satellite coordinate system; and converting the satellite coordinate system SBC into the geographic coordinate system GEO according to the conversion matrix.
In some embodiments, the method includes selecting characteristic frequency band waveform or power spectrum data in 2-level data products as input, resampling the variable magnetic field data to generate a time series product of revisit orbit observation data of a global scope and a Chinese area, labeling earthquake and magnetic index information, and obtaining 3-level data products, including: according to the CSES satellite orbit revisiting period, by using the current orbit number, the corresponding revisiting orbit data of the previous 30 days can be obtained; aiming at 2A-level waveform and power spectrum data of a half-orbit, grids are divided according to latitude intervals of 0.05 degrees and 0.1 degrees; calculating a median value in each grid through time sequence analysis; and according to the obtained median, drawing and generating a data processing report.
In some embodiments, the plotting and yield data processing report based on the median value comprises: distinguishing three components and three frequency bands from revisit orbit data of the current 5 days and sequences obtained by calculating the revisit orbit data of the previous 30 days, namely resampling time sequences, drawing a picture, and meanwhile, annotating the magnetic emotion index of the current orbit time in an icon; drawing seismic event sequences of which the magnitude is greater than 5, 6 and more than 7 occurring in the first 30 days of a longitude area range of plus or minus 5 units of each track; drawing the global position of the current orbit and the seismic events with the magnitude greater than 6 occurring within 7 days before the annotation; and generating a data processing report.
In some embodiments, selecting a characteristic frequency band waveform or power spectrum data in a 2-level data product as an input, and sliding and selecting a global-range and Chinese-region overhead variable magnetic field waveform and a specific frequency band power spectrum data to obtain a spatial distribution background field and a dynamic variation amplitude to generate a 4-level data product, including: dividing the globe into a plurality of grids at spatial intervals of 2.5 ° (latitude) x 5 ° (longitude); performing a spatial distribution analysis; and adding the seismic event and the magnetic index to generate a 4-level data processing report.
In some embodiments, a spatial distribution analysis is performed, comprising: according to 30-day time interval and 5 ° (longitude) multiplied by 2.5 ° (latitude) space interval, sliding to select 2-grade data of magnetic field over global and Chinese regions, and calculating median B of all tracks in each space intervalmQuartering site ofAnd) And a fractional bit difference IQR and an abnormality discrimination value Bc(ii) a And by applying the median value B in each intervalmCarrying out interpolation to obtain global and national spatial distribution background fields; calculating the median and the quartile difference of the varying magnetic field in each grid according to the varying magnetic field data of the current 5 days; and calculating the difference value of the median value in each interval of each day and the background median value of the previous 30 days, normalizing by using the maximum value in each grid, and interpolating the normalized values of the intervals to obtain the global and national daily dynamic change spatial distribution map.
In addition, the invention also provides a data processing system of the Zhang Heng I satellite induction type magnetometer, which comprises the following components: one or more processors; a storage device for storing one or more programs which, when executed by the one or more processors, cause the one or more processors to implement the method for processing data of the Zhang Heng I satellite inductive magnetometer according to the foregoing.
After adopting the technical scheme, the invention at least has the following advantages:
according to the Zhang Heng I satellite induction type magnetometer data processing method and system, a data product with data precision and format meeting requirements is generated through a series of operations such as format conversion, data calibration, physical quantity correction, coordinate system transformation and inversion of original data of the induction type magnetometer.
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 schematic diagram of a data naming rule of 1-4 levels;
FIG. 2 is a naming schematic of a scientific data processing report;
FIG. 3 is a naming schematic of a scientific data product image;
FIG. 4 is a schematic diagram of the file format of the induction type magnetometer 1-level data hdf 5;
FIG. 5a is a schematic representation of a waveform of an inductive magnetometer VLF;
FIG. 5b is a schematic representation of power spectrum data of the inductive magnetometer VLF;
FIG. 6 is a schematic diagram of the format of induction type magnetometer class 2 scientific data hdf 5;
FIG. 7a is a schematic diagram of a 2-stage data waveform of an inductive magnetometer;
FIG. 7b is a schematic diagram of power spectrum data of 2-stage data of the induction type magnetometer;
FIG. 8 is a schematic diagram of a 4-level data product;
FIG. 9 is a schematic diagram of the process flow of the induction type magnetometer grade 0-4 data processing;
FIG. 10 is a flow chart of inductive magnetometer 1 level data processing;
FIG. 11 is a flow chart of inductive magnetometer level 2 data processing;
FIG. 12 is a 3-level data processing flow diagram;
FIG. 13 is a 4-level data processing flow diagram.
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. The meanings of these abbreviations are shown in Table 1:
TABLE 1 abbreviations
1 introduction and data production of inductive magnetometer
An induction type Magnetometer (SCM) is a satellite-borne instrument based on Faraday's law of electromagnetic induction, and is mainly used for acquiring waveform and spectrum information of a variable magnetic field vector within a frequency range of 10Hz to 20kHz at the position of a CSES satellite orbit.
The induction type magnetometer can carry out the switching of the working mode according to the injection instruction of the satellite affair management subsystem. The data products produced by the induction magnetometer in different operating modes are shown in table 2.
Table 2 data output table of induction type magnetometer
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 induction type magnetometer comprises 0-level, 1-level, 2-level, 3-level and 4-level data products, and the definition of the standard data products of the induction type magnetometer is shown in Table 3.
TABLE 3 SCM data product hierarchy definition and data product List
2.2 naming rules
In order to facilitate the retrieval and query of each level of data products, the names of each level of data products should include satellite names, load names, track numbers, data start and stop times and other necessary identifications. An example of naming rules of induction magnetometer class 0-4 data products, images and processing reports is shown in FIG. 1.
Wherein:
(1) satellite name (4-bit character): denoted by CSES;
(2) satellite number (1 digit): starting from 1, increasing sequentially;
(3) payload encoding (3 characters): HPM, SCM, EFD, LAP, RPA, GRO, TBB, HEP represent 8 loads respectively;
(4) load number (1 digit): the method is used for distinguishing the condition that one load loads a plurality of probes with the same type of items, and sequentially increasing from 1, wherein 1 represents a probe with the number of 1,2 represents a probe with the number of 2, and sequentially increasing, and if the probes are not distinguished, the probes are represented by 0; except for the high-energy particle load, all other loads are 0; for the high-energy particle load, 1-4 sequentially represent a low-energy probe, a high-energy probe, an Italian load and an X ray;
(5) data scalable coding (1 bit sign 2 bit number): from left to right, the first bit is "L", the right two bits represent data levels, 0-4 levels are represented by 00, 01, 02/2A, 03, and 04, respectively;
(6) observed object code (2-digit number): and setting an observation object code according to the classification and code of the earthquake electromagnetic satellite survey items (submission). For 0-level data, the observed object is coded as 00;
(7) track number (5 digits): beginning with 00001, the data product used for organizing data files according to the track and incapable of marking track numbers is represented by 00000;
(8) lifting rail flag (1 digit): the ascending rail is 1, and the descending rail is 0;
(9) the data start time is represented by 14 digits, wherein the year (4 digits), the month (2 digits), the day (2 digits), the hour (2 digits), the minute (2 digits) and the second (2 digits);
(10) the end time of the data is in the same format (9);
(11) receiving station encoding (3-bit number): the magnetic field sensor is specially reserved for a tri-band beacon receiver, does not relate to ground receiving station information for an inductive magnetometer, and is marked as 000;
(12) file extension name: when the file extension is dat, representing a data file stored in binary format; when the file extension is h5, the representation data is stored in hdf format.
For scientific data processing report naming, a processing report is stored in an ASCII code format with txt as a file suffix, and the naming is formed by modifying partial fields and adding corresponding 'RP' suffixes on the basis of scientific data naming, as shown in FIG. 2.
Image naming for scientific data products is similar to process report file naming by modifying part of the fields and adding the corresponding "_ xx.png" suffix, an example of which is shown in fig. 3, based on scientific data naming.
(13) The expansion code for representing the observation object is composed of two bits of characters, wherein the first bit is used for distinguishing different areas: 1 represents a Chinese region, 2 represents a global, and 0 represents an undifferentiated region; the second position is used to distinguish multiple images of the same load, and the images are sequentially increased from 1.
2.3 data product introduction at level
2.3.11-level data product
And after format conversion, orthogonality correction and transfer function correction are carried out on the 0-level data, physical quantity data of three-component waveforms and frequency spectrums of the induction type magnetometer arranged according to time are generated, and auxiliary data required in the processing process are calibration correction data. The level 1 data product comprises: scientific data, image products, and data processing reports.
(1) Level 1 scientific data
● ULF: three-component waveform data in a 10-200 Hz frequency band in the working modes of detailed inspection, patrol and calibration;
● ELF: three-component waveform data in a 200-2.2 kHz frequency band under the working modes of detailed inspection, patrol and calibration;
● VLF: three-component waveform data in a 1.8 k-20 KHz frequency band under the working modes of detailed inspection, patrol and calibration; and (3) three-component power spectrum data in a frequency band of 1.8 k-20 KHz in a patrol and calibration working mode.
Table 4 induction type magnetometer 1-stage data file attribute description table
TABLE 5 1-STAGE DATA TABLE FORMAT FOR INDUCTIVE MAGNETO-METER
(2) Level 1 image product
The waveform and power spectrum physical quantity data of three components of the magnetometer recorded by each track in ULF, ELF and VLF frequency bands change along with time; the respective parameter displays are shown in fig. 5.
(3) Level 1 data processing report
The related information of the level 1 scientific data in the processing engineering is mainly recorded, and the related information comprises data processing time, program version information, data continuity description, load working state description and the like. The components are as follows:
processing software version number: ver0.1
Starting time:
yyyymmdd HH:MM:SS.ZZZ
inputting data: 0-level data product filename
Assistance data
Correcting parameters: orthogonality correction data; three-component three-frequency-band frequency response voltage-magnetic field conversion parameters; calibrating parameters of on-orbit calibration signals;
description of load operating condition
● electronics box temperature: is normal
● monitoring voltage data
● 3.3.3 v voltage: packet number-error case, 0- < threshold range 1- > threshold range 2-invalid value;
● 2.5.5 v voltage: packet number-error case, 0- < threshold range 1- > threshold range 2-invalid value;
● 1.9.9V Voltage: packet number-error case, 0- < threshold range 1- > threshold range 2-invalid value;
● 625Hz signal detection Voltage: packet number-error case, 0- < threshold range 1- > threshold range 2-invalid value;
● 10kHz signal detection voltage: packet number-error case, 0- < threshold range 1- > threshold range 2-invalid value;
data processing process
● Normal
● lack number: packet number-error entry;
● data file corruption: packet sequence number-time;
● method used in the treatment process: DFT calculation and IDFT conversion; FFT calculation and IFFT transformation; correcting orthogonality; correcting a transfer function;
seventhly, judging on-track calibration data: normal/abnormal;
● 625Hz time Vx Δ x Vy Δ y Vz Δ z;
● 10kHz time Vx Δ x Vy Δ y Vz Δ z;
(Vx, Vy, Vz are three-component calibration signal data; Deltax, Deltay, Deltaz are three-component calibration signal data error values)
(iii) end time
yyyymmdd HH:MM:SS.ZZZ
Ninthly, outputting data: level 1 scientific data filename
2.3.22-level data product
The 2-level data is obtained by transforming the 1-level data from a sensor coordinate system to a geographical coordinate system by using a coordinate transformation matrix to generate three-component waveform and power spectrum physical quantity data with geographical, geomagnetic coordinates and attitude information. The auxiliary data required in the processing process comprises a sensor-satellite coordinate system conversion matrix, a satellite-ground coordinate conversion matrix and the like. The level 2 data product comprises: scientific data, image products, and data processing reports.
(1) Level 2 scientific data
● ULF: three-component waveform data in a 10-200 Hz frequency band in the working modes of detailed inspection, patrol and calibration;
● ELF: three-component waveform data in 200-2.2 kHz frequency band under detailed inspection, patrol and calibration working modes
● VLF: three-component waveform data in a 1.8 k-20 KHz frequency band under the working modes of detailed inspection, patrol and calibration; and (3) three-component power spectrum data in a frequency band of 1.8 k-20 KHz in a patrol and calibration working mode.
Table 6 induction type magnetometer 2-stage data file attribute description table
Table 7 inductive magnetometer 2-stage data table format description table
(2) Level 2 image product
The three-component waveform and power spectrum physical quantity data of the inductive magnetometer recorded in each track change along with time and space in ULF, ELF and VLF frequency bands. Graphs are shown as shown in fig. 7a and 7b, including ULF, ELF and VLF X, Y, Z waveforms and VLF power spectra.
(3) Level 2 data processing report
The method mainly records relevant information of 2-level scientific data in processing engineering, including information such as data processing time, program version information, coordinate transformation matrix, data quality and the like. The components are as follows:
processing software version number Ver0.1
Starting time: yyymmdd HH MM: SS.ZZZ
Inputting data: the file names of the level 1 data products are listed.
③ auxiliary data information
● calibration parameters: a sensor satellite coordinate transformation matrix; converting the matrix in a star mode;
treatment process
● the treatment process is normal
● abnormal phenomenon and possible cause
● data file corruption: packet sequence number-time;
● method used in the treatment process: converting sensor-satellite coordinates; satellite-geographic coordinate conversion; converting geographic-geomagnetic coordinates;
whether the input function is corrected: 0-no;
1-is: correcting the starting time and the transfer function correction coefficient;
sixthly, yyymmdd HH as the end time, and SS.ZZZ as the end time
And (c) outputting data, namely listing the file names of 2-level data products.
2.3.33-level data product
The 3-level data of the induction type magnetometer is a time sequence product which is used for resampling the variable magnetic field data on the basis of the 2-level data to generate revisit orbit observation data of a global scope and a Chinese area, and is labeled with earthquake and magnetic index information. The 3-level data product comprises: scientific data, image production, and process reports.
(1) Level 2 scientific data
The 3-level scientific data file content comprises:
the median, mean, standard deviation and other statistics of the current orbit in the grid area
Revisiting statistics such as median, upper and lower quartile values and difference (i.e. quartile), standard deviation of data in grid region in a certain period
Identification exceeding upper and lower limits (0 is normal; difference between current observation value and upper and lower limits is marked if exceeding the boundary)
The format of these data is shown in tables 8 and 9.
TABLE 8 Induction magnetometer 3-level data file Attribute description
Serial number
Attribute name
Attribute content
Remarks for note
1
PAYLOADID
Instrument code
2
ORBITNUM
Track number
3
ORBITFLAG
Lifting rail mark
Lifting rail 1 and lowering rail 0
4
PSD_tRES
Time resolution of VLF power spectra
2 seconds
5
PSD_fRES
Frequency resolution of VLF power spectrum
12.5Hz
6
SOFTVERSION
Program version number
0.1
TABLE 9 3-STAGE DATA TABLE FORMAT FOR INDUCTIVE MAGNETO METER (FOUR-DIVIDE ALGORITHM is used as an example)
(2) 3-level image product
(3) Level 3 data processing report
The induction type magnetometer 3-level data processing report comprises the following components:
processing software version number: ver0.1
Starting time: yyymmdd HH MM: SS.ZZZ
Inputting data: (a) current track data, listing the file name; (b) 3-level data of the first 5 re-visited orbits (taking 30-day detection time span of the satellite as an example), and listing file names of the data; a missing data condition.
(iv) auxiliary data:
(a) seismic catalog: yyymmdd. earthquales. txt (b) magnetic susceptibility index:
(b) current 5 days of KP (3 hours), DST (1 hour), F107(1 day), AE (1 minute sampling at 0.1 ° intervals along latitude), and then sampling day-to-day sliding four-digit difference method;
fifthly, the treatment method
End time: yyymmdd HH MM: SS.ZZZ
And seventhly, outputting a product: listing file names of 3-level data products
Results of processing [ + ]
(a) The treatment process is normal
(b) The lack of the number of each grid in the treatment process; there are problems in inversion;
description of the Announcement of an abnormal situation
(a) 3 points are continuously arranged, 34 minutes are from 3 months, 23 days and 01 hours in 2020 to 36 minutes are from 3 months, 23 days and 01 hours in 2020, and the latitude and longitude range of the sub-satellite points (23-33 degrees N, 110-112 degrees E) is out of range. (b) suggesting further analysis.
2.3.44-level data product
The 4-level data of the induction type magnetometer are selected by sliding the waveform of the magnetic field changing in the air over the global range and the Chinese area and the power spectrum data of a specific frequency band according to a plurality of revisit cycle time spans and certain longitude and latitude space intervals to obtain a space distribution background field; the dynamic variation amplitude is obtained by subtracting the spatial distribution background field in a certain period of time (hereinafter, 30 days as an example) from the spatial distribution of the variable magnetic field waveform and the power spectrum of the specific frequency band observed in the current 5 days.
(1) 4-level data product
TABLE 10 Induction magnetometer 4-level data file Attribute description
Serial number
Attribute name
Attribute content
Remarks for note
1
PAYLOADID
Instrument code
2
SOFTVERSION
Program version number
0.1
3
ORBITTYPE
Track classification
TABLE 11 4-STAGE DATA TABLE FORMAT FOR INDUCTIVE MAGNETO METER
The file name of an image of global magnetic field spatial distribution, which is a 4-level data product, and subgraphs contained in the image are shown in fig. 8, and the specific contents are described as follows.
FIG. 1: the x component is the current 5-day and 625Hz power spectral density amplitude spatial distribution map;
and (2) in a subfigure: the power spectral density amplitude spatial distribution map of 625Hz 30 days before the x component;
FIG. 3 is a drawing: subtracting the difference value of the sub-image 2 from the sub-image 1, and marking the seismic events which are more than 6 levels within 30 days;
and (5) in a subfiguration of 4: the current 5-day and 625Hz power spectral density amplitude spatial distribution map of the y component;
subfigure 5: the power spectral density amplitude spatial distribution map of 625Hz in 30 days before the y component;
subfigure 6: subtracting the difference value of the subgraph 5 from the subgraph 4, and marking the seismic events which are more than 6 levels within 30 days;
subfigure 7: the current 5-day and 625Hz power spectral density amplitude spatial distribution map of the z component;
subfigure 8: the power spectral density amplitude spatial distribution map of 625Hz in 30 days before the z component;
subfigure 9: the difference after sub-graph 8 is subtracted from sub-graph 7 and seismic events greater than 6 levels or more within 30 days are noted.
The image for the chinese region is similar to the global image content. The waveform amplitude global distribution diagram after the ULF and ELF resampling is similar to the power spectral density global distribution diagram of a specific frequency point.
(3) Level 4 data processing report
The induction type magnetometer 4-level data processing report comprises the following components:
processing software version number: ver0.1
② the initial time of yyyymmdd HH to MM to SS.ZZZ
Input data
(a) All level 2 data files within the first 35 days
(b) Missing data case
(iv) auxiliary data
Seismic catalog: txt is yyymmdd _ earthquales
Magnetic susceptibility index: kp (3 hours), DST (1 hour), F107(1 day), AE (1 minute) for the current 5 days.
Fifthly, the treatment method
Sixthly, yyymmdd HH as the end time, and SS.ZZZ as the end time
Output data 4-level data file name
Results of processing [ + ]
(a) The treatment process is normal;
(b) data redundancy overflow occurs in the processing process, the program is interrupted due to the occurrence of redundancy error when the program is stopped, and related pictures and analysis results are not generated
Ninthly, the abnormal situation shows that for example, 34 minutes when 3, 23 and 01 in 2020 to 36 minutes when 3, 23 and 01 in 2020, the latitude and longitude range of the subsatellite point is (23-33 degrees N, 110-112 degrees E) and the out-of-range value is as follows: . In this case, Dst is-60 nT, kp is 4, and F107 maximum is 150.
3 data processing flow
The method comprises the following steps that (1) data format conversion is carried out on 0-level data of the induction type magnetometer, if the waveform data need to be subjected to Fourier transform to be converted into a frequency domain for calibration parameter correction and then subjected to Fourier inverse transformation to generate corrected waveform data, the frequency spectrum data are directly corrected in the frequency domain to generate correction, then power spectrum calculation is carried out to generate power spectrum data, meanwhile, on-track calibration data are extracted for error analysis to judge the working state of the induction type magnetometer, and finally a 1-level data product is generated; according to the position of the satellite point, the satellite attitude and other information, coordinate transformation is utilized to generate three-component waveform and power spectrum physical quantity data with geographic and geomagnetic coordinate information, and auxiliary information such as earthquake, magnetic index and the like is labeled to generate a 2-level data product. And selecting the characteristic frequency band waveform or power spectrum data as input to respectively obtain corresponding 3-level and 4-level data products.
4 data processing method
4.10-level data Generation of level 1 data
Fig. 10 is a processing flow of generating 1-level data from 0-level data of the induction magnetometer, specifically: the method comprises the steps of carrying out format conversion on binary 0-level data, judging whether Fourier transform is needed according to the type of scientific data, carrying out orthogonality correction on spectral lines one by one in a frequency domain range, carrying out correction and physical quantity conversion by using a transfer function of each frequency band of a magnetic sensor, if the waveform data is the waveform data, converting the waveform data into a time domain range through inverse Fourier transform, converting frequency spectrum data into power spectrum data, and finally generating corrected three-component waveform and power spectrum quantity data, corresponding processing reports and fast-view image products. In the correction process, on-orbit calibration data and calibration signal calibration parameters obtained by a ground test are extracted to carry out error analysis and judge the working state of the induction type magnetometer.
Binary conversion of 4.1.1 to decimal data
Binary scientific data (waveform data of ULF, ELF, waveform and spectral phase data of VLF) generated by the level 0 processing and self-collected engineering parameters are converted into decimal data.
(1) Scientific data format conversion of induction type magnetometer
Scientific data of the induction type magnetometer is collected by an AD chip, the data is stored as a 16-bit signed binary number, the highest bit is a sign bit, 0 represents a positive number, 1 represents a negative number, and the induction type magnetometer can represent- (2)n-1-1)…(2n-1-1) decimal values in the range (i.e. -32767 to 32767) where n is 16. And if the measuring voltage range set by the AD acquisition chip is +/-10V, the electric quantity measuring precision represented by the binary numerical value is 10V/32767. Since the range of the physical quantity collected by the AD chip has been set to ± 10V, in order to correspond the binary value collected by it (hereinafter, the binary value is written using 16 system for simplicity of expression) to the voltage value of 10 system, the voltage value corresponding to the original code 0000 of 16 system is 0V, and the voltage value corresponding to the 7FFF of 16 system is + 10V. The negative voltage is represented by complement, i.e. 8001 in 16 th order corresponds to-10V. According to the corresponding relationship, the binary value HH is converted into decimal voltage V according to the formula (4-1) and the formula (4-2).
The positive conversion formula is:
the unit is volts (V).
The negative conversion formula is:
the unit is volts (V).
Scientific data of the induction type magnetometer are classified into three types: waveform data, spectral amplitude data, and phase data. The data format conversion method is as follows:
● waveform data
The waveform data is a 16-bit signed integer, the most significant bit is the sign bit: "0" means positive; "1" means negative. The positive conversion formula is the formula (4-1) and the negative conversion formula is the formula (4-2).
● spectral magnitude data
The amplitude data is a 16-bit signed integer, the most significant bit is a sign bit, and the data is used for indicating whether the amplitude data is amplified inside the instrument: "0" means positive, no amplification; "1" means negative, magnification-100 times. The positive number is not amplified, and the conversion formula is (4-1); the negative number is data represented by complement after being amplified by the instrument to-100, and the conversion formula is an expression (4-3).
The unit is volts (V). If the data is small, the program will make 100 times magnification (multiple data bits) to improve the data accuracy when saving the data, and mark the highest bit with "1" because 16 is signed data, and convert it to decimal data value divided by-100.
● spectral phase data
The phase data is 16 bits signed data, the most significant bit is the sign bit: "0" means positive and "1" means negative, and the value is amplified 1000 times inside the apparatus.
The positive conversion formula is:
units are degrees (°).
The negative conversion formula is:
units are degrees (°). Wherein, the circumference ratio is 3.14159. Because the calculated phase data is relatively small, the program performs 1000 times of amplification (multiple data bit number) at the time of data saving to improve data accuracy.
4.1.2 Fourier transform
The waveform data is mainly corrected in a frequency domain, so that the waveform data needs to be subjected to Fourier transform into the frequency domain, amplitude correction and physical quantity conversion are carried out one by one on spectral lines, and finally, signals are converted from the frequency domain to the time domain by means of inverse Fourier transform, so that the corrected waveform data is generated.
(1) Discrete Fourier Transform (DFT) of ULF and ELF waveform data
Wherein, x (n) is discrete signal data of ULF and ELF time domains; n is the number of continuous data, and according to the number of scientific data in the data packet, the ULF waveform data takes N as 82, and the ELF waveform data takes N as 820;
x (k) is the frequency domain data obtained by the DFT computation.
(2) Fast Fourier Transform (FFT) of VLF waveform data
The signal is first transformed from the time domain to the frequency domain using an FFT. The FFT is a linear integral transform, and its principle is to decompose N-point X data into two N/2-point X1(N) and X2(N) data, and perform N/2-point DFT fourier transform on X1(N) and X2 (N). The concrete formula is shown in formula (4-7) and formula (4-8).
Wherein x (n) is a time-domain discrete signal; n is the number of continuous data, and N is 4096 for VLF waveform data; x1(k) and X2(k) are N/2-point discrete Fourier transforms of X1(N) and X2(N), respectively, see formulas (4-6); x1(n) and x2(n) are arrays obtained by parity splitting.
4.1.3 orthogonality correction
Due to errors caused by machining and assembling processes, the magnetic axis and the geometric axis of the three components of the induction type magnetometer sensor cannot be completely orthogonal. The relationship between the sensor magnetic orthogonal coordinate system and the structure orthogonal coordinate system is obtained through a magnetic orthogonality calibration test and a geometric form and position measurement test, and orthogonality correction is carried out on the measured space magnetic field change data. Therefore, the auxiliary data required for the orthogonality correction process is the orthogonality correction matrix MorthThe orthogonality correction formula is as follows:
wherein the content of the first and second substances,three-component spectral amplitude data;amplitude data after three-component orthogonal correction; morthIs an orthogonality correction matrix.
4.1.4 transfer function correction
And performing amplitude correction and physical quantity conversion on the scientific data of the voltage quantity in a frequency domain one by one through electromagnetic conversion factors of the magnetic sensor to generate corrected three-component magnetic field physical quantity data. The auxiliary data required in this process are temperature data in the sensor envelope and three-component voltage-magnetic field conversion parameters with temperature information (i.e. transfer function of the magnetic sensor):
● sensors monitor temperature data;
● ULF voltage-field switching frequency response parameter;
● ELF voltage-field switching frequency response parameter;
● VLF voltage-field switching frequency response parameter;
the temperature data of the magnetic sensor is obtained from self-acquisition temperature engineering parameters of the induction type magnetometer; and obtaining the frequency band response experiment of the magnetic sensor under the conditions of different temperatures in a laboratory by using the three-component voltage-magnetic field conversion frequency response parameters. According to the real-time temperature monitoring data of the magnetic sensor, voltage-magnetic field conversion factor data in a corresponding temperature interval are selected, the voltage quantity of a corresponding frequency point of a spectral line is converted into a magnetic field B1, and the conversion formula is as follows:
4.1.5 inverse Fourier transform
After amplitude correction is completed on the scientific data in the frequency domain, waveform data are converted into time domain by utilizing inverse Fourier transform to generate waveform physical quantity data.
(1) Inverse Discrete Fourier Transform (IDFT) of ULF and ELF
Transforming the signal from the frequency domain to the time domain by IDFT:
where X (k) is the frequency domain amplitude value of the ULF and ELF waveform data; x (n) is waveform data after IDFT conversion; n is the number of data involved in the calculation, and N is 82 for ULF and 820 for ELF.
(2) VLF inverse fast Fourier transform IFFT
The signal is transformed from a frequency domain to a time domain through IFFT, and N-point inverse Fourier transform (IDFT) is decomposed into N/2 two-point IDFT for calculation, and transformation formulas are shown as an expression (4-12) and an expression (4-13).
Wherein x1(n) and x2(n) are time domain data obtained by performing IDFT conversion on frequency domain data obtained by parity splitting using the formula (4-8); x (n) is a discrete signal in the time domain; n is the number of consecutive data VLF, and N is 4096.
4.1.6 Power spectral Density calculation
Calculating the frequency spectrum data of the VLF frequency band into Power Spectrum Density (PSD) data, wherein the power spectrum density calculation adopts the following formula:
wherein, PSDCH1(kΔf)、PSDCH2(kΔf)、PSDCH3(k Δ f) is the power spectral density of each component of the VLF in nT · Hz-1/2; k is a serial number of discrete frequency domain data, and k is 0,1,2, N/2; fVLF is the sampling frequency of the VLF band, constant, and fVLF is 51.2 kHz. N is the FFT calculation sample size of the VLF frequency band, and the value is 4096; Δ f is the frequency interval of the VLF spectrum data, and Δ f ═ fband312.5/N in Hz; b1x(k·Δf)、B1y(k·Δf)、B1z(k · Δ f) is an amplitude value of VLF component spectrum data, and is expressed in nT.
4.1.7 on-track calibration data extraction
In a calibration working mode, the induction type magnetometer has the function of generating 625Hz and 10kHz calibration signals, the calibration signals take a magnetic flux negative feedback coil as a calibration coil, a calibration magnetic field signal is generated in the magnetic sensor, and the calibration magnetic field induced electromotive force is acquired by the magnetic sensor. And extracting data in an in-orbit calibration mode according to the working state of the induction type magnetometer, extracting voltage data of 625Hz and 10kHz frequency points in a VLF frequency band in a frequency domain range, carrying out error analysis on the voltage data and calibration parameters of in-orbit calibration signals measured on the ground, and judging the working state of the induction type magnetometer.
4.21-level data Generation 2-level data
And (3) performing space-time matching on information such as the satellite affair data and the magnetic condition index with each orbit data, converting the 1-level data into a satellite coordinate system and a satellite-ground coordinate system, and labeling the geographic coordinates, the geomagnetic coordinates and the attitude and orbit control data. The process flow is shown in FIG. 11.
4.2.1 conversion of sensor coordinate System to satellite coordinate System
And (3) detecting the relation between a coordinate system (a sensor coordinate system) and a satellite structure reference cubic mirror coordinate system (a satellite coordinate system) by using the sensor geometric coordinate system, and converting the three-component magnetic field data into the satellite coordinate system. The required auxiliary parameters are a sensor-to-satellite coordinate system transformation matrix M1, which is a 3 x 3 coordinate transformation matrix, measured in the sensor assembly and satellite final assembly processes and provided by loads and satellite manufacturers;
assuming that three axes of a satellite platform coordinate system are respectively X, Y and Z, and three axes of an inductive magnetometer sensor are respectively X1，Y1And Z1(X of the sensor)1The axis is parallel to the X axis of the platform), the output signals of the sensors are respectively B1x，B1yAnd B1zThen, the conversion relationship of the observation data from the sensor coordinate system (the detection coordinate system) to the satellite coordinate system is:
wherein, B2x，B2yAnd B2zThree-component magnetic field data under a satellite coordinate system, unit, nT; m1 is provided by the load vendor after the satellite is generally assembled.
4.2.2 satellite-to-ground coordinate transformation
Let three euler angles between the satellite coordinate system with respect to the geographic coordinate system (geocentric earth-fixed system) be phi, psi, theta due to changes in satellite attitude. Given the transformation matrix from the satellite coordinate system SBC to the geographic coordinate system GEO as M2, then M2 can be expressed as:
for three-component magnetic field data in satellite coordinate system (B)2x，B2yAnd B2z) Indicating that the alternating magnetic field in the geographic coordinate system is (B)3x，B3yAnd B3z) This means that there are:
the conversion relationship between the geographic coordinate system GEO and the geomagnetic coordinate system MAG is set as follows:
wherein, theta,Respectively the geographical latitude, longitude of the dipole axis.
4.32-level data Generation of 3-level data
Based on the 2 data product of the induction type magnetometer, 3-level data product and a processing report are obtained according to the following steps, and the flow is shown in fig. 12. The method for obtaining the background value is described by taking the median of the 30-day revisit orbit data as an example.
The first step is as follows: a revisit track is searched.
The specific algorithm for searching the revisit orbit is as follows:
according to the CSES satellite orbit revisiting period, by using the current orbit number, the revising orbit data corresponding to the previous 30 days can be obtained, and the specific formula is as follows:
S＝L-n×M,n＝1,…,6 (4-19)
wherein, the number is the number of the revisiting track; is the current track number; track fix differences are revisited for CSES.
The second step is that: and (5) dividing the grids.
Aiming at 2A-level waveform and power spectrum data of a half-orbit, grids are divided according to latitude intervals of 0.05 degrees (Chinese area:) and 0.1 degrees (other areas).
The third step: and (5) performing time sequence analysis to obtain a median value in each grid.
The calculation method of the median is described by the following algorithm. Suppose that there are n magnetic field data B arranged from small to large in each grid1,…,Bn. Median value B when n +1 is an integer multiple of 4mAnd the first and third quartilesThe calculation formula is as follows:
when not an integer multiple of 4, the quartile position calculated by the above formula carries a decimal number, and at this time, the decimal number is rounded.
The method for calculating the background value comprises the following steps: the median and quartile values are found using equations (4-20) for the data falling within each grid, i.e., the waveform and/or spectral data of the last 30 days of the current track. This median value is defined as the background value of the grid where the current track is located.
The resampling method comprises the following steps: the data (waveform and/or spectral data) for the current track falling within each grid is averaged using equation (4-20) and taken as the resample value within the grid area.
A method of calculating a threshold for determining an abnormality, see equation (4-21): using the sum, a Quartile difference IQR (Inter-Quartile Range, equation (4-22)) is calculated.
Bc＝Bm±IQR (4-21)
The fourth step: mapping and output data processing reports.
(1) Calculating the revisit track data (ascending track and descending track) of the current 5 days and B obtained by the revisit track data of the previous 30 daysm、BcThe sequence, i.e. the resampled time sequence, distinguishes three components from three frequency bands and is plotted. And meanwhile, annotating the magnetic emotion index of the current track time on the icon.
(2) And (3) drawing the seismic event sequences of which the seismic levels are more than 5, 6 and more than 7 in the first 30 days of the longitude area range of plus or minus 5 units of each orbit (the seismic event selection principle of a specific icon refers to the drawing mode of the seismic event sequences in the DEMETER satellite quick view).
(3) And drawing the global position of the current orbit and marking the seismic events with the magnitude greater than 6 occurring within 7 days before the annotation.
(4) And generating a data processing report.
4.42-level data Generation of 4-level data
Based on the 2-level data product of the induction magnetometer, 4-level data products and processing reports are obtained according to the following steps, and the flow is shown in fig. 13. The method for obtaining the background value is described by taking the median of the 30-day revisit orbit data as an example.
The green network area in fig. 13 illustrates each grid size of the global area: 2.55. the blue line represents the satellite down-track trajectory and the red line represents the satellite up-track trajectory. Within a 5 day regression period, the distance between two adjacent ascending (descending) tracks is about 530 km.
The first step is as follows: and (5) dividing the grids.
The globe is divided into grids at spatial intervals of 2.5 ° (latitude) x 5 ° (longitude).
The second step is that: and (5) analyzing the spatial distribution.
(1) Calculating the background value
According to 30-day time intervals and 5 degrees (longitude) multiplied by 2.5 degrees (latitude) space intervals, 2-level data of magnetic fields above the global and Chinese areas (0-60 degrees latitude; 60-140 degrees longitude) are selected in a sliding mode, and the median, quartile () and the fractional difference of all tracks in each space interval and an abnormal judgment value are calculated; interpolation is carried out on the median values in all the intervals to obtain global and national spatial distribution background fields;
(2) current 5 days changing magnetic field data
For the current 5-day varying magnetic field data, the median and quartile difference of the varying magnetic field in each grid are calculated with reference to equation (4-20).
Calculating and observing the median of all the tracks of the magnetic field in each space interval according to the current 5-day time interval and the space interval of 5 degrees (longitude) multiplied by 2.5 degrees (latitude); interpolating the median values in each interval to obtain global and national spatial distribution maps, and updating the global and national spatial distribution maps every day in a sliding manner;
(3) difference between current 5 days and background field, and normalizing
And calculating the difference value of the median value in each interval of each day and the background median value of the previous 30 days, normalizing by using the maximum value in each grid, and interpolating the normalized values of the intervals to obtain the global and national daily dynamic change spatial distribution map.
The third step: and adding the seismic event and the magnetic index to generate a 4-level data processing report.
And marking the record of the seismic event with the grade more than 7 within 30 days on the basis of the drawing result of the second step. And writing the data errors and other exceptions in the processing process into a data processing report to generate a 4-level data product processing report.
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.