A water surface elevation inversion method based on TDCP difference of Beidou
By constructing a multi-source synchronous observation dataset and handling cycle slip breakpoints and attitude compensation, the observation error problem in water surface elevation inversion on a mobile platform was solved, and continuous and accurate output of real-time water surface elevation was achieved.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- HUBEI HAIPAI MARINE TECH DEV CO LTD
- Filing Date
- 2026-05-22
- Publication Date
- 2026-07-24
AI Technical Summary
In mobile platform scenarios, carrier phase cycle slips, severe fluctuations in platform attitude, and structural flexible deflection and radio frequency multipath coupling affect the continuous observation chain, causing water surface elevation inversion interruption errors and distorted attitudes to enter the time delay differential observation equation, affecting the continuous output of real-time water surface elevation.
By synchronously acquiring direct carrier phase observation data, reflected carrier phase observation data, and inertial high-frequency attitude data, a multi-source synchronous observation dataset is constructed. Cycle slip detection is performed, and the ambiguity at the breakpoint is retained as a new unknown parameter in the observation equation. The vertical displacement feedforward compensation is generated by combining the structural resonance deformation energy equivalent and the non-line-of-sight multipath phase divergence rate. A time-delay differential observation equation is constructed to solve for the real-time water surface elevation.
In carrier phase observation data processing, cycle slip breaks no longer interrupt the time delay differential solution chain, and attitude compensation processing reduces structural flexibility deflection and radio frequency multipath interference, thus realizing continuous output of carrier phase state observation equations and accurate solution of real-time water surface elevation.
Smart Images

