A deep water MBES seabed echo signal simulation method and system

Abstract

The application provides a deepwater MBES seabed echo signal simulation method and system, the method comprises the following steps: obtaining a local seabed terrain grid of each ping based on a platform motion model; simulating or generating sensor data and a sound velocity profile based on historical data; obtaining the coordinates and amplitude of the scattering points of the nth ping through beam forming; calculating the echo signals of the transmitted signals reflected by the scattering points to each element of the receiving array, and adding all the calculation results of the scattering points to obtain the original signals of all the receiving elements of the ping; and calculating the element original signals of all the pings according to the above two steps. The advantage of the application is that the influence factors of the echo signals considered are more comprehensive than those of the prior art, and the platform motion model effectively simplifies the dynamic process of the platform, can effectively avoid redundant calculation, improve the calculation efficiency of the receiving element original signals, and is beneficial to the verification of deepwater multi-beam sounding signal processing algorithms and the development of equipment.

Description

Technical Field

[0001] This application belongs to the fields of marine acoustic seabed topography detection technology, sonar echo signal simulation, and sonar manufacturing, specifically relating to a deep-sea MBES seabed echo signal simulation method and system. Background Technology

[0002] When a multibeam bathysound sonar operates, it forms a transmit beam strip with a narrow beam opening angle along the keel direction and a wider beam opening angle perpendicular to the keel direction. Through signal processing algorithms, it forms several receive beam strips that are wider along the keel direction and narrower perpendicular to the keel direction. By calculating the arrival time of the echoes from the overlapping areas of the transmit and receive beam strips (i.e., beam footprints), and combining this with the speed of sound and the direction of the beam footprints, the depth value at the center of each area can be calculated. Multiple bathysound results obtained from a single transmission form a bathysound strip perpendicular to the track direction. Continuously transmitting strips along the track direction enables full-coverage seabed exploration.

[0003] Simulation of seabed echo signals from multibeam bathymetry sonar plays a crucial role in the manufacturing of multibeam bathymetry systems and the research and verification of signal processing algorithms. Comprehensive signal simulation methods can accelerate algorithm research and verification, shortening the sonar development cycle. In 2012, Yoann Ladroit et al. proposed the extended scatterer model. In 2015, Zhang Wei et al. proposed an amplitude-weighted model based on the extended scatterer model, improving computation time. However, none of these models considered the impact of factors such as sound velocity field and real-time platform motion on the echo signal simulation. Summary of the Invention

[0004] The purpose of this application is to overcome the shortcomings of existing simulation methods for deep-sea MBES (Multibeam Echo Sounder) seabed echo signals under platform motion conditions, which are too simple in model and do not take into account the sound velocity field and the real-time motion of the platform.

[0005] To achieve the above objectives, this application proposes a deep-water MBES seabed echo signal simulation method, comprising:

[0006] Step 1: Obtain the local seabed topography mesh for each ping based on the platform motion model;

[0007] Step 2: Simulate or generate sensor data and sound velocity profiles based on historical data;

[0008] Step 3: First, calculate the spatial position of the transmitting array based on the sensor data at the nth ping transmission time. Second, calculate the propagation time from the local seabed topographic grid points to each element of the transmitting array through ray tracing. Then, calculate the time delay of each element of the transmitting array based on three-dimensional beam stabilization technology. Finally, obtain the coordinates and amplitude of the grid points within each transmitting beam footprint in the nth ping through beamforming. The grid points within each transmitting beam footprint are the scattering points.

[0009] Step 4: Based on the amplitude weighted model, the transmitting array is equivalent to a single transmitting point source. Based on the platform motion model, and considering the Doppler effect caused by the platform motion, the three-dimensional scattering model of seabed reverberation, the amplitude of the scattering point, and the propagation loss, the original signals of the transmitted signal reflected by the scattering point to each element of the receiving array are calculated. The calculation results of all grid points are added together to obtain the original signals of all receiving elements of the ping.

[0010] Step 5: Calculate all the primitive original signals of ping based on Step 3 and Step 4.

[0011] As an improvement to the above method, the platform motion model is as follows:

[0012] Platform motion is discretized by ping, with 1 ping representing 1 motion module;

[0013] The platform position and attitude are set to remain unchanged within a certain time period. The platform motion within the motion module is discretized according to the time length, with one time length being one motion unit.

[0014] The center of the carrier coordinate system coincides with the center of the base array.

[0015] As an improvement to the above method, step 1 includes:

[0016] Step 1-1: Grid the area that the transmitted beam may illuminate in the carrier coordinate system;

[0017] The local seabed is uniformly gridded, and the seabed coverage width on one side is determined based on the plane coordinates of the center of the launch array [0,0] and the vertical track direction. W Calculate the plane coordinates of local seabed topography grid points in the carrier coordinate system based on the width of seabed cover on one side along the track direction and the grid step size;

[0018] Step 1-2: Calculate the platform position at each ping transmission time;

[0019] Simulate and generate seabed topography data or import historical seabed topography data, set the start and end positions of the survey line, platform movement speed and launch interval, and calculate each launch time. Platform planar positions [x1,y1],[x2,y2],...,[x n ,yn ],...,[x N ,y N ];

[0020] Steps 1-3: Based on coordinate transformation and interpolation fitting techniques, calculate the vertical coordinates corresponding to the planar coordinates of the local seabed topography grid points in the nth ping carrier coordinate system, according to the seabed topography data.

[0021]

[0022] f(x j ,y j )=z j j = 1, 2, ..., K

[0023]

[0024] Where α is the heading value, and rightward rotation is positive; x represents the coordinates of the i-th grid point in the local seabed topography within the carrier coordinate system; j ,y j ,z j L represents the coordinates of the j-th seabed point in the historical topographic data X, Y, Z; W L represents the number of grid points in the vertical track direction. L K represents the number of grid points along the course of the flight path; K represents the number of seabed points. for Transform the x and y coordinates to the historical terrain data coordinate system.

[0025] As an improvement to the above method, step 2 includes:

[0026] Sensor data includes attitude data, heave data, surface sound velocity data, and GNSS data;

[0027] Sensor data can be simulated using fixed values ​​or simple trigonometric functions, or directly using measured data.

[0028] The sound velocity profile is simulated by empirical formulas or uses actual measured sound velocity profiles.

[0029] As an improvement to the above method, step 3 includes:

[0030] Step 3-1: Calculate the transmitted beam direction vector in the carrier coordinate system

[0031]

[0032]

[0033] Where φ is the angle between the projection of the transmitted beam direction vector onto the XZ plane and the negative Z-axis; θ tr The angle between the projection of the transmitted beam direction vector onto the YZ plane and the negative Z-axis;

[0034] Step 3-2: Calculate the time delay of each element of the transmitting array based on three-dimensional beam stabilization technology;

[0035] Set the center coordinates of the array as a o =[0,0,0] T Calculate the transmitted beam direction vector in the carrier coordinate system. Angle with the array axis

[0036]

[0037]

[0038]

[0039]

[0040]

[0041]

[0042]

[0043]

[0044] in, These are the direction vectors along the x, y, and z axes of the matrix in the carrier coordinate system. Γ is the bow roll angle, positive when rotating to the right; p is the pitch angle, positive when the bow is raised; r is the roll angle, positive when the port side is raised; Γ is the coordinate rotation matrix.

[0045] Calculate the time delay of each element in the transmitting array:

[0046]

[0047]

[0048] Where, τ m c is the time delay of the m-th element of the transmission array; s Indicates the surface velocity of sound; Let P and Q be the x, y, and z coordinates of the m-th primitive in the carrier coordinate system; tr These represent the number of primitives in the transmitting array along the keel direction and perpendicular to the keel direction, respectively; T is the vector transpose.

[0049] Step 3-3: Obtain the launch time through interpolation The attitude and heave data were used to calculate the coordinates of each element of the launch array at the launch moment.

[0050]

[0051]

[0052]

[0053]

[0054]

[0055] in, for Attitude values ​​at any given time; for The original coordinates of the m-th element of the emission array at any given time; for The coordinates of the m-th element of the transmitting array after attitude transformation, draft correction, and heave compensation; draft is the draft depth. for The rise and fall values ​​at any given moment;

[0056] The propagation time t from the grid point to each element of the transmitting array is calculated using ray tracing. im i = 1, 2, ..., L W ×L L m=1,2,...,P×Q tr :

[0057]

[0058] Among them, t im αi represents the propagation time from the i-th grid point to the m-th array element; x represents the horizontal distance of a single layer; t represents the propagation time of a single layer; α1 represents the grazing angle at the start of this layer; α(z) represents the grazing angle at the end of this layer; z1 represents the depth at the start of this layer; z represents the depth at the end of this layer; g represents the gradient within the layer; c represents the speed of sound.

[0059] The propagation time is calculated using a layered iterative method. The grazing angle of the first layer is fixed, and the propagation time of each layer is calculated using the above formula. The propagation time t is then accumulated and added to the grid point layer. km Compare the difference between the accumulated horizontal distance and the actual horizontal distance to see if it is less than a set threshold. If it is less, the calculation ends. If it is greater, change the first-layer grazing angle and repeat the above calculation.

[0060] Steps 3-4: Calculate the coordinates and amplitude of the scattering point using beamforming, and select the coordinates and amplitude of grid points with amplitudes greater than a set threshold as the coordinates and amplitude of the scattering point.

[0061]

[0062]

[0063] in, w is the magnitude of the k-th grid point; m f represents the weight of the m-th element of the transmitting array; c The signal center frequency is denoted by j; j is the complex unit. Amplitude Converted to decibel values.

[0064] As an improvement to the above method, step 4 includes:

[0065] Step 4-1: Calculate the shortest propagation time t for transmission and reception. min , in t min arrive Within the range, calculate the flight time for all receiving locations of the ping based on the pulse width τ. M represents the number of receiving positions; the coordinates of the receiving array elements at each receiving position in the nth ping are calculated by combining draft, heave, attitude, and GNSS data.

[0066]

[0067]

[0068] Where, x n y n for Platform planar position coordinates at launch time; for Heading angle at any moment; for Constantly receive the coordinates of the center plane of the array; for Receiving array center relative The relative position at the moment of launch; for Attitude values ​​at any given time; for Receive the original coordinates of the m′-th element of the array at any given time; for The coordinates of the m′ element of the array are received at all times after attitude transformation, draft correction, and heave compensation. for The rise and fall value at time; Q re The number of receiver array elements;

[0069] Step 4-2: Treat the transmitting array as a single emission source located at the center of the array, and calculate the propagation time from the emission source to the scattering point by ray tracing. and glancing angle Propagation time from scattering point to center of receiving array and glancing angle Combined with the launch delay τ tr Calculate the total propagation time of the scattering point:

[0070]

[0071] Among them, K s Indicates the number of scattering points;

[0072] Step 4-3: [Regarding...] Sort, and use binary search based on thresholds. The lateral angle η is calculated based on the plane coordinates of the matching scattering point and the plane coordinates of the receiving array center; the grazing angle β is calculated based on the grazing angle. ts and lateral angle β sr Calculate the scattering coefficient f(β) ts ,β sr ,η):

[0073]

[0074] Where μ is the backscattering intensity; ν is the side scattering intensity; σ is the side scattering deviation; and ΔΩ represents the deviation of the scattered sound rays from the mirror reflection direction.

[0075] Step 4-4: Based on the transmitted beam angle θ tr Platform velocity calculation: The radial velocity v of the platform relative to the transmitted beam at the moment of launch. tr According to the receiving beam angle θ re Platform velocity calculation: The radial velocity v of the platform relative to the receiving beam at the receiving moment. re According to v tr and v re Calculate the Doppler influence factor κ:

[0076]

[0077] Where c is the speed of sound;

[0078] Steps 4-5: When Where the signal pulse width is T s With a bandwidth of B, and assuming that the Doppler effect does not affect the envelope shape of the signal; calculate... The original signal B of the primitive ρ at time moment ρ (t):

[0079]

[0080] Among them, K match The number of matching scattering points; D is the normalized intensity at the k-th scattering point; rec To ensure the beam directivity of the receiving array along the flight path; S is the incident signal at the k-th scattering point; k a is the area of ​​the k-th scattering point; k Let φ be the Rayleigh distribution amplitude at the k-th scattering point; k The phase is uniformly distributed; t kρ Let be the propagation time from the k-th scattering point to the elementary element ρ;

[0081] Steps 4-6: Repeat steps 4-2 to 4-5 to obtain the original signals of all primitives for the nth ping.

[0082] This application also provides a deep-water MBES seabed echo signal simulation system, implemented based on the above method, the system comprising:

[0083] The Local Seabed Topography Mesh Module is used to obtain the local seabed topography mesh for each ping based on the platform motion model.

[0084] A module for generating sensor data and sound velocity profiles is used for simulation or to generate sensor data and sound velocity profiles based on historical data.

[0085] The coordinate and amplitude calculation module is used to first calculate the spatial position of the transmitting array based on sensor data at the nth ping transmission time, then calculate the propagation time from local seabed topographic grid points to each element of the transmitting array through ray tracking, then calculate the time delay of each element of the transmitting array based on three-dimensional beam stabilization technology, and finally obtain the coordinates and amplitudes of the grid points within each transmitting beam footprint in the nth ping through beamforming; the grid points within each transmitting beam footprint are the scattering points;

[0086] The module for calculating the original signal of the receiving element is based on the amplitude weighting model. The transmitting array is equivalent to a single transmitting point source. Based on the platform motion model, and considering the Doppler effect caused by the platform motion, the three-dimensional scattering model of seabed reverberation, the amplitude of the scattering point, and the propagation loss, the original signal of the transmitted signal reflected by the grid points to each element of the receiving array is calculated. The calculation results of all scattering points are added together to obtain the original signal of all receiving elements of the ping.

[0087] The primitive signal calculation module is used to calculate the primitive signals of all pings based on the coordinate and amplitude calculation module and the received primitive signal calculation module.

[0088] Compared with existing technologies, the advantages of this application are:

[0089] The method of this invention considers more comprehensive influencing factors of echo signals than existing methods, and the proposed platform motion model effectively simplifies the dynamic process of the platform, effectively avoiding redundant calculations and improving the calculation efficiency of the original primitive signals. In summary, this method is beneficial for the verification of deep-water multibeam bathymetry signal processing algorithms and equipment development. Attached Figure Description

[0090] Figure 1 The diagram shown is a schematic representation of the structure of the method of the present invention;

[0091] Figure 2 The image shown is a schematic diagram illustrating the calculation of local seabed topography data.

[0092] Figure 3 The diagram shown is a schematic of coordinate rotation;

[0093] Figure 4 The diagram shown is a schematic of signal transmission and reception.

[0094] Figure 5 The diagram shown is a schematic of seabed scattering. Detailed Implementation

[0095] The technical solution of this application will be described in detail below with reference to the accompanying drawings.

[0096] This invention provides a method and system for simulating deep-water MBES (Multibeam Echo Sounder) seabed echo signals. The aim is to address the shortcomings of existing simulation methods for deep-water MBES seabed echo signals under platform motion conditions, which often employ overly simplistic models and fail to consider the effects of sound velocity fields and real-time platform motion. This invention aims to provide a more comprehensive simulation method for deep-water MBES seabed echo signals. To achieve this objective, this invention proposes a platform motion model and the following technical solutions:

[0097] 1. Platform Motion Model

[0098] The platform motion model is mainly divided into three parts: ① Platform motion is discretized by ping, with 1 ping considered as 1 motion module, and the local seabed topography grid data of the motion module in the carrier coordinate system is calculated; ② Assuming that the platform position and attitude do not change within a certain time length (generally assumed to be one pulse width), the platform motion within the motion module is discretized by time length, with 1 time length considered as 1 motion unit, and the position of each motion unit is calculated, thereby calculating the total motion time of each motion unit and the unit motion time relative to the first transmission beam transmission time of that ping; ③ Assuming that the center of the carrier coordinate system and the center of the array are concurrent, the sensor data of each motion unit is obtained by interpolation based on the total motion time of each motion unit, and the positions of each element of the receiving array are calculated based on coordinate rotation and translation. Ignoring multipath effects and secondary seabed reflections, the propagation time of the transmitted signal after reflection from the scattering point to the center of the array is calculated based on the sound velocity field using ray tracking technology, and the matching unit motion time is searched. If a match is found, the signal from that scattering point to each element of the receiving array is calculated.

[0099] 2. Technical Solution

[0100] ① Simulate and generate seabed topography data or import historical seabed topography data, and obtain local seabed topography grid data for each ping based on the platform motion model; ② Generate sensor data and sound velocity profiles, including attitude data, heave data, surface sound velocity data, GNSS data, etc.; ③ First, based on the nth... The process involves: 1) Calculating the spatial location of the transmitting array using sensor data at the ping transmission time; 2) Calculating the propagation time from local seabed topographic grid points to each element of the transmitting array using ray tracking; 3) Calculating the time delay of each element of the transmitting array using three-dimensional beam stabilization technology; 4) Obtaining the coordinates and amplitude of grid points (scattering points) within the footprints of each transmitting beam in the nth ping using beamforming; 5) Calculating the echo signals of the transmitted signal reflected by the scattering points to each element of the receiving array based on the amplitude-weighted model, considering the platform motion model, the Doppler effect caused by the platform motion, the three-dimensional scattering model of seabed reverberation, the amplitude of the scattering points, and propagation loss; 6) Adding the calculation results of all scattering points to obtain the original signals of all receiving elements of the ping; 7) Repeating steps 3 and 4 to obtain the original signals of all elements of the ping.

[0101] Example 1

[0102] like Figure 1 As shown, the present invention provides a deep-water MBES seabed echo signal simulation method, comprising:

[0103] Step 1: Obtain the local seabed topography mesh for each ping based on the platform motion model.

[0104] Step 1-1: Grid the area that the transmitted beam may illuminate in the carrier coordinate system.

[0105] The local seabed is uniformly gridded, and the seabed coverage width on one side is determined based on the plane coordinates of the center of the launch array [0,0] and the vertical track direction. W , width of seabed coverage on one side along the direction of the navigation track (b) L And the grid step size d calculates the plane coordinates (X, φ) of the local seabed topography grid points in the carrier coordinate system. a ,Y a Appendix Figure 2 A schematic diagram of local seabed topographic data calculation is provided.

[0106] Step 1-2: Calculate the platform position at each ping transmission time.

[0107] Simulate and generate seabed topographic data or import historical seabed topographic data X, Y, Z, and set the start and end positions of the survey line [x s ,y s ],[x e ,y e The platform's velocity v and launch interval T are used to calculate each launch time. Platform planar positions [x1,y1],[x2,y2],...,[x n ,y n ],...,[x N ,y N ].

[0108] Steps 1-3: Based on coordinate transformation and interpolation fitting techniques, calculate the vertical coordinates corresponding to the planar coordinates of local seabed topography grid points in the coordinate system of each ping carrier, according to the seabed topography data.

[0109]

[0110]

[0111] Where α is the heading value, and rightward rotation is positive; x represents the coordinates of the i-th grid point in the local seabed topography within the carrier coordinate system; i ,y i ,z i L represents the coordinates of the i-th seabed point in the historical topographic data X, Y, Z; W L represents the number of grid points in the vertical track direction. L K represents the number of grid points along the track direction; K represents the number of seabed points.

[0112] Step 2: Simulate or generate sensor data and sound velocity profiles based on historical data.

[0113] Sensor data includes attitude data, heave data, surface sound velocity data, GNSS data, etc. Sensor data can be simulated using fixed values, simple trigonometric functions, or directly using measured data; sound velocity profiles can be simulated using empirical formulas or using actual measured sound velocity profiles.

[0114] The following steps make two assumptions: Assumption 1: The platform's position and attitude do not change within one pulse width; Assumption 2: Multipath effects and secondary seabed reflections are ignored.

[0115] Step 3: Calculate the coordinates and amplitude of the grid points (i.e., scattering points) within the footprint of each ping transmitted beam.

[0116] Step 3-1: Calculate the transmitted beam direction vector in the carrier coordinate system:

[0117]

[0118] Where φ is the angle between the projection of the transmitted beam direction vector onto the XZ plane and the negative Z-axis; θ tr The angle between the projection of the transmitted beam direction vector onto the YZ plane and the negative Z-axis.

[0119] Step 3-2: Calculate the time delay of each element of the transmitting array based on three-dimensional beam stabilization technology.

[0120] Calculate the transmitted beam direction vector in the carrier coordinate system. The angle between the center of the array and the axis of the array. Let the coordinates of the center of the array be a. o =[0,0,0] T In the carrier coordinate system, the direction vectors along the x, y, and z axes of the matrix are... The array axis affected by attitude is Appendix Figure 3 A schematic diagram of coordinate rotation is provided.

[0121]

[0122]

[0123]

[0124]

[0125]

[0126] in, Γ is the bow roll angle, positive when rotating to the right; p is the pitch angle, positive when the bow is raised; r is the roll angle, positive when the port side is raised; Γ is the coordinate rotation matrix. Calculate the time delay of each element of the transmitting array:

[0127]

[0128]

[0129] Among them, c s Indicates the surface velocity of sound; Let τ be the x, y, z coordinates of the i-th primitive in the carrier coordinate system; i For the delay of transmitting the i-th element of the array, P and Q tr These represent the number of elements in the transmitting array along the keel direction and perpendicular to the keel direction, respectively.

[0130] Step 3-3: Obtain the launch time through interpolation The attitude and heave data were used to calculate the coordinates of each element of the launch array at the launch moment. For example, time:

[0131]

[0132]

[0133] in, for Attitude values ​​at any given time; The original coordinates of the i-th element of the transmitting array; is the coordinate of the i-th element of the transmitting array after attitude transformation, draft correction, and heave compensation; draft is the draft depth; for The rise and fall value at any given moment.

[0134] The propagation time t from the grid point to each element of the transmitting array is calculated using ray tracing. ki k = 1, 2, ..., L W ×L L i = 1, 2, ..., P × Q tr :

[0135]

[0136] Among them, t ki α is the propagation time from the k-th grid point to the i-th element; α1 is the grazing angle at the beginning of this layer; α(z) is the grazing angle at the end of this layer; z1 is the depth corresponding to the beginning of this layer; z is the depth corresponding to the end of this layer; g is the gradient within the layer; x is the horizontal distance of a single layer; t is the propagation time of a single layer; c is the speed of sound.

[0137] Equation (13) is the formula for calculating the propagation time under the condition of sound speed profile. The sound ray tracking adopts layered iterative calculation. The grazing angle of the first layer (emission array element) is fixed, and the propagation time and horizontal distance of each layer are calculated by Equation (13). The propagation time and horizontal distance of the grid point layer are accumulated. The difference between the accumulated horizontal distance and the actual horizontal distance is compared with the threshold. If it is less than the threshold, the calculation ends. If it is greater than the threshold, the grazing angle of the first layer is changed and the above steps are repeated.

[0138] Steps 3-4: Calculate the coordinates and amplitudes of the grid points (i.e., scattering points) within the footprint of the transmitted beam through beamforming, and select the coordinates and amplitudes of the grid points (i.e., scattering points) with values ​​greater than -ndB as the coordinates and amplitudes of the grid points (i.e., scattering points) within the footprint of the transmitted beam.

[0139]

[0140] in, w is the amplitude at the k-th scattering point; i f represents the weight of the i-th element of the transmitting array; c The center frequency of the signal; Amplitude Converted to decibel values.

[0141] Step 4: Calculate the echo signal of the nth ping primitive j.

[0142] Appendix Figure 4 A schematic diagram of signal transmission and reception is provided.

[0143] Step 4-1: Calculate the shortest propagation time t for transmission and reception. min , in t min arrive Within the range, calculate the flight time for all receiving locations of the ping based on the pulse width τ. M represents the number of receiving positions. The coordinates of the receiving array elements at each receiving position in the nth ping are calculated by combining draft, heave, attitude, and GNSS data. For example, time:

[0144]

[0145]

[0146] Where, x n y n for Platform planar position coordinates at launch time; for Heading angle at any moment; for Constantly receive the coordinates of the center plane of the array; for Receiving array center relative The relative position at the moment of launch; for Attitude values ​​at any given time; To receive the original coordinates of the i-th element of the array; To receive the coordinates of the i-th element of the array after attitude transformation, draft correction, and heave compensation; for The rise and fall value at time; Q re The number of elements in the receiving array.

[0147] The following is based on Take the receiving position at any given time as an example.

[0148] Step 4-2: Treat the transmitting array as a single emission source located at the center of the array, and calculate the propagation time from the emission source to the scattering point by acoustic ray tracing. and glancing angle Propagation time from scattering point to center of receiving array and glancing angle K s This refers to the number of scattering points. Combined with the emission delay τ... tr Calculate the total propagation time of the scattering point:

[0149]

[0150] Step 4-3: [Regarding...] Sort, and use binary search based on thresholds. For the matching scattering point, the lateral angle η is calculated based on the plane coordinates of the scattering point and the plane coordinates of the receiving array center. (Appendix) Figure 5 A schematic diagram of seabed scattering is provided. The scattering coefficient is calculated based on the grazing angle and the lateral angle.

[0151]

[0152] Where μ is the backscattering intensity; ν is the sidescattering intensity; σ is the sidescattering deviation; and ΔΩ represents the deviation of the scattered sound rays from the direction of mirror reflection.

[0153] Step 4-4: Based on the transmitted beam angle θ tr The platform velocity v is calculated as the radial velocity v of the platform relative to the transmitted beam at the time of launch. tr According to the receiving beam angle θ re The platform velocity v is calculated as the radial velocity v of the platform relative to the receiving beam at the receiving moment. re According to v tr and v re Calculate the Doppler influence factor κ:

[0154]

[0155] Steps 4-5: When Where the signal pulse width is T s Given a bandwidth of B, and assuming the Doppler effect does not affect the envelope shape of the signal, calculate... Echo signal of time element j:

[0156]

[0157] Among them, K match The number of matching scattering points; D is the normalized intensity at the k-th scattering point; rec To ensure the beam directivity of the receiving array along the flight path; S is the incident signal at the k-th scattering point; k a is the area of ​​the k-th scattering point; k Let φ be the Rayleigh distribution amplitude at the k-th scattering point; k The phase is uniformly distributed. The propagation time t from each scattering point to the primitive is different, and different B(t) are calculated.

[0158] Steps 4-6: Repeat steps 4-2 to 4-5 to obtain the original signals of all primitives for the nth ping.

[0159] An echo signal refers to the signal form from a single scattering point to a receiving element. The original signal refers to the superposition of echo signals arriving at the same receiving element.

[0160] Step 5: Calculate all the primitive original signals of ping based on Step 3 and Step 4.

[0161] Example 2

[0162] This application also provides a deep-water MBES seabed echo signal simulation system, implemented based on the above method, the system comprising:

[0163] The Local Seabed Topography Mesh Module is used to obtain the local seabed topography mesh for each ping based on the platform motion model.

[0164] A module for generating sensor data and sound velocity profiles is used for simulation or to generate sensor data and sound velocity profiles based on historical data.

[0165] The coordinate and amplitude calculation module is used to first calculate the spatial position of the transmitting array based on the sensor data at the nth ping transmission time, then calculate the propagation time from the local seabed topographic grid points to each element of the transmitting array through ray tracking, then calculate the time delay of each element of the transmitting array based on three-dimensional beam stabilization technology, and finally obtain the coordinates and amplitude of the grid points (i.e. scattering points) within the footprints of each transmitting beam in the nth ping through beamforming.

[0166] The module for calculating the original signal of the receiving element is based on the amplitude weighting model. The transmitting array is equivalent to a single transmitting point source. Based on the platform motion model, and considering the Doppler effect caused by the platform motion, the three-dimensional scattering model of seabed reverberation, the amplitude of the scattering point, and the propagation loss, the module calculates the echo signal of the transmitted signal reflected by the scattering point to each element of the receiving array. The calculation results of all scattering points are added together to obtain the original signal of all receiving elements of the ping.

[0167] The primitive signal calculation module is used to calculate the primitive signals of all pings based on the coordinate and amplitude calculation module and the received primitive signal calculation module.

[0168] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of this application and are not intended to limit it. Although this application has been described in detail with reference to the embodiments, those skilled in the art should understand that modifications or equivalent substitutions to the technical solutions of this application do not depart from the spirit and scope of the technical solutions of this application, and should all be covered within the scope of the claims of this application.

Claims

1. A method for simulating deep-sea MBES seabed echo signals, comprising: Step 1: Obtain the local seabed topography mesh for each ping based on the platform motion model; Step 2: Simulate or generate sensor data and sound velocity profiles based on historical data; Step 3: First, according to the... The system calculates the spatial location of the transmitting array using sensor data at the ping transmission time. Next, it calculates the propagation time from local seabed topographic grid points to each element of the transmitting array using ray tracking. Then, it calculates the time delay of each element based on three-dimensional beam stabilization technology. Finally, it obtains the final position using beamforming. Ping the coordinates and amplitude of the grid points within each transmitted beam footprint; the grid points within each transmitted beam footprint are the scattering points; Step 4: Based on the amplitude weighted model, the transmitting array is equivalent to a single transmitting point source. Based on the platform motion model, and considering the Doppler effect caused by the platform motion, the three-dimensional scattering model of seabed reverberation, the amplitude of the scattering point, and the propagation loss, the echo signals of the transmitted signal reflected by the scattering point to each element of the receiving array are calculated. The calculation results of all scattering points are added together to obtain the original signal of all receiving elements of the ping. Step 5: Calculate all the original primitive signals of ping based on steps 3 and 4; The platform motion model is as follows: Platform motion is discretized by ping, with 1 ping representing 1 motion module; The platform position and attitude are set to remain unchanged within a certain time period. The platform motion within the motion module is discretized according to the time length, with one time length being one motion unit. The center of the carrier coordinate system coincides with the center of the base array.

2. The deep-water MBES seabed echo signal simulation method according to claim 1, characterized in that, Step 1 includes: Step 1-1: Grid the area that the transmitted beam may illuminate in the carrier coordinate system; The local seabed is uniformly gridded according to the plane coordinates of the center of the launch array. , width of seabed coverage on one side in the vertical track direction Calculate the plane coordinates of local seabed topography grid points in the carrier coordinate system based on the width of seabed cover on one side along the track direction and the grid step size; Step 1-2: Calculate the platform position at each ping transmission time; Simulate and generate seabed topography data or import historical seabed topography data, set the start and end positions of the survey line, platform movement speed and launch interval, and calculate each launch time. Platform Plane Position ; Steps 1-3: Based on coordinate transformation and interpolation fitting techniques, calculate the first... n Vertical coordinates corresponding to the planar coordinates of local seabed topographic grid points in the ping carrier coordinate system : in, This is the heading value; turning right is positive. For the local seabed topography in the carrier coordinate system, the first Coordinates of grid points; Historical topographic data The Coordinates of a point on the seabed; This represents the number of grid points in the vertical track direction. The number of grid points along the flight path; This represents the number of points on the seabed. , for Transformed to the coordinate system of historical terrain data x, y coordinate.

3. The deep-water MBES seabed echo signal simulation method according to claim 2, characterized in that, Step 2 includes: Sensor data includes attitude data, heave data, surface sound velocity data, and GNSS data; Sensor data can be simulated using fixed values ​​or simple trigonometric functions, or directly using measured data. The sound velocity profile is simulated by empirical formulas or uses actual measured sound velocity profiles.

4. The deep-water MBES seabed echo signal simulation method according to claim 3, characterized in that, Step 3 includes: Step 3-1: Calculate the transmitted beam direction vector in the carrier coordinate system : in, For the transmitted beam direction vector in Plane projection and negative The included angle of the axis; For the transmitted beam direction vector in Plane projection and negative The included angle of the axis; Step 3-2: Calculate the time delay of each element of the transmitting array based on three-dimensional beam stabilization technology; Set the center coordinates of the array Calculate the transmitted beam direction vector in the carrier coordinate system. Angle with the array axis : in, , , These are the direction vectors along the x, y, and z axes of the matrix in the carrier coordinate system. The bow roll angle is positive when rotated to the right. The pitch angle is positive when the bow is raised. The roll angle is positive when the port side is raised. The coordinate rotation matrix; Calculate the time delay of each element in the transmitting array: in, For the launch array The latency of each primitive; Indicates the surface velocity of sound; , , For the first The basic unit in the carrier coordinate system coordinate; and These represent the number of elements in the transmitting array along the keel direction and perpendicular to the keel direction, respectively. T Transpose of a vector; Step 3-3: Obtain the launch time through interpolation The attitude and heave data were used to calculate the coordinates of each element of the launch array at the launch moment. in, for Attitude values ​​at any given time; for Launch array at any time Primitive coordinates; for Launch array at any time The coordinates of the primitive element after attitude transformation, draft correction, and heave compensation; Draft depth; for The rise and fall values ​​at any given moment; The propagation time from the grid point to each element of the transmitting array is calculated using ray tracing. : in, For the first i From the grid point to the first m The propagation time of each array element; The horizontal distance of a single layer; The propagation time for a single layer; This is the glancing angle starting from this layer; This is the glancing angle at the end of this layer; This is the depth corresponding to the start of this layer; This is the depth corresponding to the end of this layer; Intra-layer gradient; Speed ​​of sound; The propagation time is calculated using a layered iterative method. The grazing angle of the first layer is fixed, and the propagation time of each layer is calculated using the above formula. The propagation time is then accumulated and added to the grid point layer. Compare the difference between the accumulated horizontal distance and the actual horizontal distance to see if it is less than a set threshold. If it is less, the calculation ends. If it is greater, change the first-layer grazing angle and repeat the above calculation. Steps 3-4: Calculate the coordinates and amplitude of the scattering point using beamforming, and select the coordinates and amplitude of grid points with amplitudes greater than a set threshold as the coordinates and amplitude of the scattering point. in, For the first The amplitude of each scattering point; For the launch array Primitive weights; The center frequency of the signal; j For complex units; Amplitude Converted to decibel values.

5. The deep-water MBES seabed echo signal simulation method according to claim 4, characterized in that, Step 4 includes: Step 4-1: Calculate the shortest propagation time for transmission and reception. ,exist arrive Within the range, based on pulse width Calculate the flight time for all receiving locations of this ping. ; M Indicates the number of receiving positions; calculates the number by combining draft, heave, attitude data and GNSS data. Ping the coordinates of the receiver array elements at each receiving location: in, for Platform planar position coordinates at launch time; for Heading angle at any moment; for Constantly receive the coordinates of the center plane of the array; for Receiving array center relative The relative position at the moment of launch; for Attitude values ​​at any given time; for Time receiving array Primitive coordinates; for Time receiving array The coordinates of the primitive element after attitude transformation, draft correction, and heave compensation; for The rise and fall values ​​at any given moment; The number of receiver array elements; Step 4-2: Treat the transmitting array as a single emission source located at the center of the array, and calculate the propagation time from the emission source to the scattering point by acoustic ray tracing. and glancing angle The propagation time from the scattering point to the center of the receiving array and glancing angle Combined with launch delay Calculate the total propagation time of the scattering point: in, Indicates the number of scattering points; Step 4-3: [Regarding...] Sort, and use binary search based on thresholds. The lateral angle is calculated based on the plane coordinates of the matching scattering point and the plane coordinates of the receiving array center. According to the glancing angle and lateral angle Calculate the scattering coefficient : in, This represents the backscattering intensity. This represents the intensity of lateral scattering. This is due to lateral scattering deviation; This indicates the deviation of the scattered sound rays from the direction of mirror reflection; Step 4-4: Based on the transmitted beam angle Platform velocity calculation: Radial velocity of the platform relative to the transmitted beam at the moment of launch. According to the receiving beam angle Platform velocity calculation: the radial velocity of the platform relative to the receiving beam at the receiving moment. ,according to and Calculate the Doppler influence factor : in, Speed ​​of sound; Steps 4-5: When The signal pulse width is bandwidth is Assume that the Doppler effect will not affect the envelope shape of the signal; calculate Time primitive The original signal : in, The number of matching scattering points; For the first Normalized intensity of each scattering point; To ensure the beam directivity of the receiving array along the flight path; For the first The incident signal at each scattering point; For the first Area of ​​each scattering point; For the first Rayleigh distribution amplitude at each scattering point; The phase is uniformly distributed; For the first k scattering points to the primitive The duration of transmission; Steps 4-2 to 4-5: Repeat steps 4-2 to 4-5 to obtain the result. Ping the original signals of all channel primitives.

6. A deep-water MBES seabed echo signal simulation system, implemented based on the method described in any one of claims 1-5, characterized in that, The system includes: The Local Seabed Topography Mesh Module is used to obtain the local seabed topography mesh for each ping based on the platform motion model. A module for generating sensor data and sound velocity profiles is used for simulation or to generate sensor data and sound velocity profiles based on historical data. The coordinate and magnitude calculation module is used to first calculate the coordinates and magnitude based on the first... The system calculates the spatial location of the transmitting array using sensor data at the ping transmission time. Next, it calculates the propagation time from local seabed topographic grid points to each element of the transmitting array using ray tracking. Then, it calculates the time delay of each element based on three-dimensional beam stabilization technology. Finally, it obtains the final position using beamforming. Ping the coordinates and amplitude of the grid points within each transmitted beam footprint; the grid points within each transmitted beam footprint are the scattering points; The module for calculating the original signal of the receiving element is based on the amplitude weighting model. The transmitting array is equivalent to a single transmitting point source. Based on the platform motion model, and considering the Doppler effect caused by the platform motion, the three-dimensional scattering model of seabed reverberation, the amplitude of the scattering point, and the propagation loss, the module calculates the echo signal of the transmitted signal reflected by the scattering point to each element of the receiving array. The calculation results of all scattering points are added together to obtain the original signal of all receiving elements of the ping. The primitive signal calculation module is used to calculate the primitive signals of all pings based on the coordinate and amplitude calculation module and the received primitive signal calculation module.

Images