An underwater acoustic image reconstruction method based on angle error compensation

Abstract

The present invention relates to an underwater acoustic image reconstruction method based on angle error compensation, comprising the following steps: (1) For each subsystem composed of a receiving array element and a transmitting array element, calculate the spatial coordinate information of the transmitting array element and the receiving array element and the two-way slant range history; (2) Calculate the two-way slant range history after phase center approximation; (3) Calculate the system function in the range-Doppler domain; (4) Calculate the Doppler phase error; (5) Divide the data into non-overlapping data blocks and compensate for the micro range migration error caused by the range error; (6) Obtain the signal after range delay error compensation; (7) Perform differential range migration correction processing, range direction pulse compression, secondary range pulse compression, and consistent range migration correction processing; (8) Perform azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing; (9) Perform coherent fusion processing and inverse Fourier transform to obtain the final high-resolution image. The present invention can quickly reconstruct images and meet the requirements of real-time processing.

Description

Technical Field

[0001] The present invention relates to the technical field of intelligent manufacturing, and particularly relates to an underwater acoustic image reconstruction method based on angle error compensation. Background Art

[0002] Sound waves are the only energy carriers that can propagate over long distances underwater discovered by humans so far, and underwater acoustic transducers based on sound waves are commonly used tools for underwater detection and mapping. After applying the synthetic aperture technology to the underwater acoustic field, it becomes an underwater high-resolution imaging device. Its basic principle is to use small-aperture array elements and rely on signal processing to synthesize a virtual large-aperture array to improve the resolution in the azimuth direction. The resolution of underwater synthetic aperture detection equipment has been improved by about an order of magnitude compared with traditional side-scan equipment and has been widely used in both military and civilian fields. In the military, it is mainly used for detecting sunken or suspended explosives, matching navigation, high-resolution reconnaissance in important sea areas, etc. In the civilian field, it is mainly used for underwater topographic and geomorphic mapping, seabed exploration, archaeological search and rescue, salvage of sunken objects, channel dredging, etc.

[0003] The image reconstruction based on synthetic aperture technology is based on the acoustic path differences obtained by solving geometric relationships. When the motion error exceeds one-eighth of the wavelength, it will affect the quality of image reconstruction. Therefore, it is necessary to accurately estimate and compensate the change of the acoustic path during the motion process. In the practical application of synthetic aperture technology, due to factors such as wind and waves and the maneuverability of ships, the underwater acoustic equipment platform will deviate from the ideal track. In addition, factors such as the undulation and multi-path of the propagation medium will also affect the quality of the image, resulting in a poor focusing effect of the underwater acoustic image. Therefore, high-resolution synthetic aperture underwater acoustic equipment must take motion compensation measures to reduce or even eliminate the influence of motion errors.

[0004] Generally speaking, navigation equipment can provide accurate motion information such as yaw and pitch (such as ring laser gyroscopes), and can also provide robust acceleration information in the roll and heading directions (such as accelerometers). If the motion attitude information is not used and the method published by authors such as Ma Mengbo in the "Journal of Huazhong University of Science and Technology (Natural Science Edition)" titled "CZT Imaging Algorithm for Multi-Receiver Array Synthetic Aperture Sonar" is directly used for image reconstruction processing, the performance of the image reconstruction result will be greatly reduced. If the motion attitude information is used and the method published by authors such as Ma Mengbo in the "Journal of Huazhong University of Science and Technology (Natural Science Edition)" titled "Motion Compensation of Multi-Receiver Subarray SAS Based on Inertial Navigation System" is used for image reconstruction processing, a relatively fine image reconstruction result can be obtained. However, the computational complexity of this point-by-point image reconstruction algorithm combined with attitude compensation is very large, making this method unable to meet the system timeliness processing requirements. Summary of the Invention

[0005] The technical problem to be solved by the present invention is to provide an underwater acoustic image reconstruction method based on angle error compensation that enables fast image reconstruction and meets the requirements of real-time processing.

[0006] To solve the above problems, an underwater acoustic image reconstruction method based on angle error compensation according to the present invention includes the following steps:

[0007] ⑴ For each subsystem composed of a receiving array element and a transmitting array element, calculate the spatial coordinate information of the transmitting and receiving array elements in the presence of angle errors, and calculate the two-way slant range history corresponding to any transmitting and receiving array element subsystem to obtain the two-way slant range history between the subsystem composed of the i-th receiving array element and the transmitting array element and the point target;

[0008] ⑵ According to the phase center approximation method, calculate the two-way slant range history after phase center approximation;

[0009] ⑶ According to the two-way slant range history after phase center approximation, calculate the system function of each transmitting and receiving array element subsystem in the range-Doppler domain;

[0010] ⑷ For the Doppler phase error caused by the range delay error during the phase center approximation process, calculate the Doppler phase error caused by compensating the range error in the two-dimensional time domain;

[0011] ⑸ In the range direction, divide the data into N non-overlapping data blocks, and compensate the micro range migration error caused by the range error in the range direction frequency domain for any one data block;

[0012] ⑹ For all data blocks after micro range migration compensation, arrange them in order in the range direction time domain, that is, obtain the signal after range delay error compensation;

[0013] ⑺ For the signal of each receiving array element after range delay error compensation, perform differential range migration correction processing in the range-Doppler domain in sequence, and complete range direction pulse compression, second range pulse compression, and consistent range migration correction processing while performing filtering operation in the two-dimensional frequency domain;

[0014] ⑻ For the data of each receiving array element after range direction pulse compression, second range pulse compression, and consistent range migration correction processing, perform azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing while performing filtering operation in the range-Doppler domain;

[0015] ⑼ For the data of all receiving array elements after azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing, perform coherent fusion processing in the range-Doppler domain, and perform inverse Fourier transform in the azimuth direction to obtain the final high-resolution image.

[0016] The three-dimensional space coordinates of the transmitting array element in step (1) are:

[0017]

[0018] Where: v represents the platform movement speed, in meters per second; t represents the slow time in the azimuth direction, in seconds; x t represents the azimuth coordinate of the transmitting array element at time t, in meters; t represents the coordinate of the transmitting array element at time t in meters; z t It represents the coordinate of the transmitting array element at the time t in meters;

[0019] The three-dimensional space coordinates of the i-th receiving array element are:

[0020]

[0021] Where: i (1≤i≤I) represents the index of the receiving array element; I represents the total number of receiving array elements; r represents the closest distance between the target and the motion trajectory of the underwater acoustic transducer on the slant range plane, in meters; u i represents the distance between the i-th receiving array element and the transmitting array element, in meters; represents the rotation vector, θ pitch Indicates the pitch angle error, in degrees; θ heading Indicates the yaw angle error in degrees; x i represents the coordinate of the i-th receiving array element along the azimuth at time t, in meters; i represents the coordinate of the i-th receiving array element along the ground distance at time t, in meters; z i represents the coordinate of the i-th receiving array element along the height at time t, in meters; It represents the distance that the ith receiving element moves along the azimuth direction from the time when the transmitting element transmits the signal to the time when the ith receiving element receives the target echo signal, in meters; c represents the speed of sound in water, in meters per second;

[0022] The two-way slant range from the subsystem composed of the i-th receiving array element and the transmitting array element to the point target is:

[0023]

