Method, system and equipment for joint inversion of vector flow field using space-time wave current of high-frequency ground wave radar

Abstract

The method, system and equipment for the joint inversion of space-time wave and current vector flow field by high-frequency ground wave radar belong to the technical field of high-frequency ground wave radar ocean monitoring. In order to solve the problem that the existing vector velocity inversion method has large errors due to the obvious changes in ocean currents. The present invention first adopts the spatial correlation method, and uses the radial velocity of the ocean current of the adjacent beam sea element to estimate the vector velocity of the current beam sea element, and then judges the out-of-tolerance position sea element where the ocean current obviously changes according to the law of ocean current movement, and derives the relationship between the change of ocean current and wave height in time and space through the ocean dynamic equation, calculates the tangential velocity of the ocean current, and synthesizes it with the radial velocity vector inverted by the radar to obtain the vector velocity at the current moment. The present invention is used for wave and current joint vector field inversion.

Description

Technical Field

[0001] The invention belongs to the technical field of high-frequency ground wave radar ocean monitoring, and in particular relates to a space-time joint vector field inversion method for ocean currents by a single-station radar. Background Art

[0002] High-frequency ground-wave radar can monitor areas below the horizon through the diffraction characteristics of electromagnetic waves, and has the advantages of wide range and all-weather real-time monitoring. Ocean echo is one of the main clutter interferences for high-frequency ground-wave radar target detection, but it also contains a large amount of sea state information and is the main detection object for sea state inversion.

[0003] A large number of observation results show that the generation mechanism of ocean echoes can be explained by the Bragg resonance principle, that is, when high-frequency electromagnetic waves are irradiated on the rough sea surface, the electromagnetic waves will interact with the sea surface and cause strong scattering, thereby generating first-order ocean echoes. In theory, the first-order ocean echo will form a pair of energy peaks symmetrical about zero frequency on the frequency axis of the radar echo RD spectrum, namely the Bragg peak. By detecting the frequency corresponding to the actual ocean echo spectrum peak and calculating the difference with the Bragg frequency, the radial velocity of the ocean current can be obtained, and then the vector velocity can be inverted through the radial velocity.

[0004] Usually, a more accurate vector flow velocity can be obtained by using a dual-station radar to invert the flow field. However, the construction cost of a dual-station radar is high, the location selection requirements are high, and the two radars are required to cover the same sea area at the same time. Therefore, it is more urgent and necessary to study the method of inverting the vector flow field with a single-station radar.

[0005] The traditional vector velocity inversion method mainly uses the adjacent vector synthesis method. The adjacent vector synthesis method is based on the fact that the vector velocity of two adjacent beam positions within a certain range is the same. At this time, the radial velocity corresponding to the two adjacent sea elements is the projection of the vector velocity in the two beam directions. Therefore, the vector velocity can be solved by the radial velocity of each beam, the angle between beams, and the angle of each beam. However, there may be obvious changes in ocean currents in the actual monitoring area. In this case, the traditional method of calculating the vector velocity will result in large errors. Summary of the invention

[0006] The present invention aims to solve the problem that the existing vector velocity inversion method has large errors due to the obvious changes in ocean currents.

[0007] The method for joint inversion of vector flow field by high-frequency ground wave radar space-time wave current includes the following steps:

[0008] Step 1: Use the matching peak method to detect the frequency f corresponding to the first-order spectrum peak of the actual ocean echo from the received radar echo spectrum d ;

[0009] Step 2: Based on fd and the theoretical Bragg frequency f B The difference Δf gives the radial velocity V of the ocean current. r ;

[0010] Step 3: Invert vector velocity using spatial correlation method:

[0011] Estimate the vector flow velocity of the current beam based on the radial velocity of the current beam and the adjacent left and right beams:

[0012]

[0013] Among them, u and v represent the vector velocity components of the ocean current in the x direction and y direction respectively, Vr i-1 、Vr i 、Vr i+1 Represents the radial velocities of the left beam, current beam, and right beam on the current range unit, respectively. The positive direction is toward the radar, θ i-1 ,θ i ,θ i+1 Respectively represent the angles of the left beam, current beam, and right beam starting from the positive x direction in the two-dimensional Cartesian coordinate system;

[0014] Convert formula (3) into matrix form and get:

[0015] Vr=G·X (8)

[0016]

[0017] According to the least squares theory, when the beam angle θ is known, i and radial velocity Vr i Under the condition of , we can calculate the matrix X corresponding to u and v, and get the velocity V of ocean current according to u and v. f and flow to α;

[0018] Step 4: Flow velocity V based on ocean current movement f and flow direction α, to determine the extreme position of the ocean element where the ocean current changes significantly:

[0019] Assume that the current moment is n and the next moment is n+1. According to the movement law of the ocean current, calculate the distance gate I corresponding to the next position of the ocean current of each sea element. n+1 , beam J n+1 and the time T when it flows to the next position n+1 ;

[0020] For the ocean current within the radar range, calculate the velocity difference ΔV between it and the corresponding ocean current at the next position f and the angle difference Δα, and compare it with the pre-set speed difference threshold β and angle difference threshold γ; if there is ΔV f>β or Δα>γ, or T n+1 For ocean currents that are not within the time series of radar observation data, proceed to step 5;

[0021] Step 5: Use the space-time wave-current joint method to invert the vector flow field. The specific process includes the following steps:

[0022] S51. According to the continuity equation and motion equation in the control equation group of the ocean parameter dynamic model, the discretized forms of u and v are obtained:

[0023]

[0024] Where ζ is the wave height, subscripts j and k are the unit numbers in the y and x directions in the two-dimensional Cartesian coordinate system, n represents the time, Δt represents the time step, Δx and Δy represent the grid intervals in the x and y directions in the two-dimensional Cartesian coordinate system, respectively; Ω is the value of the earth's angular velocity, φ is the local latitude; g represents the acceleration of gravity;

[0025] S52. There is a vector relationship between the tangential velocity Vs and u and v:

[0026] -Vs=usin(π-θ)+(-v)cos(π-θ) (23)

[0027] Vs=-usinθ-vcosθ (24)

[0028] Based on the form of formula (24), we can obtain the following derivative with respect to time:

[0029]

[0030] Substituting formula (22) into formula (25) yields:

[0031]

[0032] The relationship between the ocean current tangential velocity and wave height is obtained from formula (26):

[0033]

[0034] in, and Represent the tangential velocity of the ocean current on the (k, j) unit at the current moment and the previous moment, respectively. Represents the wave height on the (k,j) unit at the previous moment, and They represent the velocity components of the ocean current velocity in the x direction and the y direction in the cell (k, j) at the previous moment, and θ is the angle of the radar beam from the x direction to the y direction in the rectangular coordinate system;

[0035] Taking the movement of ocean currents in the x direction as an example, according to the CFL condition:

[0036]

[0037] Among them, σ is called the Courant number;

[0038] Get the value range of the time step Δt:

[0039]

[0040] S53, correct the ocean current at the out-of-tolerance position determined in step 4, and calculate the tangential velocity at the current moment using formula (27) according to and Get vector flow velocity and flow direction in, Represents the radial velocity of the ocean current in the cell (k, j) at time n.

[0041] Furthermore, step 2 is based on f d and the theoretical Bragg frequency f B The difference Δf gives the radial velocity V of the ocean current. r as follows:

[0042]

[0043] Among them, λ is the wavelength of the electromagnetic wave emitted by the radar, V r is the radial velocity of the ocean current, which is positive when the current moves toward the radar and negative when it moves away from the radar.

[0044] Furthermore, the theoretical Bragg frequency f B as follows:

[0045]

[0046] Where g is the acceleration due to gravity.

[0047] Furthermore, the velocity of the ocean current obtained in step 3 according to u and v is

[0048] Furthermore, the flow direction α=π-arctan(v / u) obtained in step 3 according to u and v.

[0049] Furthermore, according to the continuity equation and the motion equation in the control equation group of the ocean parameter dynamic model, the process of obtaining the discretized form of u and v includes the following steps:

[0050] The continuity equation is:

[0051]

[0052] Wherein, x, y, z represent the three axes of the three-dimensional Cartesian coordinate system; w is the velocity component of the fluid in the vertical direction, and ρ is the density of seawater;

[0053] The equation of motion is:

[0054]

[0055]

[0056]

[0057] Among them, f x 、f y 、f z is the Coriolis force component, A μ is the fluid viscosity coefficient, Γ is the tidal potential of the celestial body, g represents the gravitational acceleration, and P represents the seawater pressure;

[0058] Then, based on the continuity equation and the motion equation, the partial differentials of u and v with respect to time t are obtained:

[0059]

[0060] After gridding the data, formula (20) is transformed into a discretized form:

[0061]

[0062] Furthermore, the process of obtaining the partial differentials of u and v with respect to time t based on the continuity equation and the motion equation includes the following steps:

[0063] For the equation of motion, ignoring the Coriolis force, tidal forces of celestial bodies, and viscosity terms, when the fluid is stationary, the upward pressure gradient force is balanced by the downward force, and we get:

[0064]

[0065] The hydrostatic pressure of the ocean includes the pressure caused by sea level changes, the pressure caused by ocean density stratification, and atmospheric pressure; the vertical pressure is expressed as follows:

[0066] P=ρ0gζ (17)

[0067] Among them, ρ0 is the average density, ζ is the wave height;

[0068] Combining equation (17), the horizontal component of unit mass is expressed as follows:

[0069]

[0070]

[0071] Then we get:

[0072]

[0073] The system for joint inversion of vector flow field by high-frequency ground wave radar space-time wave current includes:

[0074] The radial velocity calculation unit of the ocean current: using the matching peak method, the frequency f corresponding to the first-order spectrum peak of the actual ocean echo is detected from the received radar echo spectrum d ; Based on f d and the theoretical Bragg frequency f B The difference Δf gives the radial velocity V of the ocean current. r ;

[0075] Spatial correlation method inversion unit: The spatial correlation method is used to invert the vector velocity. The specific inversion process includes the following steps:

[0076] Estimate the vector flow velocity of the current beam based on the radial velocity of the current beam and the adjacent left and right beams:

[0077]

[0078] Among them, u and v represent the vector velocity components of the ocean current in the x direction and y direction respectively, Vr i-1 、Vr i 、Vr i+1 Represents the radial velocities of the left beam, current beam, and right beam on the current range unit, respectively. The positive direction is toward the radar, θ i-1 ,θ i ,θ i+1 Respectively represent the angles of the left beam, current beam, and right beam starting from the positive x direction in the two-dimensional Cartesian coordinate system;

[0079] Convert formula (3) into matrix form and get:

[0080] Vr=G·X (8)

[0081]

[0082] According to the least squares theory, when the beam angle θ is known, i and radial velocity Vr i Under the condition of , we can calculate the matrix X corresponding to u and v, and get the velocity V of ocean current according to u and v. f and flow to α;

[0083] Determination unit of the sea element of the over-error position: velocity V based on the movement of ocean currents f and flow direction α, to determine the extreme position sea element where the ocean current changes significantly. The specific process includes the following steps:

[0084] Assume that the current moment is n and the next moment is n+1. According to the movement law of the ocean current, calculate the distance gate I corresponding to the next position of the ocean current of each sea element. n+1 , beam J n+1 and the time T when it flows to the next position n+1 ;

[0085] For the ocean current within the radar range, calculate the velocity difference ΔV between it and the corresponding ocean current at the next position f and the angle difference Δα, and compare it with the pre-set speed difference threshold β and angle difference threshold γ; if there is ΔV f >β or Δα>γ, or T n+1 Ocean currents that are not within the time series of radar observation data;

[0086] Space-time wave-current joint method inversion unit: The space-time wave-current joint method is used to invert the vector flow field. The specific process includes the following steps:

[0087] S51. According to the continuity equation and motion equation in the control equation group of the ocean parameter dynamic model, the discretized forms of u and v are obtained:

[0088]

[0089] Where ζ is the wave height, subscripts j and k are the unit numbers in the y and x directions in the two-dimensional Cartesian coordinate system, n represents the time, Δt represents the time step, Δx and Δy represent the grid intervals in the x and y directions in the two-dimensional Cartesian coordinate system, respectively; Ω is the value of the earth's angular velocity, φ is the local latitude; g represents the acceleration of gravity;

[0090] S52. There is a vector relationship between the tangential velocity Vs and u and v:

[0091] -Vs=usin(π-θ)+(-v)cos(π-θ) (23)

[0092] Vs=-usinθ-vcosθ (24)

[0093] Based on the form of formula (24), we can obtain the following derivative with respect to time:

[0094]

[0095] Substituting formula (22) into formula (25) yields:

[0096]

[0097] The relationship between the ocean current tangential velocity and wave height is obtained from formula (26):

[0098]

[0099] in, and Represent the tangential velocity of the ocean current on the (k, j) unit at the current moment and the previous moment, respectively. Represents the wave height on the (k,j) unit at the previous moment, and They represent the velocity components of the ocean current velocity in the x direction and the y direction in the cell (k, j) at the previous moment, and θ is the angle of the radar beam from the x direction to the y direction in the rectangular coordinate system;

[0100] Taking the movement of ocean currents in the x direction as an example, according to the CFL condition:

[0101]

[0102] Among them, σ is called the Courant number;

[0103] Get the value range of the time step Δt:

[0104]

[0105] S53, correct the ocean current at the out-of-tolerance position determined in step 4, and calculate the tangential velocity at the current moment using formula (27) according to and Get vector flow velocity and flow direction in, Represents the radial velocity of the ocean current in the cell (k, j) at time n.

[0106] A computer storage medium stores at least one instruction, which is loaded and executed by a processor to implement the high-frequency ground wave radar space-time wave and current joint inversion vector flow field method.

[0107] A high-frequency ground wave radar space-time wave current joint inversion vector flow field device, the device comprises a processor and a memory, the memory stores at least one instruction, the at least one instruction is loaded and executed by the processor to implement the high-frequency ground wave radar space-time wave current joint inversion vector flow field method.

[0108] Beneficial effects:

[0109] The present invention firstly adopts the spatial correlation method, uses the radial velocity of the ocean current of the adjacent beam sea element to estimate the vector velocity of the current beam sea element, and then judges the out-of-tolerance position sea element where the ocean current changes significantly according to the law of ocean current movement, derives the temporal and spatial relationship between the ocean current and the wave height through the ocean dynamic equation, calculates the ocean current tangential velocity, and synthesizes it with the radial velocity vector inverted by the radar to obtain the vector velocity at the current moment. The space-time wave and current joint method avoids the problem of large errors in inversion when the ocean current changes significantly in the traditional method, and optimizes the inversion effect. BRIEF DESCRIPTION OF THE DRAWINGS

[0110] Figure 1 This is the relationship diagram between the radial velocity of each beam and the vector flow velocity.

[0111] Figure 2 This is the relationship diagram between vector flow velocity and tangential velocity.

[0112] Figure 3 This is the overall flow chart of the space-time joint inversion method.

[0113] Figure 4 This is the vector flow field diagram simulated in Example 1.

[0114] Figure 5 This is the inversion result before correction of Example 1 (the inversion result with a larger flow velocity in a certain area).

[0115] Figure 6 This is the inversion result of the space-time joint method in Example 1.

[0116] Figure 7 This is the vector flow field simulated in Example 2.

[0117] Figure 8 This is the inversion result before correction of Example 2 (the inversion result with a larger flow velocity in a certain area).

[0118] Fig. 9 This is the inversion result of the time-space joint method in Example 2.

[0119] Fig.10 It is the satellite data of Copernicus Marine Service in Example 2. DETAILED DESCRIPTION

[0120] The present invention is based on a single-station radar and adopts a space-time joint method to invert the vector velocity. When the ocean current in a certain sea area is approximately unchanged or changes evenly and slowly, the space-correlation method is adopted to invert the vector velocity. When the ocean current changes unevenly or changes significantly, the space-time wave-current joint method is adopted to invert the vector velocity. The relationship between the changes of the flow field and the wave field in space and time is derived based on an ocean dynamic model. The tangential velocity at the current moment is calculated through the changes of the tangential velocity of the ocean current and the waves at the previous moment, and the vector velocity is obtained by synthesizing it with the radial velocity inverted by the radar.

[0121] Specific implementation method 1: Combination Figure 3 To explain this embodiment,

[0122] The method for joint inversion of space-time wave and current vector flow field by high-frequency ground wave radar described in this embodiment specifically includes the following steps:

[0123] Step 1: Use the matching peak method to detect the frequency f corresponding to the first-order spectrum peak of the actual ocean echo from the received radar echo spectrum d .

[0124] Step 2: Calculate f d and the theoretical Bragg frequency f B The difference Δf gives the radial velocity of the ocean current:

[0125]

[0126] Among them, λ is the wavelength of the electromagnetic wave emitted by the radar, V r is the radial velocity of the ocean current, which is positive when the current moves toward the radar and negative when it moves away from the radar.

[0127] Theoretical Bragg frequency f B The calculation of is as follows:

[0128]

[0129] Where g is the acceleration due to gravity.

[0130] Step 3: Invert vector velocity using spatial correlation method:

[0131] The vector velocity of the current beam is estimated based on the radial velocity of the current beam and the adjacent left and right beams. The vector relationship is as follows: Figure 1 As shown, the relationship is as follows:

[0132]

[0133] Among them, u and v represent the vector velocity components of the ocean current in the x direction and y direction respectively, Vr i-1 、Vr i 、Vr i+1 Respectively represent the radial velocity of the ocean current corresponding to the left beam, current beam, and right beam on the current range unit, with the positive direction toward the radar, θ i-1 ,θ i ,θ i+1 They respectively represent the angles of the left beam, the current beam, and the right beam starting from the positive x direction in the two-dimensional Cartesian coordinate system. The directions described below are all in the Cartesian coordinate system.

[0134] Convert formula (3) into matrix form:

[0135]

[0136] make

[0137]

[0138]

[0139]

[0140] get:

[0141] Vr=G·X (8)

[0142] According to the least squares theory, when the beam angle θ is known, i and radial velocity Vr i Calculate the matrix X corresponding to u and v in the case of:

[0143] X=(G T G) -1 G T Vr (9)

[0144] according to

[0145]

[0146] α=π-arctan(v / u) (11)

[0147] Calculate the velocity V of the ocean current f and flows to α.

[0148] Step 4: Flow velocity V based on ocean current movement f and flow direction α, to determine the extreme position of the ocean element where the ocean current changes significantly:

[0149] The movement of ocean currents can be carried out in a continuous and coherent manner, showing relatively stable flow characteristics without obvious changes in speed or direction. In this case, the ocean current can be regarded as maintaining a relatively consistent state of motion over a relatively long time scale, flowing in the ocean at a similar speed and direction. Let the current moment be n and the next moment be n+1. According to the movement law of the ocean current, calculate the distance gate I corresponding to the next position of the ocean current of each sea element n+1 , beam J n+1 and the time T when it flows to the next position n+1 .

[0150] For the ocean current within the radar range, calculate the velocity difference ΔV between it and the corresponding ocean current at the next position fCompared with the angle difference Δα, the pre-set speed difference threshold β (such as 0.1m-0.3m) and the angle difference threshold γ (such as 10°-20°):

[0151] If ΔV f ≤β and Δα≤γ, indicating that the error of the current ocean current inversion result is within the allowable range, and the spatial correlation method is feasible and effective for inverting ocean currents;

[0152] If ΔV f >β or Δα>γ, which means that there are large errors in the current ocean current inversion results. The spatial correlation method is not suitable for the current changes in sea conditions. The space-time wave and current joint method is used for subsequent inversion and correction.

[0153] For some special cases, such as T n+1 The ocean current is not within the time series of the radar observation data, or the ocean current is in the outermost beam and points outward, so it is impossible to find the next moment T n+1 , and then the space-time wave-current joint method is used for inversion and correction;

[0154] Step 5: Use the space-time wave-current joint method to invert the vector flow field. The specific process includes the following steps:

[0155] S51. According to the continuity equation and motion equation in the control equation group of the ocean parameter dynamic model, the continuity equation is:

[0156]

[0157] In this step, all vectors are in the three-dimensional Cartesian coordinate system, x, y, and z represent the three axes of the three-dimensional Cartesian coordinate system; u and v have the same meaning as in the previous steps, w is the velocity component of the fluid in the vertical direction, and ρ is the density of seawater.

[0158] The equation of motion is:

[0159]

[0160]

[0161]

[0162] Among them, f x 、f y 、f z is the Coriolis force component, A μ is the fluid viscosity coefficient, Γ is the tidal potential of the celestial body, g represents the gravitational acceleration, and P represents the sea water pressure.

[0163] For the equation of motion, it is assumed that the wave period is very small compared to the inertial period, and the Coriolis force, tidal force of celestial bodies, viscosity terms, etc. are ignored. Considering the case of hydrostatic equilibrium, that is, when the fluid is stationary, the upward pressure gradient force and the downward force are balanced, it can be obtained from equation (15):

[0164]

[0165] The hydrostatic pressure of the ocean includes the pressure caused by sea level changes, the pressure caused by ocean density stratification, and atmospheric pressure. Since the sea surface will adjust instantly when the atmospheric pressure changes, so that the pressure under the sea surface can be kept roughly the same, the sea surface has little effect on the horizontal gradient of the fluid, and assuming that the seawater density is evenly distributed, the vertical pressure can be expressed as follows:

[0166] P=ρ0gζ (17)

[0167] Among them, ρ0 is the average density and ζ is the wave height.

[0168] The force that generates the fluid is mainly the horizontal pressure gradient force. Combining equation (17), the horizontal component of its unit mass can be expressed as follows:

[0169]

[0170]

[0171] From formulas (13) and (14), ignoring the high-order terms, celestial tidal forces, and viscosity terms that have little effect on the fluid, we can obtain:

[0172]

[0173] According to the knowledge of ocean dynamics, the Coriolis force components in the x-direction and y-direction can be expressed as formula (21):

[0174]

[0175] Among them, u and v have the same meanings as in the previous steps, Ω is the value of the earth's angular velocity, and φ is the local latitude.

[0176] After gridding the data, formula (20) is transformed into a discretized form:

[0177]

[0178] Where j and k are the unit numbers in the y and x directions in the two-dimensional Cartesian coordinate system, n represents the time, Δt represents the time step, Δx and Δy represent the grid intervals in the x and y directions in the two-dimensional Cartesian coordinate system, respectively.

[0179] S52, such as Figure 2As shown, there is a vector relationship between the tangential velocity Vs and u and v as shown in the figure:

[0180] -Vs=usin(π-θ)+(-v)cos(π-θ) (23)

[0181] Vs=-usinθ-vcosθ (24)

[0182] The time derivative of formula (24) can be obtained as follows:

[0183]

[0184] Substituting formula (22) into formula (25) we can obtain:

[0185]

[0186] The relationship between the ocean current tangential velocity and wave height is obtained from formula (26):

[0187]

[0188] in, and Represent the tangential velocity of the ocean current on the (k, j) unit at the current moment and the previous moment, respectively. Represents the wave height on the (k,j) unit at the previous moment, and They represent the velocity components of the ocean current velocity in the x and y directions at the previous moment (k, j), j and k are the unit numbers in the y and x directions in the two-dimensional Cartesian coordinate system, n represents the moment, Δt represents the time step, Δx and Δy represent the grid spacing in the x and y directions in the two-dimensional Cartesian coordinate system, and θ is the angle of the radar beam from the x direction to the y direction in the rectangular coordinate system.

[0189] Taking the movement of ocean currents in the x direction as an example, according to the CFL condition:

[0190]

[0191] The choice of Δt needs to satisfy the requirement that the fluid movement distance uΔt is less than the grid length Δx, where σ is called the Courant number.

[0192] Get the value range of the time step Δt:

[0193]

[0194] S53, correct the ocean current at the out-of-tolerance position determined in step 4, and calculate the tangential velocity at the current moment using formula (27) According to formula (30), the vector flow velocity is and flow direction

[0195]

[0196] in, Represents the radial velocity of the ocean current in the cell (k, j) at time n.

[0197] The present invention firstly adopts the spatial correlation method, uses the radial velocity of the ocean current of the adjacent beam sea element to estimate the vector velocity of the current beam sea element, and then judges the out-of-tolerance position sea element where the ocean current changes significantly according to the law of ocean current movement, derives the temporal and spatial variation relationship between the ocean current and the wave height through the ocean dynamic equation, calculates the ocean current tangential velocity, and synthesizes it with the radial velocity vector inverted by the radar to obtain the vector velocity at the current moment. The present invention not only adopts the spatial correlation method, but also further introduces the time series on the basis of the spatial correlation method, and proposes a new method for inverting the vector flow field of space-time wave and current joint, which can significantly improve the inversion effect of the vector flow field when the ocean current changes significantly in the sea area. Specific implementation method 2:

[0199] The system for joint inversion of vector flow field by space-time wave and current of high-frequency ground wave radar described in this embodiment includes:

[0200] The radial velocity calculation unit of the ocean current: using the matching peak method, the frequency f corresponding to the first-order spectrum peak of the actual ocean echo is detected from the received radar echo spectrum d ; Based on f d and the theoretical Bragg frequency f B The difference Δf gives the radial velocity V of the ocean current. r ;

[0201] Spatial correlation method inversion unit: The spatial correlation method is used to invert the vector velocity. The specific inversion process includes the following steps:

[0202] Estimate the vector flow velocity of the current beam based on the radial velocity of the current beam and the adjacent left and right beams:

[0203]

[0204] Among them, u and v represent the vector velocity components of the ocean current in the x direction and y direction respectively, Vr i-1 、Vr i 、Vr i+1 Represents the radial velocities of the left beam, current beam, and right beam on the current range unit, respectively. The positive direction is toward the radar, θ i-1 ,θ i ,θ i+1 Respectively represent the angles of the left beam, current beam, and right beam starting from the positive x direction in the two-dimensional Cartesian coordinate system;

[0205] Convert formula (3) into matrix form and get:

[0206] Vr=G·X (8)

[0207]

[0208] According to the least squares theory, when the beam angle θ is known, i and radial velocity Vr i Under the condition of , we can calculate the matrix X corresponding to u and v, and get the velocity V of ocean current according to u and v. f and flow to α;

[0209] Determination unit of the sea element of the over-error position: velocity V based on the movement of ocean currents f and flow direction α, to determine the extreme position sea element where the ocean current changes significantly. The specific process includes the following steps:

[0210] Assume that the current moment is n and the next moment is n+1. According to the movement law of the ocean current, calculate the distance gate I corresponding to the next position of the ocean current of each sea element. n+1 , beam J n+1 and the time T when it flows to the next position n+1 ;

[0211] For the ocean current within the radar range, calculate the velocity difference ΔV between it and the corresponding ocean current at the next position f and the angle difference Δα, and compare it with the pre-set speed difference threshold β and angle difference threshold γ; if there is ΔV f >β or Δα>γ, or T n+1 Ocean currents that are not within the time series of radar observation data;

[0212] Space-time wave-current joint method inversion unit: The space-time wave-current joint method is used to invert the vector flow field. The specific process includes the following steps:

[0213] S51. According to the continuity equation and motion equation in the control equation group of the ocean parameter dynamic model, the discretized forms of u and v are obtained:

[0214]

[0215] Where ζ is the wave height, subscripts j and k are the unit numbers in the y and x directions in the two-dimensional Cartesian coordinate system, n represents the time, Δt represents the time step, Δx and Δy represent the grid intervals in the x and y directions in the two-dimensional Cartesian coordinate system, respectively; Ω is the value of the earth's angular velocity, φ is the local latitude; g represents the acceleration of gravity;

[0216] S52. There is a vector relationship between the tangential velocity Vs and u and v:

[0217] -Vs=usin(π-θ)+(-v)cos(π-θ) (23)

[0218] Vs=-usinθ-vcosθ (24)

[0219] Based on the form of formula (24), we can obtain the following derivative with respect to time:

[0220]

[0221] Substituting formula (22) into formula (25) yields:

[0222]

[0223] The relationship between the ocean current tangential velocity and wave height is obtained from formula (26):

[0224]

[0225] in, and Represent the tangential velocity of the ocean current on the (k, j) unit at the current moment and the previous moment, respectively. Represents the wave height on the (k,j) unit at the previous moment, and They represent the velocity components of the ocean current velocity in the x direction and the y direction in the cell (k, j) at the previous moment, and θ is the angle of the radar beam from the x direction to the y direction in the rectangular coordinate system;

[0226] Taking the movement of ocean currents in the x direction as an example, according to the CFL condition:

[0227]

[0228] Among them, σ is called the Courant number;

[0229] Get the value range of the time step Δt:

[0230]

[0231] S53, correct the ocean current at the out-of-tolerance position determined in step 4, and calculate the tangential velocity at the current moment using formula (27) according to and Get vector flow velocity and flow direction in, Represents the radial velocity of the ocean current in the cell (k, j) at time n. Specific implementation method three:

[0233] This embodiment is a computer storage medium, in which at least one instruction is stored. The at least one instruction is loaded and executed by a processor to implement the high-frequency ground wave radar space-time wave and current joint inversion vector flow field method.

[0234] It should be understood that the instructions include computer program products, software or computerized methods corresponding to any method described in the present invention; the instructions can be used to program a computer system, or other electronic device. Computer storage media may include readable media on which instructions are stored, which may include but are not limited to magnetic storage media, optical storage media; magneto-optical storage media include read-only memory ROM, random access memory RAM, erasable programmable memory (e.g., EPROM and EEPROM) and flash memory layers, or other types of media suitable for storing electronic instructions. Specific implementation method four:

[0236] This embodiment is a high-frequency ground wave radar space-time wave current joint inversion vector flow field device, the device includes a processor and a memory, it should be understood that including any device including a processor and a memory described in the present invention, the device may also include other units and modules that perform display, interaction, processing, control, etc. and other functions through signals or instructions;

[0237] At least one instruction is stored in the memory, and the at least one instruction is loaded and executed by the processor to implement the above.

[0238] Example

[0239] The embodiment of the present invention uses measured data from Weihai, Shandong Province. The sea conditions in the sea area monitored by the radar station in Weihai, Shandong Province are stable. In actual monitoring, the ocean currents usually change slowly. In the experiment, the ocean currents in a certain area are simulated to have a large difference from the nearby ocean currents. Therefore, taking a small range of drastic changes on a beam as an example, the vector velocity of some sea elements is increased to simulate the uneven ocean currents. The corresponding wave height inversion results can reflect the actual changes in the ocean currents.

[0240] Example 1

[0241] The data used in this embodiment comes from the measured data of Weihai Radar Station on October 16, 2022. In the inversion process, firstly, the inversion result of the actual vector flow field is obtained according to steps 1 to 2. On the basis of the actual vector flow field, the vector flow velocity corresponding to the 20th to 50th distance unit of the third beam is increased to obtain the simulated flow field as shown in FIG. Figure 4 shown.

[0242] According to the simulated ocean current radial velocity, execute step 3 and obtain Figure 5From the inversion results shown, it can be seen that due to the large change in unilateral radial velocity, there are large errors in the ocean current inversion at the corresponding positions of the 2nd and 4th beams. At this time, the inversion results are not accurate enough and it is necessary to make subsequent corrections.

[0243] Execute step 4. The ocean currents of the 2nd and 4th beams do not meet the uniformity condition of the ocean currents at adjacent positions and have large errors. Therefore, they are determined as out-of-tolerance positions that need to be corrected.

[0244] Execute step 5 for the second and fourth beams, solve the tangential velocity according to the ocean dynamic equation and calculate the vector velocity, and obtain the inversion results of this method as follows: Figure 6 As shown in the figure, it can be seen that the corrected ocean current is consistent with the overall trend of the ocean current in the region. Figure 5 Compared with the inversion effect, it is obviously better.

[0245] After calculation, the inversion results of the second and fourth beam current velocities are compared with the simulated vector flow field. The root mean square error of the inversion results before correction using this method is 0.93m / s, while the root mean square error of the inversion results after correction is 0.04m / s. It can be seen that this method can avoid the estimation error caused by the drastic changes in ocean currents.

[0246] Comparing the forecast data with the inversion data, the flow directions of the two are similar, and the inversion method has a better inversion effect.

[0247] Example 2

[0248] The data used in this embodiment is derived from the measured data of Weihai Radar Station on May 29, 2021. The inversion process is the same as that in Example 1. First, the inversion result of the actual vector flow field is obtained according to steps 1 to 2. On the basis of the actual vector flow field, the vector flow velocity corresponding to the 20th to 50th distance unit of the third beam is increased to obtain the simulated flow field as shown in Figure 7 shown.

[0249] According to the simulated ocean current radial velocity, execute step 3 and obtain Figure 8 From the inversion results shown, it can be seen that due to the large change in unilateral radial velocity, there are large errors in the ocean current inversion at the corresponding positions of the 2nd and 4th beams. At this time, the inversion results are not accurate enough.

[0250] Execute step 4. The ocean currents of the 2nd and 4th beams do not meet the uniformity condition of the ocean currents at adjacent positions and have large errors. Therefore, they are determined as out-of-tolerance positions that need to be corrected.

[0251] Execute step 5 for the second and fourth beams, solve the tangential velocity according to the ocean dynamic equation and calculate the vector velocity, and obtain the inversion results of this method as follows: Fig. 9As shown in the figure, it can be seen that the corrected ocean current is consistent with the overall trend of the ocean current in the region. Figure 8 Compared with the inversion effect, it is obviously better.

[0252] In the first embodiment, the inversion result is compared with the China National Marine Forecast Data, and it can be seen that the inversion result is relatively accurate. In this example, the satellite data of Copernicus Marine Service is used to compare with the inversion result to increase the credibility of the data. Fig.10 The ocean current data downloaded from Copernicus Marine Service. It can be seen that the results obtained by the inversion method are consistent with the downloaded ocean current data, which verifies the feasibility and accuracy of the method.

[0253] The above calculation examples of the present invention are only used to explain the calculation model and calculation process of the present invention in detail, and are not intended to limit the implementation methods of the present invention. For ordinary technicians in the relevant field, other different forms of changes or modifications can be made based on the above description. It is impossible to list all the implementation methods here. All obvious changes or modifications derived from the technical solution of the present invention are still within the scope of protection of the present invention.

Claims

1. A method for joint inversion of vector flow field using space-time wave currents by high-frequency ground wave radar, characterized in that: The following steps are involved: Step 1: Use the matching peak method to detect the frequency f corresponding to the first-order spectrum peak of the actual ocean echo from the received radar echo spectrum d ; Step 2: Based on f d and the theoretical Bragg frequency f B The difference Δf gives the radial velocity V of the ocean current. r ; Step 3: Invert vector velocity using spatial correlation method: Estimate the vector flow velocity of the current beam based on the radial velocity of the current beam and the adjacent left and right beams: Among them, u and v represent the vector velocity components of the ocean current in the x direction and y direction respectively, Vr i-1 、Vr i 、Vr i+1 Represents the radial velocities of the left beam, current beam, and right beam on the current range unit, respectively. The positive direction is toward the radar, θ i-1 ,θ i ,θ i+1 Respectively represent the angles of the left beam, current beam, and right beam starting from the positive x direction in the two-dimensional Cartesian coordinate system; Convert formula (3) into matrix form and get: Vr=G·X (8) According to the least squares theory, when the beam angle θ is known, i and radial velocity Vr i Under the condition of , we can calculate the matrix X corresponding to u and v, and get the velocity V of ocean current according to u and v. f and flow to α; Step 4: Flow velocity V based on ocean current movement f and flow direction α, to determine the extreme position of the ocean element where the ocean current changes significantly: Assume that the current moment is n and the next moment is n+1. According to the movement law of the ocean current, calculate the distance gate I corresponding to the next position of the ocean current of each sea element. n+1 , beam J n+1 and the time T when it flows to the next position n+1 ; For the ocean current within the radar range, calculate the velocity difference ΔV between it and the corresponding ocean current at the next position f and the angle difference Δα, and compare it with the pre-set speed difference threshold β and angle difference threshold γ; if there is ΔV f >β or Δα>γ, or T n+1 For ocean currents that are not within the time series of radar observation data, proceed to step 5; Step 5: Use the space-time wave-current joint method to invert the vector flow field. The specific process includes the following steps: S51. According to the continuity equation and motion equation in the control equation group of the ocean parameter dynamic model, the discretized forms of u and v are obtained: Where ζ is the wave height, subscripts j and k are the unit numbers in the y and x directions in the two-dimensional Cartesian coordinate system, n represents the time, Δt represents the time step, Δx and Δy represent the grid intervals in the x and y directions in the two-dimensional Cartesian coordinate system, respectively; Ω is the value of the earth's angular velocity, φ is the local latitude; g represents the acceleration of gravity; S52. There is a vector relationship between the tangential velocity Vs and u and v: -Vs=usin(π-θ)+(-v)cos(π-θ) (23) Vs=-usinθ-vcosθ (24) Based on the form of formula (24), we can obtain the following derivative with respect to time: Substituting formula (22) into formula (25) yields: The relationship between the ocean current tangential velocity and wave height is obtained from formula (26): in, and Represent the tangential velocity of the ocean current on the (k, j) unit at the current moment and the previous moment, respectively. Represents the wave height on the (k,j) unit at the previous moment, and They represent the velocity components of the ocean current velocity in the x direction and the y direction in the cell (k, j) at the previous moment, and θ is the angle of the radar beam from the x direction to the y direction in the rectangular coordinate system; Taking the movement of ocean currents in the x direction as an example, according to the CFL condition: Among them, σ is called the Courant number; Get the value range of the time step Δt: S53, correct the ocean current at the out-of-tolerance position determined in step 4, and calculate the tangential velocity at the current moment using formula (27) according to and Get vector flow velocity and flow direction in, Represents the radial velocity of the ocean current in the cell (k, j) at time n.

2. The high-frequency ground wave radar space-time wave current joint inversion vector flow field method according to claim 1 is characterized in that: Step 2 is based on f d and the theoretical Bragg frequency f B The difference Δf gives the radial velocity V of the ocean current. r as follows: Among them, λ is the wavelength of the electromagnetic wave emitted by the radar, V r is the radial velocity of the ocean current, which is positive when the current moves toward the radar and negative when it moves away from the radar.

3. The high-frequency ground wave radar space-time wave current joint inversion vector flow field method according to claim 2 is characterized in that: The theoretical Bragg frequency f B as follows: Where g is the acceleration due to gravity.

4. The high-frequency ground wave radar space-time wave current joint inversion vector flow field method according to claim 3 is characterized in that: The velocity of the ocean current obtained from u and v in step 3 5. The high-frequency ground wave radar space-time wave current joint inversion vector flow field method according to claim 4 is characterized in that: The flow direction α=π-arctan(v / u) obtained in step 3 based on u and v.

