Real-time inversion method and system of three-dimensional ocean wave spectrum based on single-frequency GNSS buoy

Abstract

The present application belongs to the field of three-dimensional ocean wave spectrum calculation, and specifically relates to a real-time inversion method and system for three-dimensional ocean wave spectrum based on a single-frequency GNSS buoy. The method adopts the time difference of carrier phase (TDCP) technology, and only uses the GNSS data of the L1 frequency band to realize the real-time calculation of the three-dimensional ocean wave spectrum. By using real-time broadcast ephemeris instead of post-precision ephemeris, the present invention greatly shortens the calculation delay and can realize real-time calculation without waiting. In addition, since the TDCP technology has a low demand for computing resources, the present invention can integrate the calculation process into the micro-control system inside the GNSS buoy, thereby realizing real-time observation of the three-dimensional spectrum of ocean waves. This innovative technical solution not only reduces the system cost, but also improves the real-time and efficiency of ocean observation.

Description

Technical Field

[0001] The present application belongs to the field of three-dimensional ocean wave spectrum calculation, and specifically relates to a real-time inversion method and system for three-dimensional ocean wave spectrum based on a single-frequency GNSS buoy. Background Art

[0002] At present, using GNSS buoys to observe sea surface height changes and then invert information such as the directional spectrum of waves has become an important means of ocean observation. However, existing technologies mainly rely on the precise point positioning (PPP) solution method of post-precision ephemeris to obtain high-precision buoy displacement information. This method usually requires waiting for 1 to 2 days to obtain precise ephemeris, resulting in delays in data processing. In addition, the PPP solution process has high requirements for computing resources and can usually only be completed on high-performance computers, making it difficult to achieve real-time calculations inside the buoy.

[0003] The PPP method usually requires ionospheric elimination calculations, which requires the use of geodetic multi-frequency GNSS systems. The hardware cost of such systems is high, which not only increases the cost of observation equipment, but also limits the widespread application of this technology. Therefore, the existing GNSS buoy technology has obvious deficiencies in cost and real-time performance, and cannot meet the needs of real-time calculation of wave directional spectra.

[0004] In addition, some current technical solutions can also obtain three-dimensional wave spectra using inertial measurement units (IMUs). However, the angle measurement of IMUs often drifts significantly during long-term use, which can lead to a decrease in measurement accuracy. Therefore, in order to ensure the accuracy of the measurement, IMU systems usually need to be calibrated regularly, which not only increases the difficulty of maintenance, but also may cause downtime during the observation process, further affecting the real-time and reliability of the system. Therefore, the IMU combination system also has certain limitations in practical applications. Summary of the invention

[0005] In order to overcome these technical deficiencies, the present invention proposes a real-time inversion method for three-dimensional wave spectrum based on a low-cost single-frequency GNSS module. The method adopts the time difference of carrier phase (TDCP) technology, and can realize the real-time calculation of the three-dimensional wave spectrum using only the GNSS data of the L1 frequency band. By using real-time broadcast ephemeris instead of post-precision ephemeris, the present invention greatly shortens the calculation delay and can realize real-time calculation without waiting. In addition, since the TDCP technology has a low demand for computing resources, the present invention can integrate the calculation process into the micro-control system inside the GNSS buoy, thereby realizing real-time observation of the three-dimensional spectrum of waves. This innovative technical solution not only reduces the system cost, but also improves the real-time and efficiency of ocean observation. Its technical solution is:

[0006] A real-time inversion method for three-dimensional ocean wave spectrum based on a single-frequency GNSS buoy comprises the following steps:

[0007] S1. Establish a carrier phase measurement model;

[0008] S2. Establish a difference-based observation equation to solve the receiver position change and clock bias, and construct the error equation of carrier phase time difference TDCP;

[0009] S3.TDCP cycle slip detection and repair;

[0010] S4. Calculate the three-dimensional wave spectrum using the three-dimensional velocity information.

[0011] Preferably, in step S1, the carrier phase measurement model is as follows:

[0012] In step S1, the carrier phase measurement model is as follows:

[0013] ;

[0014] yes The wavelength of the frequency signal; Indicates the band code of the satellite signal. This application uses a single frequency band , the frequency is 1575.42 MHz;

[0015] It is a receiver. It is a receiver. With satellite Between Carrier phase observations of frequency;

[0016] is the geometric distance between the receiver and the satellite;

[0017] It is the speed of light;

[0018] and are the clock errors of the receiver and satellite respectively;

[0019] is the integer ambiguity of the carrier phase;