[0024] in: It represents the slant range of the i-th receiving array element in the case of angular error, with the unit of meter; t represents the slow time in the azimuth direction, with the unit of second; r represents the closest distance between the target and the motion trajectory of the underwater acoustic transducer in the slant range plane, with the unit of meter; h represents the height of the underwater acoustic transducer from the bottom of the water, with the unit of meter; θ represents the slant angle of the underwater acoustic transducer, with the unit of degree; the three-dimensional coordinates of the point target in the azimuth direction, ground distance direction, and height direction are 0, rsinθ, and h respectively, and the unit of all is meter; c represents the speed of sound in water, with the unit of meter per second; v represents the speed of the underwater acoustic transducer platform, with the unit of meter per second; u i It represents the distance between the i-th receiving array element and the transmitting array element, with the unit of meter; θ pitch It represents the pitch angle error, with the unit of degree; θ heading It represents the yaw error, with the unit of degree.

[0025] The calculation formula for the two-way slant range history after phase center approximation in step (2) is:

[0026]

[0027] Among them: u i It represents the distance between the i-th (1≤i≤I) receiving array element and the transmitting array element, with the unit of meter; I represents the total number of receiving array elements; r represents the closest distance between the target and the motion trajectory of the underwater acoustic transducer in the slant range plane, with the unit of meter; ΔR i It represents the distance delay error after phase center approximation, with the unit of meter, ΔR i The specific expression of

[0028]

[0029] Among them: It represents the slant range of the i-th receiving array element in the case of angular error, with the unit of meter; r represents the closest distance between the target and the motion trajectory of the underwater acoustic transducer in the slant range plane, with the unit of meter; h represents the height of the underwater acoustic transducer from the bottom of the water, with the unit of meter; θ pitch It represents the pitch angle error, with the unit of degree; θ heading It represents the yaw error, with the unit of degree.

[0030] The calculation formula for the system function of each transceiver array element subsystem in the range-Doppler domain in step (3) is:

[0031]

[0032] Among them: sS i (τ,f a ; r) represents the system function of the i-th receiving array element in the range-Doppler domain, It represents the window function of the transmitted signal, Wa (f a -f dc ) represents the combined directivity function of the receiving and transmitting array elements; j 2 = -1 represents the imaginary unit; τ represents the fast time in the range direction, with the unit of seconds; c represents the speed of sound in water, with the unit of meters per second; ΔR i represents the range delay error after phase center approximation, with the unit of meters; f a represents the azimuth Doppler frequency, with the unit of Hertz; f dc represents the Doppler center frequency caused by the squint of the underwater acoustic transducer, with the unit of Hertz; r represents the closest distance between the target and the receiving array in the slant range plane, with the unit of meters; θ i (τ, f a ; r) represents the phase of the system function, and the specific expression is:

[0033]

[0034] Where: u i represents the distance between the i-th (1 ≤ i ≤ I) receiving array element and the transmitting array element, with the unit of meters; I represents the total number of receiving array elements; f0 represents the center frequency of the transmitted chirp signal, with the unit of Hertz; v represents the moving speed of the underwater acoustic transducer platform, with the unit of meters per second; K m represents the range chirp rate considering the coupling between azimuth and range, with the unit of Hertz per second, and the specific expression is:

[0035]

[0036] Where: K r represents the chirp rate of the transmitted chirp signal, with the unit of Hertz per second; D(f a ) The physical meaning is that there is coupling between azimuth and range, which is called the range curvature factor, and the specific expression is:

[0037]

[0038] In step (4), the Doppler phase error caused by compensating the range error in the two-dimensional time domain is calculated according to the following formula:

[0039]

[0040] Where: represents the Doppler phase error, f0 represents the center frequency of the transmitted chirp signal, with the unit of Hertz; ΔR i represents the range delay error after phase center approximation, with the unit of meters; c represents the speed of sound in water, with the unit of meters per second.

[0041] In step (5), the width Δr of each data block is determined by the following formula:

[0042]

[0043] where: abs represents the absolute value operation, denotes the value at the center distance r = r n,c of the nth data block; denotes the value at the edge distance r = r n,c - 0.5Δr of the nth data block; r n,c represents the center distance of the nth data block, in meters; the expression (r n,c - 0.5Δr) represents the edge distance of the nth data block, in meters; Δr represents the width of the data block, in meters;

[0044] The total number N of data blocks is:

[0045]

[0046] where: r max represents the maximum distance that the underwater acoustic transducer can detect, in meters; r min represents the minimum distance detected by the underwater acoustic transducer, in meters; ceil represents the ceiling operation;

[0047] The compensation function of the nth data block is:

[0048]

[0049] where: Ω i (t, f r , n) represents the compensation function of the micro-range migration error in the nth data block, f r represents the range instantaneous frequency, in hertz; c represents the speed of sound in water, in meters per second; ΔR i,n represents the range delay error ΔR after phase center approximation i at the center distance r = r n,c of the nth data block, in meters, and the specific expression is:

[0050]

[0051] where: u i represents the distance between the i-th receiving array element and the transmitting array element, in meters; denotes the value of the i-th receiving array element at the center distance of the nth data block under the condition of angle error, in meters; r n,crepresents the value at the center distance of the nth data block, with the unit of meter; h represents the height of the underwater acoustic transducer from the bottom of the water, with the unit of meter; θ pitch represents the pitch angle error, with the unit of degree; θ heading represents the yaw error, with the unit of degree.

[0052] The method of arranging in sequence according to the order of data blocks in the range-time domain in step (6) is as follows:

[0053]

[0054] where: represents the data after micro-range migration compensation for the first data block of the ith receiving array element, and so on. represents the data after micro-range migration compensation for the nth data block of the ith receiving array element. The subscript n (1≤n≤N) represents the index of the data block of the ith receiving array element. represents the signal after arranging all N data blocks of the ith receiving array element in sequence.

[0055] The filtering function during the differential range migration correction processing in step (7) in the range-Doppler domain is:

[0056]

[0057] where: Γ i (τ, f a ) represents the filtering function for performing differential range migration correction processing on the ith receiving array element in the range-Doppler domain. j 2 =-1 represents the imaginary unit. K m represents the chirp rate in the case of fusing azimuth and range coupling, with the unit of Hertz / second; D(f a ) represents the range curvature factor. D(f dc ) represents the range curvature factor D(f a ) at the Doppler center f a =f dc The value at this point. c represents the speed of sound in water, with the unit of meter / second; r ref represents the reference distance, with the unit of meter; τ represents the fast time in the range direction, with the unit of second.

[0058] The filtering function when completing range-direction pulse compression, second-range pulse compression, and consistent range migration correction processing is:

[0059]

[0060] where: H i (f r , f a ; r ref) represents the filtering function for performing range - direction pulse compression, second - order range - pulse compression, and consistent range migration correction processing on the i - th receiving array element in the two - dimensional frequency domain, D(f a ) represents the range curvature factor, D(f dc ) represents the value of the range curvature factor D(f a ) at the Doppler center f a = f dc ; K m represents the chirp rate in the case of fusing azimuth and range coupling, with the unit of Hertz per second; f r represents the instantaneous frequency in the range direction, with the unit of Hertz; r ref represents the reference range, with the unit of meter; c represents the speed of sound in water, with the unit of meter per second.

[0061] The compensation function during azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing in step (8) is:

[0062]

[0063] where: Λ i (τ, f a ) represents the compensation function for performing azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing in the range - Doppler domain, D(f a ) represents the range curvature factor, D(f dc ) represents the value of the range curvature factor D(f a ) at the Doppler center f a = f dc ; f0 represents the center frequency of the transmitted linear frequency - modulated signal, with the unit of Hertz; K m represents the chirp rate in the case of fusing azimuth and range coupling, with the unit of Hertz per second; f a represents the azimuth Doppler frequency, with the unit of Hertz; r ref represents the reference range, with the unit of meter; r represents the shortest distance between the target and the movement trajectory of the underwater acoustic transducer in the slant - range plane, with the unit of meter; c represents the speed of sound in water, with the unit of meter per second.