6. The method for joint inversion of vector flow field by high frequency ground wave radar space-time wave current according to any one of claims 1 to 5, characterized in that: According to the continuity equation and motion equation in the control equation group of the ocean parameter dynamic model, the process of obtaining the discretized form of u and v includes the following steps: The continuity equation is: Wherein, x, y, z represent the three axes of the three-dimensional Cartesian coordinate system; w is the velocity component of the fluid in the vertical direction, and ρ is the density of seawater; The equation of motion is: Among them, f x 、f y 、f z is the Coriolis force component, A μ is the fluid viscosity coefficient, Γ is the tidal potential of the celestial body, g represents the gravitational acceleration, and P represents the seawater pressure; Then, based on the continuity equation and the motion equation, the partial differentials of u and v with respect to time t are obtained: After gridding the data, formula (20) is transformed into a discretized form:

7. The high-frequency ground wave radar space-time wave current joint inversion vector flow field method according to claim 6 is characterized in that: The process of obtaining the partial differentials of u and v with respect to time t based on the continuity equation and the motion equation includes the following steps: For the equation of motion, ignoring the Coriolis force, tidal forces of celestial bodies, and viscosity terms, when the fluid is stationary, the upward pressure gradient force is balanced by the downward force, and we get: The hydrostatic pressure of the ocean includes the pressure caused by sea level changes, the pressure caused by ocean density stratification, and atmospheric pressure; the vertical pressure is expressed as follows: P=ρ0gζ (17) Among them, ρ0 is the average density, ζ is the wave height; Combining equation (17), the horizontal component of unit mass is expressed as follows: Then we get:

8. A system for joint inversion of vector flow field by high-frequency ground wave radar space-time wave current, characterized in that: include: The radial velocity calculation unit of the ocean current: using the matching peak method, the frequency f corresponding to the first-order spectrum peak of the actual ocean echo is detected from the received radar echo spectrum d ; Based on f d and the theoretical Bragg frequency f B The difference Δf gives the radial velocity V of the ocean current. r ; Spatial correlation method inversion unit: The spatial correlation method is used to invert the vector velocity. The specific inversion process includes the following steps: Estimate the vector flow velocity of the current beam based on the radial velocity of the current beam and the adjacent left and right beams: Among them, u and v represent the vector velocity components of the ocean current in the x direction and y direction respectively, Vr i-1 、Vr i 、Vr i+1 Represents the radial velocities of the left beam, current beam, and right beam on the current range unit, respectively. The positive direction is toward the radar, θ i-1 ,θ i ,θ i+1 Respectively represent the angles of the left beam, current beam, and right beam starting from the positive x direction in the two-dimensional Cartesian coordinate system; Convert formula (3) into matrix form and get: Vr=G·X (8) According to the least squares theory, when the beam angle θ is known, i and radial velocity Vr i Under the condition of , we can calculate the matrix X corresponding to u and v, and get the velocity V of ocean current according to u and v. f and flow to α; Determination unit of the sea element of the over-error position: velocity V based on the movement of ocean currents f and flow direction α, to determine the extreme position sea element where the ocean current changes significantly. The specific process includes the following steps: Assume that the current moment is n and the next moment is n+1. According to the movement law of the ocean current, calculate the distance gate I corresponding to the next position of the ocean current of each sea element. n+1 , beam J n+1 and the time T when it flows to the next position n+1 ; For the ocean current within the radar range, calculate the velocity difference ΔV between it and the corresponding ocean current at the next position f and the angle difference Δα, and compare it with the pre-set speed difference threshold β and angle difference threshold γ; if there is ΔV f >β or Δα>γ, or T n+1 Ocean currents that are not within the time series of radar observation data; Space-time wave-current joint method inversion unit: The space-time wave-current joint method is used to invert the vector flow field. The specific process includes the following steps: S51. According to the continuity equation and motion equation in the control equation group of the ocean parameter dynamic model, the discretized forms of u and v are obtained: Where ζ is the wave height, subscripts j and k are the unit numbers in the y and x directions in the two-dimensional Cartesian coordinate system, n represents the time, Δt represents the time step, Δx and Δy represent the grid intervals in the x and y directions in the two-dimensional Cartesian coordinate system, respectively; Ω is the value of the earth's angular velocity, φ is the local latitude; g represents the acceleration of gravity; S52. There is a vector relationship between the tangential velocity Vs and u and v: -Vs=usin(π-θ)+(-v)cos(π-θ) (23) Vs=-usinθ-vcosθ (24) Based on the form of formula (24), we can obtain the following derivative with respect to time: Substituting formula (22) into formula (25) yields: The relationship between the ocean current tangential velocity and wave height is obtained from formula (26): in, and Represent the tangential velocity of the ocean current on the (k, j) unit at the current moment and the previous moment, respectively. Represents the wave height on the (k,j) unit at the previous moment, and They represent the velocity components of the ocean current velocity in the x direction and the y direction in the cell (k, j) at the previous moment, and θ is the angle of the radar beam from the x direction to the y direction in the rectangular coordinate system; Taking the movement of ocean currents in the x direction as an example, according to the CFL condition: Among them, σ is called the Courant number; Get the value range of the time step Δt: S53, correct the ocean current at the out-of-tolerance position determined in step 4, and calculate the tangential velocity at the current moment using formula (27) according to and Get vector flow velocity and flow direction in, Represents the radial velocity of the ocean current in the cell (k, j) at time n.

9. A computer storage medium, characterized in that The storage medium stores at least one instruction, and the at least one instruction is loaded and executed by the processor to implement the high-frequency ground wave radar space-time wave and current joint inversion vector flow field method as described in any one of claims 1 to 7.

10. A high-frequency ground wave radar space-time wave current joint inversion vector flow field device, characterized in that: The device includes a processor and a memory, wherein the memory stores at least one instruction, and the at least one instruction is loaded and executed by the processor to implement the high-frequency ground wave radar space-time wave and current joint inversion vector flow field method as described in any one of claims 1 to 7.

Images