[0020] and are the delay errors in the ionosphere and troposphere, respectively;

[0021] Random measurement noise including multipath and other residual errors.

[0022] Preferably, in step S2, the carrier phase difference calculation formula is:

[0023] ;

[0024] The symbol represents the amount of change between consecutive observed data points.

[0025] Preferably, the error equation of TDCP is:

[0026] ;

[0027] Represents the residual between the observed and calculated values;

[0028] is the observation vector, which contains the differential carrier phase observation value and the differential values ​​of various errors;

[0029] ;

[0030] The unknowns include the change in the receiver's position between successive epochs. and receiver clock bias changes , use the least squares method to solve for the unknowns.

[0031] Preferably, cycle slips must be detected and repaired before velocity estimation, the steps are as follows:

[0032] S31. Detection: Detect cycle slips using a variety of algorithms, including the Melbourne-Wübbena combination, the TurboEdit algorithm, or consistency checks based on carrier phase and pseudorange observations;

[0033] S32. Repair: Once a cycle slip is detected, the phase data needs to be repaired or an estimation method needs to be used to compensate for the jump in this data;

[0034] S33. Verification: After repair, it is necessary to verify whether the cycle slip has been successfully corrected, which can be achieved by analyzing the residual error of the repaired data;

[0035] S34. Data integration: After confirming the data quality, the repaired carrier phase data is used for subsequent velocity estimation.

[0036] Preferably, the buoy velocity is calculated using the least square method, and the steps are as follows:

[0037] a. Design matrix construction: Based on the TDCP observation equation, construct the design matrix , each row in the design matrix represents the coefficients corresponding to a differential observation equation, which are related to the change of the receiver position and clock difference changes related;

[0038] b. Normal equation establishment: Use the design matrix and differential carrier phase observations to construct the normal equation:

[0039] ;

[0040] in Indicates A The transposed matrix of is the weight matrix, which represents the accuracy of each observation; is the observed quantity vector; use matrix inversion to solve for the unknown quantity , i.e., the displacement between epochs and the receiver clock change;

[0041] c. Residual calculation and evaluation: Calculate the residual To assess how well the observations fit the model:

[0042] ;

[0043] The size and distribution of residuals can be used to test data quality and model suitability;

[0044] d. Receiver velocity solution: using the position change and the corresponding time interval , solve for the receiver's velocity ,Right now .

[0045] Preferably, the three-dimensional wave spectrum is calculated using the three-dimensional velocity information, and the steps are as follows:

[0046] S41. Set the data sampling to 2 Hz, which corresponds to the sampling frequency of GNSS data;

[0047] S42. Set the sea depth of the GNSS buoy anchor position ;

[0048] S43. Set the output result of the single-frequency GNSS buoy TDCP algorithm to construct the three-dimensional data moment (east, north, vertical);

[0049] S44. To estimate the wave directional spectrum, set the calculation frequency range to 0.005 to 0.7 Hz, the step size to 0.001 Hz, and set the discretization angle range of the wave propagation direction to 0 to 360 degrees, with an angle resolution of 1 degree;

[0050] S45. Calculate the directional spectrum of the wave using the discrete Fourier transform method (DFTM). In this step, the signal is transformed into the frequency domain using a 256-point fast Fourier transform and iterated 20 times to optimize the estimated result of the directional spectrum; in addition, the estimated directional spectrum is smoothed using a smoothing algorithm to reduce the impact of high-frequency noise on the result;

[0051] When implementing DFTM calculation, firstly, the three-dimensional displacement signals (vertical displacement, eastward displacement and northward displacement) collected by the GNSS buoy are detrended. Assume that the displacement data collected by the buoy is , detrend it to ensure data concentration:

[0052] ;

[0053] in, For time, is the mean value of the displacement signal.

[0054] When implementing DFTM calculation, it is necessary to calculate the cross power spectrum. and The displacement signal is subjected to fast Fourier transform (FFT) to obtain the cross power spectrum in the frequency domain. , for each frequency , cross power spectrum Calculated by the following formula:

[0055] ;

[0056] in, is the direction parameter The Fourier transform of is the direction parameter The complex conjugate of the Fourier transform of ;

[0057] S46. Calculate the wave number and direction discretization. For each frequency , calculate the corresponding wave number :

[0058] ;

[0059] in, is the acceleration due to gravity, It is the depth of water;

[0060] By wave number and the relative geometric layout between sensors, calculating different angles The direction discretization uses the following angle range:

[0061] ;

[0062] in, is the direction resolution, usually 32;

[0063] S47. Directional response function calculation:

[0064] For each direction parameter ,Depending on the type (such as east-west displacement, north-south displacement, or vertical displacement), the transfer function is used to calculate its response to waves in different directions and frequencies. The transfer function is , the function is calculated by the following formula:

[0065] For the east-west displacement sensor, its transfer function is:

[0066] ;

[0067] For the north-south displacement sensor, its transfer function is:

[0068] ;

[0069] For the vertical displacement sensor, the transfer function is assumed to be constant:

[0070] ;

[0071] in, , represents the imaginary unit, depth response Calculated by the following formula:

[0072] ;

[0073] To improve the calculation stability, if , then set its value to 0.1; is the vertical position of the sensor;

[0074] S48. Directional spectrum estimation: The discrete Fourier transform method (DFTM) is used to calculate the directional spectrum of the waves. By combining the cross power spectrum with the transfer function matrix, the directional spectrum of different frequencies can be obtained. and angle Wave energy distribution under :

[0075] ;

[0076] in: Indicates the number of directions, the value is 3;

[0077] is the direction parameter and The cross power spectrum of;

[0078] is the direction parameter Transfer function, is the direction parameter The complex conjugate of the transfer function.

[0079] The final directional spectrum matrix is ​​obtained by summing different directions and frequencies through the above formula , this matrix represents the distribution of wave energy in various directions at different frequencies;

[0080] S49. Post-processing and smoothing of results:

[0081] In order to reduce the impact of noise on the results, the directional spectrum results are smoothed; the smoothing operation uses a two-dimensional convolution filter and is processed through the following kernel matrix:

[0082] ;

[0083] The resulting matrix after smoothing , indicating that at each frequency and angle Wave energy distribution under the condition of wave motion. After smoothing, the result is more stable.

[0084] At this point, the calculation of the wave directional spectrum has been completed, and the three-dimensional spectrum of waves observed using a single-frequency GNSS buoy has been realized.

[0085] A three-dimensional ocean wave spectrum real-time inversion system based on a single-frequency GNSS buoy, used for a three-dimensional ocean wave spectrum real-time inversion method based on a single-frequency GNSS buoy, comprising a single-frequency GNSS receiver, an information processing module, a power supply module and a communication module, wherein the single-frequency GNSS receiver is responsible for collecting GNSS signals, the information processing module is responsible for system control and analytical calculation of GNSS signals, the power supply module is responsible for system power supply, and the communication module is responsible for sending the results to a server;

[0086] The single-frequency GNSS receiver communicates with the microcontroller through the communication module;

[0087] The information processing module includes a single-chip microcomputer and an ARM A7 processor. The single-chip microcomputer is responsible for basic equipment management tasks, including data collection of GNSS receivers, power management of equipment, and control of the ARM A7 processor. The ARM A7 processor also uses the constructed carrier phase measurement model to calculate the three-dimensional wave spectrum.

[0088] Power supply module: Solar panels and lithium batteries are used to power the entire system.

[0089] Preferably, when the ARM A7 processor is in a sleep state, it is activated by a wake-up signal sent by the single-chip microcomputer only after the single-chip microcomputer completes 20 minutes of GNSS data collection; at this time, the ARM A7 starts and undertakes advanced processing of the data, including executing standardized processing of GNSS data and carrier phase time difference TDCP algorithm to calculate and generate three-dimensional spectral data of waves; after the data processing is completed, the ARM A7 transmits the results to the single-chip microcomputer, and automatically returns to the sleep state according to the single-chip microcomputer instructions to reduce energy consumption.

[0090] Compared with the prior art, the present invention has the following beneficial effects:

[0091] Compared with the prior art, this patent has the following significant beneficial effects:

[0092] 1. Low cost: The use of single-frequency GNSS receivers and simplified system design significantly reduces the manufacturing and maintenance costs of the equipment.

[0093] 2. No need for regular direction calibration: Since the measurement principle of the GNSS system does not involve the integration of internal states and does not drift over time like IMU devices, there is no need for regular direction calibration.

[0094] 3. Real-time observation: Using the TDCP algorithm of a single-frequency GNSS buoy, the present invention can achieve real-time monitoring of the marine environment. BRIEF DESCRIPTION OF THE DRAWINGS

[0095] Figure 1 This is the flow chart of this application.

[0096] Figure 2 A schematic diagram of the system.

[0097] Figure 3 is the three-dimensional ocean wave spectrum calculated by single-frequency GNSS.