[0064] The coherent fusion processing method in step (9) is:

[0065]

[0066] where: represents the result after coherent fusion of the data of each receiving array element in the range - Doppler domain, It represents the result after azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing on the data of the i-th receiving array element in the range-Doppler domain. The subscript i (1 ≤ i ≤ I) represents the index of the receiving array element, and I represents the total number of receiving array elements;

[0067] The final high-resolution image is obtained as follows:

[0068]

[0069] Among them: IFFT represents the inverse Fourier transform, and ff(τ,t) represents the final high-resolution image.

[0070] The present invention has the following advantages compared with the prior art:

[0071] Based on the angular error of the navigation data, the present invention compensates for the adverse effects of the deviation of the underwater acoustic equipment system from non-uniform linear motion on the image by fusing the correction of the angular error in each processing step of the fast image reconstruction method, thereby improving the performance of the image reconstruction method and the efficiency of image reconstruction at the same time, and can provide great support for the real-time processing of high-performance underwater acoustic equipment. Description of the Drawings

[0072] The following further details the specific embodiments of the present invention with reference to the drawings.

[0073] Figure 1 It is the flow chart of the image reconstruction method of the present invention. Among them: Range FFT represents the Fourier transform in the range direction; Range IFFT represents the inverse Fourier transform in the range direction; Azimuth FFT represents the Fourier transform in the azimuth direction; Azimuth IFFT represents the inverse Fourier transform in the azimuth direction.

[0074] Figure 2 It is the two-dimensional image reconstruction geometry of a multi-subarray synthetic aperture underwater acoustic equipment with angular error in the present invention.

[0075] Figure 3 It is the coordinate system in which the transmitting and receiving array elements form a linear array in the present invention.

[0076] Figure 4 It is the system approximation error under different angular error conditions in the present invention. Among them: (a) Pitch angle is 0.5°, yaw angle is 1°; (b) Pitch angle is 1°, yaw angle is 2°; (c) Pitch angle is 3°, yaw angle is 4°; (d) Pitch angle is 5°, yaw angle is 6°.

[0077] Figure 5This is the system approximation error of the underwater acoustic transducer under different beam widths when the pitch angle is 2° and the yaw angle is 3° in the present invention. Wherein: (a) The beam width is 3°; (b) The beam width is 6°; (c) The beam width is 9°; (d) The beam width is 12°.

[0078] Figure 6 This is the angle error recorded by the navigation device in the present invention. Wherein: (a) Yaw angle error; (b) Pitch angle error.

[0079] Figure 7 This is the processing result of the measured data at sea in the present invention. Wherein: (a) Angle error compensation is not performed; (b) Only pitch error is compensated; (c) Only yaw error is compensated; (d) Both yaw and pitch errors are compensated. Detailed implementation manners

[0080] As Figure 1 shown, an underwater acoustic image reconstruction method based on angle error compensation includes the following steps:

[0081] ⑴ For each subsystem composed of a receiving element and a transmitting element, calculate the spatial coordinate information of the transmitting element and the receiving element in the presence of angle errors, and calculate the two-way slant range history corresponding to any transmitting-receiving element subsystem to obtain the two-way slant range history between the subsystem composed of the i-th receiving element and the transmitting element and the point target.

[0082] Wherein: The three-dimensional spatial coordinates of the transmitting element are:

[0083]

[0084] Wherein: v represents the platform movement speed, in meters per second; t represents the slow time in the azimuth direction, in seconds; x t represents the coordinate of the transmitting element in the azimuth direction at time t, in meters; y t represents the coordinate of the transmitting element in the ground range at time t, in meters; z t represents the coordinate of the transmitting element in the height at time t, in meters.

[0085] The three-dimensional spatial coordinates of the i-th receiving element are:

[0086]

[0087] Wherein: i (1≤i≤I) represents the index of the receiving element; I represents the total number of receiving elements; r represents the shortest distance between the target and the movement trajectory of the underwater acoustic transducer in the slant range plane, in meters; u i represents the distance between the i-th receiving element and the transmitting element, in meters; represents the rotation vector, θpitch Indicates the pitch angle error, in degrees; θ heading Indicates the yaw angle error in degrees; x i represents the coordinate of the i-th receiving array element along the azimuth at time t, in meters; i represents the coordinate of the i-th receiving array element along the ground distance at time t, in meters; z i represents the coordinate of the i-th receiving array element along the height at time t, in meters; It represents the distance that the ith receiving element moves along the azimuth direction from the time when the transmitting element transmits the signal to the time when the ith receiving element receives the target echo signal, in meters; c represents the speed of sound in water, in meters per second;

[0088] The two-way slant range from the subsystem composed of the i-th receiving array element and the transmitting array element to the point target is:

[0089]

[0090] in: represents the slant range of the i-th receiving array element in the case of angle error, in meters; t represents the slow time in azimuth, in seconds; r represents the closest distance between the target and the motion trajectory of the hydroacoustic transducer on the slant range plane, in meters; h represents the height of the hydroacoustic transducer from the bottom of the water, in meters; θ represents the slant angle of the hydroacoustic transducer, in degrees; the three-dimensional coordinates of the point target in azimuth, ground distance, and height are 0, rsinθ, h, all in meters; c represents the speed of sound in water, in meters / second; v represents the platform speed of the hydroacoustic transducer, in meters / second; u i represents the distance between the i-th receiving array element and the transmitting array element, in meters; θ pitch Indicates the pitch angle error, in degrees; θ heading Indicates the yaw error in degrees.

[0091] (2) According to the phase center approximation method, the two-way slant range history after phase center approximation is calculated as follows:

[0092]

[0093] Where: u i represents the distance between the i-th (1≤i≤I) receiving array element and the transmitting array element, in meters; I represents the total number of receiving array elements; r represents the closest distance between the target and the motion trajectory of the underwater acoustic transducer on the slant range plane, in meters; ΔR i Represents the distance delay error after phase center approximation, in meters, ΔR i The specific expression is:

[0094]

[0095] Wherein: represents the slant range of the i-th receiving array element in the case of angular error, with the unit of meter; r represents the closest distance between the target and the movement trajectory of the underwater acoustic transducer in the slant range plane, with the unit of meter; h represents the height of the underwater acoustic transducer from the bottom of the water, with the unit of meter; θ pitch represents the pitch angle error, with the unit of degree; θ heading represents the yaw error, with the unit of degree.

[0096] (3) Calculate the system function of each transceiver array element subsystem in the range-Doppler domain according to the two-way slant range history after phase center approximation. The calculation formula is:

[0097]

[0098] Wherein: sS i (τ, f a ; r) represents the system function of the i-th receiving array element in the range-Doppler domain, represents the window function of the transmitted signal, W a (f a - f dc ) represents the combined directivity function of the receiving and transmitting array elements; j 2 = -1 represents the imaginary unit; τ represents the fast time in the range direction, with the unit of second; c represents the speed of sound in water, with the unit of meter per second; ΔR i represents the range delay error after phase center approximation, with the unit of meter; f a represents the azimuth Doppler frequency, with the unit of hertz; f dc represents the Doppler center frequency caused by the slant of the underwater acoustic transducer, with the unit of hertz; r represents the closest distance between the target and the receiving array in the slant range plane, with the unit of meter; θ i (τ, f a ; r) represents the phase of the system function, and the specific expression is:

[0099]

[0100] Wherein: u i represents the distance between the i-th (1 ≤ i ≤ I) receiving array element and the transmitting array element, with the unit of meter; I represents the total number of receiving array elements; f0 represents the center frequency of the transmitted chirp signal, with the unit of hertz; v represents the movement speed of the underwater acoustic transducer platform, with the unit of meter per second; K m represents the range chirp rate considering the coupling of azimuth and range, with the unit of hertz per second, and the specific expression is:

