CN114839354A - Beidou/GPS soil humidity measurement method based on sliding algorithm and weighting strategy - Google Patents
Beidou/GPS soil humidity measurement method based on sliding algorithm and weighting strategy Download PDFInfo
- Publication number
- CN114839354A CN114839354A CN202210770775.1A CN202210770775A CN114839354A CN 114839354 A CN114839354 A CN 114839354A CN 202210770775 A CN202210770775 A CN 202210770775A CN 114839354 A CN114839354 A CN 114839354A
- Authority
- CN
- China
- Prior art keywords
- signal
- noise ratio
- soil humidity
- beidou
- gps
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 239000002689 soil Substances 0.000 title claims abstract description 117
- 238000004422 calculation algorithm Methods 0.000 title claims abstract description 57
- 238000000691 measurement method Methods 0.000 title claims abstract description 20
- 238000005259 measurement Methods 0.000 claims abstract description 75
- 238000000034 method Methods 0.000 claims abstract description 60
- 238000012935 Averaging Methods 0.000 claims abstract description 6
- 230000009977 dual effect Effects 0.000 claims abstract description 5
- 238000004364 calculation method Methods 0.000 claims description 36
- 238000001914 filtration Methods 0.000 claims description 19
- 238000012545 processing Methods 0.000 claims description 19
- 230000008569 process Effects 0.000 claims description 15
- 230000003044 adaptive effect Effects 0.000 claims description 11
- 238000005070 sampling Methods 0.000 claims description 11
- 230000009467 reduction Effects 0.000 claims description 9
- 238000010183 spectrum analysis Methods 0.000 claims description 7
- 238000012886 linear function Methods 0.000 claims description 6
- 238000007781 pre-processing Methods 0.000 claims description 6
- 238000001228 spectrum Methods 0.000 claims description 4
- 238000006243 chemical reaction Methods 0.000 claims description 3
- 238000013527 convolutional neural network Methods 0.000 claims description 3
- 230000008030 elimination Effects 0.000 claims description 3
- 238000003379 elimination reaction Methods 0.000 claims description 3
- 230000000737 periodic effect Effects 0.000 claims description 3
- 238000010187 selection method Methods 0.000 claims description 3
- 239000011159 matrix material Substances 0.000 description 16
- 238000009795 derivation Methods 0.000 description 8
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 4
- 238000007796 conventional method Methods 0.000 description 2
- 238000012937 correction Methods 0.000 description 2
- 238000001514 detection method Methods 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 239000000523 sample Substances 0.000 description 2
- 230000009897 systematic effect Effects 0.000 description 2
- 230000004075 alteration Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000012272 crop production Methods 0.000 description 1
- 238000001035 drying Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 239000003673 groundwater Substances 0.000 description 1
- 230000002262 irrigation Effects 0.000 description 1
- 238000003973 irrigation Methods 0.000 description 1
- 238000007726 management method Methods 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000012544 monitoring process Methods 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N33/00—Investigating or analysing materials by specific methods not covered by groups G01N1/00 - G01N31/00
- G01N33/24—Earth materials
- G01N33/246—Earth materials for water content
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
- G01S19/01—Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
- G01S19/13—Receivers
- G01S19/33—Multimode operation in different systems which transmit time stamped messages, e.g. GPS/GLONASS
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/02—Preprocessing
- G06F2218/04—Denoising
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A40/00—Adaptation technologies in agriculture, forestry, livestock or agroalimentary production
- Y02A40/10—Adaptation technologies in agriculture, forestry, livestock or agroalimentary production in agriculture
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Remote Sensing (AREA)
- General Physics & Mathematics (AREA)
- Radar, Positioning & Navigation (AREA)
- Life Sciences & Earth Sciences (AREA)
- Chemical & Material Sciences (AREA)
- Theoretical Computer Science (AREA)
- Mathematical Physics (AREA)
- Data Mining & Analysis (AREA)
- Health & Medical Sciences (AREA)
- General Life Sciences & Earth Sciences (AREA)
- Mathematical Analysis (AREA)
- Medicinal Chemistry (AREA)
- Analytical Chemistry (AREA)
- Biochemistry (AREA)
- General Health & Medical Sciences (AREA)
- Immunology (AREA)
- Pathology (AREA)
- Algebra (AREA)
- Computational Mathematics (AREA)
- Geology (AREA)
- Food Science & Technology (AREA)
- Mathematical Optimization (AREA)
- Environmental & Geological Engineering (AREA)
- Pure & Applied Mathematics (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- General Engineering & Computer Science (AREA)
- Computer Networks & Wireless Communication (AREA)
- Geophysics And Detection Of Objects (AREA)
- Position Fixing By Use Of Radio Waves (AREA)
Abstract
The invention discloses a Beidou/GPS soil humidity measurement method based on a sliding algorithm and a weighting strategy, which comprises the following steps: s1, establishing a reflected signal-to-noise ratio phase reference value database; s2, calculating the soil humidity measured by the Beidou system by using the soil humidity measurement model and combining the signal-to-noise ratio phase reference value in the signal-to-noise ratio phase reference value database of the Beidou; s3, calculating the soil humidity measured by the GPS system by using the soil humidity measurement model and combining the signal-to-noise ratio phase reference value in the signal-to-noise ratio phase reference value database of the GPS; and S4, combining the soil humidity measured by the Beidou and GPS dual system in a multi-frequency mode, and averaging the results. The method considers the influence of the signal-to-noise ratio selection length, the effective reflection high fluctuation of the satellite and the frequency characteristic difference of the satellite on the soil humidity measurement precision, effectively solves the problems of low measurement precision, unstable performance and the like of a single system, can meet the requirements of soil humidity measurement time and spatial resolution, and can meet the requirements of soil humidity measurement with high precision, high real-time performance and high stability.
Description
Technical Field
The invention relates to the technical field of navigation satellite remote sensing inversion, in particular to a Beidou/GPS soil humidity measuring method based on a sliding algorithm and a weighting strategy.
Background
The soil humidity, also called soil water content, can be used for the research of global atmospheric water circulation, is closely related to the material transfer and energy exchange among the global water circulation, the atmospheric circulation and the biological ecological circulation, and is an important component of the global ecological system. Meanwhile, the soil humidity not only influences regional rainfall, groundwater supply and weather change, but also can be used as an important index for agricultural drought monitoring, and information is provided for irrigation of crops. Therefore, the method and the device for real-time and effective soil humidity measurement are provided and established, which can provide effective information for water vapor ecological research, crop production, management decision and the like, and are one of the problems to be solved in the current soil humidity measurement.
Currently, there are three main types of soil moisture measurement methods. The first method is a traditional measurement method, mainly including a drying method and a probe measurement. Although the traditional method can obtain higher measurement precision, the method needs field operation and processing, has lower time and space convenience rate, cannot meet the requirement of large-range real-time measurement, and has higher labor cost. The second type is a soil humidity measuring method of a professional remote sensing satellite, the method utilizes an environment remote sensing satellite to measure, a large area range can be measured, time and spatial resolution are high, an optical lens and a sensor which are depended by the remote sensing satellite are easily influenced by the environment, and measurement failure or accuracy reduction based on the remote sensing satellite method can be caused by severe weather such as cloud, fog, haze and the like. The third type is soil humidity measurement based on a GNSS (Global Navigation Satellite System) Navigation Satellite, and the method utilizes the phase fluctuation of the signal-to-noise ratio of the reflection signal of the GNSS Satellite to invert the soil humidity, thereby acquiring the soil humidity measurement result in a large range in real time.
However, the existing method does not consider the influence of the signal-to-noise ratio selection length on the effective reflection of the computed satellite and the subsequent phase detection amount accuracy based on the signal-to-noise ratio, but adopts the fixed signal-to-noise ratio of 5 to 30 degrees to select the length, so that the final soil humidity measurement accuracy is reduced. Meanwhile, the influence of the high fluctuation of the effective reflection of the satellite along with the time is not considered in the conventional method, so that the accuracy of the signal-to-noise ratio phase of subsequent calculation is reduced. In addition, the existing soil humidity inversion usually adopts single frequency for inversion, and the influence of the difference of signal-to-noise ratios of different satellite frequency bands is ignored, so that the measurement precision is reduced.
Disclosure of Invention
The invention provides a Beidou/GPS soil humidity measurement method based on a sliding algorithm and a weighting strategy, which is used for researching the influence of the selection length of a signal-to-noise ratio, the effective reflection high fluctuation of a satellite and the frequency characteristic difference of the satellite on the soil humidity measurement precision, effectively solving the problems of low measurement precision, unstable performance and the like of a single system, meeting the requirements of time and space resolution of soil humidity measurement, meeting the requirements of measurement precision, instantaneity, stability and the like, and providing powerful support for soil humidity measurement.
In order to solve the technical problems, the technical scheme of the invention is as follows:
a Beidou/GPS soil humidity measurement method based on a sliding algorithm and a weighting strategy comprises the following steps:
s1, establishing a database of the phase reference value of the signal-to-noise ratio of the reflected signal
Connecting and collecting original observation values of the Beidou and the GPS on rainy days, wherein the original observation values comprise multi-frequency signal-to-noise ratios, pseudo ranges and carrier phases, calculating elevation angles of corresponding epoch moments, processing data of the original observation values, calculating to obtain signal-to-noise ratio phase reference values of the Beidou and the GPS, and storing the signal-to-noise ratio phase reference values as a database;
s2, calculating the current soil humidity of each signal-to-noise phase reference value by using a soil humidity measurement model and combining the signal-to-noise phase reference values in a signal-to-noise phase reference value database of the Beidou, and finally calculating the average value of the soil humidity obtained by calculating all frequencies of all the Beidou satellites to obtain the soil humidity measured by the Beidou system;
s3, calculating the current soil humidity of each signal-to-noise phase reference value by using the soil humidity measurement model and combining the signal-to-noise phase reference values in the signal-to-noise phase reference value database of the GPS, and finally, calculating the average value of the soil humidity obtained by calculating all frequencies of the GPS dual-frequency satellite to obtain the soil humidity measured by the GPS system;
and S4, combining the multi-frequency measurement of the Beidou and GPS dual systems, averaging the results, and outputting the measurement result.
Preferably, the method for obtaining the signal-to-noise ratio phase reference value in step S1 includes:
s1-1, processing the original data
S1-1-1, filtering the influence of the signal-to-noise ratio direct signal of the original observation value through a high-order Chebyshev polynomial fitting algorithm, and reserving the signal-to-noise ratio signal only containing the reflected signal and random noise;
s1-1-2, carrying out noise reduction processing on the signal-to-noise ratio signal processed in the step S1-1-1 through a Gaussian low-pass filtering noise reduction algorithm to obtain a signal-to-noise ratio signal only retaining a reflection signal;
s1-1-3, calculating the effective reflection height of the satellite through an adaptive sliding window and a weighting constraint strategy;
s1-2, calculating a signal-to-noise phase reference value through an extended Kalman filtering algorithm according to the satellite effective reflection height obtained in the step S1-1-3 and by combining the signal-to-noise ratio signal of the reflection signal obtained in the step S1-1-2.
Preferably, the principle of the high-order chebyshev polynomial fitting algorithm in the step S1-1-1 is as follows:
in a time period t 0 , t 0 +Δt]The signal-to-noise ratio signal is fitted with a Chebyshev polynomial, where t 0 Is the starting epoch of the snr, and at is the snr selected length,
to normalize the signal-to-noise ratio, the time is converted as follows:
at this time, the parameter τ ϵ [ -1, 1], and thus, the signal-to-noise ratio can be expressed as a Chebyshev polynomial:
in the formula,the signal-to-noise ratio of the fitted direct signal is represented,nfor the order of the chebyshev polynomial,C i chebyshev polynomial coefficients, chebyshev polynomials, being signal-to-noise ratio componentsT i The recurrence formula is:
after fitting, the direct signal can be eliminated to obtain a signal only retaining the signal-to-noise ratio of the reflected signal and the high-frequency random noise, and the direct signal-to-noise ratio elimination formula is expressed as follows:
in the formula:S r_n for signals that contain only the signal-to-noise ratio and noise of the reflected signal,S o representing the signal-to-noise ratio that was originally just received,S d the signal-to-noise ratio of the direct signal of the corresponding point is obtained after Chebyshev polynomial fitting.
Preferably, the gaussian low-pass filtering noise reduction algorithm processing in the step S1-1-2 is expressed as follows:
assuming that the signal to be processed containing the signal-to-noise ratio and noise of the reflected signal is represented asS r_n (x) And the signal after filtering which only contains the signal-to-noise ratio of the reflected signal is expressed asS r The calculation procedure is as follows:
in the formula,G(x) Is a gaussian function, the symbol denotes the convolution operation,xrepresenting the SNR sequence, the high-frequency random noise can be effectively filtered through the processing of the process, the signal only retaining the SNR of the reflected signal is obtained,
wherein G (x) is a Gaussian kernel function having a width ofIt is decided that,can be calculated from the standard deviation of the signal-to-noise ratio, as a distribution parameter of the gaussian function,eandπare natural constants respectively.
Preferably, in step S1-1-3, the method for acquiring the satellite effective reflection is as follows:
(1) determination of signal-to-noise ratio signal by adaptive sliding algorithmS r The data set is selected from a set of data,
the selection method comprises the following steps: the starting angle of the altitude angle is still 5 degrees, but the ending angle is selected in a sliding way between 25 and 35 degrees, a fixed value is not adopted, but the data groups are selected in a sliding way in a self-adaptive way to be respectively calculated in the middle, the selection rule is determined by the data sampling rate,
if the data sampling rate is 30s, the signal-to-noise ratio between 25 degrees and 35 degrees needs to be calculated in groups, when 20 epochs exist between 25 degrees and 35 altitudes, the signals are divided into 20 groups, the first group of data consists of epochs between 5 degrees and 25 degrees of altitude plus the first epoch between 25 degrees and 35 degrees, the second group of data consists of epochs between 5 degrees and 25 degrees of altitude plus the first two epochs between 25 degrees and 35 degrees, and the like, so that 20 groups of data to be processed are formed;
if the data sampling rate is 15s, the step length is determined to be 2, namely, a group of data is selected every other epoch for calculation;
if the data sampling rate is 5s, the step length is 3, that is, a group of data is selected for calculation every three epochs. If the time is 1s, the step length is 5, a group of data is selected for calculation every five epochs,
after the adaptive sliding algorithm processing, the signal-to-noise ratio signal to be processed is divided intomGroup, can be represented as: That is to say havemThe signals with different signal-to-noise ratio lengths can be used for calculating the effective reflection height of the satellite;
(2) the Lomb-Scargle spectrum analysis method is adopted to solve the satellite effective reflection height for each group of signal-to-noise ratio signals,
wherein, Lomb-Sacragle spectrum analysis method is used for solving the power spectrum of the signal, and the expression is as follows:
in the formula,P x ( f ) Is at a frequency off The power of the periodic signal of (a);S r (x) In order to contain the signal-to-noise ratio of the reflected signal,Nin order to be the length of the sequence,τthe time shift invariant can be calculated by the following formula:
can find outP x ( f ) At the time of peakf The value is calculated according to the relation between the frequency and the effective reflection heighth e The conversion relation is:
in the formula,λis the wavelength of the signal, and is a known quantity, therefore, the effective reflection height of the satellite can be obtained by the above formulah e ,
Repeating the above calculation processmRespectively resolving the group signal-to-noise ratio signals to obtainmThe effective reflection is high, expressed as:the satellite to be calculated isThe effective reflection height is averaged to finally obtain the effective reflection height of the satellite in the dayWhereiniTo show the day
(3) Repeating the steps (1) and (2), resolving data of ten days before the measurement day, and obtaining the satellite effective reflection height of nearly ten days, which is expressed asWeighting the satellite effective reflection height of the ten days by adopting a weighting constraint strategy, wherein a specific weighting model is expressed as follows:
in the formula,w i which represents the weighting coefficient(s) of the,idata representing the number of days, ten days in total, measured on the day of the dayi=10, and so on, the weights decrease gradually.
Preferably, in step S1-2, the signal-to-noise ratio phase reference value obtaining method includes:
high effective reflection through satelliteCombining the signal-to-noise ratio signal containing the reflected signalThe phase value of the signal-to-noise ratio can be obtained by the relationship between the two, and the specific formula is as follows:
in the formula,Aandthe amplitude and phase representing the signal-to-noise ratio are of the formulaAn unknown quantity.EIs the satellite altitude, a known quantity.
The solution is performed by using an extended kalman filter algorithm, as follows:
first, for a nonlinear system, the state equation and observation equation of its state space can be expressed as:
wherein the first term is a state equation, the second term is an observation equation,x k andz k the actual state vector and the observation vector, respectively. Wherein,x k is formed byAAndthe unknown state quantity of the unknown quantity composition can be expressed asx k (A,),z k Is formed byThe observed state quantity, i.e. the output quantity of the system, of known quantity composition is expressed asz k (),u k As a function of the system, as input quantities to the equation of state, from known quantitiesAndEthe components of the composition are as follows,kwhich is indicative of a sequence of signals,f to relate tok-1 andka non-linear function of the state of the system at the moment,his a state quantityx k And observed quantityz k The non-linear function of interest is,w k andv k respectively process noise and observationsNoise, and is white Gaussian noise with mutually independent mean value of zero,
in practical applications, it cannot be determinedw k Andv k the specific value at each time instant, but the estimate of the state vector and the observed quantity for which the noise value comes at time k can be ignored, and therefore the above equation can be approximately written as:
the extended kalman filter algorithm transforms the nonlinear system into a linear system by performing a taylor first order expansion of the nonlinear system near the optimal estimation point of its state. In this case, the above formula can be expressed as:
therefore, the signal-to-noise ratio phase reference value can be obtained by solvingSum amplitude valueA. In this case only phase values 。
Preferably, the flow of the taylor first-order expansion algorithm is as follows:
1) inputting initial conditions, inputtingAnd, is thatk-1 a posteriori state estimates of the system at time instance,is thatk-1 time instant systematic error covariance;
2) parameter prediction, the state prediction equation is:
whereinis thatkThe posterior state estimation of the time system, the error covariance prediction equation is: , is thatkState of the momentThe covariance of the prediction is determined by the prediction,Ais thatfFunction pairxThe jacobian matrix after the partial derivation is solved,Wis thatfTo pairwThe jacobian matrix after the partial derivation is solved,Qis the process noise covariance matrix.
A correction module, the system gain equation expressed as:
whereinHis thathTo pairxThe jacobian matrix after the partial derivation is solved,Vis thathTo pairvThe jacobian matrix after the partial derivation is solved,Ris an observed noise covariance matrix, where the filter equation is expressed as:
in the formula, the meaning of the parameters is introduced before, and the error covariance update equation is expressed as:
Preferably, the specific method of step S2 is as follows:
s2-1, reading the signal-to-noise ratio phase reference value of the Beidou system from the signal-to-noise ratio phase reference value database, and performing combined calculation by combining the signal-to-noise ratio phase value of the day to obtain the soil humidity of the measurement day, wherein the calculation formula is as follows:
in the formula,the difference between the signal-to-noise ratio phase value and the signal-to-noise ratio phase reference value on the measurement day is obtained by preprocessing the original observation value obtained by the Beidou system on the measurement day through VWC, VWC represents the soil humidity required, and the signal-to-noise ratio phase value on the measurement day is obtained through step S1;
s2-2, respectively calculating corresponding soil humidity by using data of three frequencies of all Beidou observation satellites, and then calculating an average value to obtain the final soil humidity measured by the Beidou system, wherein the final soil humidity is expressed as VWC B 。
Further, three frequencies of the Beidou satellite are respectively as follows: b1 is 1575.42MHz, B2 is 1176.45MHz, and B3 is 1268.52MHz, which belongs to the common knowledge in the field.
Preferably, the specific method of step S3 is as follows:
s3-1, reading the signal-to-noise ratio phase reference value of the GPS system from the signal-to-noise ratio phase reference value database, and performing combined calculation by combining the signal-to-noise ratio phase value of the day to obtain the soil humidity of the measurement day, wherein the calculation formula is as follows:
in the formula,the difference between the signal-to-noise ratio phase value and the signal-to-noise ratio phase reference value on the measurement day is obtained by preprocessing the original observation value acquired by the GPS on the measurement day through the step S1, wherein the VWC represents the soil humidity required by the measurement day;
s3-2, respectively calculating corresponding soil humidity by using dual-frequency data of all observation satellites of the GPS, and then calculating an average value to obtain the final soil humidity measured by the GPS, wherein the final soil humidity is expressed as VWC G 。
It should be noted that, currently, all satellites of the GPS can transmit dual-frequency signals, and only some satellites transmit triple-frequency signals. Therefore, to maintain uniformity, measurements are taken here using both frequency signals. The dual-frequency signal frequencies of the GPS satellite are respectively as follows: l1 is 1575.42MHz, L2 is 1227.60MHz, and is common knowledge in the field.
Preferably, the implementation method of step S4 is:
and (4) solving the soil humidity of the final measurement day by using the Beidou and GPS soil humidity obtained by calculation in the step S3 and the step S4, namely, taking the average value of the Beidou and GPS soil humidity, and expressing as follows:
in the formula, VWC B Is the soil humidity, VWC, measured by the Beidou system G Is the soil humidity, VWC, measured by the GPS system Final (a Chinese character of 'gan') I.e. the final measured soil moisture.
The invention has the following characteristics and beneficial effects:
by adopting the technical scheme, the existing Beidou/GPS navigation system satellite is utilized to measure the soil humidity, and compared with the existing method based on a remote sensing satellite measurement method and a probe measurement method, the method can solve the problems of high measurement cost, high environmental influence and the like, and can solve the problems of low measurement time and spatial resolution, incapability of meeting large-scale measurement and the like. Compared with the existing method based on the signal-to-noise ratio of the GPS signal, the method can solve the problems that the high precision of effective reflection of the satellite is influenced by the selection length of the signal-to-noise ratio, and the like, and can also solve the problems that the effective reflection of the satellite is high in fluctuation, the measurement precision is low and the measurement is unstable due to neglecting the frequency characteristic difference. Meanwhile, the Beidou and the GPS are combined to carry out dual-system multi-frequency combined measurement, so that the accuracy and the stability of the measurement can be guaranteed, and the soil humidity measurement with high precision, high real-time performance and high stability is met.
In addition, the signal-to-noise ratio of the direct signal is eliminated by using a high-order Chebyshev polynomial fitting method, denoising is carried out by using a Gaussian low-pass filtering denoising algorithm, and high-frequency random noise is filtered out, so that a signal only containing the signal-to-noise ratio of the reflected signal is obtained. Because the signal-to-noise ratio selection length directly influences the precision of the effective reflection height of the follow-up satellite, the effective reflection height of each frequency of each satellite is independently calculated by adopting a self-adaptive sliding algorithm, and the effective reflection heights of the satellite for continuous ten days are integrated by adopting a weighting constraint strategy, so that the precision of the effective reflection height of the satellite is effectively improved. And the problem of nonlinear solution is solved by adopting an extended Kalman filtering algorithm, and the solution accuracy of the phase value of the signal-to-noise ratio of the reflection signal is improved. The difference between different frequency characteristics of different satellites is considered, and the accuracy of soil humidity measurement is guaranteed by adopting a method of independently calculating all frequencies of all satellites and then calculating the mean value. And the Beidou and GPS navigation systems are combined to carry out dual-system multi-frequency combined detection, so that the precision, the real-time performance and the stability of soil humidity measurement are further ensured.
Drawings
In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings used in the description of the embodiments or the prior art will be briefly described below, and it is obvious that the drawings in the following description are only some embodiments of the present invention, and for those skilled in the art, other drawings can be obtained according to these drawings without creative efforts.
Fig. 1 is a flow chart for establishing a snr-phase reference value database based on a sliding algorithm and a weighting strategy according to the present embodiment.
Fig. 2 is a flowchart of the Beidou/GPS soil humidity measurement method based on the sliding algorithm and the weighting strategy in the embodiment.
Detailed Description
It should be noted that the embodiments and features of the embodiments may be combined with each other without conflict.
In the description of the present invention, it is to be understood that the terms "center", "longitudinal", "lateral", "up", "down", "front", "back", "left", "right", "vertical", "horizontal", "top", "bottom", "inner", "outer", and the like, indicate orientations or positional relationships based on those shown in the drawings, and are used only for convenience in describing the present invention and for simplicity in description, and do not indicate or imply that the referenced devices or elements must have a particular orientation, be constructed and operated in a particular orientation, and thus, are not to be construed as limiting the present invention. Furthermore, the terms "first", "second", etc. are used for descriptive purposes only and are not to be construed as indicating or implying relative importance or implicitly indicating the number of technical features indicated. Thus, a feature defined as "first," "second," etc. may explicitly or implicitly include one or more of that feature. In the description of the present invention, "a plurality" means two or more unless otherwise specified.
In the description of the present invention, it should be noted that, unless otherwise explicitly specified or limited, the terms "mounted," "connected," and "connected" are to be construed broadly, e.g., as meaning either a fixed connection, a removable connection, or an integral connection; can be mechanically or electrically connected; they may be connected directly or indirectly through intervening media, or they may be interconnected between two elements. The specific meaning of the above terms in the present invention can be understood by those of ordinary skill in the art through specific situations.
The invention provides a Beidou/GPS soil humidity measurement method based on a sliding algorithm and a weighting strategy, which comprises the following steps as shown in figure 1:
s1, establishing a phase reference value database of the signal-to-noise ratio of the reflected signal
The method comprises the steps of connecting and collecting original observation values of the Beidou and the GPS in rainy days, wherein the original observation values comprise multi-frequency signal-to-noise ratios, pseudo ranges and carrier phases, calculating elevation angles of corresponding epoch moments, processing data of the original observation values, calculating to obtain signal-to-noise ratio phase reference values of the Beidou and the GPS, and storing the signal-to-noise ratio phase reference values as a database.
It should be noted that the signal-to-noise ratio may be directly extracted from the observation file, and the altitude information corresponding to the epoch may be calculated by the pseudorange, the carrier phase, and the satellite ephemeris. This process is common knowledge in the field of positioning and the calculation process is not listed here in detail. For the convenience of subsequent presentation, the raw SNR extracted here is used in subsequent stepsS o Indicating, for calculated altitude angleEAnd (4) showing.
Specifically, as shown in fig. 2, the method includes the following sub-steps:
s1-1, processing the original data
S1-1-1, filtering the influence of the signal-to-noise ratio direct signal of the original observation value through a high-order Chebyshev polynomial fitting algorithm, and reserving the signal-to-noise ratio signal only containing the reflected signal and random noise. The chebyshev polynomial fitting is based on the principle of the best approximation of a function, which is formed by fitting a chebyshev polynomial as a basis function, using the function as a minimum of the sum of the variances between the function value at a given point and the known function value at that point, and the principle of the chebyshev polynomial fitting algorithm is expressed as follows:
in the time periodt 0 , t 0 +Δt]The signal-to-noise ratio signal is fitted with a Chebyshev polynomial, wheret 0 Is the starting epoch of the signal-to-noise ratio, and ΔtThe length is chosen for the signal-to-noise ratio,
to normalize the signal-to-noise ratio, the time is converted as follows:
at this time, the parameter τ ϵ [ -1, 1], and thus, the signal-to-noise ratio can be expressed as a Chebyshev polynomial:
in the formula,the signal-to-noise ratio of the fitted direct signal is represented,nin order to improve the fitting accuracy of the signal-to-noise ratio, in this embodiment, a polynomial of order 5 is selected for fitting, that is, n =5, so as to improve the fitting accuracy of the signal-to-noise ratio.C i Chebyshev polynomial coefficients, chebyshev polynomials, being signal-to-noise ratio componentsT i The recurrence formula is:
after fitting, the direct signal can be eliminated to obtain a signal only retaining the signal-to-noise ratio of the reflected signal and the high-frequency random noise, and the direct signal-to-noise ratio elimination formula is expressed as follows:
in the formula:S r_n for signals that contain only the signal-to-noise ratio and noise of the reflected signal,S o representing the signal-to-noise ratio that was originally just received,S d the signal-to-noise ratio of the direct signal of the corresponding point is obtained after Chebyshev polynomial fitting.
S1-1-2, noise reduction processing is carried out on the signal-to-noise ratio signal processed in the step S1-1-1 through a Gaussian low-pass filtering noise reduction algorithm, and a signal-to-noise ratio signal only retaining the reflection signal is obtained. Since the high-frequency random noise included in the signal-to-noise ratio obeys normal distribution, the gaussian low-pass filter is used for denoising in the embodiment, so that the denoising effect is optimal. The gaussian low-pass filtering algorithm is essentially a linear smooth filtering algorithm, and mainly selects a weight value by using the shape of a gaussian function, so that not only can noise be accurately filtered, but also useful signals can be completely reserved. Since the snr signal in this embodiment is one-dimensional, a one-dimensional gaussian function with zero mean is used for filtering, which is expressed as follows:
assuming that the signal to be processed containing the signal-to-noise ratio and noise of the reflected signal is represented asS r_n (x) The signal after filtering which only contains the signal-to-noise ratio of the reflected signal is represented asS r The calculation procedure is as follows:
in the formula,G(x) Is a gaussian function, the symbol denotes the convolution operation,xrepresenting the SNR sequence, the high-frequency random noise can be effectively filtered through the processing of the process, the signal only retaining the SNR of the reflected signal is obtained,
in the formula,G(x) Is a Gaussian kernel function, the width of the Gaussian function isIt is decided that,can be calculated from the standard deviation of the signal-to-noise ratio, as a distribution parameter of the gaussian function,eandπrespectively, are natural constants.
By the technical scheme, high-frequency random noise can be effectively filtered, a signal only retaining the signal-to-noise ratio of the reflection signal is obtained, and the signal is expressed in the subsequent steps.
And S1-1-3, calculating the effective reflection height of the satellite through an adaptive sliding window and a weighted constraint strategy.
The process can be divided into three steps, wherein the first step is to determine the signal-to-noise ratio selection length, and the second step is to solve the satellite effective reflection height through a Lomb-Scargle spectrum analysis method. And thirdly, obtaining the final effective reflection height of the satellite by adopting a weighting constraint strategy. The specific process is as follows:
(1) determination of signal-to-noise ratio signal by adaptive sliding algorithmS r The selected data set.
It should be noted that, in the conventional method, the snr is selected to solve at an altitude angle of 5 to 30 degrees. However, the snr characteristics vary greatly from satellite to satellite and from frequency to frequency, resulting in a difference in effective reflection for satellites calculated at different snr lengths. Therefore, an adaptive sliding algorithm is used for the calculation. The specific process is as follows: the starting angle of the elevation angle is still 5 degrees, but the ending angle is between 25 and 35 degrees, and is not calculated with a fixed value, but with an adaptively chosen value between them. The selection rule is determined by the data sampling rate.
The specific selection method comprises the following steps: the starting angle of the elevation angle is still 5 degrees, but the ending angle is chosen in a sliding manner between 25 and 35 degrees, and is not calculated by using a fixed value, but by using adaptive sliding data sets. The selection rule is determined by the data sampling rate. Such as: if the data sampling rate is 30s, the snr of 25-35 degrees is calculated in groups (i.e. assuming 20 epochs between 25-35 altitudes, the signal is divided into 20 groups, the first group of data is composed of epochs at 5-25 degrees in altitude plus the first epoch at 25-35 degrees, the second group of data is composed of epochs at 5-25 degrees in altitude plus the first two epochs at 25-35 degrees, and so on to form 20 groups of data to be processed). If the time is 15s, the step length is determined to be 2, i.e. a group of data is selected for calculation every other epoch. If it is 5s, the step size is 3, i.e. a group of data is selected for calculation every three epochs. If the time is 1s, the step length is 5, and a group of data is selected for calculation every five epochs. After the adaptive sliding algorithm processing, it is assumed that the signal-to-noise ratio signal to be processed is divided intomGroups, which can be represented as:that is to say havemA set of signals of different signal-to-noise ratio lengths can be used to calculate the satellite effective reflection height.
(2) The effective reflection height of the satellite is obtained by adopting a Lomb-Scargle spectrum analysis method for each group of signal-to-noise ratio signals, the power spectrum of the signals is firstly calculated, and the corresponding frequency when the energy is highest is found outf. Then pass throughfAnd solving the corresponding satellite effective reflection height according to the relation between the value and the satellite effective reflection height.
Wherein, Lomb-Sacragle spectrum analysis method is used for solving the power spectrum of the signal, and the expression is as follows:
in the formula,P x ( f ) Is at a frequency off The power of the periodic signal of (a);S r (x) In order to contain the signal-to-noise ratio of the reflected signal,Nin order to be the length of the sequence,τfor time shift invariants, this can be calculated by:
can find P x (f) F value when reaching peak value, and then calculating the effective reflection height of the corresponding satellite according to the relationship between the frequency and the effective reflection heightThe conversion relation is:
in the formula, λ is the wavelength of the signal and is a known quantity, and therefore, the effective reflection height of the satellite can be obtained by the above formula,
It should be noted that λ is known quantity and can be directly obtained, and the wavelength of the GPS L1 frequency band is: 19.0cm, and the wavelength of the Beidou B1 frequency band is 19.2 cm. Belonging to the common knowledge in the field. Therefore, the effective reflection height of the satellite can be obtained by the above equation.
Repeating the above calculation processmRespectively resolving the group signal-to-noise ratio signals to obtainmThe effective reflection is high, expressed as:averaging the calculated satellite effective reflection heights to obtain the satellite effective reflection height of the dayWhereiniRepresents the day of the day;
(3) repeating the steps (1) and (2), resolving data of ten days before the measurement day, and obtaining the satellite effective reflection height of nearly ten days, which is expressed asWeighting the satellite effective reflection height of the ten days by adopting a weighting constraint strategy, wherein a specific weighting model is expressed as follows:
in the formula,w i which represents the weighting coefficient(s) of the,idata representing the number of days, ten days in total, measured on the day of the dayi=10, and so on, the weights decrease gradually.
It will be appreciated that the weighting patterns are weighted according to the proximity of time, i.e. the closer to the time of the measurement day, the higher the weight, the further away the time the lower the weight.
S1-2, calculating a signal-to-noise phase reference value through an extended Kalman filtering algorithm according to the satellite effective reflection height obtained in the step S1-1-3 and by combining the signal-to-noise ratio signal of the reflection signal obtained in the step S1-1-2.
High effective reflection through satelliteCombining the signal-to-noise ratio signal containing the reflection signal, the phase value of the signal-to-noise ratio can be obtained through the relation between the signal-to-noise ratio signal and the reflection signal, and the specific formula is as follows:
in the formula, A andthe amplitude and phase, which represent the signal-to-noise ratio, are unknowns in the formula.
Furthermore, considering that the signal-to-noise ratio and the altitude angle are in a nonlinear relationship at the moment, in order to improve the precision of parameter calculation, an extended kalman filter algorithm is adopted for the calculation. The extended kalman filter algorithm is called by an existing function module in the Matlab platform, and here we only make a simple algorithm description, as shown below:
first, for a nonlinear system, the state equation and observation equation of its state space can be expressed as:
wherein the first term is a state equation, the second term is an observation equation,x k andz k the actual state vector and the observation vector, respectively. Wherein,x k is formed byAAndthe unknown state quantity of the unknown quantity composition can be expressed asx k (A,);z k Is formed byThe observed state quantity, i.e. the output quantity of the system, of known quantity composition is expressed asz k ();u k As a function of the system, as input quantities to the equation of state, from known quantitiesAndEand (4) forming.kWhich is indicative of a sequence of signals,f to relate tok-1 andka non-linear function of the state of the system at the moment,his a state quantityx k And observed quantityz k The non-linear function of interest is,w k andv k respectively process noise and observation noise, and is Gaussian white noise with mutually independent mean value of zero,
in practical applications, w cannot be determined k And v k The specific value at each time instant, but the estimate of the state vector and the observed quantity for which the noise value comes at time k can be ignored, and therefore the above equation can be approximately written as:
the extended kalman filter algorithm transforms the nonlinear system into a linear system by performing a taylor first order expansion of the nonlinear system near the optimal estimation point of its state. In this case, the above formula can be expressed as:
Wherein, the flow of the Taylor first-order expansion algorithm is as follows:
1) inputting initial conditions, inputtingAnd, is thatk-1 a posteriori state estimates of the system at time instance,is thatk-1 time instant systematic error covariance;
2) parameter prediction, the state prediction equation is:
wherein, in the process,is thatkThe posterior state estimation of the time system, the error covariance prediction equation is:
,is thatkState of the momentThe covariance of the prediction is determined by the prediction,Ais thatfFunction pairxThe jacobian matrix after the partial derivation is solved,Wis thatfFor is towThe jacobian matrix after the partial derivation is solved,Qis the process noise covariance matrix;
a correction module, the system gain equation expressed as:
whereinHis thathTo pairxAfter the derivation of the partial derivatives, the jacobian matrix is obtained,Vis thathTo pairvThe jacobian matrix after the partial derivation is solved,Ris an observed noise covariance matrix, where the filter equation is expressed as:
in the formula, the meaning of the parameters is introduced before, and the error covariance update equation is expressed as:
It should be noted that, in this embodiment, the signal-to-noise ratio phase reference value is obtained through extended kalman filter calculation, and may be obtained through a least square algorithm.
It can be understood that, in the above technical solution, the signal-to-noise ratio phase reference value database is established. By continuously acquiring the navigation data of the Beidou and GPS systems in non-rainy days, and by utilizing the process of the step S1, the sub-system respectively calculates the corresponding signal-to-noise ratio phase values according to the satellite sub-frequency, establishes a database according to the system type, the satellite number and the frequency number and stores the database on a local computer in a text form. For example, the phase value of the B1 frequency band of the Beidou B01 satelliteThe phase value of B2 frequency band is. By analogy, a signal-to-noise ratio phase reference value database of multi-frequency of the Beidou and GPS dual system can be established.
S2, calculating the current soil humidity of each signal-to-noise ratio phase reference value by using the soil humidity measurement model and combining the signal-to-noise ratio phase reference value database of the Beidou, and finally calculating the average value of the soil humidity obtained by calculating all frequencies of all satellites of the Beidou to obtain the soil humidity measured by the Beidou system.
S2-1, reading the signal-to-noise ratio phase reference value of the Beidou system from the signal-to-noise ratio phase reference value database, and performing combined calculation by combining the signal-to-noise ratio phase value of the day to obtain the soil humidity of the measurement day, wherein the calculation formula is as follows:
in the formula,the difference between the signal-to-noise ratio phase value and the signal-to-noise ratio phase reference value on the measurement day is obtained by processing an original observation value acquired through a Beidou system on the measurement day through step S1, wherein VWC represents the soil humidity required by the measurement day;
it should be noted that the formula is an empirical formula for solving the soil moisture by the signal-to-noise ratio, and belongs to the common knowledge in the field. Therefore, this embodiment will not be described in detail.
S2-2, respectively calculating corresponding soil humidity by using data of three frequencies of all Beidou observation satellites, and then solving the mean value to obtain the final soil humidity measured by the Beidou system, wherein the final soil humidity is expressed as VWC B 。
Further, three frequencies of the Beidou satellite are respectively as follows: b1 is 1575.42MHz, B2 is 1176.45MHz, and B3 is 1268.52MHz, which belongs to the common knowledge in the field.
S3, calculating the current soil humidity of each signal-to-noise phase reference value by using a soil humidity measurement model and combining the signal-to-noise phase reference value database of the GPS, and finally, calculating the average value of the soil humidity obtained by calculating all frequencies of the GPS dual-frequency satellite to obtain the soil humidity measured by the GPS system;
s3-1, reading the signal-to-noise ratio phase reference value of the GPS system from the signal-to-noise ratio phase reference value database, and performing combined calculation by combining the signal-to-noise ratio phase value of the day to obtain the soil humidity of the measurement day, wherein the calculation formula is as follows:
in the formula,the VWC represents the soil humidity required, and the signal-to-noise ratio phase value on the measurement day is a value obtained by preprocessing the original observation value acquired by the GPS on the measurement day through step S1.
S3-2, respectively calculating corresponding soil humidity by using dual-frequency data of all observation satellites of the GPS, and then calculating an average value to obtain the final soil humidity measured by the GPS, wherein the final soil humidity is expressed as VWC G 。
It should be noted that, currently, all satellites of the GPS transmit dual-frequency signals, and only some satellites transmit triple-frequency signals. Therefore, to maintain uniformity, measurements are taken here using both frequency signals. The dual-frequency signal frequencies of the GPS satellite are respectively as follows: l1 is 1575.42MHz, L2 is 1227.60MHz, and is common knowledge in the field.
And S4, combining the multi-frequency measurement of the Beidou and GPS dual systems, averaging the results, and outputting the measurement result.
And (4) solving the soil humidity of the final measurement day by using the Beidou and GPS soil humidity obtained by calculation in the step S3 and the step S4, namely, taking the average value of the Beidou and GPS soil humidity, and expressing as follows:
in the formula, VWC B Is the soil humidity, VWC, measured by the Beidou system G Is the soil humidity, VWC, measured by the GPS system Final (a Chinese character of 'gan') I.e. the final measured soil moisture.
The embodiments of the present invention have been described in detail with reference to the accompanying drawings, but the present invention is not limited to the described embodiments. It will be apparent to those skilled in the art that various changes, modifications, substitutions and alterations can be made in these embodiments, including the components, without departing from the principles and spirit of the invention, and still fall within the scope of the invention.
Claims (9)
1. A Beidou/GPS soil humidity measurement method based on a sliding algorithm and a weighting strategy is characterized by comprising the following steps:
s1, establishing a phase reference value database of the signal-to-noise ratio of the reflected signal
Connecting and collecting original observation values of the Beidou and the GPS on rainy days, wherein the original observation values comprise multi-frequency signal-to-noise ratios, pseudo ranges and carrier phases, calculating elevation angles of corresponding epoch moments, processing data of the original observation values, calculating to obtain signal-to-noise ratio phase reference values of the Beidou and the GPS, and storing the signal-to-noise ratio phase reference values as a database;
s2, calculating the current soil humidity of each signal-to-noise phase reference value by using a soil humidity measurement model and combining the signal-to-noise phase reference values in a signal-to-noise phase reference value database of the Beidou, and finally calculating the average value of the soil humidity obtained by calculating all frequencies of all the Beidou satellites to obtain the soil humidity measured by the Beidou system;
s3, calculating the current soil humidity of each signal-to-noise phase reference value by using the soil humidity measurement model and combining the signal-to-noise phase reference values in the signal-to-noise phase reference value database of the GPS, and finally, calculating the average value of the soil humidity obtained by calculating all frequencies of the GPS dual-frequency satellite to obtain the soil humidity measured by the GPS system;
and S4, combining the multi-frequency measurement of the Beidou and GPS dual systems, averaging the results, and outputting the measurement result.
2. The Beidou/GPS soil humidity measurement method based on sliding algorithm and weighting strategy as claimed in claim 1, wherein the acquisition method of signal-to-noise ratio phase reference value in step S1 is:
s1-1, processing the original data
S1-1-1, filtering the influence of a signal-to-noise ratio direct signal in an original observation value through a high-order Chebyshev polynomial fitting algorithm, and reserving the signal-to-noise ratio signal only containing a reflected signal and random noise;
s1-1-2, carrying out noise reduction processing on the signal-to-noise ratio signal processed in the step S1-1-1 through a Gaussian low-pass filtering noise reduction algorithm to obtain a signal-to-noise ratio signal only retaining a reflection signal;
s1-1-3, calculating the effective reflection height of the satellite through an adaptive sliding window algorithm and a weighting constraint strategy;
s1-2, calculating a signal-to-noise phase reference value through an extended Kalman filtering algorithm according to the satellite effective reflection height obtained in the step S1-1-3 and by combining the signal-to-noise ratio signal of the reflection signal obtained in the step S1-1-2.
3. The Beidou/GPS soil moisture measurement method based on sliding algorithm and weighting strategy according to claim 2, characterized in that the principle of the high-order Chebyshev polynomial fitting algorithm in the step S1-1-1 is expressed as follows:
in the time periodt 0 , t 0 +Δt]The signal-to-noise ratio signal is fitted with a Chebyshev polynomial, wheret 0 Is the starting epoch of the signal-to-noise ratio, and ΔtThe length is chosen for the signal-to-noise ratio,
to normalize the signal-to-noise ratio, the time is converted as follows:
at this time, the parameter τ ϵ [ -1, 1], and thus, the signal-to-noise ratio can be expressed as a Chebyshev polynomial:
in the formula,the signal-to-noise ratio of the fitted direct signal is shown,nfor the order of the chebyshev polynomial,C i chebyshev polynomial coefficients, chebyshev polynomials, being signal-to-noise ratio componentsT i The recurrence formula is:
after fitting, the direct signal can be eliminated to obtain a signal only retaining the signal-to-noise ratio of the reflected signal and the high-frequency random noise, and the direct signal-to-noise ratio elimination formula is expressed as follows:
in the formula:S r_n for signals that contain only the signal-to-noise ratio and noise of the reflected signal,S o representing the signal-to-noise ratio that was originally just received,S d the signal-to-noise ratio of the direct signal of the corresponding point is obtained after Chebyshev polynomial fitting.
4. The Beidou/GPS soil moisture measurement method based on sliding algorithm and weighting strategy according to claim 3, characterized in that the Gaussian low-pass filtering noise reduction algorithm processing in the step S1-1-2 is expressed as follows:
assuming that the signal to be processed containing the signal-to-noise ratio and noise of the reflected signal is represented asS r_n (x) And the signal after filtering which only contains the signal-to-noise ratio of the reflected signal is expressed asS r The calculation procedure is as follows:
in the formula,G(x) Is a gaussian function, the symbol denotes the convolution operation,xrepresenting the SNR sequence, the high-frequency random noise can be effectively filtered through the processing of the process, the signal only retaining the SNR of the reflected signal is obtained,
5. The Beidou/GPS soil moisture measurement method based on sliding algorithm and weighting strategy as claimed in claim 4, wherein in the step S1-1-3, the satellite effective reflection high obtaining method comprises:
(1) determination of signal-to-noise ratio signal by adaptive sliding algorithmS r The data set is selected from a set of data,
the selection method comprises the following steps: the starting angle of the altitude angle is still 5 degrees, but the ending angle is selected in a sliding way between 25 and 35 degrees, a fixed value is not adopted, but the data groups are selected in a sliding way in a self-adaptive way to be respectively calculated in the middle, the selection rule is determined by the data sampling rate,
if the data sampling rate is 30s, the signal-to-noise ratio between 25 degrees and 35 degrees needs to be calculated in groups, when 20 epochs exist between 25 degrees and 35 altitudes, the signals are divided into 20 groups, the first group of data consists of epochs between 5 degrees and 25 degrees of altitude plus the first epoch between 25 degrees and 35 degrees, the second group of data consists of epochs between 5 degrees and 25 degrees of altitude plus the first two epochs between 25 degrees and 35 degrees, and the like, so that 20 groups of data to be processed are formed;
if the data sampling rate is 15s, the step length is determined to be 2, namely, a group of data is selected for calculation every other epoch;
if the data sampling rate is 5s, the step length is 3, namely, a group of data is selected for calculation every three epochs,
if the time is 1s, the step length is 5, a group of data is selected for calculation every five epochs,
after the adaptive sliding algorithm processing, the signal-to-noise ratio signal to be processed is divided intomGroups, which can be represented as:that is to say havemThe signals with different signal-to-noise ratio lengths can be used for calculating the effective reflection height of the satellite;
(2) the effective reflection height of the satellite is obtained by adopting a Lomb-Scargle spectrum analysis method for each group of signal-to-noise ratio signals,
wherein, Lomb-Sacragle spectrum analysis method is used for solving the power spectrum of the signal, and the expression is as follows:
in the formula,P x ( f ) Is at a frequency off The power of the periodic signal of (a);S r (x) In order to contain the signal-to-noise ratio of the reflected signal,Nin order to be the length of the sequence,τthe time shift invariant can be calculated by the following formula:
can find outP x ( f ) At the time of peakfValue of,Then the corresponding satellite effective reflection height is obtained through the relation between the frequency and the effective reflection height ,Has a conversion relation of:
In the formula,λis the wavelength of the signal, and is a known quantity, therefore, the effective reflection height of the satellite can be obtained by the above formulah e ,
Repeating the above calculation processmRespectively resolving the group signal-to-noise ratio signals to obtainmThe effective reflection is high, expressed as:averaging the calculated satellite effective reflection heights to obtain the satellite effective reflection height of the dayWhereiniRepresents the day of the day;
(3) repeating the steps (1) and (2), resolving data of ten days before the measurement day, and obtaining the satellite effective reflection height of nearly ten days, which is expressed asWeighting the satellite effective reflection height of the ten days by adopting a weighting constraint strategy, wherein a specific weighting model is expressed as follows:
in the formula,w i which represents the weighting coefficient(s) of the,idata representing the number of days, ten days in total, measured on the day of the dayi=10, and so on, the weights decrease gradually.
6. The Beidou/GPS soil moisture measurement method based on sliding algorithm and weighting strategy as claimed in claim 5, wherein in the step S1-2, the method for obtaining the SNR phase reference value comprises:
high effective reflection through satelliteCombining signal-to-noise ratio signals containing the reflected signalThe phase value of the signal-to-noise ratio can be obtained by the relationship between the two, and the specific formula is as follows:
in the formula, A andthe amplitude and phase representing the signal-to-noise ratio are unknown quantities in the formula, and E is the altitude angle at the epoch time and is a known quantity, which is obtained by preprocessing in step S1Signal-to-noise ratio signal obtained in step S1-1-2 and including the reflected signalAlso known quantity, from step S1-1-3, and therefore, can be solved for A and A by the extended Kalman Filter AlgorithmThe amount of the unknown quantity is,
the solution is performed by using an extended kalman filter algorithm, as follows:
first, for a nonlinear system, the state equation and observation equation of its state space can be expressed as:
wherein the first term is a state equation and the second term is an observation equation,x k Andz k respectively, the actual state vector and the observation vector, wherein,x k is formed byAAndthe unknown state quantity of the unknown quantity composition can be expressed asx k (A, );z k Is formed byThe observed state quantity, i.e. the output quantity of the system, of known quantity composition is expressed asz k ();u k As a system function, as input quantities of a state equation, from a known quantity andEthe components of the composition are as follows,kwhich is indicative of a sequence of signals,
f is related tok-1 andka non-linear function of the state of the system at the moment,his a state quantityx k And observed quantityz k The non-linear function of interest is,w k andv k respectively process noise and observation noise, and is Gaussian white noise with mutually independent mean value of zero,
in practical applications, it cannot be determinedw k Andv k specific values at each moment, but negligible noise values comekThe state vector at that time and the estimate of the observed quantity, and therefore the above equation can be approximated by:
the extended kalman filter algorithm transforms the nonlinear system into a linear system by performing taylor first order expansions of the nonlinear system near the optimal estimation point of its state, which can be expressed as:
7. The Beidou/GPS soil moisture measurement method based on the sliding algorithm and the weighting strategy as claimed in any one of claims 1 to 6, wherein the specific method of the step S2 is as follows:
s2-1, reading the signal-to-noise ratio phase reference value of the Beidou system from the signal-to-noise ratio phase reference value database, and performing combined calculation by combining the signal-to-noise ratio phase value of the day to obtain the soil humidity of the measurement day, wherein the calculation formula is as follows:
in the formula,the method comprises the steps that the difference between a signal-to-noise ratio phase value and a signal-to-noise ratio phase reference value on the day of measurement is obtained by preprocessing an original observation value obtained through a Beidou system on the day of measurement in step S1, wherein VWC represents the soil humidity required by the day of measurement;
s2-2, respectively calculating corresponding soil humidity by using data of three frequencies of all Beidou observation satellites, and then calculating an average value to obtain the final soil humidity measured by the Beidou system, wherein the final soil humidity is expressed as VWC B Three frequencies of the Beidou satelliteRespectively as follows: b1 is 1575.42MHz, B2 is 1176.45MHz, and B3 is 1268.52MHz, which belongs to the common knowledge in the field.
8. The Beidou/GPS soil moisture measurement method based on sliding algorithm and weighting strategy as claimed in claim 7, wherein the specific method of the step S3 is as follows:
s3-1, reading the signal-to-noise ratio phase reference value of the GPS system from the signal-to-noise ratio phase reference value database, and performing combined calculation by combining the signal-to-noise ratio phase value of the day to obtain the soil humidity of the measurement day, wherein the calculation formula is as follows:
in the formula,the difference between the signal-to-noise ratio phase value and the signal-to-noise ratio phase reference value on the measurement day is obtained by preprocessing the original observation value acquired by the GPS on the measurement day through the step S1, wherein the VWC represents the soil humidity required by the measurement day;
s3-2, respectively calculating corresponding soil humidity by using dual-frequency data of all observation satellites of the GPS, and then calculating an average value to obtain the final soil humidity measured by the GPS, wherein the final soil humidity is expressed as VWC G The dual-frequency signal frequencies of the GPS satellite are respectively as follows: l1 is 1575.42MHz, L2 is 1227.60MHz, and is common knowledge in the field.
9. The Beidou/GPS soil moisture measurement method based on sliding algorithm and weighting strategy as claimed in claim 8, wherein the implementation method of the step S4 is:
and (4) solving the soil humidity of the final measurement day by using the Beidou and GPS soil humidity obtained by calculation in the step S3 and the step S4, namely, taking the average value of the Beidou and GPS soil humidity, and expressing as follows:
in the formula, VWC B Is the soil humidity, VWC, measured by the Beidou system G Is the soil humidity, VWC, measured by the GPS system Final (a Chinese character of 'gan') I.e. the final measured soil moisture.
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210770775.1A CN114839354B (en) | 2022-07-02 | 2022-07-02 | Beidou and GPS soil humidity measurement method based on sliding algorithm and weighting strategy |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210770775.1A CN114839354B (en) | 2022-07-02 | 2022-07-02 | Beidou and GPS soil humidity measurement method based on sliding algorithm and weighting strategy |
Publications (2)
Publication Number | Publication Date |
---|---|
CN114839354A true CN114839354A (en) | 2022-08-02 |
CN114839354B CN114839354B (en) | 2022-11-18 |
Family
ID=82573543
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202210770775.1A Active CN114839354B (en) | 2022-07-02 | 2022-07-02 | Beidou and GPS soil humidity measurement method based on sliding algorithm and weighting strategy |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN114839354B (en) |
Citations (21)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101900692A (en) * | 2010-06-18 | 2010-12-01 | 武汉大学 | Method for measuring large-area soil humidity |
CN102968552A (en) * | 2012-10-26 | 2013-03-13 | 郑州威科姆科技股份有限公司 | Satellite orbit data estimation and correction method |
CN104020180A (en) * | 2014-06-19 | 2014-09-03 | 武汉大学 | Soil humidity inversion method based on low elevation signal received by Beidou base station |
CN104224181A (en) * | 2014-09-26 | 2014-12-24 | 中国科学院生物物理研究所 | SAR (Specific Absorption Rate) real-time monitoring system and method of multi-channel magnetic resonance imaging equipment |
CN105352979A (en) * | 2015-12-08 | 2016-02-24 | 武汉大学 | Soil humidity estimation method based on Beidou GEO satellite signals |
CN105787474A (en) * | 2016-03-29 | 2016-07-20 | 武汉瑞能电力设备有限责任公司 | Processing method of bridge vibration monitoring data |
CN106125106A (en) * | 2016-08-10 | 2016-11-16 | 清华大学 | The method measuring soil moisture based on the ground Big Dipper/GPS dual-mode survey station |
CN106290408A (en) * | 2016-07-21 | 2017-01-04 | 清华大学 | Based on the soil moisture measurement method running GNSS station signal-to-noise ratio data continuously |
CN107796484A (en) * | 2017-01-11 | 2018-03-13 | 中南大学 | One kind is based on BDStar navigation system signal-to-noise ratio data observed stage changing method |
CN109932735A (en) * | 2019-03-25 | 2019-06-25 | 中国铁路设计集团有限公司 | The localization method of the short baseline single-frequency simple epoch solution of Beidou |
CN111399015A (en) * | 2020-04-07 | 2020-07-10 | 中船重工鹏力(南京)大气海洋信息系统有限公司 | Dual-mode satellite fusion positioning method suitable for ship traffic management system |
CN111474209A (en) * | 2020-03-13 | 2020-07-31 | 山东航向电子科技有限公司 | Real-time soil humidity measuring method based on intelligent terminal equipment |
CN212692782U (en) * | 2020-09-10 | 2021-03-12 | 安徽阡陌网络科技有限公司 | High-precision lightweight handheld land area measuring device |
CN112505068A (en) * | 2020-11-03 | 2021-03-16 | 桂林理工大学 | Surface soil humidity multi-satellite combined inversion method based on GNSS-IR |
CN112782689A (en) * | 2020-12-29 | 2021-05-11 | 西南交通大学 | Multi-satellite data fusion GNSS-IR soil humidity monitoring method |
CN113049777A (en) * | 2021-03-12 | 2021-06-29 | 北京航空航天大学 | Device for measuring soil humidity through GNSS direct reflection signal carrier interference |
WO2022007211A1 (en) * | 2020-07-10 | 2022-01-13 | 自然资源部第一海洋研究所 | Gnss-based real-time high-precision wave measurement method and apparatus |
CN113959329A (en) * | 2021-12-17 | 2022-01-21 | 西南交通大学 | Snow depth inversion method based on multi-satellite data fusion |
CN215953863U (en) * | 2020-12-25 | 2022-03-04 | 武汉大学 | Beidou/GNSS water vapor real-time inversion device |
CN114355411A (en) * | 2021-12-22 | 2022-04-15 | 杭州电子科技大学 | Flood detection method based on Beidou or GPS carrier-to-noise ratio observation value |
CN114397425A (en) * | 2021-12-22 | 2022-04-26 | 杭州电子科技大学 | GNSS-IR soil humidity inversion method based on generalized continuation approximation |
-
2022
- 2022-07-02 CN CN202210770775.1A patent/CN114839354B/en active Active
Patent Citations (21)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101900692A (en) * | 2010-06-18 | 2010-12-01 | 武汉大学 | Method for measuring large-area soil humidity |
CN102968552A (en) * | 2012-10-26 | 2013-03-13 | 郑州威科姆科技股份有限公司 | Satellite orbit data estimation and correction method |
CN104020180A (en) * | 2014-06-19 | 2014-09-03 | 武汉大学 | Soil humidity inversion method based on low elevation signal received by Beidou base station |
CN104224181A (en) * | 2014-09-26 | 2014-12-24 | 中国科学院生物物理研究所 | SAR (Specific Absorption Rate) real-time monitoring system and method of multi-channel magnetic resonance imaging equipment |
CN105352979A (en) * | 2015-12-08 | 2016-02-24 | 武汉大学 | Soil humidity estimation method based on Beidou GEO satellite signals |
CN105787474A (en) * | 2016-03-29 | 2016-07-20 | 武汉瑞能电力设备有限责任公司 | Processing method of bridge vibration monitoring data |
CN106290408A (en) * | 2016-07-21 | 2017-01-04 | 清华大学 | Based on the soil moisture measurement method running GNSS station signal-to-noise ratio data continuously |
CN106125106A (en) * | 2016-08-10 | 2016-11-16 | 清华大学 | The method measuring soil moisture based on the ground Big Dipper/GPS dual-mode survey station |
CN107796484A (en) * | 2017-01-11 | 2018-03-13 | 中南大学 | One kind is based on BDStar navigation system signal-to-noise ratio data observed stage changing method |
CN109932735A (en) * | 2019-03-25 | 2019-06-25 | 中国铁路设计集团有限公司 | The localization method of the short baseline single-frequency simple epoch solution of Beidou |
CN111474209A (en) * | 2020-03-13 | 2020-07-31 | 山东航向电子科技有限公司 | Real-time soil humidity measuring method based on intelligent terminal equipment |
CN111399015A (en) * | 2020-04-07 | 2020-07-10 | 中船重工鹏力(南京)大气海洋信息系统有限公司 | Dual-mode satellite fusion positioning method suitable for ship traffic management system |
WO2022007211A1 (en) * | 2020-07-10 | 2022-01-13 | 自然资源部第一海洋研究所 | Gnss-based real-time high-precision wave measurement method and apparatus |
CN212692782U (en) * | 2020-09-10 | 2021-03-12 | 安徽阡陌网络科技有限公司 | High-precision lightweight handheld land area measuring device |
CN112505068A (en) * | 2020-11-03 | 2021-03-16 | 桂林理工大学 | Surface soil humidity multi-satellite combined inversion method based on GNSS-IR |
CN215953863U (en) * | 2020-12-25 | 2022-03-04 | 武汉大学 | Beidou/GNSS water vapor real-time inversion device |
CN112782689A (en) * | 2020-12-29 | 2021-05-11 | 西南交通大学 | Multi-satellite data fusion GNSS-IR soil humidity monitoring method |
CN113049777A (en) * | 2021-03-12 | 2021-06-29 | 北京航空航天大学 | Device for measuring soil humidity through GNSS direct reflection signal carrier interference |
CN113959329A (en) * | 2021-12-17 | 2022-01-21 | 西南交通大学 | Snow depth inversion method based on multi-satellite data fusion |
CN114355411A (en) * | 2021-12-22 | 2022-04-15 | 杭州电子科技大学 | Flood detection method based on Beidou or GPS carrier-to-noise ratio observation value |
CN114397425A (en) * | 2021-12-22 | 2022-04-26 | 杭州电子科技大学 | GNSS-IR soil humidity inversion method based on generalized continuation approximation |
Non-Patent Citations (3)
Title |
---|
MINGKUN SU等: "BeiDou system satellite-induced pseudorange multipath bias mitigation based on different orbital characteristic for static applications", 《IET RADAR SONAR NAVIG.》 * |
SIBO ZHANG等: "Use of GNSS SNR data to retrieve soil moisture and vegetation variables over a wheat crop", 《HYDROL. EARTH SYST. SCI. DISCUSS.》 * |
张双成等: "BDS/GPS多卫星解译土壤湿度变化研究", 《测绘科学》 * |
Also Published As
Publication number | Publication date |
---|---|
CN114839354B (en) | 2022-11-18 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN114518586B (en) | GNSS precise single-point positioning method based on spherical harmonic expansion | |
Vousdoukas et al. | Performance of intertidal topography video monitoring of a meso-tidal reflective beach in South Portugal | |
CN109001845B (en) | rainfall forecasting method | |
CN110530901B (en) | Medium and small scale soil water monitoring system and method integrating cosmic ray neutron method and unmanned aerial vehicle remote sensing | |
CN111751286B (en) | Soil moisture extraction method based on change detection algorithm | |
CN116519913B (en) | GNSS-R data soil moisture monitoring method based on fusion of satellite-borne and foundation platform | |
CN115980751A (en) | Power law model InSAR troposphere delay correction method | |
CN115047499A (en) | Inversion method and system for satellite-borne GNSS-R soil temperature and humidity | |
CN114924241A (en) | Frequency correction method and system for satellite-borne rainfall measurement radar and ground-based weather radar | |
CN112782701A (en) | Visibility perception method, system and equipment based on radar | |
Yu et al. | Generating daily 100 m resolution land surface temperature estimates continentally using an unbiased spatiotemporal fusion approach | |
CN107273797B (en) | Rice sub-pixel identification method based on water body index variation coefficient | |
Levin et al. | Observation impacts on the Mid-Atlantic Bight front and cross-shelf transport in 4D-Var ocean state estimates: Part I—Multiplatform analysis | |
Wang et al. | Soil moisture retrieval from sentinel-1 and sentinel-2 data using ensemble learning over vegetated fields | |
CN114821360A (en) | Method and device for intelligently extracting crop leaf area index and electronic equipment | |
CN114611699A (en) | Soil moisture downscaling method and device, electronic equipment and storage medium | |
CN115980317B (en) | Foundation GNSS-R data soil moisture estimation method based on corrected phase | |
CN114839354B (en) | Beidou and GPS soil humidity measurement method based on sliding algorithm and weighting strategy | |
CN112733445A (en) | Large-area-scale soil moisture inversion method based on evapotranspiration vegetation index spatial characteristics | |
CN116029162B (en) | Flood disaster inundation range monitoring method and system by using satellite-borne GNSS-R data | |
CN114252875B (en) | High-precision meshing method for imaging altitude data | |
CN113625307A (en) | Landslide monitoring system and method based on GNSS | |
CN111044489B (en) | Method for obtaining atmosphere refractive index height distribution profile based on multi-wavelength measurement | |
CN110543835A (en) | Satellite sea surface salinity remote sensing product precision evaluation method based on triple matching theory | |
Farah | Assessment of UNB3M neutral atmosphere model and EGNOS model for near-equatorial-tropospheric delay correction |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |