Real-time signal processing method for three-dimensional imaging sonar based on CPU+GPU heterogeneous hardware platform

Abstract

This invention discloses a real-time signal processing method for 3D imaging sonar based on a CPU+GPU heterogeneous hardware platform, comprising the following steps: acquiring distance information from slice data for near-field and far-field discrimination, and performing compensation processing on near-field data; optimizing the near-field and far-field beamforming algorithms of the planar array based on Fast Fourier Transform (FFT) focusing; taking the modulus of the planar array beamforming results; performing sidelobe suppression processing on the beamforming results of each slice; and performing filtering interpolation processing according to the selected filtering method. This application fully leverages the synergistic advantages of CPU logic control and GPU parallel computing, optimizing and parallelizing the core step of planar array beamforming for 3D imaging sonar signals using the FFT algorithm, and specifically optimizing the data read / write method during GPU computation. This effectively reduces the computational load of the system and improves memory utilization, enabling efficient and reliable real-time 3D imaging of underwater targets.

Description

Technical Field

[0001] This invention belongs to the research field of Underwater Acoustic Signal Processing (UASP), and more specifically, this invention relates to the design and implementation of a three-dimensional imaging sonar signal processing algorithm based on a CPU+GPU heterogeneous hardware platform. Background Technology

[0002] Three-dimensional imaging sonar systems transmit pulse signals to perform acoustic imaging of underwater targets and scenes, enabling tasks such as underwater target imaging, target detection, classification and recognition, and terrain mapping. Beamforming algorithms play a crucial role in sonar signal processing, providing spatial filtering and localization imaging. In recent years, three-dimensional imaging sonar utilizing array beamforming methods has received widespread attention both domestically and internationally. However, three-dimensional imaging sonar beamforming algorithms are based on two-dimensional beamforming methods using planar arrays, requiring the computation of a large number of beam signals. Direct implementation in hardware systems demands substantial computational resources, making it difficult to meet real-time requirements. Summary of the Invention

[0003] To address the aforementioned issues, there is an urgent need to design an efficient and reliable real-time 3D imaging sonar signal processing method that can reduce the computational and memory requirements of signal processing, improve computational efficiency, and enable the 3D imaging sonar system to efficiently and reliably complete real-time 3D imaging processing of underwater targets and acquire 3D point cloud acoustic image information of underwater scenes.

[0004] This invention discloses a real-time signal processing method for three-dimensional imaging sonar based on a CPU+GPU heterogeneous hardware platform. The method includes the following steps:

[0005] Step 1: Obtain the distance information of the slice data of the point cloud image of the three-dimensional sonar, perform near-field discrimination, and perform compensation processing on the near-field data;

[0006] Step 2: Process the three-dimensional imaging sonar signal based on the Fast Fourier Transform, and replace the phase shift factor in the frequency domain beamforming with the transform factor of the Fourier Transform to complete the planar array beamforming.

[0007] Step 3: Perform mode extraction on the planar array beam to generate beam mode matrix data for multiple slices;

[0008] Step 4: Perform sidelobe suppression processing on the beam mode matrix data of each slice to generate suppressed beam domain data;

[0009] Step 5: Filter and interpolate the suppressed beam domain data to obtain three-dimensional point cloud data.

[0010] Furthermore, step 1 also includes the following steps:

[0011] Step 11: Calculate the boundary threshold based on the sonar array element aperture size and signal wavelength. Based on the distance values ​​corresponding to the input slice data, it is divided into far-field data and near-field data;

[0012] The threshold value is:

[0013] ;

[0014] Where D is the maximum aperture size of the sonar planar array. The operating wavelength of the sonar;

[0015] When the target distance corresponding to the slice data When this time, it is determined to be near-field data;

[0016] Step 12: Perform distance compensation on the near-field data, and store the compensated data continuously in a horizontal-then-vertical manner.

[0017] For the receiving array is The expression for frequency domain beamforming of a planar array is:

[0018] ;