Figure CN122260348B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of satellite navigation and positioning, and more specifically, to a method for water surface elevation inversion based on BeiDou TDCP differential. Background Technology
[0002] Global Navigation Satellite System (GNSS) reflectance measurement utilizes the observational characteristics of L-band signals emitted by satellites after reflection from the water surface to invert water level. Existing technologies have developed a processing approach based on carrier phase and dual-frequency combination for water level inversion. For example, patent CN101846746A describes a carrier phase altimetry device and method based on GNSS-R technology, which uses open-loop tracking of reflected signals and time-difference phase inversion to invert water level height; patent CN117310770A describes a sea level height inversion method based on BeiDou ionospheric-free combination precise single-point positioning, which completes sea level height inversion through dual-frequency ionospheric-free combination, parameter estimation, and spectral analysis.
[0003] The above scheme illustrates that existing technologies can perform water surface elevation inversion based on carrier phase observation data, dual-frequency combined observations, and differential processing. However, in mobile platform scenarios, carrier phase cycle slips, severe platform attitude fluctuations, and structural flexible deflection combined with radio frequency multipath coupling still collectively affect the continuous observation chain. Especially when solving consecutively between adjacent observation epochs, if the integer ambiguity repair approach is still used, or inertial attitude compensation is directly employed, discontinuity errors and distorted attitude can easily enter subsequent observation equations, thus affecting the construction of the time-delay differential observation equations and the continuous output of real-time water surface elevation.
[0004] To address the aforementioned problems, a technical solution is provided. Summary of the Invention
[0005] To overcome the aforementioned deficiencies in the prior art, embodiments of the present invention provide a water surface elevation inversion method based on BeiDou TDCP differential. This method involves synchronously acquiring direct carrier phase observation data, reflected carrier phase observation data, and inertial high-frequency attitude data to construct a multi-source synchronous observation dataset. Then, cycle slip detection is performed on the carrier phase observation data, and the ambiguity at the breakpoints is retained as a new unknown parameter in the observation equation. Subsequently, the structural resonance deformation energy equivalent, non-line-of-sight multipath phase divergence rate, and rigid-flexible physical fidelity coefficient are jointly generated to obtain the vertical displacement feedforward compensation. Finally, a time-delay differential observation equation is constructed, and the real-time water surface elevation is solved to address the problems mentioned in the background art.
[0006] To achieve the above objectives, the present invention provides the following technical solution: Simultaneously receive direct carrier phase observation data and reflected carrier phase observation data from BeiDou satellites, and acquire inertial high-frequency attitude data of the mobile platform. Strictly align the data with timestamps to construct a multi-source synchronous observation dataset. Cycle slip detection is performed on carrier phase observation data using a combination of wide-lane and ionospheric residuals. When a cycle slip breakpoint is identified, integer ambiguity repair is suspended, and the ambiguity at the corresponding breakpoint is retained as a new unknown parameter in the observation equation. The equivalent energy of structural resonant deformation is calculated based on inertial high-frequency attitude data, and the phase divergence rate of non-line-of-sight multipath is calculated based on carrier phase observation data. Rigid-flexible physical fidelity coefficients are generated, and the vertical displacement feedforward compensation is obtained accordingly. The vertical displacement feedforward compensation is injected into the observation model to construct the time delay difference observation equation between adjacent observation epochs. The new unknown parameters are eliminated by algebraic subtraction of epochs, and the real-time water surface elevation is solved according to the functional relationship between the time delay difference observation equation and the geometric path of water surface reflection.
[0007] Furthermore, a multi-source synchronous observation dataset is constructed, including: using the observation epoch of carrier phase observation data as the main index, pairing direct observation records and reflected observation records formed by the same BeiDou satellite at the same observation epoch to form a satellite communication link; attaching the inertial high-frequency attitude data and its associated time domain segments that are strictly aligned with the observation epoch to the corresponding satellite communication link; and recording the adjacency relationship between the current observation epoch and the immediately preceding observation epoch.
[0008] Furthermore, it synchronously receives direct carrier phase observation data and reflected carrier phase observation data from BeiDou satellites, and strictly aligns the data with timestamps, including: receiving direct signals through a right-hand circularly polarized antenna facing the zenith, and receiving reflected signals through a left-hand circularly polarized antenna facing the water surface; and the main control chip establishes a unified time reference based on second pulse interrupts and standard time messages, interpolating the inertial high-frequency attitude data surrounding the observation epoch to the corresponding observation epoch time.
[0009] Furthermore, cycle slip detection is performed using wide-lane combinations and ionospheric residual combinations, including: constructing wide-lane combinations and ionospheric residual combinations along each satellite communication link, and comparing the corresponding detection quantities according to adjacent observation epochs; when any detection quantity satisfies the breakpoint discrimination boundary, the current observation epoch is written into the cycle slip breakpoint marker, and the corresponding satellite communication link is divided into consecutive observation segments according to the cycle slip breakpoint marker, and the corresponding consecutive observation segment number is written in.
[0010] Furthermore, the ambiguity at the corresponding breakpoint is retained as a new unknown parameter in the observation equation, including: assigning a new unknown parameter to each continuous observation segment, and making each observation epoch in the same continuous observation segment refer to the corresponding new unknown parameter in the carrier phase state observation equation; when the current observation epoch and the immediately preceding observation epoch belong to the same continuous observation segment, a differential call flag is written.
[0011] Furthermore, the equivalent energy of structural resonant deformation is calculated based on inertial high-frequency attitude data, including: extracting vertical acceleration from inertial high-frequency attitude data and integrating it to obtain a vertical velocity sequence; performing frequency domain decomposition on the vertical velocity sequence; extracting the high-frequency power spectral density within the structural resonant frequency band after removing the low-frequency wave envelope; and combining it with the pre-calibrated rigid body mass of the moving platform to obtain the equivalent energy of structural resonant deformation.
[0012] Furthermore, the non-line-of-sight multipath phase divergence rate is calculated based on carrier phase observation data, including: combining the dual-frequency pseudorange corresponding to the carrier phase observation data with the dual-frequency carrier observation values to obtain the micro-multipath high-frequency composite residual, and performing zero-mean line crossing statistics on the micro-multipath high-frequency composite residual within the same continuous observation segment to form the non-line-of-sight multipath phase divergence rate organized according to the satellite communication link and observation epoch.
[0013] Furthermore, rigid-flexible physical fidelity coefficients are generated, and vertical displacement feedforward compensation is obtained accordingly. This includes: inputting the structural resonant deformation energy equivalent and the non-line-of-sight multipath phase divergence rate into a prefabricated cross-modal autoencoder to obtain reconstruction results; generating rigid-flexible physical fidelity coefficients based on the reconstruction results; dynamically scaling the state update covariance of the Kalman filter based on the rigid-flexible physical fidelity coefficients; and mapping the three-dimensional attitude fluctuation components of the platform into vertical displacement feedforward compensation under the constraint of the state update covariance.
[0014] Furthermore, the vertical displacement feedforward compensation is injected into the observation model to construct the time delay differential observation equation between adjacent observation epochs. This includes: for satellite communication links with differential call markers, extracting the direct carrier phase observation data, reflected carrier phase observation data, and vertical displacement feedforward compensation of the current observation epoch and the immediately preceding observation epoch, constructing the corresponding carrier phase state observation equation, and performing algebraic subtraction on the carrier phase state observation equation.
[0015] Furthermore, the real-time water surface elevation is solved based on the functional relationship between the time delay differential observation equation and the water surface reflection geometric path. This includes substituting the phase difference residual after algebraic subtraction of epochs, along with the antenna height, satellite elevation angle, and the real-time water surface elevation of the next adjacent epoch, into the water surface reflection geometric path, and using a least squares iterative algorithm to output the real-time water surface elevation of the current epoch.
[0016] The technical effects and advantages of the water surface elevation inversion method based on BeiDou TDCP differential in this invention are as follows: In carrier phase observation data processing, this invention no longer treats cycle slip breaks as outliers that must be repaired before continued use. Instead, it first retains the ambiguity at the breakpoint as a new unknown parameter in the observation equation, and then eliminates this parameter during the algebraic subtraction process between adjacent observation epochs. Therefore, cycle slip breaks do not directly disrupt the subsequent time delay differential solution chain, and the carrier phase state observation equation can maintain a clear organizational relationship within continuous observation segments, making it suitable for continuously outputting water surface elevation results around the observation epoch.
[0017] In attitude compensation processing, this invention does not directly provide the compensation amount based on inertial high-frequency attitude data. Instead, it first calculates the equivalent energy of structural resonant deformation and the phase divergence rate of non-line-of-sight multipath propagation, then generates rigid-flexible physical fidelity coefficients, and uses these coefficients to constrain the state update covariance to form the vertical displacement feedforward compensation amount. Thus, the compensation amount entering the observation model is not simply determined by attitude change, but is based on joint discrimination by the inertial and carrier sides, which helps reduce the direct interference to the compensation chain when structural flexible deflection and radio frequency multipath propagation coexist. Attached Figure Description
[0018] Figure 1 This is a flowchart illustrating the overall process of a water surface elevation inversion method based on BeiDou TDCP differential positioning according to the present invention. Figure 2 This is a diagram of the observation architecture of a shipborne GNSS-R dual-antenna and inertial measurement unit in a water surface elevation inversion method based on BeiDou TDCP differential according to the present invention. Figure 3 This is a schematic diagram of the data organization of a multi-source synchronous observation dataset in a water surface elevation inversion method based on BeiDou TDCP differential in this invention; Figure 4 This is a schematic diagram of cycle slip detection and ambiguity variable segmentation isolation in a water surface elevation inversion method based on BeiDou TDCP differential in this invention; Figure 5 This is a schematic diagram of the generation of rigid-flexible physical fidelity coefficients and the output of vertical displacement feedforward compensation in a water surface elevation inversion method based on BeiDou TDCP differential in this invention. Figure 6 This is a schematic diagram of the single-epoch observation model and time delay difference solution in the water surface elevation inversion method based on BeiDou TDCP differential in this invention. Detailed Implementation
[0019] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.
[0020] Please see Figures 1-6 This invention provides a water surface elevation inversion method based on BeiDou TDCP differential, comprising: S1: Simultaneously receive direct carrier phase observation data and reflected carrier phase observation data from BeiDou satellites, and acquire inertial high-frequency attitude data of the mobile platform. Strictly align the timestamps of the data to construct a multi-source synchronous observation dataset. S2: Cycle slip detection is performed on the carrier phase observation data using a wide-lane combination and an ionospheric residual combination. When a cycle slip breakpoint is identified, integer ambiguity repair is suspended, and the ambiguity at the corresponding breakpoint is retained as a new unknown parameter in the observation equation. S3: Calculate the structural resonant deformation energy equivalent based on the inertial high-frequency attitude data, calculate the non-line-of-sight multipath phase divergence rate based on the carrier phase observation data, generate rigid-flexible physical fidelity coefficients, and obtain the vertical displacement feedforward compensation amount accordingly. S4: Inject the vertical displacement feedforward compensation into the observation model to construct the time delay difference observation equation between adjacent observation epochs. Eliminate the new unknown parameters by algebraic subtraction of epochs, and solve the real-time water surface elevation according to the functional relationship between the time delay difference observation equation and the water surface reflection geometric path.
[0021] Specifically, this invention focuses on processing direct carrier phase observation data, reflected carrier phase observation data, and inertial high-frequency attitude data from mobile platforms of the BeiDou satellite. First, a multi-source synchronous observation dataset is constructed using a unified time reference, enabling each observation epoch to form a continuously accessible satellite communication link within the same data organization framework. Then, cycle slip detection is performed on the carrier phase observation data, and after identifying cycle slip breakpoints, the ambiguity at the breakpoints is preserved with entirely new unknown parameters, eliminating the need for integer ambiguity repair in subsequent processing. Based on this, the equivalent energy of structural resonance deformation and the non-line-of-sight multipath phase divergence rate are extracted from the inertial and carrier sides, respectively, to form rigid-flexible physical fidelity coefficients and further obtain the vertical displacement feedforward compensation. Finally, the vertical displacement feedforward compensation is injected into the observation model to construct the time delay differential observation equation between adjacent observation epochs, and the real-time water surface elevation is output by combining the functional relationship of the water surface reflection geometric path.
[0022] In the BeiDou L-band water surface elevation inversion chain, the direct carrier phase observation data, reflected carrier phase observation data, and the mobile platform's inertial high-frequency attitude data have different sources, sampling frequencies, and time bases. If directly processed by carrier phase observation, the satellite communication link and attitude state at the same observation epoch may lose their correspondence. Therefore, before proceeding with cycle slip detection and attitude compensation calculations, a unified data organization foundation needs to be established around the observation epochs to ensure that various observations can be paired, linked, and recorded under the same time base. Thus, step S1's task is to collect heterogeneous sampling sources, unify timing, and organize them into a directly callable multi-source synchronous observation dataset.
[0023] S101. Trigger dual-channel observation and acquisition and establish a unified time reference.
[0024] Upon receiving the system acquisition command, the main control chip first simultaneously activates the data acquisition channels of the GNSS-R dual-antenna system and the inertial measurement unit using the same hardware time base, ensuring that the carrier phase observation data and the inertial high-frequency attitude data are aligned from the start of acquisition. The GNSS-R dual-antenna system used here consists of a right-hand circularly polarized antenna facing the zenith and a left-hand circularly polarized antenna facing the water surface. The right-hand circularly polarized antenna is only responsible for receiving direct signals from the BeiDou satellite, while the left-hand circularly polarized antenna is only responsible for receiving echoes of the same source signal reflected by the water surface, thus distinguishing the direct observation link from the reflection observation link at the source. The satellite communication link refers to the pairing result of one direct observation record and one reflection observation record formed by the same BeiDou satellite at the same observation epoch. This pairing result serves as the smallest processing unit for cycle slip detection, continuous observation segment division, and time delay differential solution in subsequent steps.
[0025] The main control chip demodulates and encapsulates the BeiDou B1C and B2a carrier phase observation data output from the two antennas, forming direct carrier phase observation sequences and reflected carrier phase observation sequences arranged recursively according to the observation time. Simultaneously, the inertial measurement unit continuously outputs three-dimensional attitude angle, angular velocity, and acceleration measurements, which the main control chip retains as a high-frequency inertial attitude data sequence. To ensure that both types of data enter the same time reference frame, the main control chip uses a second pulse interrupt as an integer second anchor point and the satellite time provided by the standard time message as an absolute time tag, binding each carrier phase observation epoch and each inertial sampling point to the same continuous time axis.
[0026] In the shipborne embodiment, S101 to S103 below are all executed by the main control chip. When the ship enters the survey line and issues a data acquisition command, the right-hand circularly polarized antenna begins to record the direct carrier phase of a certain Beidou satellite, and the left-hand circularly polarized antenna simultaneously records the carrier phase of the satellite after reflection from the water surface. The inertial measurement unit outputs the ship's pitch, roll, heading, and corresponding high-frequency motion in parallel. The main control chip does not perform any prior deletion on any of the data links, but first completes the unified time stamp to retain complete input basis for subsequent epoch-based pairing.
[0027] S102. Use carrier phase observation epochs as the main index to match inertial high-frequency attitude data.
[0028] After unifying the time stamp, the main control chip uses the low-frequency epoch of the carrier phase observation data as the primary index to perform same-scale pairing on the direct carrier phase observation sequence, the reflected carrier phase observation sequence, and the inertial high-frequency attitude data sequence. Specifically, the main control chip first pairs the direct carrier phase observation sequence and the reflected carrier phase observation sequence with the same epoch according to the satellite identifier and frequency point, retaining only the satellite communication link where both direct and reflected observation values exist at the same observation time, thereby forming epoch-based observation units for carrier phase observation data; each epoch-based observation unit contains a clear observation time, satellite communication link, and the direct and reflected carrier phases under that link.
[0029] After obtaining the epoch-ratio observation unit, the main control chip uses this observation time as a reference to retrieve adjacent sampling points surrounding that time from the inertial high-frequency attitude data sequence and performs time interpolation alignment. For the inertial high-frequency attitude data vector... When the carrier phase observation epoch is And satisfy At that time, the main control chip generates the epoch according to the following formula. Strictly aligned inertial high-frequency attitude data: ; In the formula, and These represent the inertial high-frequency attitude data vectors located on either side of the carrier phase observation epoch. and Indicates the corresponding sampling time. This represents the carrier phase observation epoch time to be aligned. The purpose of this formula is to project the inertial high-frequency attitude data onto a unified time scale where the carrier phase observation epoch is located, so that the subsequent step S3 can directly call the attitude quantity at the same scale as the current carrier observation when calculating the structural resonance deformation energy equivalent, while retaining the original high-frequency sequence for frequency domain processing.
[0030] In the shipborne embodiment, if a certain carrier phase observation epoch corresponds to a phase of rapid hull undulation, the main control chip does not only capture a single attitude value, but also generates... At the same time, A continuous segment of nearby inertial high-frequency attitude data is attached to the epoch-level observation unit as an auxiliary time-domain segment. In this way, when extracting vertical acceleration, integrating to obtain the vertical velocity sequence, and calculating the high-frequency power spectral density in step S3, data can be directly read from this auxiliary time-domain segment without the need for cross-source backtracking retrieval.
[0031] S103. Construct a multi-source synchronous observation dataset and output the subsequent calling relationships.
[0032] After completing the epochal matching between the observation unit and the inertial high-frequency attitude data, the main control chip constructs a multi-source synchronous observation dataset according to the two-layer organization rule of the observation epoch-satellite communication link. Each record in this dataset consists of the direct carrier phase, the reflected carrier phase, the inertial high-frequency attitude data strictly aligned with that time, and the corresponding inertial high-frequency attitude data's associated time domain segment at the same observation time. The observation epoch is used for subsequent step S4 to perform adjacent epoch difference, the satellite communication link is used for subsequent step S2 to perform link-by-link cycle slip detection and ambiguity parameter isolation at breakpoints, and the inertial high-frequency attitude data and its associated time domain segment are used for subsequent step S3 to calculate the structural resonance deformation energy equivalent and complete the state update of the attitude compensation link, respectively.
[0033] To maintain a closed loop in subsequent calls, the main control chip records the adjacency relationships of each observation epoch while writing the multi-source synchronous observation dataset. This allows each current observation epoch to directly index the corresponding record on the same satellite communication link of the immediately preceding observation epoch. This adjacency relationship is fixed and retained as a data organization rule. Its function is as follows: step S2 can perform wide-lane combination and ionospheric residual combination operations on the carrier phase observation data of adjacent epochs along the same satellite communication link; step S4 can directly extract the observation equations of the current observation epoch and the immediately preceding observation epoch and perform algebraic subtraction.
[0034] In the shipborne embodiment, when a BeiDou satellite is simultaneously observed by both a right-hand circularly polarized antenna and a left-hand circularly polarized antenna within two consecutive observation epochs, the main control chip will continuously write the direct carrier phase, reflected carrier phase, and their respective matched inertial high-frequency attitude data of the satellite's communication link into the multi-source synchronous observation dataset. When proceeding to step S2, it can directly determine whether a carrier phase jump breakpoint occurs along the link. When proceeding to step S4, it can directly extract the records of these two adjacent epochs to establish the time delay differential observation equation without reorganizing the original data.
[0035] After processing in step S1, the main control chip has formed a multi-source synchronous observation dataset organized around the observation epoch-satellite communication link. This dataset includes direct carrier phase and reflected carrier phase observation data available for step S2, inertial high-frequency attitude data and associated time-domain segments available for step S3, and retains the traceable adjacency relationships between adjacent epochs. At this point, the heterogeneous sampling sources have been compressed into the same time reference and the same calling framework, and subsequent steps can directly proceed to cycle slip detection, physical fidelity determination, and time delay differential solution without requiring secondary time synchronization or pairing.
[0036] Corresponding to step S2, after the multi-source synchronous observation dataset has been formed, although the carrier phase observation data has the conditions to be processed according to the observation epoch and satellite communication link, the phase-locked instability of the mobile platform under high sea states will still cause discontinuous jumps in the carrier phase sequence. Once such jumps are directly introduced into the state observation equation, they will not only disrupt the continuity within the same satellite communication link, but also break the difference relationship between the ambiguity term and adjacent observation epochs. Therefore, it is necessary to first identify the breakpoints and clarify the ambiguity retention method at the breakpoints. Thus, step S2 continues to perform cycle slip detection, breakpoint segmentation, and ambiguity variable isolation, outputting the carrier phase state observation equation organization result that can be directly called in step S4.
[0037] S201. Construct a link-by-link dual-frequency cycle slip detection sequence.
[0038] Step S2 starts from the multi-source synchronous observation dataset output in step S1. All subsequent steps are executed by the main control chip. The main control chip first extracts the direct carrier phase observation data and reflected carrier phase observation data of adjacent observation epochs sequentially from each satellite communication link according to the two-layer index of the observation epoch-satellite communication link. Then, under each observation type, it reads the B1C carrier phase observation value and the B2a carrier phase observation value, forming a dual-frequency detection input sequence for the same satellite communication link. Since step S1 has preserved the adjacency relationship of the immediately preceding observation epoch, the main control chip can directly and synchronously obtain the corresponding record of the previous observation epoch when reading the current observation epoch, thereby establishing a link-by-link recursive timing detection window.
[0039] After obtaining the dual-frequency detection input sequence, the main control chip performs a difference operation on the B1C carrier phase observations and B2a carrier phase observations under the same satellite communication link and the same observation type, after normalizing the wavelength to the cycle unit, to generate a wide-lane ambiguity combination sequence. At the same time, the dual-frequency carrier phase observations are first differentially divided according to adjacent observation epochs, and then inversely combined according to the frequency square relationship between the two frequency points to generate an ionospheric residual combination sequence.
[0040] Specifically, for the first The satellite communication link is in the observation epoch. The wide-lane ambiguity combination at the location is first converted by the main control chip into cycle units for both the B1C and B2a carrier phase observations, denoted as follows: and Press again Construct the wide-lane ambiguity combination value, and according to The variation in width of the channel is constructed as a large cycle slip detection quantity.
[0041] Specifically, the main control chip multiplies the B1C carrier phase observation value and the B2a carrier phase observation value by the corresponding carrier wavelength to obtain the carrier length difference between adjacent observation epochs. and ,in , Press again Construct the ionospheric residual combination values, and according to Construct a small cycle slip detection quantity.
[0042] Wide-lane ambiguity combination sequences are used to amplify the abrupt change characteristics of large cycle slips in the dual-frequency difference, while ionospheric residual combination sequences are used to preserve the residual abrupt change characteristics of small cycle slips between adjacent observation epochs. Both types of sequences are continuously written to the current processing cache using satellite communication links and observation epochs as indices.
[0043] In the shipborne embodiment, when the ship experiences a short-term heave and causes a sudden jump in the B2a carrier phase observation value of the reflection channel, the main control chip first extracts the current observation epoch and the B1C carrier phase observation value and the B2a carrier phase observation value immediately adjacent to the previous observation epoch from the reflection carrier phase observation data corresponding to the satellite communication link. Then, it sequentially generates the wide-lane ambiguity combination value and the ionospheric residual combination value of the link, so that the same jump has a comparable temporal position in both types of detection sequences, providing a unified entry point for subsequent breakpoint determination.
[0044] S202. Determining cycle slip breakpoints based on wide lane combination and ionospheric residual combination.
[0045] After obtaining the link-by-link dual-frequency cycle slip detection sequence, the main control chip calculates the changes in the wide-lane ambiguity combination sequence and the ionospheric residual combination sequence for adjacent observation epochs according to the satellite communication link. For any current observation epoch, the main control chip first calculates the absolute difference between the current wide-lane ambiguity combination value and the wide-lane ambiguity combination value of the immediately preceding observation epoch, and then calculates the absolute difference between the current ionospheric residual combination value and the ionospheric residual combination value of the immediately preceding observation epoch, thereby forming the large cycle slip detection quantity and the small cycle slip detection quantity at the same epoch position.
[0046] The historical sequence of changes in the width of the alley is composed of the 20 most recent valid adjacent observation epochs within the same continuous observation segment. The composition is as follows: when there are fewer than 20 valid samples within a continuous observation segment, the historical sequence of wide-lane changes recorded during the factory static water calibration phase of the satellite communication link is used as a substitute, and the 90th percentile of the substitute sequence is taken as the boundary for wide-lane discontinuity discrimination. The historical sequence of ionospheric residual changes is composed of the most recent 20 valid adjacent observation epochs within the same continuous observation segment. When there are fewer than 20 valid samples in a continuous observation period, the historical sequence of ionospheric residual changes recorded during the factory static water calibration stage of the satellite communication link is used as the replacement, and the 90th percentile of the replacement sequence is taken as the ionospheric residual breakpoint discrimination boundary.
[0047] The main control chip first performs a coarse judgment on large cycle slip detections using the wide-lane breakpoint discrimination boundary. If the current absolute difference exceeds the wide-lane breakpoint discrimination boundary, the current observation epoch is directly marked as a cycle slip breakpoint. If the current absolute difference does not exceed the wide-lane breakpoint discrimination boundary, the chip continues to perform a fine judgment on small cycle slip detections using the ionospheric residual breakpoint discrimination boundary. If the current absolute difference exceeds the ionospheric residual breakpoint discrimination boundary, the current observation epoch is also marked as a cycle slip breakpoint. The main control chip writes the above judgment results, marked as cycle slip breakpoints, into the satellite communication link record corresponding to the multi-source synchronous observation dataset, and uses this mark as a boundary to divide the satellite communication link into consecutive observation segments.
[0048] Once a cycle slip breakpoint is written to an observation epoch, the main control chip immediately updates the continuous observation segment number along the satellite communication link. This ensures that the continuous observation records before the breakpoint retain their original continuous observation segment numbers, while the observation epoch at the breakpoint and subsequent continuous observation records enter a new continuous observation segment. The resulting continuous observation segment is not an auxiliary description, but rather a direct index for subsequent ambiguity variable assignments.
[0049] In the shipborne embodiment, if a reflective satellite communication link experiences a significant abrupt change in wide-lane ambiguity combination at the current observation epoch after the bow strikes the water, the main control chip immediately writes a cycle slip breakpoint marker at that epoch and cuts that epoch out of the original continuous observation segment. If another direct-satellite communication link does not experience a significant abrupt change in wide-lane ambiguity combination, but the ionospheric residual combination exhibits a brief jump between adjacent observation epochs, the main control chip writes the same cycle slip breakpoint marker at that epoch using a small cycle slip path. Both paths ultimately output a unified cycle slip breakpoint marker and continuous observation segment number, and subsequent processing no longer distinguishes between detection sources.
[0050] S203. Assign entirely new unknown parameters to the observation equations and form subsequent call tags.
[0051] After the cycle slip breakpoint markers and continuous observation segment numbers are generated, the main control chip no longer uses the ambiguity representation method shared across segments of the same satellite communication link. Instead, it rewrites the ambiguity terms of the carrier phase state observation equation according to the continuous observation segments. Specifically, the main control chip assigns a completely new unknown parameter to each continuous observation segment of each satellite communication link and establishes a unique index for it using the satellite communication link number and the continuous observation segment number. All observation epochs within the same continuous observation segment refer to the same completely new unknown parameter when writing the carrier phase state observation equation, while different continuous observation segments refer to their own independent completely new unknown parameters. For the newly assigned unknown parameter, the main control chip uses the compensated observation value of the current observation epoch minus the geometric path value obtained by geometric extrapolation according to the immediate previous observation epoch as the initial value, and sets its initial covariance to a preset multiple of the upper limit of the factory static calibration ambiguity covariance to maintain the independent degrees of freedom of this variable relative to the previous continuous observation segment within the new continuous observation segment.
[0052] Correspondingly, for any satellite communication link within a certain continuous observation segment, the main control chip first converts the reflected carrier phase observation value and the direct carrier phase observation value at the current observation epoch into length domain observation values and calculates the difference to establish the carrier phase state observation equation for the satellite communication link at that observation epoch. That is, the length domain phase difference observation value consists of a geometric path term, an ambiguity constant term, a common error term, and a measurement noise term. The geometric path term represents the propagation path change of the satellite communication link at the current observation epoch; the common error term represents the environmental and equipment errors that change slowly with the observation epoch over a short time scale; and the measurement noise term represents the random disturbance at the current observation epoch. When no cycle slip breakpoint occurs, the same satellite communication link shares a single ambiguity constant term across epochs. The rewriting refers to the fact that after step S202 outputs the cycle slip breakpoint marker and the continuous observation segment number, the main control chip no longer allows all observation epochs of the same satellite communication link to share the same ambiguity constant term, but instead replaces the ambiguity constant term with an ambiguity variable indexed by the continuous observation segment number. If the current observation epoch belongs to the consecutive observation segment numbered as Then, the ambiguity variable is uniformly referenced for each observation epoch within the continuous observation segment. Once a new cycle slip breakpoint is identified and a new continuous observation segment number is generated, the observation epoch at which the breakpoint occurs and the subsequent continuous observation records are changed to reference the ambiguity variable corresponding to the next continuous observation segment. Thus, the ambiguity variable remains unchanged within the same continuous observation segment, while the ambiguity variables between different continuous observation segments are independent of each other. The location of the cycle slip is thus preserved in the carrier phase state observation equation in a parameter-isolated manner.
[0053] To enable step S4 to directly extract adjacent observation epoch pairs that can be subtracted, the main control chip, while writing the ambiguity variable, further compares the consecutive observation segment numbers of the current observation epoch with those of the immediately preceding observation epoch. When the consecutive observation segment numbers are the same, the main control chip writes a difference call flag, indicating that the adjacent observation epoch pair shares the same entirely new unknown parameter, and can then directly enter the time-delay differential observation equation. When the consecutive observation segment numbers switch at the cycle slip breakpoint, the main control chip does not write a difference call flag, but only retains the index of the new ambiguity variable corresponding to the current observation epoch. The adjacent observation epoch pair that can be called by step S4 is then formed after the next observation epoch enters the same consecutive observation segment.
[0054] In the shipborne embodiment, when the reflective satellite communication link of a Beidou satellite has stably recorded multiple observation epochs in the first continuous observation segment, the main control chip will cause the carrier phase state observation equations of these observation epochs to jointly reference the same new unknown parameter. When a subsequent observation epoch is determined to be a cycle slip breakpoint in step S202, the main control chip will start the second continuous observation segment at that epoch and write another new unknown parameter for that segment. In this way, adjacent observation epochs in the first continuous observation segment can be directly subtracted in step S4, and the second continuous observation segment will form a new differentially callable link starting from the second observation epoch within it.
[0055] After processing in step S2, each satellite communication link in the multi-source synchronous observation dataset has obtained cycle slip breakpoint markers, continuous observation segment numbers, new unknown parameters indexed by continuous observation segments, and differential call markers. The ambiguity constant term in the carrier phase state observation equation has also been segmented and isolated. Therefore, subsequent step S3 can still directly call the carrier phase observation data that has not been renamed, and step S4 can extract adjacent observation epoch pairs that share the same new unknown parameters based on the differential call markers, and the constant term will be naturally eliminated when the epochs are subtracted.
[0056] After the segmentation relationship of the discontinuities in the carrier phase observation data has been clarified, the impact of the mobile platform's attitude fluctuations on the observation model cannot still be directly compensated using the inertial high-frequency attitude data with equal weights. This is because structural flexible deflection and radio frequency multipath disturbances under high sea states may occur simultaneously within the same observation epoch, causing deviations between the mechanical side measurement results and the carrier side residual performance. At this point, step S3 needs to extract discriminant quantities that can characterize mechanical deformation and radio frequency divergence from the inertial high-frequency attitude data and the carrier phase observation data, respectively, to determine the constraint method when attitude compensation is incorporated into the observation model.
[0057] S301. Calculate the structural resonant deformation energy equivalent based on inertial high-frequency attitude data.
[0058] Step S3 begins by starting with the multi-source synchronous observation dataset output from step S1. The main control chip takes the inertial high-frequency attitude data corresponding to the current observation epoch and its associated time-domain segment as input. First, according to the pitch, roll, and yaw angles already aligned for the current observation epoch, it rotates the body coordinate system acceleration components output by the inertial measurement unit to a local horizontal coordinate system. Then, it extracts the vertical acceleration sequence along the vertical direction of the local horizontal coordinate system. Since the associated time-domain segment shares the same time reference as the current observation epoch, the main control chip can directly integrate this vertical acceleration sequence in chronological order without cross-source resampling to obtain the vertical velocity sequence corresponding to the current observation epoch.
[0059] After obtaining the vertical velocity sequence, the main control chip performs a Fast Fourier Transform on the sequence and constructs a power spectral density curve with frequency as the horizontal axis and the squared vertical velocity spectral density as the vertical axis. The main control chip first locates the low-frequency wave envelope peak on this power spectral density curve, then searches for the first energy valley point from this peak towards the high-frequency side, determining the frequency corresponding to this valley point as the lower boundary of the structural resonant frequency band. Subsequently, the Nyquist frequency sampled by the inertial measurement unit is used as the upper boundary of the structural resonant frequency band, retaining only the high-frequency power spectral density within this band to characterize the structural flutter intensity between the rigid hull base and the antenna mounting location. The main control chip further combines this with the pre-calibrated rigid body mass of the mobile platform to calculate the equivalent structural resonant deformation dynamic energy using the following formula: ; In the formula, This represents the structural resonance deformation energy equivalent of the current observation epoch. This refers to the pre-calibrated rigid body mass of the mobile platform. This represents the vertical velocity power spectral density obtained after removing the low-frequency wave envelope. This indicates the lower boundary of the structural resonant frequency band. This represents the upper boundary of the structural resonant frequency band. The calculation result is written to the processing buffer with the current observation epoch as the index, serving as the first input to the subsequent cross-modal autoencoder.
[0060] In the shipborne embodiment, S301 to S303 described below are all executed by the main control chip. When the hull passes through a wave ridge and experiences short-term high-frequency flutter, the main control chip extracts the vertical acceleration sequence from the time-domain segment attached to the inertial high-frequency attitude data of the current observation epoch. After attitude rotation, time integration, and spectral decomposition, an independent high-frequency energy band can be identified to the right of the main peak of the low-frequency wave envelope. This high-frequency energy band is integrated into the structural resonant deformation kinetic energy equivalent, thereby converting the visible hull flutter into a physical quantity that can be used for subsequent calculations.
[0061] S302. Calculate the non-line-of-sight multipath phase divergence rate based on carrier phase observation data.
[0062] After obtaining the structural resonance deformation energy equivalent of the current observation epoch, the main control chip continues to read the carrier phase observation data corresponding to the current observation epoch and which has been assigned continuous observation segment numbers by step S2 from the multi-source synchronous observation dataset. The dual-frequency pseudorange introduced here is the auxiliary ranging observation value output by the BeiDou receiving channel at the same observation epoch as the carrier phase observation data. The main control chip first reads the dual-frequency pseudorange and dual-frequency carrier observation values of the B1C frequency point and the B2a frequency point at the current observation epoch according to the satellite communication link. Then, according to the correspondence of the same observation epoch, the same satellite communication link, and the same frequency point, the above four observations are organized into a set of radio frequency residual calculation inputs. Among them, the dual-frequency pseudorange and dual-frequency carrier observation value of the B1C frequency point constitute the first frequency point input, and the dual-frequency pseudorange and dual-frequency carrier observation value of the B2a frequency point constitute the second frequency point input. In subsequent processing, the main control chip performs dual-frequency non-ionospheric combination on the two sets of inputs respectively, and constructs the combined residual at the same observation epoch with the combined pseudorange result and the combined carrier result. In order to avoid bringing the breakpoint error at the cycle slip breakpoint into the residual statistics, the main control chip only organizes the observation sequence that traces back from the current observation epoch within the same continuous observation segment. Any position that crosses the cycle slip breakpoint mark is used as the cutoff boundary of the current statistical window.
[0063] Within the same continuous observation segment, the main control chip first constructs a dual-frequency pseudorange combination and a dual-frequency carrier combination using the dual-frequency pseudorange and carrier combination coefficients, respectively, based on the dual-frequency pseudorange and carrier observations. Then, the difference between the dual-frequency pseudorange combination and the dual-frequency carrier combination forms the combined residual for the current observation epoch. Subsequently, the main control chip calculates the moving average of the combined residuals for the most recent 10 valid observation epochs within the current statistical window and uses this moving average as the low-frequency trend term for the current observation epoch. When there are fewer than 10 valid samples in the current continuous observation segment, the average of the combined residuals of all available samples in the continuous observation segment is used as the low-frequency trend term. The low-frequency trend term specifically corresponds to the average offset and cross-epoch gradual offset that persist between multiple adjacent observation epochs within the continuous observation segment; its rate of change is significantly lower than the high-frequency fluctuations caused by micro-multipath propagation. The main control chip subtracts the low-frequency trend term from the combined residual for the current observation epoch to obtain the micro-multipath high-frequency composite residual.
[0064] For the dual-frequency pseudorange and dual-frequency carrier observations at frequency points 1 and 2, the main control chip respectively... and A dual-frequency ionospheric pseudorange combination and a dual-frequency ionospheric carrier combination are constructed, and the difference between the two is used as the basic quantity of the combined residual. The microscopic multipath high-frequency composite residual is calculated as follows: Construction, in which For the 10 most recent valid observation epochs within the same continuous observation segment Moving mean; when there are fewer than 10 valid samples in the current continuous observation segment, the mean of all available samples in that continuous observation segment is used instead. In statistics At that time, the main control chip only records a valid zero-mean line crossing event when the zero-mean residuals of two adjacent sampling points have opposite signs and both of their absolute values are greater than the median value of still water noise in the current continuous observation segment; when the current statistical window contains only one valid sample, the main control chip records the satellite communication link as a valid event. Record it as 0.
[0065] For the current satellite communication link within the current statistical window, the [number]th [item / item / etc.] For each sampling point, the main control chip first obtains the zero-mean sequence of the microscopic multipath high-frequency composite residual, then counts the number of times the zero-mean sequence crosses the zero-mean line, and converts it into the non-line-of-sight multipath phase divergence rate according to the window duration: ; In the formula, Indicates the epoch of the current observation. Non-line-of-sight multipath phase divergence rate of a satellite communication link This indicates the number of times the microscopic multipath high-frequency composite residual of the satellite communication link crosses the zero mean line within the current statistical window. This indicates the duration covered by the statistical window. The larger the value, the denser the high-frequency multipath phase flipping in the radio frequency environment of the satellite communication link.
[0066] Since the current observation epoch typically corresponds to multiple available satellite communication links, after completing the link-by-link calculation, the main control chip summarizes the non-line-of-sight multipath phase divergence rates of all satellite communication links marked with differential call flags for the same epoch, and takes the median as the non-line-of-sight multipath phase divergence rate for the current observation epoch. The purpose of median aggregation is to maintain the ability of the current observation epoch to characterize the common radio frequency environment of most satellite communication links, while preventing individual abnormal links from dominating subsequent fidelity determination.
[0067] In the shipborne embodiment, when the ship rolls and causes changes in the near-water surface obstruction of the reflecting antenna, the microscopic multipath high-frequency composite residuals of certain reflecting satellite communication links will repeatedly and densely cross the zero-mean line within continuous observation segments. The main control chip counts the number of these crossings along each continuous observation segment, and then summarizes the statistical results of multiple available satellite communication links at the current observation epoch into a non-line-of-sight multipath phase divergence rate, so that the degree of radio frequency disturbance enters the next processing stage with a unified epoch value.
[0068] S303. Generates rigid-flexible physical fidelity coefficients based on a cross-modal autoencoder and outputs vertical displacement feedforward compensation.
[0069] The structural resonance deformation energy equivalent at the current observation epoch is obtained. Non-line-of-sight multipath phase divergence rate Then, the main control chip unifies the dimensions of both according to the normalized scale fixed at the factory, and constructs the cross-modal autoencoder input vector. The cross-modal autoencoder is an unsupervised neural network model pre-trained using massive amounts of homomorphic data from calm water areas at the time of device delivery. During execution, it only performs forward reconstruction calculations; the main control chip processes the input vector... After being fed into the encoder and decoder, the reconstructed vector is obtained. The Euclidean distance between the input vector and the reconstructed vector is used as the reconstruction error for the current observation epoch, and then rigid-flexible physical fidelity coefficients are generated according to the negative exponential normalization rule. ; In the formula, This represents the rigid-flexible physics fidelity coefficient for the current observation epoch. This represents the Euclidean norm. The closer it is to 1, the closer the structural resonant morphological change energy equivalent and the non-line-of-sight multipath phase divergence rate are to the homomorphic relationship in the training samples of the steady water area. The smaller the value, the more pronounced the physical disconnect between mechanical deformation and radio frequency divergence.
[0070] After obtaining the rigid-flexible physical fidelity coefficients, the main control chip extracts the platform's three-dimensional attitude fluctuation components from the inertial high-frequency attitude data that is strictly aligned with the current observation epoch. The platform's three-dimensional attitude fluctuation components are defined as the deviations of the pitch, roll, and yaw angles at the current observation epoch from the low-frequency attitude baseline within the same subordinate time domain segment. The main control chip organizes these components into an attitude fluctuation vector. Then, the rigid-flexible physical fidelity coefficients are applied to the state update covariance of the Kalman filter to obtain the dynamic state update covariance of the current observation epoch: ; In the formula, This represents the baseline state update covariance at the current observation epoch without the introduction of rigid-flexible physical fidelity coefficients. This represents the dynamically scaled state update covariance. The main control chip uses... The Kalman filter update gain calculation is used to automatically reduce the attitude fluctuation update weights when the rigid-flexible physical fidelity coefficient is low, and the output is constrained by the state update covariance to form the platform's three-dimensional attitude fluctuation components. .
[0071] The low-frequency attitude baseline is formed by taking the pitch, roll, and yaw angle sequences from the corresponding time-domain segment of the current observation epoch, and passing them through cutoff frequencies no higher than the lower boundary of the structural resonant frequency band. The low-pass filtering is used to obtain the platform's three-dimensional attitude fluctuation components, which are formed by subtracting the corresponding low-frequency attitude baseline from the original attitude sequence and then taking the values at the current observation epoch. Reference state update covariance. The noise matrix, composed of the pitch, roll, and yaw measurement noise variances obtained from the factory calibration of the inertial measurement unit, is added to the attitude propagation covariance of the current observation epoch. Its diagonal elements correspond to the reference update variances of the pitch, roll, and yaw states, respectively.
[0072] After obtaining the three-dimensional attitude fluctuation components of the platform constrained by the state update covariance, the main control chip introduces the antenna phase center lever arm vector. The antenna phase center lever arm vector is a fixed spatial vector determined during equipment installation and calibration, pointing from the origin of the inertial measurement unit coordinate system to the antenna phase center. The main control chip denotes it as... Simultaneously, based on the local horizontal coordinate system of the current observation epoch, the unit vector of the bottom reflection water surface normal is denoted as... Under small-angle attitude fluctuation conditions, the main control chip performs a cross product between the platform's three-dimensional attitude fluctuation components constrained by the state update covariance and the antenna phase center lever arm vector to obtain the instantaneous displacement vector of the antenna phase center. This instantaneous displacement vector is then projected along the unit vector normal to the bottom-level reflected water surface to obtain the vertical displacement feedforward compensation. ; In the formula, This represents the vertical displacement feedforward compensation amount for the current observation epoch. The main control chip writes this vertical displacement feedforward compensation amount into the carrier phase state observation equation call cache according to the observation epoch, so that it can be directly injected into the observation model in step S4 when constructing the time delay differential observation equation.
[0073] In one embodiment, the cross-modal autoencoder can employ a symmetrical fully connected autoencoder network structure. The input layer consists of two neurons, corresponding to the structural resonant deformation energy equivalent and the non-line-of-sight multipath phase divergence rate, respectively. The encoder is configured with 8 neurons, 4 neurons, and 2 neurons sequentially, while the decoder is configured with 4 neurons, 8 neurons, and 2 neurons sequentially. Except for the output layer, all hidden layers employ modified linear activation functions, and the output layer employs a linear activation function to maintain consistency between the reconstruction result and the input dimensions. The training dataset can be continuously collected in calm waters before the device leaves the factory, for example, collecting 1.2 million sets of synchronous samples. Each set of samples consists of the structural resonant deformation energy equivalent and the non-line-of-sight multipath phase divergence rate at the same observation epoch, and is paired with the satellite communication link according to the timestamp before being written into the training set. Of these, 800,000 sets are used for training, 200,000 sets for validation, and 200,000 sets for testing. Before training, the structural resonant deformation energy equivalent and the non-line-of-sight multipath phase divergence rate are normalized. The normalization parameters can be obtained statistically from the training set; for example, the mean of the structural resonant deformation energy equivalent is 0.42 and the standard deviation is 0.18, and the mean of the non-line-of-sight multipath phase divergence rate is 1.36 and the standard deviation is 0.41. These mean and standard deviations are then stored as pre-training parameters in the main control chip. The training method can be unsupervised reconstruction training. The loss function is the mean square error between the input vector and the reconstruction vector. The optimizer uses an adaptive moment estimation algorithm with a learning rate of 0.001, a batch size of 256, and 200 training epochs. Training stops when the validation set loss does not decrease for 20 consecutive epochs, and the set of network weights and biases with the minimum validation set loss is saved as pre-training parameters. During device operation, only the stored network weights, biases, and the normalized parameters are called to perform forward reconstruction calculations on the structural resonant deformation energy equivalent and the non-line-of-sight multipath phase divergence rate input at the current observation epoch.
[0074] In the shipborne embodiment, when the deck of the ship exhibits significant elastic torsion on a certain segment of the track and the non-line-of-sight multipath phase divergence rate of multiple satellite communication links increases synchronously, the reconstruction error output by the cross-modal autoencoder will be significantly amplified. The main control chip then obtains a smaller rigid-flexible physical fidelity coefficient and expands the state update covariance of the current observation epoch accordingly. Subsequently, the main control chip performs constrained updates on the pitch, roll, and heading angle fluctuations of the same observation epoch, and then converts the attitude fluctuations into the vertical displacement feedforward compensation of the antenna phase center relative to the bottom water surface through the antenna phase center lever arm vector, thereby completing the amplitude constraint of the attitude information contaminated by mechanical flexibility before entering the observation equation. Corresponding to step S3, this step calculates the structural resonant deformation energy equivalent based on inertial high-frequency attitude data and calculates the non-line-of-sight multipath phase divergence rate based on carrier phase observation data. The two are then input into a prefabricated cross-modal autoencoder to generate rigid-flexible physical fidelity coefficients. On this basis, the state update covariance is dynamically scaled according to the rigid-flexible physical fidelity coefficients, and the platform's three-dimensional attitude fluctuation components are mapped to the vertical displacement feedforward compensation amount of the antenna phase center relative to the bottom water surface, thus giving the attitude quantities clear credibility constraints before entering the observation model.
[0075] Corresponding to step S4, after the continuous observation segments of the carrier phase observation data have been divided, the new unknown parameters have been retained in the observation equations, and the vertical displacement feedforward compensation has been formed, the current processing chain has the foundation to construct the differential relationship between adjacent observation epochs. At this point, step S4 needs to incorporate the observation-side compensation, the segmented carrier phase state observation equations, and the water surface reflection geometry into the same solution framework, so that the ambiguity terms that remain constant within the same continuous observation segment are eliminated through epoch algebraic subtraction, and the remaining observations can directly correspond to the water surface elevation.
[0076] S401. Construct a single-epoch observation model for injecting vertical displacement feedforward compensation.
[0077] Step S4 starts simultaneously from the multi-source synchronous observation dataset, differential call markers, and vertical displacement feedforward compensation. The main control chip first scans the current observation epoch records one by one according to the satellite communication links, retaining only the records with differential call markers. Then, for each retained satellite communication link, it synchronously extracts the direct carrier phase observation data, reflected carrier phase observation data, continuous observation segment number, and the new unknown parameter index allocated in step S2 for the current observation epoch. Since the differential call marker is only written when the current observation epoch and the immediately preceding observation epoch belong to the same continuous observation segment, the new unknown parameters obtained by the main control chip at this stage naturally satisfy the condition of remaining unchanged across two observation epochs, and can directly proceed to subsequent epochs for difference calculation.
[0078] At the single-epoch level, the main control chip first converts the reflected carrier phase observation data and the direct carrier phase observation data of the same satellite communication link into length domain observation values, and then subtracts the vertical displacement feedforward compensation amount corresponding to the current observation epoch from the observation side to form the single-epoch observation model after injection compensation: ; In the formula, Indicates the first The carrier wavelength currently used in this satellite communication link Indicates the current observation epoch. Phase observations of the reflected carrier of the satellite communication link. Indicates the current observation epoch. Direct carrier phase observations for each satellite communication link. This represents the output of step S3 and the vertical displacement feedforward compensation amount corresponding to the current observation epoch. Indicates the first The geometric path function of water surface reflection for the satellite communication link at the current observation epoch. This indicates the real-time water surface elevation to be determined at the current observation epoch. This indicates that step S2 is in the continuous observation segment The internally isolated new and unknown parameters are preserved. This represents the common error term that changes slowly over a short timescale at the current observation epoch. This represents the measurement noise term. The main control chip writes the single-epoch observation model into the current solution cache according to the observation epoch and satellite communication link.
[0079] In the shipborne embodiment, the following steps are all executed by the main control chip. When multiple satellite communication links simultaneously have differential call markers at a certain observation epoch, the main control chip reads the direct carrier phase observation data, reflected carrier phase observation data, and vertical displacement feedforward compensation amount of these links respectively, and forms a single-epoch observation model after injection compensation for each link; if a satellite communication link is located at the beginning observation epoch of a new continuous observation segment, the record does not yet have a differential call marker, and the main control chip only caches the single-epoch observation model without entering the differential calculation in this round.
[0080] S402. Perform algebraic subtraction on adjacent observation epochs to form a time-delay difference observation equation.
[0081] After the single-epoch observation model with compensation is established, the main control chip, according to the observation epoch adjacency relationship retained in step S1, synchronously extracts the two single-epoch observation models of the current observation epoch and the immediately preceding observation epoch for each satellite communication link marked with differential call flags, and directly performs algebraic subtraction while keeping the same continuous observation segment number and the same new unknown parameter index unchanged. The resulting time-delay differential observation equation is: ; in, ; In the formula, Indicates the first The compensated time delay differential observation values between the current observation epoch and the immediately preceding observation epoch of the satellite communication link. This represents the real-time water surface elevation immediately preceding the previous observation epoch. This is because the two single-epoch observation models share the same entirely new unknown parameter. The constant term is naturally eliminated when subtracting.
[0082] After performing algebraic subtraction, the main control chip further performs a short-time stationarity determination on the common error term based on whether the time interval between the current observation epoch and the immediately preceding observation epoch are adjacent and extremely short timescales. When the time interval meets the current sampling period and both epochs are within the same continuous observation segment, the main control chip will... Treating the low dynamic common error residuals as approximately zero, the time delay differential observation equation simplifies to a pure phase difference residual equation containing only the geometric path difference of water surface reflection and the differential noise term: ; In the formula, This represents the noise term after epoch difference. The main control chip uses this simplified time-delay difference observation equation as the formal solution equation for the current observation epoch, and excludes records that do not meet the difference call conditions from the solution set for this round.
[0083] In the shipborne embodiment, when the same satellite communication link remains locked for two consecutive observation epochs and no new cycle slip breakpoint marker is triggered, the main control chip directly performs a difference operation on the single-epoch observation models of these two observation epochs. After the difference operation, the equation no longer contains the new unknown parameters that were isolated in step S2, nor does it explicitly retain the slowly changing common error terms. If a cycle slip breakpoint marker occurs in the intermediate epoch, the link will not generate a difference call marker at the breakpoint. The main control chip will not forcibly perform a difference operation across the breakpoint, but will wait for a new pair of adjacent observation epochs that can be differenced to be re-formed within the new continuous observation segment.
[0084] S403. Solve by performing least squares iteration based on the water surface reflection geometric path function.
[0085] After generating all usable time-delay differential observation equations for the current observation epoch, the main control chip begins to organize the water surface reflection geometric path function. For the next observation epoch... A satellite communication link is established. The main control chip reads the antenna altitude and satellite elevation angle of the current observation epoch, and establishes the water surface reflection geometric path function according to the two-way trigonometric projection relationship. ; In the formula, Indicates the antenna height at the current observation epoch. Indicates the current observation epoch. The satellite elevation angle of the satellite communication link. Correspondingly, the geometric path function for water surface reflection in the immediate preceding observation epoch is... The main control chip uses the real-time water surface elevation obtained immediately preceding the previous observation epoch as the initial value for the current observation epoch, denoted as... Then, the time delay differential observation equations of all available satellite communication links are stacked in link order to form the equation set for the current observation epoch.
[0086] in, Take the installation calibration height of the antenna phase center relative to the equipment measurement datum plane; when the vertical displacement feedforward compensation has been injected on the observation side, the geometric path function... Maintain the installation height as specified and avoid repeated stacking. Satellite elevation angle The angle between the satellite line-of-sight vector and the local horizontal coordinate system is determined after the navigation message and ephemeris parameters obtained from the demodulation of the direct-view channel at the current observation epoch are calculated.
[0087] The main control chip then uses a least-squares iterative algorithm to calculate the real-time water surface elevation for the current observation epoch. For the ... In the next iteration, the main control chip first estimates the current water surface elevation. Linearize the surface reflection geometric path function of each satellite communication link to form the observation residual vector. and design vector , of which The design vector elements corresponding to each satellite communication link are: ; The observed residual vector is the first The item is: ; Based on this, the main control chip calculates the current iteration increment: ; ; Update the real-time water surface elevation estimate for the current observation epoch. When When the value is not greater than the smaller of the median absolute value of the differential residuals of all available satellite communication links at the current observation epoch and the upper limit of the single-epoch delay differential noise obtained from the factory static water calibration, the main control chip stops iterating and outputs a value. As the real-time water surface elevation for the current observation epoch; if only one satellite communication link satisfies the differential call flag for the current observation epoch, the main control chip will still solve the problem according to the same iterative process, but will retain the differential residual corresponding to that link for continuous verification in subsequent epochs.
[0088] In the shipborne embodiment, when multiple subtractable satellite communication links from different satellites exist simultaneously at the current observation epoch, the main control chip stacks the compensated time delay differential observations of these links sequentially and uniformly substitutes them into the least squares iterative algorithm. After completing several iterations, the real-time water surface elevation of the current observation epoch is written back to the observation epoch record corresponding to the multi-source synchronous observation dataset and continues to participate in the subsequent recursive solution as the solution value of the next observation epoch that is immediately adjacent to the previous observation epoch.
[0089] After processing in step S4, the vertical displacement feedforward compensation is injected into the carrier phase state observation equation, constructing the time delay differential observation equation between adjacent observation epochs. Then, based on the differential call mark, the observation records within the same continuous observation segment are extracted and epochal algebraic subtraction is performed to eliminate the new unknown parameters during the differential process. Subsequently, the differential residual is solved by combining the functional relationship between antenna height, satellite elevation angle and water surface reflection geometry path, forming the real-time water surface elevation continuously output according to the observation epoch.
[0090] The above embodiments are merely preferred embodiments of the present invention and are not intended to limit the scope of protection of the present invention. Any modifications, equivalent substitutions, or improvements made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.
Claims
1. A method for water surface elevation inversion based on BeiDou TDCP differential, characterized in that, Including the following steps: Simultaneously receive direct carrier phase observation data and reflected carrier phase observation data from BeiDou satellites, and acquire inertial high-frequency attitude data of the mobile platform. Strictly align the data with timestamps to construct a multi-source synchronous observation dataset. Cycle slip detection is performed on carrier phase observation data using a combination of wide-lane and ionospheric residuals. When a cycle slip breakpoint is identified, integer ambiguity repair is suspended, and the ambiguity at the corresponding breakpoint is retained as a new unknown parameter in the observation equation. The equivalent energy of structural resonant deformation is calculated based on inertial high-frequency attitude data, and the phase divergence rate of non-line-of-sight multipath is calculated based on carrier phase observation data. Rigid-flexible physical fidelity coefficients are generated, and the vertical displacement feedforward compensation is obtained accordingly. The vertical displacement feedforward compensation is injected into the observation model to construct the time delay difference observation equation between adjacent observation epochs. The new unknown parameters are eliminated by algebraic subtraction of epochs, and the real-time water surface elevation is solved according to the functional relationship between the time delay difference observation equation and the geometric path of water surface reflection.
2. The water surface elevation inversion method based on BeiDou TDCP differential as described in claim 1, characterized in that, Construct a multi-source synchronous observation dataset, including: using the observation epoch of carrier phase observation data as the main index, pairing direct observation records and reflected observation records formed by the same BeiDou satellite at the same observation epoch to form a satellite communication link; attaching the inertial high-frequency attitude data and its associated time domain segments that are strictly aligned with the observation epoch to the corresponding satellite communication link; and recording the adjacency relationship between the current observation epoch and the immediately preceding observation epoch.
3. The water surface elevation inversion method based on BeiDou TDCP differential as described in claim 2, characterized in that, The system synchronously receives direct carrier phase observation data and reflected carrier phase observation data from BeiDou satellites, and performs strict timestamp alignment on the data. This includes receiving direct signals through a right-hand circularly polarized antenna facing the zenith and receiving reflected signals through a left-hand circularly polarized antenna facing the water surface. The main control chip establishes a unified time reference based on second pulse interrupts and standard time messages, and interpolates the inertial high-frequency attitude data surrounding the observation epoch to the corresponding observation epoch time.
4. The water surface elevation inversion method based on BeiDou TDCP differential as described in claim 1, characterized in that, Cycle slip detection is performed using wide-lane combination and ionospheric residual combination, including: constructing wide-lane combination and ionospheric residual combination along each satellite communication link, and comparing the corresponding detection quantities according to adjacent observation epochs; when any detection quantity meets the breakpoint discrimination boundary, the current observation epoch is written into the cycle slip breakpoint mark, and the corresponding satellite communication link is divided into consecutive observation segments according to the cycle slip breakpoint mark, and the corresponding consecutive observation segment number is written.
5. The water surface elevation inversion method based on BeiDou TDCP differential as described in claim 4, characterized in that, The ambiguity at the corresponding breakpoint is retained as a new unknown parameter in the observation equation, including: assigning a new unknown parameter to each continuous observation segment, and making each observation epoch in the same continuous observation segment refer to the corresponding new unknown parameter in the carrier phase state observation equation; when the current observation epoch and the immediately preceding observation epoch belong to the same continuous observation segment, a differential call flag is written.
6. The water surface elevation inversion method based on BeiDou TDCP differential as described in claim 1, characterized in that, The calculation of structural resonant deformation energy equivalent based on inertial high-frequency attitude data includes: extracting vertical acceleration from inertial high-frequency attitude data and integrating it to obtain a vertical velocity sequence; performing frequency domain decomposition on the vertical velocity sequence; extracting the high-frequency power spectral density within the structural resonant frequency band after removing the low-frequency wave envelope; and combining it with the pre-calibrated rigid body mass of the moving platform to obtain the structural resonant deformation energy equivalent.
7. The water surface elevation inversion method based on BeiDou TDCP differential as described in claim 1, characterized in that, The calculation of non-line-of-sight multipath phase divergence rate based on carrier phase observation data includes: combining the dual-frequency pseudorange corresponding to the carrier phase observation data with the dual-frequency carrier observation values to obtain the micro-multipath high-frequency composite residual, and performing zero-mean line crossing statistics on the micro-multipath high-frequency composite residual within the same continuous observation segment to form the non-line-of-sight multipath phase divergence rate organized according to the satellite communication link and observation epoch.
8. The water surface elevation inversion method based on BeiDou TDCP differential as described in claim 1, characterized in that, The process of generating rigid-flexible physical fidelity coefficients and obtaining vertical displacement feedforward compensation amounts includes: inputting the structural resonant deformation energy equivalent and the non-line-of-sight multipath phase divergence rate into a prefabricated cross-modal autoencoder to obtain reconstruction results; generating rigid-flexible physical fidelity coefficients based on the reconstruction results; dynamically scaling the state update covariance of the Kalman filter based on the rigid-flexible physical fidelity coefficients; and mapping the three-dimensional attitude fluctuation components of the platform into vertical displacement feedforward compensation amounts under the constraints of the state update covariance.
9. The water surface elevation inversion method based on BeiDou TDCP differential as described in claim 1, characterized in that, The vertical displacement feedforward compensation is injected into the observation model to construct the time delay differential observation equation between adjacent observation epochs. This includes: for satellite communication links with differential call markers, extracting the direct carrier phase observation data, reflected carrier phase observation data, and vertical displacement feedforward compensation of the current observation epoch and the immediately preceding observation epoch, constructing the corresponding carrier phase state observation equation, and performing algebraic subtraction on the carrier phase state observation equation.
10. A water surface elevation inversion method based on BeiDou TDCP differential as described in claim 9, characterized in that, The real-time water surface elevation is solved based on the functional relationship between the time-delay differential observation equation and the water surface reflection geometry path. This includes substituting the phase difference residual after algebraic subtraction of epochs, along with the antenna height, satellite elevation angle, and the real-time water surface elevation of the current observation epoch, into the water surface reflection geometry path, and then using a least squares iterative algorithm to output the real-time water surface elevation of the current observation epoch.