[0098] Figure 4 To show the comparative verification results of Wave Knight. DETAILED DESCRIPTION

[0099] In order to facilitate understanding of the present application, the present application will be described more fully below with reference to the relevant drawings. Embodiments of the present application are provided in the drawings. However, the present application can be implemented in many different forms and is not limited to the embodiments described herein. On the contrary, the purpose of providing these embodiments is to make the disclosure of the present application more thorough and comprehensive.

[0100] The present invention proposes a real-time inversion method for three-dimensional ocean wave spectrum based on single-frequency GNSS buoy, and its technical solution covers the design of hardware system and the inversion calculation process of three-dimensional ocean wave spectrum. The design goal of the whole system is to achieve low-cost, high-precision, real-time ocean wave observation to solve the problems of calculation delay, high system cost and dependence on precise ephemeris in the prior art.

[0101] Figure 2 As shown, a system for real-time inversion of three-dimensional wave spectrum based on single-frequency GNSS buoy includes a single-frequency GNSS receiver, an information processing module, a power supply module and a communication module;

[0102] The single-frequency GNSS receiver communicates with the microcontroller through the communication module;

[0103] The information processing module includes a single-chip microcomputer and an ARM A7 processor. The single-chip microcomputer is responsible for basic equipment management tasks, including data collection of GNSS receivers, power management of equipment, and control of the ARM A7 processor. The ARM A7 processor also uses the constructed carrier phase measurement model to calculate the three-dimensional wave spectrum.

[0104] Power supply module: Solar panels and high-performance lithium batteries are used to power the entire system.

[0105] 1. The hardware system design is as follows:

[0106] 1.1 Hardware core system:

[0107] The hardware system of the present invention focuses on ensuring the accuracy of three-dimensional wave spectrum measurement while significantly reducing costs. The core hardware components of the system include a single-frequency GNSS receiver, an economical GNSS antenna, a microcontroller, a power supply system, and a communication module. These components meet the dual requirements of cost-effectiveness and operational performance.

[0108] The single-frequency GNSS receiver is a key component of the system, which provides a cost-effective solution. Compared with the traditional multi-frequency geodetic GNSS receiver, the single-frequency receiver has a greater cost advantage, and through the advanced signal processing and data analysis algorithms of the present invention, it can still maintain good measurement accuracy in complex marine environments. The use of this receiver greatly reduces the cost of the entire system, while meeting the basic accuracy requirements of three-dimensional wave spectrum measurement.

[0109] Economical GNSS antennas are used to further reduce the cost of the system. Although these economical antennas may have certain disadvantages in signal quality and anti-interference ability compared with professional measurement choke ring antennas, the time differential carrier phase (TDCP) technology used in this system focuses on measuring speed changes rather than absolute position, making the signal quality requirements relatively low.

[0110] The key advantage of TDCP technology is that it can provide extremely high velocity measurement accuracy, reaching the millimeter per second (mm / s) level, which is accurate enough to meet the needs of real-time inversion of three-dimensional wave spectra. This algorithm estimates the velocity of the buoy by analyzing the slight changes in the received GNSS signal, thereby indirectly reflecting the dynamic characteristics of the waves. Even economical antennas can accurately capture key velocity changes with the support of this technology, ensuring data accuracy and system reliability.

[0111] The microcontroller in the system combines a low-power MCU and a high-performance ARM A7 processor to optimize the efficiency of data processing and system instruction execution. The low-power MCU microcontroller selected has the characteristics of high performance and energy saving. It can run continuously at a power consumption of only 0.1W, which is suitable for long-term marine environment monitoring. The high-performance ARM A7 processor is mainly responsible for the execution of complex algorithms, including the standardized conversion of GNSS data and the calculation of the carrier phase time difference (TDCP) algorithm. Its large memory and powerful computing power ensure the efficient completion of these tasks.

[0112] The design of the power supply system fully considers the special requirements of the marine environment. The system uses solar panels and high-performance lithium batteries, which not only ensures long-term continuous power supply, but also improves the system's self-sufficiency in an unmanned environment. In addition, the power supply system also includes intelligent power management technology, which can automatically adjust energy consumption according to the operating status of the equipment and optimize energy efficiency.

[0113] Finally, the communication module ensures real-time data transmission and remote control capabilities of the system. The communication technology used supports long-distance, low-power data transmission, and ensures communication stability and reliability even in harsh marine environments. The two-way communication function of the communication module allows operators to monitor the system status from a remote location and make necessary adjustments, enhancing the operational flexibility of the system.

[0114] 1.2 Coupling design of MCU and ARM A7:

[0115] The integrated design of the MCU and ARM A7 optimizes the system's power management to ensure extended operating time for the buoy. In this design, the low-power MCU is responsible for basic device management tasks, including data acquisition from the GNSS receiver, power management of the device, and control of the ARM A7 processor. The MCU remains in a continuous operation state to monitor the system status in real time and execute basic instructions.

[0116] Specifically, the ARM A7 processor is in sleep mode to save energy and is activated only by a wake-up signal from the MCU after the MCU completes 20 minutes of GNSS data collection. At this point, the ARM A7 starts and undertakes advanced processing of the data, including performing normalization of GNSS data and the time difference carrier phase (TDCP) algorithm to calculate and generate three-dimensional spectral data of the waves. After the data processing is completed, the ARM A7 transmits the results to the MCU and automatically returns to sleep mode according to the MCU instructions to minimize energy consumption.

[0117] The entire operation process is automated and follows a preset operating cycle, ensuring efficient operation of the system while reducing energy consumption. This coupled design not only enhances the energy efficiency of the system, but also improves the equipment's ability to continue working in an unmanned monitoring environment, making it an ideal marine monitoring solution.

[0118] In the present invention, the signal processing algorithm is the key to achieve accurate three-dimensional wave spectrum measurement. The microcontroller system controls the GNSS receiver to collect carrier phase observation data (L1 frequency) according to the set parameters (frequency is 2Hz). The collected data is then processed by the time difference carrier phase (TDCP) algorithm to solve the velocity information of the buoy.

[0119] Figure 1 As shown, the real-time inversion method of three-dimensional wave spectrum based on single-frequency GNSS buoy includes the following steps:

[0120] S1. Establish a carrier phase measurement model;

[0121] S2. Establish a difference-based observation equation to solve the receiver position change and clock bias, and construct the TDCP error equation;

[0122] S3.TDCP cycle slip detection and repair;

[0123] S4. Calculate the three-dimensional wave spectrum using the three-dimensional velocity information.

[0124] 2.1 Carrier phase measurement model establishment:

[0125] This step is used to calculate the distance between the receiver and the satellite. The model contains carrier phase observations that incorporate geometric distance, clock error, atmospheric delay, and integer ambiguity.

[0126] ;

[0127] The meanings of the parameters are as follows:

[0128] yes The wavelength of the frequency signal; Indicates the band code of the satellite signal. This application uses a single frequency band , the frequency is 1575.42 MHz;

[0129] It is a receiver. It is a receiver. With satellite Between Carrier phase observations of frequency;

[0130] is the geometric distance between the receiver and the satellite;

[0131] It is the speed of light;

[0132] and are the clock errors of the receiver and satellite respectively;

[0133] is the integer ambiguity of the carrier phase;

[0134] and are the delay errors in the ionosphere and troposphere, respectively;

[0135] Random measurement noise including multipath and other residual errors.

[0136] 2.2 Carrier phase difference:

[0137] Carrier phase differential processing is an important step in extracting velocity information from serial carrier phase data. By processing the differences between consecutive observations, it effectively eliminates a variety of error sources that are not directly related to velocity estimation, such as the clock offsets of the receiver and satellite and constant integer ambiguity, and significantly reduces the impact of atmospheric delay.

[0138] The carrier phase difference calculation formula performed in this process is:

[0139] ;

[0140] Here The symbol represents the change between consecutive observation data points. By difference, the distance change and other related differential items can reveal the position change of the receiver and its velocity information. This differential process helps to accurately eliminate the errors caused by the consistency or slow changes of the satellite and receiver time bias, while effectively mitigating the effects of ionospheric and tropospheric delays.

[0141] 2.3 Error equation: A difference-based observation equation is established to solve the receiver position change and clock bias.

[0142] The error equation is established to minimize these errors by adjusting the model parameters, thereby obtaining a more accurate velocity estimate. Using the aforementioned differential satellite distance equation, we insert it into the carrier phase differential observation equation to construct the TDCP error equation:

[0143] ;

[0144] in represents the residual between the observed and calculated values, and is the observation vector, which contains the differential carrier phase observation value and the differential values ​​of various errors:

[0145] ;

[0146] In this model, the unknowns include the change in the receiver position between successive epochs. and receiver clock bias changes When enough satellites (more than four) are observed, these unknowns can be solved using the least squares method.

[0147] 2.4 Cycle slip detection and repair:

[0148] In the TDCP method, a cycle slip refers to the interruption of carrier phase continuity due to receiver lock loss or other factors. A cycle slip can cause a sudden change in the carrier phase observation value, thus affecting the accuracy of speed estimation. Therefore, possible cycle slips must be detected and repaired before speed estimation. The steps are as follows:

[0149] a. Detection: Cycle slips are detected using a variety of algorithms. Common methods include the Melbourne-Wübbena combination, the TurboEdit algorithm, or consistency checks based on carrier phase and pseudorange observations.

[0150] b. Repair: Once a cycle slip is detected, the phase data needs to be repaired. This may involve discarding the affected observations or using estimation methods to compensate for the slips in these observations.

[0151] c. Verification: After repair, it is necessary to verify whether the cycle slip has been corrected successfully, which is usually achieved by analyzing the residual data after repair.

[0152] d. Data integration: After confirming the data quality, the repaired carrier phase data can be used for subsequent velocity estimation.

[0153] 2.5 Least squares solution for buoy velocity:

[0154] Using observation data from more than four satellites, the receiver position change is solved by the least squares method. When using the TDCP (Time-Differenced Carrier Phase) method for velocity estimation, the least squares solution step can be divided into the following sub-steps:

[0155] a. Design matrix construction: Based on the TDCP observation equation, construct the design matrix Each row in the design matrix represents the coefficients of a differential observation equation, which are related to the change of the receiver position. and clock difference changes related.

[0156] b. Normal equation establishment: Use the design matrix and differential carrier phase observations to construct the normal equation ,in is the weight matrix, which represents the accuracy of each observation; is the observation vector, i.e., the difference observation value. It is usually solved using matrix inversion or other numerical methods such as Cholesky decomposition.

[0157] Residual Calculation and Evaluation: Calculating Residuals , to assess the fit between the observed values ​​and the model. The size and distribution of the residuals can be used to test the quality of the data and the suitability of the model.

[0158] d. Receiver velocity solution: using the position change and the corresponding time interval , solve for the receiver's velocity ,Right now .

[0159] 3. Calculate the 3D wave spectrum using 3D velocity information:

[0160] S41. Set the sampling frequency of GNSS data;

[0161] S42. Set the sea depth of the GNSS buoy anchor position;

[0162] S43. Set the output result of the single-frequency GNSS buoy TDCP algorithm to construct three-dimensional data moments in the east, north and vertical directions;

[0163] S44. Setting the calculation frequency range, and at the same time setting the discretization angle range of the wave propagation direction to 0 to 360 degrees;

[0164] S45. Calculate the directional spectrum of waves using the discrete Fourier transform method DFTM;

[0165] The three-dimensional displacement signal collected by the GNSS buoy is detrended. Suppose the displacement data collected by the buoy is , detrending it:

[0166] ;

[0167] in, For time, is the mean value of the displacement signal;

[0168] By the direction parameter and The displacement signal is subjected to fast Fourier transform FFT to obtain the cross power spectrum in the frequency domain ; For each frequency , cross power spectrum Calculated by the following formula:

[0169] ;

[0170] in, is the direction parameter The Fourier transform of is the direction parameter The complex conjugate of the Fourier transform of ;

[0171] S46. Calculate the wave number and direction discretization; for each frequency, calculate the corresponding wave number, the wave number is calculated by the following formula:

[0172] ;

[0173] in, is the acceleration due to gravity, It is the depth of water;

[0174] By wave number As well as the relative geometric layout between sensors, the phase difference of wave propagation in different directions is calculated, and the direction discretization adopts the following angle range :

[0175]

[0176] in, is the directional resolution;

[0177] S47. Direction transfer function calculation: The transfer function is , the function is calculated by the following formula:

[0178] For the east-west displacement sensor, its transfer function is:

[0179] ;

[0180] For the north-south displacement sensor, its transfer function is:

[0181] ;

[0182] For the vertical displacement sensor, the transfer function is assumed to be constant:

[0183] ;

[0184] in, Represents imaginary unit, depth response Calculated by the following formula:

[0185] ;

[0186] like , then set its value to 0.1; is the vertical position of the sensor;

[0187] S48. Directional spectrum estimation:

[0188] The discrete Fourier transform method DFTM is used to calculate the directional spectrum of the waves; by combining the cross power spectrum with the transfer function matrix, the directional spectrum of different frequencies is obtained. and angle The wave energy distribution under , is calculated as:

[0189] Directional Spectrum It can be expressed as:

[0190] ;

[0191] in is the direction parameter and The cross power spectrum is expressed as:

[0192] ;

[0193] is the direction parameter The transfer function, is the direction parameter The complex conjugate of the transfer function.

[0194] In order to reduce the influence of noise on the results, a two-dimensional convolution filter is used to filter the directional spectrum. Smoothing process, wave energy distribution after smoothing , the results are more stable.

[0195] At this point, the calculation of the wave directional spectrum has been completed, and the three-dimensional spectrum of waves observed using a single-frequency GNSS buoy has been realized.

[0196] Figure 3 The wave directional spectrum calculated using this method is shown in Figure 2. Figure 4 The wave directional spectrum calculated using this method is shown in Figure 2. Figure 3 and Figure 4 The three-dimensional wave spectrum structures of the two are generally the same. The consistency shows that the method of the present invention has the same accuracy as the wave rider device in capturing the main wave direction, and the main wave direction is evenly distributed around 280°.

[0197] The above description of the disclosed embodiments enables one skilled in the art to implement or use the present invention. Various modifications to these embodiments will be apparent to one skilled in the art, and the general principles defined herein may be implemented in other embodiments without departing from the spirit or scope of the present invention. Therefore, the present invention will not be limited to the embodiments shown herein, but rather to the widest scope consistent with the principles and novel features disclosed herein.

Claims

1. A real-time inversion method for three-dimensional ocean wave spectrum based on a single-frequency GNSS buoy, characterized in that: The following steps are involved: S1. Establish a carrier phase measurement model; S2. Establish a difference-based observation equation to solve the receiver position change and clock bias, and construct the error equation of carrier phase time difference TDCP; S3. Carrier phase time difference TDCP cycle slip detection and repair; S4. Calculate the three-dimensional wave spectrum using the three-dimensional velocity information; The directional spectrum of waves is calculated using the discrete Fourier transform method DFTM: The three-dimensional displacement signal collected by the GNSS buoy is detrended. Suppose the displacement data collected by the buoy is X(t), and detrended: x(t)=x(t)-mean(x(t)); Where t is time, mean(x(t)) is the mean of the displacement signal; By performing fast Fourier transform FFT on the displacement signals of the directional parameters m and n, the cross power spectrum P in the frequency domain is obtained. mn (f), for each frequency f, the cross power spectrum P mn (f) is calculated by the following formula: Among them, X m (f) is the Fourier transform of the direction parameter m, is the conjugate complex number of the Fourier transform of the direction parameter n; Compute the wave number and direction discretization: For each frequency f, calculate the corresponding wavenumber k(f): Among them, g is the acceleration due to gravity, h is the water depth, π is the circumference of a circle, The wave propagation phase difference at different angles θ is calculated by using the wave number k(f) and the relative geometric layout between the sensors. The direction discretization uses the following angle range: Among them, D d is the directional resolution; Directional transfer function calculation: The transfer function of the directional parameter m is H m (f, θ), the function is calculated in three directions as follows: For the east-west displacement sensor, its transfer function is: H dspx (f,θ)=ii·K z ·cos(θ); For the north-south displacement sensor, its transfer function is: H dspy (f,θ)=ii·K z ·sin(θ); For the vertical displacement sensor, the transfer function is assumed to be constant: H elev (f,θ)=1; in, Represents imaginary unit, depth response K z Calculated by the following formula: If K z <0.1, then set its value to 0.1; z elev is the vertical position of the sensor; Directional spectrum estimation: The discrete Fourier transform method DFTM is used to calculate the directional spectrum of the waves. By combining the cross power spectrum with the transfer function matrix, the wave energy distribution S(f, θ) at different frequencies f and angles θ is obtained. The calculation formula is: Where: M represents the number of directions, P mn (f) is the cross power spectrum of directional parameters m and n; H n (f, θ) is the transfer function of the directional parameter n, is the complex conjugate of the transfer function with directional parameter m; By summing different directions and frequencies through the above formula, the final directional spectrum matrix S(f, θ) is obtained, which represents the distribution of wave energy in various directions at different frequencies; Result post-processing and smoothing: The smoothing operation uses a two-dimensional convolution filter with the following kernel matrix: The smoothed result matrix S′(f, θ)=S(f, θ)*Knu represents the wave energy distribution at each frequency f and angle θ.