[0019] Indicates the planar array receiving signal The frequency domain representation, Indicates the row number of the planar array. Indicates the planar array number, Represents the distance variable. Indicates the beam frequency. Indicates a horizontal angle. Indicates the vertical angle. This represents the phase shift parameter related to the angle and array index number;

[0020] Add near-field distance compensation item The frequency domain data expression of the near-field data received signal is obtained as follows:

[0021] ;

[0022] Indicates distance The delay parameter at that location, The index number represents the frequency; the remaining parameters are the same as those in the far-field condition.

[0023] Furthermore, step 2 also includes the following steps:

[0024] Step 21: Weight the compensated near-field and far-field data according to the window coefficient to generate windowed data;

[0025] Step 22: Interpolate the sparse array elements of the windowed data into full array data;

[0026] Step 23: Perform FFT calculations on the full array data in both horizontal and vertical directions to obtain two-dimensional beam domain data.

[0027] Furthermore, in step 23, under far-field conditions, the expression for frequency domain beamforming is as follows:

[0028] ;

[0029] in, Indicates azimuth angle And array element index number The relevant phase shift parameters, Indicates pitch angle And array element index number The relevant phase shift parameters;

[0030] Perform a two-dimensional Fourier transform on the frequency domain beam:

[0031] ;

[0032] express The two-dimensional Fourier transform result, Indicates the azimuth frequency index. Indicates pitch frequency index. The length of the azimuth-to-Fourier transform. The pitch-to-Fourier transform length;

[0033] Based on the expressions for frequency domain beamforming and the two-dimensional Fourier transform of frequency domain signals, we obtain... and The correspondence between them;

[0034] ;

[0035] ;

[0036] Indicates the operating wavelength of the sonar. Indicates the spacing between azimuth array elements. This indicates the pitch array element spacing.

[0037] Furthermore, in step 3, the planar array beam includes a horizontal beam and a vertical beam. The horizontal and vertical beams are decomposed using Euler's formula, and the exponential terms are expanded:

[0038] ;

[0039] ;

[0040] Multiply the two and rearrange:

[0041] ;

[0042] The real part is:

[0043] ;

[0044] The imaginary part is:

[0045] ;

[0046] Then take the modulus of the beamforming result:

[0047] .

[0048] Furthermore, step 4 also includes the following steps:

[0049] Step 41: Calculate the maximum value of all slices in each beam direction. ;

[0050] Step 42, select all wavelengths of the beam less than 0.3. The value is set to zero.

[0051] Furthermore, step 5 also includes the following steps:

[0052] Step 51: Perform maximum value filtering.

[0053] Defining the neighborhood and input / output: The beamforming result after modulus extraction is a two-dimensional matrix. The neighborhood of the maximum value filter is The filtered output matrix is Each pixel From the neighborhood The maximum value of the input pixels under the coverage is determined;

[0054] Neighborhood Coverage and Pixel Indexing: For the first pixel of the output matrix... Pixels, neighborhood The range of input pixels covered is:

[0055] ;

[0056] in, Indicates the size of the output matrix. Indicates the row index of the output matrix. Indicates the column index of the output matrix. This indicates rounding down; an index must be present. ;

[0057] Maximum value calculation: Output pixels It is the maximum value of all input pixels in the neighborhood, that is:

[0058] ;

[0059] This represents the beam result after modulus extraction. Only the maximum value is retained in this beam direction, and the rest are set to zero.

[0060] Furthermore, step 5 also includes the following steps:

[0061] Step 52: Threshold filtering is performed.

[0062] Assume the beamforming result after mode extraction is a two-dimensional matrix. Each of its elements Indicates the beam direction amplitude response ;

[0063] threshold Including proportions based on maximum value, proportions based on mean + standard deviation, and fixed thresholds:

[0064] Based on the maximum value ratio: a fixed percentage of the maximum beam amplitude is taken.

[0065] ;

[0066] in This is the proportionality coefficient. It is the global maximum value of the beamforming result;

[0067] Based on mean + standard deviation: Suitable for scenarios with relatively uniform noise distribution.

[0068] ;

[0069] in It is the global average of the beam amplitude. It is the standard deviation. For adjustment coefficients;

[0070] Fixed threshold: A constant is set directly based on system sensitivity or detection requirements. ;

[0071] Element-by-element thresholding and filtering rules: for a matrix Each element in Processed according to the "threshold comparison - amplitude replacement" rule, the output matrix is... elements Defined as:

[0072] ;

[0073] The "threshold comparison - amplitude replacement" rule is "low-pass threshold filtering". If it is necessary to retain signals below the threshold, the inequality is reversed.

[0074] Furthermore, step 5 also includes the following steps:

[0075] Step 53: Perform beam domain interpolation on the maximum value filtering or threshold filtering results. While performing maximum value filtering in each beam direction, use parabolic interpolation to interpolate the beam data, with one-point parabolic interpolation between adjacent data points; beam densification is achieved using parabolic interpolation.

[0076]

[0077] Represents the ordinate of the parabola. Represents the x-coordinate of the parabola. Represents the quadratic parameter. Represents the parameter of the first-order term. Represents a constant term;

[0078] After beam interpolation of the 108×108 maximum values, three-dimensional point cloud data is output. Finally, the three-dimensional acoustic image information of the underwater scene can be obtained through color mapping.

[0079] The beneficial effects achieved by this invention are:

[0080] The signal processing method for three-dimensional imaging sonar planar array based on CPU+GPU heterogeneous hardware platform provided in this application effectively reduces the computational load and memory requirements of the system, and can effectively and reliably complete the real-time three-dimensional imaging task of underwater targets. Attached Figure Description

[0081] Appendix Figure 1 This is the main implementation flow of the real-time signal processing method for three-dimensional imaging sonar based on a CPU+GPU heterogeneous hardware platform proposed in this invention.

[0082] Appendix Figure 2 This is a detailed data processing flow for the parallel implementation of the algorithm based on a CPU+GPU heterogeneous hardware platform proposed in this invention.

[0083] Appendix Figure 3 This invention proposes a horizontal FFT beamforming access method.

[0084] Appendix Figure 4 This invention proposes a vertical FFT beamforming access method.

[0085] Appendix Figure 5 This invention proposes a method for storing and retrieving beam data modulation results.

[0086] Appendix Figure 6 This refers to the experimental test scenario and target in the embodiments of the present invention;

[0087] Appendix Figure 7 These are the experimental results from the embodiments of the present invention: maximum value filtering and threshold filtering imaging results of floating bridges and stone slab targets in the harbor. Detailed Implementation

[0088] The present invention will be further described below with reference to specific embodiments, and the advantages and features of the present invention will become clearer as a result. However, these embodiments are merely exemplary and do not constitute any limitation on the scope of the present invention. Those skilled in the art should understand that modifications or substitutions can be made to the details and form of the technical solutions of the present invention without departing from the spirit and scope of the present invention, but all such modifications and substitutions fall within the protection scope of the present invention.

[0089] Example conditions:

[0090] (1) Frequency domain data (slice data and distance information) was obtained by testing the floating bridge and underwater rock targets on the lake, and the basic parameters of the sonar used (sonar array aperture size, signal wavelength, etc.) were obtained. The test scenario and targets are attached. Figure 6 .

[0091] (2) The transmitter sends a 300kHz single-frequency signal for 100µs. Each frame of the image is divided into 2048 slices. There are 8 sparse array signal acquisition boards, each with 64 acquisition channels. The sampling rate of each channel is 2MSPS and the sampling precision is 12bit. In the FPGA, the 12bit precision signal is padded to 16bit. Each channel needs to acquire at least 131072 points of 12bit data. After quadrature demodulation and CIC filtering with a 64x downsampling, 2048 32bit frequency domain complex signals are obtained, representing the amplitude and phase of the current signal.

[0092] (3) The data transmission size of each frame of the sparse array is 32Mb (32 bits × 512 channels × 2048 sections). The transmission rate of Gigabit Ethernet is 1000 Mbps.

[0093] (4) The host computer hardware architecture is CPU (I5-8500, Intel) + GPU (GTX-1060, NVIDIA), with 16GB of video memory.

[0094] like Figure 1 As shown, this embodiment provides a real-time three-dimensional imaging sonar signal processing method based on a CPU+GPU architecture, including:

[0095] Step 1: Obtain the distance information of the slice data to distinguish between near and far fields, and perform compensation processing on the near field data.

[0096] Step 11: Calculate the boundary threshold based on the sonar array element aperture size and signal wavelength. Based on the distance values ​​corresponding to the input slice data, it is divided into far-field data and near-field data.

[0097] Among them, the slice data is the sampling data of the point cloud image of the three-dimensional sonar at a certain moment. The point cloud image of the three-dimensional sonar is composed of multibeam data at different distances, and the data at each distance is the slice data. The distance corresponding to the slice is calculated from the sampling time and the sound speed.

[0098] The core of judging near-field and far-field data lies in the relationship between the aperture size of the sonar array elements and the signal wavelength (Rayleigh distance). () serves as the dividing line.

[0099]

[0100] Where D is the maximum aperture size of the sonar planar array. This refers to the operating wavelength of the sonar.

[0101] When the target distance corresponding to the slice data When the sound wave propagates in the form of a spherical wave, the phase difference of the signal received by the array element changes with the distance, and the spatial correlation between the signals is strong.

[0102] Step 12: Perform distance compensation on the near-field data, and store the compensated near-field data and far-field data consecutively in a horizontal-then-vertical manner.

[0103] Assume the receiving array is The expression for frequency domain beamforming of a planar array is:

[0104]

[0105] Add near-field distance compensation item This yields the frequency domain data expression of the near-field data received signal.

[0106]

[0107] Step 2: Process the three-dimensional imaging sonar signal based on Fast Fourier Transform (FFT), and replace the phase shift factor in the frequency domain beamforming with the transform factor of the Fourier transform to complete the beamforming of the planar array.

[0108] Step 21, Data Windowing: The compensated near-field and far-field data are weighted according to the generated window coefficients. The data storage location of the processed windowed data remains unchanged.

[0109] The window type can be either Hamming window or Chebyshev window. The window coefficients have been pre-calculated and stored in memory, and can be directly called when needed. The characteristics of the Hamming window function are as follows:

[0110]

[0111] The characteristics of the Chebyshev window function are as follows:

[0112]

[0113] The product of the channel data and the window function is used as the weighted result of the data window coefficient.

[0114] Windowing is applied to the primitive data of 512 channels in each slice. This module uses a total of 2048 blocks, each block processes one slice of data; each block has a total of 512 threads, each thread processes one primitive data.

[0115] Step 22: Interpolate the sparse array elements of the windowed data into a full array by padding with zeros to form a 48×48 array data.

[0116] The 512-element sparse array element field data is padded with zeros to form a 48×48 array data. First, a full 48×48 matrix of all zeros is generated. Then, the channel data corresponding to the sparse element number is assigned a value. Elements without a sparse element number remain zero.

[0117] This module uses a total of 2048 blocks, each block processes one slice of data; each block has a total of 512 threads, each thread processes one primitive data.

[0118] Step 23: Perform FFT calculations on the full array data in both horizontal and vertical directions to obtain two-dimensional beam domain data;

[0119] Under far-field conditions, the expression for frequency domain beamforming is as follows:

[0120]

[0121] in Indicates azimuth angle And array element index number The relevant phase shift parameters, Indicates pitch angle And array element index number The relevant phase shift parameters.

[0122] A two-dimensional Fourier transform is performed on the frequency domain data of the received signal, and the result is as follows:

[0123] ;

[0124] Based on the expressions for frequency domain beamforming and the two-dimensional Fourier transform of frequency domain signals, we obtain... and The correspondence between them;

[0125]

[0126]

[0127] Therefore, the expression for two-dimensional beamforming under far-field conditions can be decomposed into the product of two one-dimensional beamforming expressions. That is, two-dimensional beamforming of a planar array can be simplified into one-dimensional linear array beamforming in the horizontal and vertical directions.

[0128] In frequency domain beamforming, the phase shift factor and the transform factor of the Fourier transform have the same expression form, so the phase shift beamforming of a planar array can be completed using the Fast Fourier Transform (FFT).

[0129] Horizontal far-field beamforming;

[0130] Perform FFT calculations on each slice of data along the horizontal direction, that is, perform FFT calculations on each row of each slice of data matrix, and then store the results column-wise, as shown in the appendix. Figure 3 To reduce the processing workload of data position adjustment and the computational workload of subsequent vertical FFT beamforming, the storage position of the data corresponding to the effective beams (the number of effective beams is calculated to be 108 based on the required field of view) is adjusted only.

[0131] This module uses a total of 2048 blocks, with 108 threads per block. Each thread adjusts only one data position at a time, and only after the first data is processed will the next data be processed, until all data has been processed.

[0132] Vertical far-field beamforming;

[0133] The horizontal FFT results for each slice are then processed using an FFT in the vertical direction. The horizontal FFT results are stored column-wise, and the output is the vertical processing result. Therefore, FFT calculations are performed on each row of the data matrix for each slice, and the data storage location remains unchanged. The module input is the horizontal FFT result, and the output is the FFT beamforming result. See the appendix for specific storage and retrieval methods. Figure 4This module uses a total of 2048 blocks, with 108 threads per block. Each thread adjusts only one data position at a time, processing the next data only after the previous one is finished, until all data has been processed. The result is a horizontal × vertical beam matrix of 2048 slices.

[0134] Step 3: Perform modulus processing on the horizontal and vertical beams, that is, take the absolute value of the complex result to obtain the real result. The module input is the two-dimensional FFT beamforming result and the number of slices, and the output is the modulus value of the beam result.

[0135] Using Euler's formula Decompose the horizontal × vertical beam matrix and expand the exponential terms:

[0136]

[0137]

[0138] Multiply the two and rearrange:

[0139]

[0140] therefore, The real part is:

[0141]

[0142] The imaginary part is:

[0143]

[0144] Then take the modulus of the beamforming result:

[0145]

[0146] This module uses 2048 blocks, with 108 threads per block. Each thread performs modulus calculation on only one data point at a time, processing the next data point only after the previous one is completed, until all data has been processed. Data with the same row number in the beam modulus matrix data of different slices is stored contiguously, while the data storage locations for all beam modulus values ​​within the same slice are not contiguous. Therefore, the data storage locations need to be adjusted, as detailed in the appendix. Figure 5 As shown.

[0147] Step 4: Perform sidelobe suppression processing on the beam data of each slice, that is, set the data that is lower than 0.3 times (-10dB) of the slice maximum value to zero.

[0148] Step 4 specifically includes:

[0149] Step 41: Calculate the maximum value of all slices in each beam direction. ;

[0150] Step 42, select all wavelengths of the beam less than 0.3. The value is set to zero.

[0151] This module uses a total of 2048 blocks, each block performs sidelobe suppression processing on the beam data of one slice; each block has a total of 512 threads, each thread processes one beam data at a time, and after the processing is completed, it will proceed to the processing of the next data until all data has been processed.

[0152] Step 5: Filter and interpolate the beam domain data. Filtering methods include maximum value filtering and threshold filtering.

[0153] Step 51: If the filtering method is maximum value filtering, the maximum value filtering selects a maximum value in each beam direction of all slices, that is, range direction filtering. The filtering result is 108×108 maximum value results. Maximum value filtering and beam interpolation are performed simultaneously with maximum value filtering in each beam direction, and beam data interpolation is performed using parabolic interpolation. One-point parabolic interpolation is performed between two adjacent data points to obtain 215×215 beam results.

[0154] The beamforming results after modulus extraction are subjected to maximum value filtering, specifically including:

[0155] Defining the neighborhood and input / output: Assume the beamforming result after modulus extraction is a two-dimensional matrix. The neighborhood of the maximum value filter is The filtered output matrix is Each pixel From the neighborhood The maximum value of the input pixels under the coverage is determined.

[0156] Neighborhood Coverage and Pixel Indexing: For the first pixel of the output matrix... Pixels, neighborhood The range of input pixels covered is:

[0157]

[0158] in This indicates rounding down (e.g., when k=3). ). Indexing needs to be ensured. (If the range is exceeded, zeros can be added or the data can be mirrored and extended).

[0159] Maximum value calculation: Output pixels It is the maximum value of all input pixels in the neighborhood, that is:

[0160]

[0161] Only the maximum value is retained in this beam direction, and the rest are set to zero.

[0162] This module uses a total of 23 blocks, each with 512 threads. Each thread performs maximum value search on the beam data of all slices in a beam direction and performs parabolic interpolation and data storage according to the processing requirements of beams in different regions.

[0163] Step 52: If the filtering method is threshold filtering. Threshold filtering filters the beam data based on the input maximum and minimum thresholds. All beam data exceeding the minimum threshold are retained, and after subtracting the minimum threshold, the difference between the maximum and minimum thresholds is normalized. Normalized results greater than 1 are assigned a value of 1. Threshold filtering and beam interpolation refer to beam interpolation based on rectangular region interpolation criteria. If the rectangular region meets the interpolation conditions, region filtering is performed first, followed by beam interpolation using 1-point parabolic interpolation. Finally, the interpolated beam data undergoes threshold filtering and normalization again, with boundary data reused between adjacent rectangular regions.

[0164] If the filtering method is threshold filtering, let the beamform result after modulus extraction be a two-dimensional matrix. (M is the number of rows, N is the number of columns), each element I(i,j) represents the amplitude response of the beam in direction (i,j). The core of threshold filtering is to perform binarization or amplitude clipping on matrix elements by setting a threshold T, specifically including:

[0165] Definition and determination of threshold T: Threshold T is the core parameter of filtering. There are three common ways to define it, and the choice should be made based on the characteristics of the actual data.

[0166] (1) Based on the maximum value ratio: Take a fixed ratio of the maximum beam amplitude:

[0167]

[0168] in This is the proportionality coefficient. It is the global maximum value of the beam result.

[0169] (2) Based on mean + standard deviation: suitable for scenarios with relatively uniform noise distribution:

[0170]

[0171] in It is the global average of the beam amplitude. It is the standard deviation. This is an adjustment factor (usually taken as 1 to 3).

[0172] (3) Fixed threshold: A constant is set directly according to the system sensitivity or detection requirements. (If the noise floor is known to be) (Then the beam with higher noise level is retained).

[0173] Element-by-element thresholding and filtering rules: For a matrix Each element in Processed according to the "threshold comparison - amplitude replacement" rule, the output matrix is... elements Defined as:

[0174]

[0175] This rule is for "low-pass threshold filtering" (retaining signals above the threshold). If it is necessary to retain signals below the threshold (such as for noise extraction), the inequality can be reversed. Saved at that time (Time setting 0).

[0176] This module uses a total of 2048 blocks, with 512 threads in each block. Each thread performs threshold filtering on the beam data of a rectangular area and performs parabolic interpolation and data storage according to the processing requirements of different areas.

[0177] Step 53: Perform beam domain interpolation on the maximum value filtering or threshold filtering results.

[0178] The directional function is parabolic in shape near the main lobe, and beam densification can be achieved by using parabolic interpolation.

[0179]

[0180] After beam interpolation of the 108×108 maximum values, three-dimensional point cloud data is output. Finally, the three-dimensional acoustic image information of the underwater scene can be obtained through color mapping.

[0181] This method has been validated through lake trials, and the results are attached. Figure 7 During the experiment, the researchers unanimously agreed that:

[0182] The signal processing method for three-dimensional imaging sonar planar array based on CPU+GPU heterogeneous hardware platform provided in this application effectively reduces the computational load and memory requirements of the system, and can effectively and reliably complete the real-time three-dimensional imaging task of underwater targets.

[0183] The above are merely preferred embodiments of the present invention and do not constitute any limitation on the scope of protection of the present invention; all technical solutions formed by equivalent transformations or equivalent substitutions fall within the scope of protection of the present invention; the parts of the present invention not described in detail are well known to those skilled in the art.

Claims

1. A real-time signal processing method for three-dimensional imaging sonar based on a CPU+GPU heterogeneous hardware platform, characterized in that, The real-time signal processing method for three-dimensional imaging sonar based on a CPU+GPU heterogeneous hardware platform includes the following steps: Step 1: Obtain the distance information of the slice data of the point cloud image of the three-dimensional sonar, perform near-field discrimination, and perform compensation processing on the near-field data; Step 2: Process the three-dimensional imaging sonar signal based on the Fast Fourier Transform, and replace the phase shift factor in the frequency domain beamforming with the transform factor of the Fourier Transform to complete the planar array beamforming. Step 3: Perform mode extraction on the planar array beam to generate beam mode matrix data for multiple slices; Step 4: Perform sidelobe suppression processing on the beam mode matrix data of each slice to generate suppressed beam domain data; Step 5: Filter and interpolate the suppressed beam domain data to obtain three-dimensional point cloud data.

2. The real-time signal processing method for three-dimensional imaging sonar based on a CPU+GPU heterogeneous hardware platform according to claim 1, characterized in that, Step 1 also includes the following steps: Step 11: Calculate the boundary threshold based on the sonar array element aperture size and signal wavelength. Based on the distance values ​​corresponding to the input slice data, it is divided into far-field data and near-field data; The threshold value is: ； Where D is the maximum aperture size of the sonar planar array. The operating wavelength of the sonar; When the target distance corresponding to the slice data When this time, it is determined to be near-field data; Step 12: Perform distance compensation on the near-field data, and store the compensated data continuously in a horizontal-then-vertical manner. For the receiving array is The expression for frequency domain beamforming of a planar array is: ； Indicates the planar array receiving signal The frequency domain representation, Indicates the row number of the planar array. Indicates the planar array number, Represents the distance variable. Indicates the beam frequency. Indicates a horizontal angle. Indicates the vertical angle. This represents the phase shift parameter related to the angle and array index number; Add near-field distance compensation item The frequency domain data expression of the near-field data received signal is obtained as follows: ； Indicates distance The delay parameter at that location, The index number represents the frequency; the remaining parameters are the same as those in the far-field condition.

3. The real-time signal processing method for three-dimensional imaging sonar based on a CPU+GPU heterogeneous hardware platform according to claim 1, characterized in that, Step 2 also includes the following steps: Step 21: Weight the compensated near-field and far-field data according to the window coefficient to generate windowed data; Step 22: Interpolate the sparse array elements of the windowed data into full array data; Step 23: Perform FFT calculations on the full array data in both horizontal and vertical directions to obtain two-dimensional beam domain data.

4. The real-time signal processing method for three-dimensional imaging sonar based on a CPU+GPU heterogeneous hardware platform according to claim 3, characterized in that, In step 23, under far-field conditions, the expression for frequency domain beamforming is as follows: ； in, Indicates azimuth angle And array element index number The relevant phase shift parameters, Indicates pitch angle And array element index number The relevant phase shift parameters; Perform a two-dimensional Fourier transform on the frequency domain beam: ； express The two-dimensional Fourier transform result, Indicates the azimuth frequency index. Indicates pitch frequency index. The length of the azimuth-to-Fourier transform. The pitch-to-Fourier transform length; Based on the expressions for frequency domain beamforming and the two-dimensional Fourier transform of frequency domain signals, we obtain... and The correspondence between them; ； ； Indicates the operating wavelength of the sonar. Indicates the spacing between azimuth array elements. This indicates the pitch array element spacing.