[0101]

[0102] Where: K r represents the chirp rate of the transmitted chirp signal, with the unit of Hertz per second; D(f a ) has the physical meaning that there is coupling in the azimuth and range directions, which is called the range curvature factor, and the specific expression is:

[0103]

[0104] ⑷ For the Doppler phase error caused by the range delay error in the process of phase center approximation, calculate the Doppler phase error caused by compensating the range error in the two-dimensional time domain according to the following formula:

[0105]

[0106] Where: represents the Doppler phase error, f0 represents the center frequency of the transmitted chirp signal, with the unit of Hertz; ΔR i represents the range delay error after phase center approximation, with the unit of meter; c represents the speed of sound in water, with the unit of meter per second.

[0107] ⑸ In the range direction, divide the data into N non-overlapping data blocks, and compensate the micro-range migration error caused by the range error in the range frequency domain for any one data block.

[0108] The width Δr of each data block is determined by the following formula:

[0109]

[0110] Where: abs represents the absolute value operation, represents the value at the center distance r = r n,c of the nth data block; represents the value at the edge distance r = r n,c -0.5Δr of the nth data block; r n,c represents the center distance of the nth data block, with the unit of meter; the expression (r n,c -0.5Δr) represents the edge distance of the nth data block, with the unit of meter; Δr represents the width of the data block, with the unit of meter;

[0111] The total number N of data blocks is:

[0112]

[0113] Where: r max represents the maximum distance that the underwater acoustic transducer can detect, with the unit of meter; r minrepresents the minimum distance detected by the underwater acoustic transducer, with the unit of meter; ceil represents the operation of rounding up to the nearest integer;

[0114] The compensation function for the nth data block is:

[0115]

[0116] where: Ω i (t, f r , n) represents the compensation function for the micro-range migration error in the nth data block, f r represents the range instantaneous frequency, with the unit of hertz; c represents the speed of sound in water, with the unit of meter per second; ΔR i,n represents the range delay error ΔR after phase center approximation i at the center distance r = r n,c of the nth data block, with the unit of meter, and the specific expression is:

[0117]

[0118] where: u i represents the distance between the ith receiving array element and the transmitting array element, with the unit of meter; represents the value of the ith receiving array element at the center distance of the nth data block based on the angle error, with the unit of meter; r n,c represents the value at the center distance of the nth data block, with the unit of meter; h represents the height of the underwater acoustic transducer from the bottom of the water, with the unit of meter; θ pitch represents the pitch angle error, with the unit of degree; θ heading represents the yaw error, with the unit of degree.

[0119] ⑹ For all data blocks after micro-range migration compensation, arrange them in sequence in the range-time domain according to the order of the data blocks, and then the signal after range delay error compensation can be obtained.

[0120] The arrangement method is as follows:

[0121]

[0122] where: represents the data after micro-range migration compensation for the first data block of the ith receiving array element, and so on, represents the data after micro-range migration compensation for the nth data block of the ith receiving array element. The subscript n (1 ≤ n ≤ N) represents the index of the data block of the ith receiving array element, represents the signal after arranging all N data blocks of the ith receiving array element in sequence.

[0123] ⑺ For the signals after distance-delay error compensation for each receiving array element, perform differential range migration correction processing in the range-Doppler domain in sequence, and complete range-direction pulse compression, second-order range pulse compression, and uniform range migration correction processing while performing filtering operations in the two-dimensional frequency domain.

[0124] Among them: The filtering function during differential range migration correction processing in the range-Doppler domain is:

[0125]

[0126] Among them: Γ i (τ,f a ) represents the filtering function for performing differential range migration correction processing on the i-th receiving array element in the range-Doppler domain, j 2 =-1 represents the imaginary unit, K m represents the chirp rate in the case of combined azimuth and range coupling, with the unit of Hertz per second; D(f a ) represents the range curvature factor, D(f dc ) represents the value of the range curvature factor D(f a ) at the Doppler center f a =f dc ; c represents the speed of sound in water, with the unit of meters per second; r ref represents the reference distance, with the unit of meters; τ represents the fast time in the range direction, with the unit of seconds;

[0127] The filtering function for completing range-direction pulse compression, second-order range pulse compression, and uniform range migration correction processing is:

[0128]

[0129] Among them: H i (f r ,f a ; r ref ) represents the filtering function for performing range-direction pulse compression, second-order range pulse compression, and uniform range migration correction processing on the i-th receiving array element in the two-dimensional frequency domain, D(f a ) represents the range curvature factor, D(f dc ) represents the value of the range curvature factor D(f a ) at the Doppler center f a =f dc ; K m represents the chirp rate in the case of combined azimuth and range coupling, with the unit of Hertz per second; f r represents the instantaneous frequency in the range direction, with the unit of Hertz; r ref represents the reference distance, with the unit of meters; c represents the speed of sound in water, with the unit of meters per second.

[0130] ⑻ For the data after range - Doppler domain processing including range - direction pulse compression, second - order range - pulse compression, and consistent range migration correction for each receiving array element, azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing are performed simultaneously while filtering in the range - Doppler domain. The compensation function is as follows:

[0131]

[0132] Where: Λ i (τ, f a ) represents the compensation function for azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing in the range - Doppler domain. D(f a ) represents the range curvature factor. D(f dc ) represents the value of the range curvature factor D(f a ) at the Doppler center f a = f dc ; f0 represents the center frequency of the transmitted chirp signal, with the unit of Hertz; K m represents the chirp rate in the case of fusing azimuth and range coupling, with the unit of Hertz per second; f a represents the azimuth Doppler frequency, with the unit of Hertz; r ref represents the reference range, with the unit of meter; r represents the closest distance between the target and the movement trajectory of the underwater acoustic transducer in the slant - range plane, with the unit of meter; c represents the speed of sound in water, with the unit of meter per second.

[0133] ⑼ For the data after azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing for all receiving array elements, coherent fusion processing is performed in the range - Doppler domain, and the final high - resolution image is obtained by performing inverse Fourier transform in the azimuth direction.

[0134] The coherent fusion processing method is as follows:

[0135]

[0136] Where: represents the result after coherent fusion of the data of each receiving array element in the range - Doppler domain. represents the result after azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing for the data of the i - th receiving array element in the range - Doppler domain. The subscript i (1 ≤ i ≤ I) represents the index of the receiving array element, and I represents the total number of receiving array elements;

[0137] The way to obtain the final high - resolution image is as follows:

[0138]

[0139] Wherein, IFFT represents the inverse Fourier transform, and ff(τ,t) represents the final high-resolution image.

[0140] Embodiment

[0141] First, considering the angular error, the present invention calculates the spatial coordinate positions of the transmitting and receiving array elements, and further calculates the two-way slant range history between the spatial target and the transmitting and receiving array elements. Then, according to the phase center approximation method, the approximated two-way slant range history is obtained. On this basis, the system function of any transmitting and receiving array element subsystem in the range-Doppler domain is calculated. For the signals after distance delay error compensation for each receiving array element, differential range migration correction processing is sequentially performed in the range-Doppler domain, range-Doppler domain pulse compression and second-order range pulse compression are completed while performing filtering operations in the two-dimensional frequency domain, and consistent range migration correction processing is performed. Azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing are realized while performing filtering operations in the range-Doppler domain. For the data after azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing for all receiving array elements, coherent fusion processing is performed in the range-Doppler domain, and the final high-resolution image is obtained by performing the inverse Fourier transform in the azimuth direction.

[0142] As Figure 2 shown is the geometry of two-dimensional image reconstruction of a multi-subarray synthetic aperture underwater acoustic equipment in the presence of angular error. In Figure 2 , the movement direction of the underwater acoustic equipment, i.e., the x-axis, is called the azimuth direction, the y-axis is called the ground range direction, and the z-axis is called the height direction. The xyz forms a three-dimensional rectangular coordinate system. In fact, the underwater acoustic synthetic aperture equipment is a two-dimensional imaging device and can only obtain a two-dimensional image. Its imaging plane is composed of the target P in Figure 2 and the x-axis. That is to say, the underwater acoustic synthetic aperture equipment projects the three-dimensional information of the target in the xyz three-dimensional rectangular coordinate system onto the two-dimensional imaging plane. In Figure 2 , r represents the shortest distance between the target and the movement trajectory of the underwater acoustic transducer in the slant range plane. There is a linear array integrating a transmitting array element and multiple receiving array elements on the x-axis. The black rectangle in the figure represents the transmitting array element, and the remaining rectangles represent the receiving array elements. For each subsystem composed of a receiving array element and a transmitting array element, the spatial coordinate information of the transmitting array element and the receiving array element is calculated. The three-dimensional spatial coordinates of the transmitting array element are:

[0143]

[0144] Wherein: v represents the platform movement speed, with the unit of meters per second; t represents the slow time in the azimuth direction, with the unit of seconds; x t represents the coordinate of the transmitting array element in the azimuth direction at time t, with the unit of meters; y tDenotes the coordinate of the transmitting array element in the ground distance at time t, with the unit of meter; z t Denotes the coordinate of the transmitting array element in the altitude at time t, with the unit of meter.

[0145] The three-dimensional space coordinates of the i-th receiving array element are:

[0146]

[0147] Where: i (1 ≤ i ≤ I) represents the index of the receiving array element, and I represents the total number of receiving array elements; r represents the shortest distance between the target and the movement trajectory of the underwater acoustic transducer in the slant range plane, with the unit of meter; u i Represents the distance between the i-th receiving array element and the transmitting array element, with the unit of meter. Represents the rotation vector. θ pitch Represents the pitch angle error, with the unit of degree, and represents the angle error caused by rotation around the z Figure 3 axis in a axis. θ heading Represents the yaw angle error, with the unit of degree, and represents the angle error caused by rotation around the x Figure 3 axis in a axis, Figure 3 is the three-dimensional direct coordinate system of the array established with the transmitting array element in the Figure 2 linear array as the coordinate origin, where O a represents the coordinate origin, x a represents the azimuth direction, y a ground distance direction, z a represents the altitude direction. In the figure, -z a is given, representing the opposite direction of the altitude direction. x i represents the coordinate of the i-th receiving array element along the azimuth direction at time t, y i represents the coordinate of the i-th receiving array element along the ground distance at time t, z i represents the coordinate of the i-th receiving array element along the altitude at time t. Represents the distance that the i-th receiving array element moves along the azimuth direction during the period from the transmitting array element emitting a signal to the i-th receiving array element receiving the target echo signal, with the unit of meter. c represents the speed of sound in water, with the unit of m / s.

[0148] For each subsystem composed of a receiving array element and a transmitting array element, calculate the two-way slant range history corresponding to any transceiver array element subsystem, and obtain the two-way slant range history from the subsystem composed of the i-th receiving array element and the transmitting array element to the point target as:

[0149]

[0150] Where: Denote the slant range of the \(i\)-th receiving array element in the presence of angular errors, with the unit of meters; \(t\) represents the slow time in the azimuth direction, with the unit of seconds; \(r\) represents the minimum distance between the target and the motion trajectory of the underwater acoustic transducer in the slant range plane, with the unit of meters; \(h\) represents the height of the underwater acoustic transducer from the bottom of the water, with the unit of meters; \(\theta\) represents the oblique angle of the underwater acoustic transducer, with the unit of degrees. Here, the three-dimensional coordinates of the point target in the azimuth direction, ground distance direction, and height direction are \(0\), \(r\sin\theta\), and \(h\), respectively, all with the unit of meters. \(c\) represents the speed of sound in water, with the unit of meters per second; \(v\) represents the speed of the underwater acoustic transducer platform, with the unit of meters per second. \(u\) i Denote the distance between the \(i\)-th receiving array element and the transmitting array element, with the unit of meters. \(\theta\) pitch Denote the pitch angle error, with the unit of degrees; \(\theta\) heading Denote the yaw error, with the unit of degrees.

[0151] For each subsystem composed of a receiving array element and a transmitting array element, calculate the two-way slant range history after phase center approximation, and its calculation formula is:[[]]

[0152]

[0153] where: \(\Delta R\) i Denote the distance delay error after phase center approximation, and its calculation formula is:[[]]

[0154]

[0155] where \(u\) i Denote the distance between the \(i\)-th (\(1\leq i\leq I\)) receiving array element and the transmitting array element; \(I\) represents the total number of receiving array elements; Denote the slant range of the \(i\)-th receiving array element in the presence of angular errors, with the unit of meters; \(r\) represents the minimum distance between the target and the motion trajectory of the underwater acoustic transducer in the slant range plane, with the unit of meters; \(h\) represents the height of the underwater acoustic transducer from the bottom of the water, with the unit of meters. \(\theta\) pitch Denote the pitch angle error, with the unit of degrees; \(\theta\) heading Denote the yaw error, with the unit of degrees.

[0156] According to the two-way slant range history from the subsystem composed of the \(i\)-th receiving array element and the transmitting array element to the point target and the two-way slant range history after phase center approximation, the system approximation error can be calculated, and the expression is:[[]]

[0157]

[0158] Figure 4 The system approximation errors under different angular error conditions. Among them Figure 4 (a) The yaw angle is \(0.5^{\circ}\) and the pitch angle is \(1^{\circ}\), Figure 4 (b) The yaw angle is \(1^{\circ}\) and the pitch angle is \(2^{\circ}\),Figure 4 (c) The yaw angle is 3° and the pitch angle is 4°. Figure 4 (d) The yaw angle is 5° and the pitch angle is 6°. It can be easily found from Figure 4 that under different angle error conditions, the system approximation error is small, all meeting the requirement that the error is less than π / 4.

[0159] Figure 5 The system approximation error of the underwater acoustic transducer under different beam widths when the yaw angle is 2° and the pitch angle is 3° is shown. Among them Figure 5 (a) The beam width is 3°. Figure 5 (b) The beam width is 6°. Figure 5 (c) The beam width is 9°. Figure 5 (d) The beam width is 12°. It can be clearly seen from the figure that under different beam width error conditions, the system approximation error is small, all meeting the requirement that the error is less than π / 4.

[0160] That is to say, considering the angle error, the two-way slant range history given by the present invention based on phase center approximation can better approximate the accurate two-way slant range history, which is also a prerequisite for obtaining accurate image reconstruction results.

[0161] According to the two-way slant range history after phase center approximation, calculate the system function of each transceiver array element subsystem in the range-Doppler domain. Its calculation formula is:

[0162]

[0163] Where: sS i (τ, f a ; r) represents the system function of the i-th receiving array element in the range-Doppler domain. represents the window function of the transmitted signal, W a (f a - f dc ) represents the combined directivity function of the transmitting and receiving array elements. τ represents the fast time in the range direction, with the unit of seconds. c represents the speed of sound in water, with the unit of meters per second; ΔR i represents the range delay error after phase center approximation, with the unit of meters. f a represents the azimuth Doppler frequency, with the unit of hertz; f dc represents the Doppler center frequency caused by the tilt of the underwater acoustic transducer, with the unit of hertz. r represents the closest distance between the target and the receiving array in the slant range plane, with the unit of meters. j 2 = -1 represents the imaginary unit. θ i (τ, f a ; r) represents the phase of the system function in the range-Doppler domain, and its expression is:

[0164]

[0165] Among them: f0 represents the center frequency of the transmitted chirp signal, with the unit of Hertz; u i represents the distance between the i-th (1 ≤ i ≤ I) receiving array element and the transmitting array element, with the unit of meter; I represents the total number of receiving array elements. K m represents the chirp rate in the case of combined azimuth and range coupling, with the unit of Hertz per second.

[0166]

[0167] Among them: K r represents the chirp rate of the transmitted signal, with the unit of Hertz per second; D(f a ) The physical meaning is that there is coupling in the azimuth and range directions, which is called the range curvature factor, and its expression is:

[0168]

[0169] Among them, v represents the moving speed of the underwater acoustic transducer platform, with the unit of meter per second.

[0170] For each subsystem composed of a receiving array element and a transmitting array element, compensate the Doppler phase error caused by the range error in the two-dimensional time domain:

[0171]

[0172] Among them represents the Doppler phase error, f0 represents the center frequency of the transmitted chirp signal, with the unit of Hertz; ΔR i represents the range delay error after phase center approximation, with the unit of meter; c represents the speed of sound in water, with the unit of meter per second.

[0173] For each subsystem composed of a receiving array element and a transmitting array element, divide the data into N non-overlapping data blocks in the range direction, and the width Δr of each data block is determined by the following formula:

[0174]

[0175] Among them, abs represents the absolute value operation, represents the value at the center distance r = r n,c of the n-th data block, represents the value at the edge distance r = r n,c -0.5Δr of the n-th data block, r n,c represents the center distance of the n-th data block, with the unit of meter; the expression (r n,c-0.5Δr) represents the edge distance of the nth data block, in meters; Δr represents the width of the data block, in meters. Thus, the total number of data blocks is:

[0176]

[0177] where r max represents the maximum distance that the underwater acoustic equipment can detect, in meters; r min represents the minimum distance detected by the underwater acoustic equipment, in meters; ceil represents the ceiling operation.

[0178] For any data block of the data of the subsystem composed of each receiving array element and transmitting array element, compensate for the micro-range migration error caused by the range error in the range frequency domain. The compensation function for the nth data block is:

[0179]

[0180] where Ω i (t, f r , n) represents the compensation function of the micro-range migration error in the nth data block, f r represents the instantaneous frequency in the range direction, in hertz; c represents the speed of sound in water, in meters per second; ΔR i,n represents the range delay error ΔR after phase center approximation i The error value at the center distance r = r n,c of the nth data block, in meters, and its calculation formula is:

[0181]

[0182] where u i represents the distance between the ith receiving array element and the transmitting array element, in meters; represents the value of the ith receiving array element at the center distance of the nth data block based on the angle error, in meters; r n,c represents the value at the center distance of the nth data block, in meters; h represents the height of the underwater acoustic transducer from the bottom of the water, in meters. θ pitch represents the pitch angle error, in degrees; θ heading represents the yaw error, in degrees.

[0183] For all data blocks after micro-range migration compensation of the subsystem composed of each receiving array element and transmitting array element, arrange them in order in the range time domain according to the data blocks. The arrangement method is as follows:

[0184]

[0185] where Denotes the data after micro-range migration compensation for the first data block of the i-th receiving array element, and so on. Denotes the data after micro-range migration compensation for the n-th data block of the i-th receiving array element. The subscript n (1 ≤ n ≤ N) represents the index of the data block of the i-th receiving array element. Denotes the signal after arranging all N data blocks of the i-th receiving array element in sequence.

[0186] For the signal after arranging all data blocks of each receiving array element in sequence, differential range migration correction processing is performed in the range-Doppler domain, and its filtering function is:

[0187]

[0188] where Γ i (τ, f a ) represents the filtering function for performing differential range migration correction processing on the i-th receiving array element in the range-Doppler domain. j 2 = -1 represents the imaginary unit. K m represents the chirp rate in the case of fusing azimuth and range coupling, with the unit of Hertz per second; D(f a ) represents the range curvature factor. D(f dc ) represents the value of the range curvature factor D(f a ) at the Doppler center f a = f dc ; c represents the speed of sound in water, with the unit of meters per second; r ref represents the reference distance, with the unit of meters per second; τ represents the fast time in the range direction, with the unit of seconds.

[0189] For the data after differential range migration correction processing of each receiving array element, filtering operations are performed in the two-dimensional frequency domain to simultaneously complete range-direction pulse compression, second-order range pulse compression, and uniform range migration correction processing. Its filtering function is:

[0190]

[0191] where H i (f r , f a ; r ref ) represents the filtering function for performing range-direction pulse compression, second-order range pulse compression, and uniform range migration correction processing on the i-th receiving array element in the two-dimensional frequency domain. D(f a ) represents the range curvature factor. D(f dc ) represents the value of the range curvature factor D(f a ) at the Doppler center f a = f dc ; K mDenotes the chirp rate in the case of combined azimuth and range coupling, with the unit of Hertz per second; f r Denotes the instantaneous range frequency, with the unit of Hertz; r ref Denotes the reference range, with the unit of meter; c denotes the speed of sound in water, with the unit of meter per second.

[0192] For the data after range pulse compression, second range pulse compression, and consistent range migration correction processing in the range direction for each receiving array element, azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing are performed in the range-Doppler domain, and its compensation function is:

[0193]

[0194] Where Λ i (τ, f a ) denotes the compensation function for azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing in the range-Doppler domain, D(f a ) denotes the range curvature factor, D(f dc ) denotes the range curvature factor D(f a ) at the Doppler center f a = f dc The value at, f0 denotes the center frequency of the transmitted chirp signal, with the unit of Hertz; K m Denotes the chirp rate in the case of combined azimuth and range coupling, with the unit of Hertz per second; f a Denotes the azimuth Doppler frequency, with the unit of Hertz; r ref Denotes the reference range, with the unit of meter; r denotes the shortest distance between the target and the motion trajectory of the underwater acoustic transducer in the slant range plane, with the unit of meter; c denotes the speed of sound in water, with the unit of meter per second.

[0195] For the data after azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing for all receiving array elements, coherent fusion processing is performed in the range-Doppler domain, and its processing method is:

[0196]

[0197] Where Denotes the result after coherent fusion of the data of each receiving array element in the range-Doppler domain, Denotes the result after azimuth focusing processing, residual phase error compensation operation, and azimuth walk correction processing for the data of the i-th receiving array element in the range-Doppler domain. The subscript i (1 ≤ i ≤ I) represents the index of the receiving array element, and I represents the total number of receiving array elements.

[0198] For all the data coherently fused by the receiving array elements, the azimuth inverse Fourier transform is performed to obtain the final high-resolution image, and the processing method is as follows:

[0199]

[0200] Among them, IFFT represents the inverse Fourier transform, and ff(τ,t) represents the final high-resolution image.

[0201] The following gives a set of typical parameters of an underwater acoustic synthetic aperture equipment system, and no weighting is used in the range direction and azimuth direction processing. The center frequency of the chirp signal is 75 kHz, the pulse width is 20 ms, the signal bandwidth is 20 kHz, the platform speed is 2.5 m / s, the pulse repetition frequency is 0.576 s, the number of receiving array elements is 36, the length of the transmitting array element along the azimuth direction is 0.16 m, the length of the receiving array element along the azimuth direction is 0.08 m, and the height of the underwater acoustic equipment from the bottom of the water is 50 m. Under the condition of using the typical synthetic aperture equipment system parameters, different yaw angle errors and pitch angle errors are adopted, and the approximate errors of the simulated system are as Figure 4 shown, where Figure 4 (a) The yaw angle is 0.5° and the pitch angle is 1°, Figure 4 (b) The yaw angle is 1° and the pitch angle is 2°, Figure 4 (c) The yaw angle is 3° and the pitch angle is 4°, Figure 4 (d) The yaw angle is 5° and the pitch angle is 6°. From Figure 4 it is not difficult to find that under different angle error conditions, the approximate error of the system is small, and all meet the requirement that the error is less than π / 4. The simulation results of the system approximate error show that the method given by the present invention can solve the problem of high-precision image reconstruction in the case of yaw angle error and pitch angle error. When the yaw angle is 2° and the pitch angle is 3°, the approximate errors of the system under different beam widths are as Figure 5 shown, where Figure 5 (a) The beam width is 3°, Figure 5 (b) The beam width is 6°, Figure 5 (c) The beam width is 9°, Figure 5 (d) The beam width is 12°. It is not difficult to see from the figure that under different beam width error conditions, the approximate error of the system is small, and all meet the requirement that the error is less than π / 4. The simulation results of the system approximate error show that the method given by the present invention can solve the problem of high-precision image reconstruction under different beam widths. For the equipment designed with the same parameters, the angle errors extracted during a certain voyage at sea are as Figure 6 shown, where Figure 6 (a) is the pitch angle error, Figure 6 (b) is the yaw angle error. If not Figure 6Compensate for the indicated angular error. The result after image reconstruction is as follows Figure 7 (a) shows that most traditional image reconstruction algorithms do not consider pose compensation. For example, the method proposed by authors such as Ma Mengbo in the article "CZT Imaging Algorithm for Multi-receiver Array Synthetic Aperture Sonar" published in the Journal of Huazhong University of Science and Technology (Natural Science Edition). Its image reconstruction performance is comparable to that of Figure 7 (a), which will greatly affect subsequent target detection and recognition. Only Figure 6 Compensate for the pitch angle error shown in Figure 7 (b). The result after image reconstruction is as follows Figure 7 (a). It can be found that the result after compensating for the pitch angle error is more accurate. Only Figure 6 Compensate for the yaw angle error shown in Figure 7 (c). The result after image reconstruction is as follows Figure 7 (a). It can be found that the improvement in the image reconstruction result after compensating for the yaw angle error is relatively obvious. At the same time, Figure 6 Compensate for the pitch angle error and yaw angle error shown in Figure 6 (b). The result after image reconstruction is as follows Figure 7 (d). Comparing the results shown in Figure 7 (a), Figure 7 (b), and Figure 7 (c), it can be found that the image reconstruction effect is more refined after compensating for the pitch angle error and yaw angle error simultaneously, further verifying the effectiveness and reliability of the method described in the present invention. It should be noted here that: Figure 7 The color bar on the right side of each figure represents the dynamic range of the signal energy amplitude of each pixel in the figure. In addition, based on Matlab software, the operation efficiency of the method of the present invention and the point-by-point image reconstruction algorithm based on angular error compensation for processing data with 4000 azimuth points and 8000 range points is compared. The image reconstruction time of the method of the present invention is about 26s, and the image reconstruction time using the method proposed by authors such as Ma Mengbo in the article "Motion Compensation of Multi-receiver Subarray SAS Based on Inertial Navigation System" published in the Journal of Huazhong University of Science and Technology (Natural Science Edition) is about 7060s. The image reconstruction efficiency of the method of the present invention is improved by about 272 times, further reflecting the high efficiency of the method of the present invention.

[0202] It should be understood that those of ordinary skill in the art can make improvements or transformations according to the above description, and all such improvements and transformations should fall within the protection scope of the appended claims of the present invention.

Claims

1. A method for underwater acoustic image reconstruction based on angle error compensation, comprising the following steps: (1) For each subsystem composed of receiving array elements and transmitting array elements, the spatial coordinate information of the transmitting array elements and the receiving array elements is calculated in the presence of angle errors, and the two-way slant range history corresponding to any transmitting and receiving array element subsystem is calculated to obtain the two-way slant range history between the subsystem composed of the i-th receiving array element and the transmitting array element and the point target; (2) According to the phase center approximation method, calculate the two-way slant range history after the phase center approximation; ⑶ According to the two-way slant range history after phase center approximation, calculate the system function of each transceiver array element subsystem in the range-Doppler domain; (4) For the Doppler phase error caused by the distance delay error in the phase center approximation processing, the Doppler phase error caused by the distance error is calculated and compensated in the two-dimensional time domain; ⑸ Divide the data into N non-overlapping data blocks, and compensate for the micro-range migration error caused by the range error in the range-to-frequency domain for any data block; (6) For all the data blocks after micro-distance migration compensation, arrange them in the order of the data blocks in the distance time domain to obtain the signal after distance delay error compensation; ⑺For each receiving array element after range delay error compensation, differential range migration correction processing is performed in the range-Doppler domain in turn, and range pulse compression, secondary range pulse compression, and consistent range migration correction processing are completed while filtering operations are performed in the two-dimensional frequency domain; ⑻For each receiving array element, after the data has been processed by range pulse compression, secondary range pulse compression, and consistent range migration correction, azimuth focusing, residual phase error compensation, and azimuth movement correction are performed while filtering in the range-Doppler domain; ⑼ For the data after azimuth focusing processing, residual phase error compensation operation and azimuth movement correction processing of all receiving array elements, coherent fusion processing is performed in the range-Doppler domain, and azimuth inverse Fourier transform is performed to obtain the final high-resolution image.

2. The method for underwater acoustic image reconstruction based on angle error compensation according to claim 1, characterized in that: The three-dimensional space coordinates of the transmitting array element in step (1) are: Where: v represents the platform movement speed, in meters per second; t represents the slow time in the azimuth direction, in seconds; x t represents the azimuth coordinate of the transmitting array element at time t, in meters; t represents the coordinate of the transmitting array element at time t in meters; z t It represents the coordinate of the transmitting array element at the time t in meters; The three-dimensional space coordinates of the i-th receiving array element are: Where: i (1≤i≤I) represents the index of the receiving array element; I represents the total number of receiving array elements; r represents the closest distance between the target and the motion trajectory of the underwater acoustic transducer on the slant range plane, in meters; u i represents the distance between the i-th receiving array element and the transmitting array element, in meters; represents the rotation vector, θ pitch Indicates the pitch angle error in degrees; θ heading Indicates the yaw angle error in degrees; x i represents the coordinate of the i-th receiving array element along the azimuth at time t, in meters; i represents the coordinate of the i-th receiving array element along the ground distance at time t, in meters; z i represents the coordinate of the i-th receiving array element along the height at time t, in meters; It represents the distance that the ith receiving element moves along the azimuth direction from the time when the transmitting element transmits the signal to the time when the ith receiving element receives the target echo signal, in meters; c represents the speed of sound in water, in meters per second; The two-way slant range from the subsystem composed of the i-th receiving array element and the transmitting array element to the point target is: in: It represents the slant range of the i-th receiving array element in the case of angle error, in meters; h represents the height of the underwater acoustic transducer from the water bottom, in meters; θ represents the slant angle of the underwater acoustic transducer, in degrees; the three-dimensional coordinates of the point target in azimuth, ground distance, and height are 0, rsinθ, h, all in meters.