2. The method for real-time inversion of three-dimensional ocean wave spectrum based on single-frequency GNSS buoy according to claim 1, characterized in that: In step S1, the carrier phase measurement model is as follows: λ i YesL i The wavelength of the frequency signal; i represents the band code of the satellite signal; The receiver is the L between the receiver r and the satellite s i Carrier phase observations of frequency; is the geometric distance between the receiver and the satellite; C is the speed of light; dt r and dT s are the clock errors of the receiver and satellite respectively; is the integer ambiguity of the carrier phase; and are the delay errors in the ionosphere and troposphere, respectively; ε φ Random measurement noise including multipath and other residual errors.

3. The method for real-time inversion of three-dimensional ocean wave spectrum based on single-frequency GNSS buoy according to claim 2 is characterized in that: In step S2, the carrier phase difference calculation formula is: The Δ symbol represents the amount of change between consecutive observation data points.

4. The method for real-time inversion of three-dimensional ocean wave spectrum based on single-frequency GNSS buoy according to claim 3 is characterized in that: The error equation of carrier phase time difference TDCP is: v represents the residual between the observed value and the calculated value; is the unit direction vector from the receiver position to the satellite position at epoch j; R(j) is the coordinate vector of the satellite at time j; ΔO r The change in the receiver position; r (j-1) is the position of the receiver at time j-1; l is the observation vector, which contains the differential carrier phase observation value and the differential values ​​of various errors; The unknowns include the change in the receiver position between successive epochs, ΔO r and receiver clock deviation change Δdt r , use the least squares method to solve for the unknowns.

5. The method for real-time inversion of three-dimensional ocean wave spectrum based on single-frequency GNSS buoy according to claim 1, characterized in that: Cycle slips must be detected and corrected before velocity estimation. The steps are as follows: S31. Detection: Cycle slip detection using a variety of algorithms, including the Melbourne-Wübbena combination, the TurboEdit algorithm, or consistency checks based on carrier phase and pseudorange observations; S32. Repair: Once a cycle slip is detected, the phase data needs to be repaired or an estimation method needs to be used to compensate for the jump in this data; S33. Verification: After repair, it is necessary to verify whether the cycle slip has been successfully corrected, which can be achieved by analyzing the residuals of the repaired data; S34. Data integration: After confirming the data quality, the repaired carrier phase data is used for subsequent velocity estimation.

6. The method for real-time inversion of three-dimensional ocean wave spectrum based on single-frequency GNSS buoy according to claim 4, characterized in that: Use the least square method to solve the buoy velocity. The steps are as follows: a. Design matrix construction: Based on the TDCP observation equation, construct the design matrix A. Each row in the design matrix represents the coefficient corresponding to a differential observation equation; b. Normal equation establishment: Use the design matrix and differential carrier phase observations to construct the normal equation: A T PAx=A T Pl; Among them, A T represents the transposed matrix of A; P is the weight matrix, which represents the accuracy of each observation; l is the observation vector, and the matrix inversion is used to solve the unknown quantity x, that is, the displacement between epochs and the receiver clock change; c. Residual calculation and evaluation: Calculate the residual v to evaluate the fit between the observation and the model: v = l - Ax; d. Receiver velocity solution: Use the position change ΔO r and the corresponding time interval t j -t j-1 , solve for the receiver's velocity Sp, that is 7. A three-dimensional ocean wave spectrum real-time inversion system based on a single-frequency GNSS buoy, used to execute the three-dimensional ocean wave spectrum real-time inversion method based on a single-frequency GNSS buoy according to any one of claims 1 to 6, characterized in that: It includes a single-frequency GNSS receiver, an information processing module, a power supply module and a communication module; The single-frequency GNSS receiver communicates with the microcontroller through the communication module; The information processing module includes a single-chip microcomputer and an ARM A7 processor. The single-chip microcomputer is responsible for basic equipment management tasks, including data collection of GNSS receivers, power management of equipment, and control of the ARM A7 processor. The ARM A7 processor also uses the constructed carrier phase measurement model to calculate the three-dimensional wave spectrum. Power supply module: Solar panels and lithium batteries are used to power the entire system.

8. The three-dimensional ocean wave spectrum real-time inversion system based on single-frequency GNSS buoy according to claim 7, characterized in that: When the ARM A7 processor is in sleep mode, it is activated by a wake-up signal from the microcontroller only after the microcontroller completes 20 minutes of GNSS data collection. At this time, the ARM A7 starts and undertakes advanced data processing, including executing the standardized processing of GNSS data and the carrier phase time difference TDCP algorithm to calculate and generate three-dimensional spectral data of ocean waves. After the data processing is completed, the ARM A7 transmits the results to the microcontroller and automatically returns to sleep mode according to the microcontroller instructions.

Images