5. The real-time signal processing method for three-dimensional imaging sonar based on a CPU+GPU heterogeneous hardware platform according to claim 1, characterized in that, In step 3, the planar array beam includes a horizontal beam and a vertical beam. Euler's formula is used to decompose the horizontal and vertical beams, expanding the exponential terms: ； ； Multiply the two and rearrange: ； The real part is: ； The imaginary part is: ； Then take the modulus of the beamforming result: 。 6. The real-time signal processing method for three-dimensional imaging sonar based on a CPU+GPU heterogeneous hardware platform according to claim 1, characterized in that, Step 4 also includes the following steps: Step 41: Calculate the maximum value of all slices in each beam direction. ; Step 42, select all wavelengths of the beam less than 0.

3. The value is set to zero.

7. The real-time signal processing method for three-dimensional imaging sonar based on a CPU+GPU heterogeneous hardware platform according to claim 1, characterized in that, Step 5 also includes the following steps: Step 51: Perform maximum value filtering. Defining the neighborhood and input / output: The beamforming result after modulus extraction is a two-dimensional matrix. The neighborhood of the maximum value filter is The filtered output matrix is Each pixel From the neighborhood The maximum value of the input pixels under the coverage is determined; Neighborhood Coverage and Pixel Indexing: For the first pixel of the output matrix... Pixels, neighborhood The range of input pixels covered is: ； in, Indicates the size of the output matrix. Indicates the row index of the output matrix. Indicates the column index of the output matrix. This indicates rounding down; an index must be present. ; Maximum value calculation: Output pixels It is the maximum value of all input pixels in the neighborhood, that is: ； This represents the beam result after modulus extraction. Only the maximum value is retained in this beam direction, and the rest are set to zero.

8. The real-time signal processing method for three-dimensional imaging sonar based on a CPU+GPU heterogeneous hardware platform according to claim 1, characterized in that, Step 5 also includes the following steps: Step 52: Threshold filtering is performed. Assume the beamforming result after mode extraction is a two-dimensional matrix. Each of its elements Indicates the beam direction amplitude response ; threshold Including proportions based on maximum value, proportions based on mean + standard deviation, and fixed thresholds: Based on the maximum value ratio: a fixed percentage of the maximum beam amplitude is taken. ； in This is the proportionality coefficient. It is the global maximum value of the beamforming result; Based on mean + standard deviation: Suitable for scenarios with relatively uniform noise distribution. ； in It is the global average of the beam amplitude. It is the standard deviation. For adjustment coefficients; Fixed threshold: A constant is set directly based on system sensitivity or detection requirements. ; Element-by-element thresholding and filtering rules: for a matrix Each element in Processed according to the "threshold comparison - amplitude replacement" rule, the output matrix is... elements Defined as: ； The "threshold comparison - amplitude replacement" rule is "low-pass threshold filtering". If you need to retain signals below the threshold, reverse the inequality.

9. The real-time signal processing method for three-dimensional imaging sonar based on a CPU+GPU heterogeneous hardware platform according to claim 1, characterized in that, Step 5 also includes the following steps: Step 53: Perform beam domain interpolation on the maximum value filtering or threshold filtering results. While performing maximum value filtering in each beam direction, use parabolic interpolation to interpolate the beam data, with one-point parabolic interpolation between adjacent data points; beam densification is achieved using parabolic interpolation. ； Represents the ordinate of the parabola. Represents the x-coordinate of the parabola. Represents the quadratic parameter. Represents the parameter of the first-order term. Represents a constant term; After beam interpolation of the 108×108 maximum values, three-dimensional point cloud data is output. Finally, the three-dimensional acoustic image information of the underwater scene can be obtained through color mapping.

Images