3. The method for underwater acoustic image reconstruction based on angle error compensation according to claim 1, characterized in that: The calculation formula of the two-way slant range history after the phase center is approximated in step (2) is: Where: v represents the platform movement speed, in meters per second; t represents the slow time in the azimuth direction, in seconds; c represents the speed of sound in water, in meters per second; u i represents the distance between the i-th (1≤i≤I) receiving array element and the transmitting array element, in meters; I represents the total number of receiving array elements; r represents the closest distance between the target and the motion trajectory of the underwater acoustic transducer on the slant range plane, in meters; ΔR i Represents the distance delay error after phase center approximation, in meters, ΔR i The specific expression is: in: represents the slant distance of the i-th receiving array element in the case of angle error, in meters; h represents the height of the underwater acoustic transducer from the bottom of the water, in meters; θ represents the slant angle of the underwater acoustic transducer, in degrees; θ pitch Indicates the pitch angle error in degrees; θ heading Indicates the yaw angle error in degrees.

4. The method for underwater acoustic image reconstruction based on angle error compensation according to claim 1, characterized in that: The calculation formula of the system function of each transceiver array element subsystem in the range-Doppler domain in step (3) is: Where: sS i (τ,f a ; r) represents the system function of the i-th receiving array element in the range-Doppler domain, Represents the window function of the transmitted signal, W a (f a -f dc ) represents the joint directivity function of the receiving and transmitting array elements; j 2 =-1 represents the imaginary unit; τ represents the time from distance to speed, in seconds; c represents the speed of sound in water, in meters per second; ΔR i It represents the distance delay error after phase center approximation, in meters; f a Indicates the azimuth Doppler frequency in Hertz; f dc represents the Doppler center frequency caused by the squint of the hydroacoustic transducer, in Hertz; r represents the closest distance between the target and the receiving array on the slant range plane, in meters; θ i (τ,f a ; r) represents the phase of the system function, and the specific expression is: Where: u i represents the distance between the i-th (1≤i≤I) receiving array element and the transmitting array element, in meters; I represents the total number of receiving array elements; f0 represents the center frequency of the transmitted linear frequency modulation signal, in Hertz; v represents the movement speed of the underwater acoustic transducer platform, in meters per second; K m It represents the range modulation frequency considering the coupling of azimuth and distance, in Hertz / second. The specific expression is: Where: K r Indicates the frequency of the transmitted linear frequency modulation signal, in Hertz / second; D(f a )The physical meaning is that there is coupling between the azimuth and the distance, which is called the distance bending factor. The specific expression is:

5. The method for underwater acoustic image reconstruction based on angle error compensation according to claim 1, characterized in that: In step (4), the Doppler phase error caused by compensating the distance error in the two-dimensional time domain is calculated as follows: in: represents the Doppler phase error, f0 represents the center frequency of the transmitted linear frequency modulation signal, in Hertz; ΔR i It represents the distance delay error after phase center approximation, in meters; c represents the speed of sound in water, in meters per second.

6. The method for underwater acoustic image reconstruction based on angle error compensation according to claim 1, characterized in that: The width Δr of each data block in step (5) is determined by the following formula: Among them: abs means taking the absolute value operation, express The distance between the center of the nth data block is r = r n,c The value at express At the edge distance of the nth data block r = r n,c The value at -0.5Δr; r n,c Represents the center distance of the nth data block in meters; the expression (r n,c -0.5Δr) represents the edge distance of the nth data block, in meters; Δr represents the width of the data block, in meters; The total number of data blocks N is: Where: r max Indicates the maximum distance that the underwater acoustic transducer can detect, in meters; r min Indicates the minimum distance detected by the underwater acoustic transducer, in meters; ceil means rounding up to an integer; The compensation function of the nth data block is: Where: Ω i (t,f r ,n) represents the compensation function of the micro-distance migration error in the nth data block, f r represents the instantaneous frequency in the distance direction, in Hertz; c represents the speed of sound in water, in meters per second; ΔR i,n Represents the distance delay error ΔR after phase center approximation i The distance between the center of the nth data block is r = r n,c The error value at , in meters, is expressed as: Where: u i represents the distance between the i-th receiving array element and the transmitting array element, in meters; Represents the value of the i-th receiving array element at the center distance of the n-th data block based on the angle error, in meters; r n,c represents the value at the center distance of the nth data block, in meters; h represents the height of the underwater acoustic transducer from the bottom of the water, in meters; θ pitch Indicates the pitch angle error in degrees; θ heading Indicates the yaw angle error in degrees.

7. The method for underwater acoustic image reconstruction based on angle error compensation according to claim 1, characterized in that: The method of arranging the data blocks in the time domain in the distance direction in step (6) is as follows: in: represents the data after micro-distance migration compensation for the first data block of the i-th receiving array element, and so on. represents the data after micro-distance migration compensation for the nth data block of the i-th receiving array element, and the subscript n (1≤n≤N) represents the index of the i-th receiving array element data block. It represents the signal after all N data blocks of the i-th receiving array element are arranged in sequence.

8. The method for underwater acoustic image reconstruction based on angle error compensation according to claim 1, characterized in that: The filtering function of step (7) when performing differential range migration correction processing in the range-Doppler domain is: Where: i (τ,f a ) represents the filter function for performing differential range motion correction processing on the i-th receiving array element in the range-Doppler domain, j 2 =-1 represents the imaginary unit, K m It represents the modulation frequency when the azimuth and distance are coupled, in Hz / s; D(f a ) represents the distance bending factor, D(f dc ) represents the distance bending factor D(f a ) at the Doppler center f a =f dc The value at the point where c represents the speed of sound in water in meters per second; r ref represents the reference distance in meters; τ represents the distance acceleration time in seconds; The filter function for completing range pulse compression, secondary range pulse compression, and consistent range migration correction is: Where: H i (f r ,f a ; r ref ) represents the filtering function for performing range pulse compression, secondary range pulse compression, and consistent range migration correction on the i-th receiving array element in the two-dimensional frequency domain, D(f a ) represents the distance bending factor, D(f dc ) represents the distance bending factor D(f a ) at the Doppler center f a =f dc The value at K m represents the modulation frequency when the azimuth and distance are coupled, in Hz / s; f r Represents the instantaneous frequency in distance, in Hertz; r ref represents the reference distance in meters; c represents the speed of sound in water in meters per second.

9. The method for underwater acoustic image reconstruction based on angle error compensation according to claim 1, characterized in that: The compensation function for performing azimuth focusing processing, residual phase error compensation operation and azimuth movement correction processing in step (8) is: Where: i (τ,f a ) represents the compensation function for azimuth focusing, residual phase error compensation and azimuth motion correction in the range-Doppler domain, D(f a ) represents the distance bending factor, D(f dc ) represents the distance bending factor D(f a ) at the Doppler center f a =f dc The value at f0 represents the center frequency of the transmitted linear frequency modulation signal, in Hertz; K m represents the modulation frequency when the azimuth and distance are coupled, in Hz / s; f a Indicates the azimuth Doppler frequency in Hertz; r ref represents the reference distance in meters; r represents the shortest distance between the target and the motion trajectory of the hydroacoustic transducer on the slant range plane in meters; c represents the speed of sound in water in meters per second.

10. The method for underwater acoustic image reconstruction based on angle error compensation according to claim 1, characterized in that: The coherent fusion processing method in step (9) is: in: It represents the result of coherent fusion of each receiving array data in the range-Doppler domain. represents the result of azimuth focusing processing, residual phase error compensation operation and azimuth movement correction processing for the i-th receiving array data in the range-Doppler domain, the subscript i (1≤i≤I) represents the index of the receiving array element, and I represents the total number of receiving array elements; The final high-resolution image is obtained as follows: Where: IFFT represents inverse Fourier transform, and ff(τ,t) represents the final high-resolution image.

Images