Underwater sound image reconstruction method based on Legendre polynomial approximation

Abstract

The invention relates to a Legendre polynomial approximation-based underwater acoustic image reconstruction method, which comprises the following steps of: S1, redefining azimuth slow time by using an intermediate variable and synthetic aperture time, and calculating a two-way slant range process of a multi-subarray synthetic aperture sonar with a double-root-sign form of a function related to the intermediate variable; s2, obtaining an expression in a power series form about azimuth slow time; s3, calculating a baseband echo signal; s4, calculating a phase expression of the two-dimensional frequency domain system function; s5, calculating a phase expression in a power series form about the instantaneous frequency in the distance direction; s6, obtaining data of each receiving array element after range pulse compression and secondary range compression processing; s7, obtaining data after range migration correction of each receiving array element; s8, acquiring data of each receiving array element after azimuth pulse compression processing; and S9, performing coherent superposition and inverse Fourier transform to obtain a reconstructed underwater acoustic high-resolution image. The reconstructed synthetic aperture sonar underwater acoustic image is more accurate.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of intelligent manufacturing, and particularly relates to a water sound image reconstruction method based on Legendre polynomial approximation. BACKGROUND

[0002] The multi-subarray synthetic aperture sonar contains a transmitting array and multiple receiving array elements, and the two-way slant range history of the transmitting and receiving array elements includes the distance between the transmitting array element and the target and the distance between the target and the receiving array element, that is, the two-way slant range history of the transmitting and receiving array elements includes two expressions with square roots; when the phase stationary point is solved by using the phase stationary principle, the solution of the monomial eighth power equation is faced, so that the analytic phase stationary point cannot be obtained, and further the two-dimensional frequency domain system function in the form of analytic expression cannot be obtained. Since the precondition of the water sound image fast reconstruction algorithm is that the system must have the two-dimensional frequency domain system function in the form of analytic expression, the water sound image fast reconstruction algorithm cannot be designed.

[0003] In order to solve this problem, the traditional processing method generally performs phase center approximation on the two-way slant range history with the double square root form, and the phase center approximation is essentially only the result based on the second-order Taylor series approximation, which brings a large approximation error. SUMMARY

[0004] The technical problem to be solved by the present application is to provide a more accurate water sound image reconstruction method based on Legendre polynomial approximation.

[0005] To solve the above problems, the water sound image reconstruction method based on Legendre polynomial approximation provided by the present application comprises the following steps:

[0006] S1, calculating the synthetic aperture time T corresponding to the point target at the position of the farthest detection distance of the multi-subarray synthetic aperture sonar s , using an intermediate variable y with a value range of [-1, 1] and the synthetic aperture time T s , redefining the azimuth slow time t, and recalculating the two-way slant range history with the double square root form of the function about the intermediate variable y based on the newly defined azimuth slow time t and the original two-way slant range history;

[0007] S2, performing fourth-order Legendre polynomial approximation on the two-way slant range history with the double square root form based on the zero-order, first-order, second-order, third-order and fourth-order basis functions derived from the Rodrigues expression of the Legendre polynomial, and arranging the two-way slant range history after the fourth-order Legendre polynomial approximation into an expression in the form of power series about the azimuth slow time t

[0008] S3, based on the expression in the form of power series about the azimuth slow time t The baseband echo signal ss of the point target is calculated based on the wideband signal p(τ) emitted by the transmitting elements in the multi-subarray synthetic aperture sonar i (τ, t; r);

[0009] S4, the baseband echo signal ss of the point target is calculated based on the phase expression of the two-dimensional frequency domain system function i (τ, t; r), the phase expression of the two-dimensional frequency domain system function is calculated based on the phase stationary principle and the series inversion method;

[0010] S5, the phase expression of the two-dimensional frequency domain system function in the form of power series with respect to the instantaneous frequency in the distance direction is further calculated by the series approximation method for the calculated phase expression of the two-dimensional frequency domain system function;

[0011] S6, the phase for the secondary range compression in the two-dimensional frequency domain is extracted based on the phase expression of the two-dimensional frequency domain system function in the form of power series with respect to the instantaneous frequency in the distance direction in step S5, and the data after the range direction pulse compression and the secondary range compression processing of each receiving element in the multi-subarray synthetic aperture sonar is obtained by performing the range direction pulse compression and the secondary range compression processing in the two-dimensional frequency domain;

[0012] S7, the distance curvature curve expression in the range-Doppler domain is extracted based on the phase expression of the two-dimensional frequency domain system function in the form of power series with respect to the instantaneous frequency in the distance direction in step S5, and the data after the range migration correction of each receiving element is obtained by performing the inverse Fourier transform in the range direction, then performing the range migration correction in the range-Doppler domain by the sinc interpolation for the data after the range direction pulse compression and the secondary range compression processing of each receiving element in step S6;

[0013] S8, the phase for the azimuth direction pulse compression in the range-Doppler domain is extracted based on the phase expression of the two-dimensional frequency domain system function in the form of power series with respect to the instantaneous frequency in the distance direction in step S5, and the data after the azimuth direction pulse compression processing of each receiving element is obtained by performing the azimuth direction pulse compression in the range-Doppler domain for the data after the range migration correction of each receiving element in step S7;

[0014] S9, the data of all receiving elements in the range-Doppler domain is coherently superimposed for the data after the azimuth direction pulse compression processing of each receiving element, and the reconstructed underwater high-resolution image is obtained after the inverse Fourier transform in the azimuth direction.

[0015] The calculation formula of the synthetic aperture time T s in step S1 is as follows:

[0016]

[0017] In the formula, rfar θ represents the coordinates of the farthest target in the distance direction, in meters; -3dB The beamwidth of a single transmitting element is expressed in radians; V represents the tow speed of the sonar, expressed in meters per second.

[0018] The formula for calculating the azimuth time t is:

[0019]

[0020] The two-way slant distance history with double square roots is as follows:

[0021]

[0022] In the formula: r represents the coordinate of any target in the range direction, in meters; C represents the speed of sound in water, in meters per second; b i R represents the distance between the i-th receiving element and the transmitting element, in meters; i (y;r) represents the two-way slant range history of the i-th receiving array element, in meters.

[0023] The two-way slant distance history after the fourth-order Legendre polynomial approximation in step S2 is as follows:

[0024] R i (y; r)=c0p0(y)+c1p1(y)+c2p2(y)+c3p3(y)+c4p4(y)

[0025] In the formula: p0(y), p1(y), p2(y), p3(y), and p4(y) are the zeroth, first, second, third, and fourth order basis functions derived from the Rodrigues expression of the Legendre polynomial, respectively; and c is the coefficient of each order after approximation by the fourth order Legendre polynomial. n The formula for calculating (n=0,1,2,3,4) is:

[0026]

[0027] in, <R i (y;r),p n (y)> represents a two-way slant distance history with double square roots and basis functions of various orders p. n The inner product of (y)(n=0,1,2,3,4) over the interval [-1,1]; n denotes the order of the Legendre polynomial;

[0028] The two-way slant range history approximated by the fourth-order Legendre polynomial is presented as a power series expression in terms of the azimuth slow time t. for:

[0029]

[0030] The expressions for the coefficients are as follows:

[0031]

[0032] Where: c0, c1, c2, c3, and c4 are the coefficients of each order after the fourth-order Legendre polynomial approximation of the accurate two-way slant distance history.

[0033] The baseband echo signal ss of the point target in step S3 i The expression for (τ,t;r) is:

[0034]

[0035] In the formula: p is the broadband signal emitted by the multi-subarray synthetic aperture sonar transmitting element; τ represents the range-direction fast time, in seconds; is the echo delay time of the i-th receiving element, in seconds; r represents the coordinates of any target in the range direction, in meters; C represents the speed of underwater acoustics, in meters per second; f c The center frequency of the transmitted signal is represented by Hertz; j is the imaginary unit.

[0036] The phase expression of the two-dimensional frequency domain system function in step S4 is:

[0037]

[0038] In the formula: f τ f t These represent the instantaneous range frequency (for the fast range time) and the Doppler frequency (for the slow azimuth time), respectively, both in Hertz; f c The signal's center frequency is represented by Hertz; C represents the speed of sound in water by meters per second; R... cen_i k 1_i k 2_i k 3_i k 4_i These are the coefficients of the expression in power series form.

[0039] The phase expression of the two-dimensional frequency domain system function in step S5, in power series form with respect to the instantaneous frequency in the range direction.

[0040]

[0041] Here, the coefficients φ in the above formula a_i (f t ), φ rcm_i (f t ;r) The expressions are as follows:

[0042]

[0043]

[0044] In the formula: f τ f t These represent the instantaneous range frequency (for the fast range time) and the Doppler frequency (for the slow azimuth time), respectively, both in Hertz; f c The signal's center frequency is represented by Hertz; C represents the speed of sound in water by meters per second; R... cen_i k 1_i k 2_i k 3_i k 4_i These are the coefficients of the expression in power series form.

[0045] In step S6, the phase that undergoes secondary range compression in the two-dimensional frequency domain is extracted, and its expression is φ. src_i (f τ ,f t The filtering function H performs range pulse compression and secondary range compression processing on the echo signal of each receiving array element in the two-dimensional frequency domain. r_i (f τ )for:

[0046]

[0047] Here, P(f) τ f represents the spectrum of the broadband signal emitted by the transmitting elements of a multi-subarray synthetic aperture sonar, which is dimensionless; τ This represents the instantaneous frequency in the range direction corresponding to the fast time in the range direction, in Hertz; * indicates the conjugate operation; r s φ represents the distance to the reference target in the survey area, in meters; r represents the coordinate of any target in the range direction, in meters; src_i (f τ ,f t The sign indicates the square of the instantaneous frequency in the distance direction in step S5. coefficient; f t This represents the Doppler frequency corresponding to the azimuth slow time, in Hertz; Represents expression P * (f τ )exp{-jφ src_i (f τ ,f t At the reference target distance r s The phase at that point.

[0048] The expression for the distance curvature curve extracted in the distance-Doppler domain in step S7 is as follows: Furthermore, the range migration ΔR within the range-Doppler domain can be obtained. i (f t The expression for r is:

[0049]

[0050] Here, r represents the coordinates of any target in the range direction, in meters; f t This represents the Doppler frequency corresponding to the azimuth slow time, in Hertz; φ rcm_i (f t The expression for r is:

[0051]

[0052] In the formula: f c The signal's center frequency is represented by Hertz; C represents the speed of sound in water by meters per second; R... cen_i k 1_i k 2_i k 3_i k 4_i These are the coefficients of the expression in power series form.

[0053] In step S8, the phase φ that undergoes azimuth pulse compression in the range-Doppler domain is extracted. a_i (f t ); Filter function H used for azimuth pulse compression a_i (f t ;r) is:

[0054] H a_i (f t ;r)=exp{-jφ a_i (f t )}

[0055] Here, r represents the coordinates of any target in the range direction, in meters; f t This represents the Doppler frequency corresponding to the azimuth slow time, in Hertz; φ a_i (f t The expression for ) is:

[0056]

[0057] In the formula: f c The signal's center frequency is represented by Hertz; C represents the speed of sound in water by meters per second; R... cen_i k 1_i k 2_i k 3_i k 4_i These are the coefficients of the expression in power series form.

[0058] The result sS(f) after coherent superposition in step S9 t ;r) is:

[0059]

[0060] In the formula: sS i (f t ;r) represents the result of azimuth pulse compression performed by the i-th receiving element in the range-Doppler domain; N is the total number of receiving elements; r represents the coordinates of any target in the range direction, in meters; f t This represents the Doppler frequency corresponding to the azimuth slow time, and the unit is Hertz.

[0061] Compared with the prior art, the present invention has the following advantages:

[0062] 1. This invention uses the fourth-order Legendre polynomial to more accurately approximate the two-way slant range history in order to reduce the influence of distance approximation error. Based on this, the phase dwell principle and series response method are used to obtain a two-dimensional frequency domain system function with analytical expression form, and a rapid reconstruction method for underwater acoustic images of multi-subarray synthetic aperture sonar is designed accordingly.

[0063] 2. Compared with the phase center approximation method, the method of the present invention can achieve accurate compensation for distance migration, obtain a more accurate approximate distance, obtain a more accurate two-dimensional frequency domain system function, and thus obtain a more accurate synthetic aperture sonar underwater acoustic image reconstruction result. Attached Figure Description

[0064] The specific embodiments of the present invention will be described in further detail below with reference to the accompanying drawings.

[0065] Figure 1 This is a flowchart of the present invention. Wherein: FT represents Fourier transform, and IFT represents inverse Fourier transform.

[0066] Figure 2 This represents the distance delay error of the traditional phase center approximation method when the spacing between the transmit and receive array elements is 0.5m. In the figure, λ represents the corresponding center frequency f. c The wavelength.

[0067] Figure 3 This represents the distance delay error of the traditional phase center approximation method when the spacing between the transmit and receive array elements is 3m. In the figure, λ represents the corresponding center frequency f. c The wavelength.

[0068] Figure 4 This represents the distance delay error when the spacing between the transceiver array elements of this invention is 0.5m. In the figure, λ represents the corresponding center frequency f. c The wavelength.

[0069] Figure 5 This represents the distance delay error when the spacing between the transceiver array elements of this invention is 3m. In the figure, λ represents the corresponding center frequency f. c The wavelength.

[0070] Figure 6 This is an azimuthal cross-sectional view of the point target reconstruction result of the present invention. Detailed Implementation

[0071] This invention primarily uses a polynomial expression to more accurately approximate the two-way slant range history of multi-subarray synthetic aperture sonar with double square root form, in order to obtain a more accurate two-dimensional frequency domain system function, thereby enabling more accurate reconstruction results of multi-subarray synthetic aperture sonar images.

[0072] like Figure 1 As shown, a method for underwater acoustic image reconstruction based on Legendre polynomial approximation includes the following steps:

[0073] S1, calculate the synthetic aperture time T for the point target at the furthest detection range location corresponding to the multi-subarray synthetic aperture sonar. s The synthesis pore size time T s The formula for calculating (in seconds) is:

[0074]

[0075] In the formula: r far θ represents the coordinates of the farthest target in the distance direction, in meters; -3dB The beamwidth of a single transmitting element is expressed in radians; V represents the tow speed of the sonar, expressed in meters per second.

[0076] Define an intermediate variable y with a value range of [-1, 1]. Then, use this intermediate variable y with a value range of [-1, 1] and the synthesis aperture time T. s Redefining the azimuth slow time t (in seconds), the result is:

[0077] The formula for calculating the azimuth time t is:

[0078]

[0079] Based on the newly defined azimuth slow time t and the original two-way slant range history The calculation of the two-way slant range history of a multi-subarray synthetic aperture sonar with a function of the intermediate variable y is as follows:

[0080]

[0081] In the formula: r represents the coordinate of any target in the range direction, in meters; C represents the speed of sound in water, in meters per second; bi R represents the distance between the i-th receiving element and the transmitting element, in meters; i (y;r) represents the two-way slant range history of the i-th receiving array element, in meters.

[0082] S2, based on the zeroth, first, second, third, and fourth order basis functions derived from the Rodrigues expression of Legendre polynomials, approximates the two-way slant range history with double radical forms using fourth-order Legendre polynomials. The approximated two-way slant range history is then presented as a power series expression about the azimuth slow time t.

[0083] The two-way slant distance history after approximation with the fourth-order Legendre polynomial is:

[0084] R i (y; r)=c0p0(y)+c1p1(y)+c2p2(y)+c3p3(y)+c4p4(y)

[0085] In the formula: p0(y), p1(y), p2(y), p3(y), and p4(y) are the zeroth, first, second, third, and fourth order basis functions derived from the Rodrigues expression of the Legendre polynomial, respectively; and c is the coefficient of each order after approximation by the fourth order Legendre polynomial. n The formula for calculating (n=0,1,2,3,4) is:

[0086]

[0087] in, <R i (y;r),p n (y)> represents a two-way slant distance history with double square roots and basis functions of various orders p. n The inner product of (y)(n=0,1,2,3,4) over the interval [-1,1]; n denotes the order of the Legendre polynomial;

[0088] The two-way slant range history approximated by the fourth-order Legendre polynomial is presented as a power series expression in terms of the azimuth slow time t. for:

[0089]

[0090] The expressions for the coefficients are as follows:

[0091]

[0092] Where: c0, c1, c2, c3, and c4 are the coefficients of each order after the fourth-order Legendre polynomial approximation of the accurate two-way slant distance history.

[0093] S3, an expression based on the power series form of the azimuth direction at a slow time t. Calculate the baseband echo signal ss of the point target using the broadband signal p(τ) emitted by the transmitting element in a multi-subarray synthetic aperture sonar. i (τ,t;r).

[0094] Assuming the broadband signal emitted by the sonar transmitting element is p(τ), then the echo delay time of the i-th receiving element corresponding to an ideal point target with azimuth coordinate 0 and range coordinate r is: This represents the transformation of the transmit / receive slant range history, approximated by a fourth-order polynomial, into a power series expression with respect to the azimuth slow time t; based on the echo delay time of the i-th receiving element. (Unit: seconds) The echo signal of the receiving array element after time-domain modulation can be obtained. for:

[0095]

[0096] Where: f c τ represents the center frequency of the transmitted signal, measured in Hertz; τ represents the distance-time, measured in seconds.

[0097] The time-domain modulated echo signal calculated above Demodulation yields the baseband echo signal ss of the i-th receiving element in the time domain. i (τ,t;r), its expression is:

[0098]

[0099] In the formula: p is the broadband signal emitted by the multi-subarray synthetic aperture sonar transmitting element; is the echo delay time of the i-th receiving array element, in seconds; j is an imaginary unit.

[0100] S4, for the baseband echo signal ss of a point target. i (τ,t;r), calculate the phase expression of the two-dimensional frequency domain system function based on the phase dwell principle and series inversion method.

[0101] The baseband echo signal ss of the point target obtained in the previous step i By performing a two-dimensional frequency domain transformation on (τ, t; r), the corresponding two-dimensional frequency domain system function SS can be obtained. i (f τ ,f t ;r) is:

[0102]

[0103] Here, P(f)τ f represents the spectrum of the broadband signal p(τ) emitted by the sonar transmitting element; τ f t These represent the instantaneous range frequency (for the fast range time) and the Doppler frequency (for the slow azimuth time), respectively, both in Hertz; f c This indicates the center frequency of the transmitted signal, measured in Hertz (Hz). This means that the two-way slant range history after approximation with the fourth-order Legendre polynomial is presented as a power series expression in terms of the azimuth slow time.

[0104] For two-dimensional frequency domain system functions SS i (f τ ,f t The phase expression in the integral term, based on the phase dwell principle, shows that the two-dimensional frequency domain system function SS i (f τ ,f t The phase expression in the integral term has a first derivative of 0 with respect to the azimuth direction and slow time; the corresponding expression is:

[0105]

[0106] The two-way slant range history, approximated by the fourth-order Legendre polynomial, is presented as a power series expression in terms of the azimuth slow time. It can be further obtained

[0107]

[0108] The phase dwell point can be obtained by using the series inversion method. for:

[0109]

[0110] Here, the first-order, second-order, and third-order coefficients are respectively

[0111] Ignoring the spectrum P(f) of the broadband signal p(τ) emitted by the sonar transmitting elements τ Based on the above two-dimensional frequency domain system function SS i (f τ ,f t ;r) and phase dwell point The phase expression of the two-dimensional frequency domain system function corresponding to the baseband echo signal of a point target can be calculated, and the result is:

[0112]

[0113] S5. For the phase expression of the calculated two-dimensional frequency domain system function, the power series phase expression of the two-dimensional frequency domain system function with respect to the instantaneous frequency in the range direction is further calculated using the series approximation method.

[0114] For three expressions containing instantaneous frequency in the range direction Regarding f τ Using the second-order approximation of 0, we can obtain the following three approximate expressions:

[0115]

[0116]

[0117] Based on the above three approximate expressions and the phase expression of the two-dimensional frequency domain system function obtained in the previous step, the phase expression of the two-dimensional frequency domain system function is further rearranged into a power series phase expression with respect to the instantaneous frequency in the range direction. The result is:

[0118]

[0119] Here, the coefficients φ in the above formula a_i (f t ), φ rcm_i (f t ;r) The expressions are as follows:

[0120]

[0121] S6, based on the phase expression of the two-dimensional frequency domain system function in step S5 in the power series form with respect to the instantaneous frequency in the range direction, extract the phase that undergoes secondary range compression in the two-dimensional frequency domain, and its expression is φ. src_i (f τ ,f t The filtering function H performs range pulse compression and secondary range compression processing on the echo signal of each receiving element in a multi-subarray synthetic aperture sonar in the two-dimensional frequency domain. r_i (f τ )for:

[0122]

[0123] Here, P(f) τ f represents the spectrum of the broadband signal emitted by the transmitting elements of a multi-subarray synthetic aperture sonar, which is dimensionless; τ This represents the instantaneous frequency in the range direction corresponding to the fast time in the range direction, in Hertz; * indicates the conjugate operation; r s φ represents the distance to the reference target in the survey area, in meters; r represents the coordinate of any target in the range direction, in meters; src_i(f τ ,f t The sign indicates the square of the instantaneous frequency in the distance direction in step S5. coefficient; f t This represents the Doppler frequency corresponding to the azimuth slow time, in Hertz; Represents expression P * (f τ )exp{-jφ src_i (f τ ,f t At the reference target distance r s The phase at that point.

[0124] Finally, the data for each receiving array element after range pulse compression and secondary range compression processing is obtained.

[0125] S7, based on the phase expression of the two-dimensional frequency domain system function in step S5 in the power series form with respect to the instantaneous frequency in the range direction, the range curvature curve expression in the range-Doppler domain is extracted as follows:

[0126] φ rcm_i (f t The expression for r is:

[0127]

[0128] For the data processed by range pulse compression and secondary range compression for each receiving element in step S6, an inverse Fourier transform is performed in the range direction, and then range migration correction is performed in the range-Doppler domain using Sink interpolation. This further yields the range migration ΔR in the range-Doppler domain. i (f t The expression for ;r (in meters) is:

[0129]

[0130] Finally, the distance migration correction data for each receiving array element is obtained.

[0131] S8, based on the phase expression of the two-dimensional frequency domain system function in step S5 in the power series form with respect to the instantaneous frequency in the range direction, extract the phase φ for azimuth pulse compression in the range-Doppler domain. a_i (f t Its expression is:

[0132]

[0133] Then, for the range migration corrected data of each receiving element in step S7, azimuth pulse compression is performed in the range-Doppler domain. The filter function H used for azimuth pulse compression is... a_i(f t ;r) is:

[0134] H a_i (f t ;r)=exp{-jφ a_i (f t )}

[0135] Finally, the data after azimuth pulse compression processing for each receiving array element is obtained.

[0136] S9, for the data after azimuth pulse compression processing of each receiving array element, coherently superimpose the data of all receiving array elements in the range-Doppler domain, and the result of the coherent superposition is sS(f t ;r) is:

[0137]

[0138] In the formula: sS i (f t ;r) represents the result of azimuth pulse compression of the i-th receiving element in the range-Doppler domain; N is the total number of receiving elements.

[0139] Then, after performing an inverse Fourier transform in the azimuth direction, the reconstructed high-resolution underwater acoustic image is obtained.

[0140] To verify the feasibility of the method of this invention, a simulation experiment was designed. The simulation parameters are as follows: the bandwidth of the transmitted linear frequency modulated signal is 20000Hz, the carrier frequency of the transmitted signal is 150000Hz, the pulse repetition frequency of the transmitted pulse signal is 0.4s, the sonar tow speed is 4m / s, the actual apertures of the transceiver array elements in the azimuth direction are 2cm and 4cm, respectively, and the total number of receiving array elements is 80. It is assumed that there exists an ideal point target in space with coordinates of 152m in the range direction and 35m in the azimuth direction.

[0141] Figure 2 , Figure 3 The distance errors of the traditional phase center approximation method are calculated for the transmit / receive array element spacing of 0.5m and 3m, respectively. Figure 4 , Figure 5 The figures show the range errors of the method of the present invention when the distance between the transmitting and receiving array elements is 0.5m and 3m, respectively. By comparison, it is easy to see that the method of the present invention has a smaller range error, verifying the conclusion that the method of the present invention can more accurately approximate the two-way slant range history of multi-subarray synthetic aperture sonar with double square root form.

[0142] After reconstructing the target, the reconstruction result is shown in the azimuth profile as follows. Figure 6 As shown, the azimuth profile of the traditional phase center approximation method has many burrs, while the azimuth profile of the method of the present invention is smoother and has no burrs.

[0143] The peak-to-side-lobe ratio, integral-to-side-lobe ratio, and azimuth resolution of the traditional phase center approximation method and the method of this invention were statistically compared to further evaluate their reconstruction performance. The peak-to-side-lobe ratio, integral-to-side-lobe ratio, and azimuth resolution of the traditional phase center approximation method were -14.45 dB, -9.98 dB, and 0.02 m, respectively. The peak-to-side-lobe ratio, integral-to-side-lobe ratio, and azimuth resolution of the method of this invention were -14.91 dB, -10.17 dB, and 0.02 m, respectively. This shows that the peak-to-side-lobe ratio and integral-to-side-lobe ratio of the method of this invention are improved by 0.46 dB and 0.19 dB, respectively, compared to the traditional phase center approximation method. The azimuth resolution of both methods is consistent. Figure 6 The azimuth profile shown and the statistically obtained peak sidelobe ratio, integral sidelobe ratio, azimuth resolution, and other indicators clearly demonstrate that the method of this invention has superior reconstruction performance compared to the traditional phase center approximation method.

[0144] It should be understood that those skilled in the art can make improvements or modifications based on the above description, and all such improvements and modifications should fall within the protection scope of the appended claims.

Claims

1. A method for underwater acoustic image reconstruction based on Legendre polynomial approximation, comprising the following steps: S1, calculate the synthetic aperture time T for the point target at the furthest detection range location corresponding to the multi-subarray synthetic aperture sonar. s Using an intermediate variable y with a value range of [-1,1] and the synthesis aperture time T s The azimuth slow time t is redefined, and the two-way slant range history of the multi-subarray synthetic aperture sonar with double radical form is recalculated based on the newly defined azimuth slow time t and the original two-way slant range history. S2, based on the zeroth, first, second, third, and fourth order basis functions derived from the Rodrigues expression of Legendre polynomials, approximates the two-way slant range history with double radical forms using fourth-order Legendre polynomials. The approximated two-way slant range history is then presented as a power series expression about the azimuth slow time t. S3, an expression based on the power series form of the azimuth direction at a slow time t. Calculate the baseband echo signal ss of the point target using the broadband signal p(τ) emitted by the transmitting element in a multi-subarray synthetic aperture sonar. i (τ,t;r); S4, for the baseband echo signal ss of a point target. i (τ,t;r), calculate the phase expression of the two-dimensional frequency domain system function based on the phase dwell principle and series inversion method; S5. For the phase expression of the calculated two-dimensional frequency domain system function, the power series phase expression of the two-dimensional frequency domain system function with respect to the instantaneous frequency in the range direction is further calculated using the series approximation method. S6. Based on the phase expression of the two-dimensional frequency domain system function in step S5 in the form of a power series with respect to the instantaneous frequency in the range direction, extract the phase that has undergone secondary range compression in the two-dimensional frequency domain. Furthermore, range pulse compression and secondary range compression processing are performed on the echo signal of each receiving element in the multi-subarray synthetic aperture sonar in the two-dimensional frequency domain to obtain the data of each receiving element after range pulse compression and secondary range compression processing. S7. Based on the phase expression of the two-dimensional frequency domain system function in step S5 in the form of a power series with respect to the instantaneous frequency in the range direction, the range curvature curve expression in the range-Doppler domain is extracted. For the data after range pulse compression and secondary range compression processing for each receiving element in step S6, an inverse Fourier transform is performed in the range direction. Then, range migration correction is performed in the range-Doppler domain using Sink interpolation to obtain the range migration corrected data for each receiving element. S8. Based on the phase expression of the two-dimensional frequency domain system function in step S5 in the form of a power series with respect to the instantaneous frequency in the range direction, extract the phase of the azimuth pulse compression in the range-Doppler domain, and perform azimuth pulse compression in the range-Doppler domain on the data after range migration correction for each receiving element in step S7 to obtain the data after azimuth pulse compression processing for each receiving element. S9: For the data after pulse compression processing in the azimuth direction of each receiving array element, the data of all receiving array elements are coherently superimposed in the range-Doppler domain, and after performing inverse Fourier transform in the azimuth direction, the reconstructed high-resolution underwater acoustic image is obtained.

2. The underwater acoustic image reconstruction method based on Legendre polynomial approximation as described in claim 1, characterized in that: The pore size synthesis time T in step S1 s The calculation formula is: In the formula: r far θ represents the coordinates of the farthest target in the distance direction, in meters; 3dB The beamwidth of a single transmitting element is expressed in radians; V represents the tow speed of the sonar, expressed in meters per second. The formula for calculating the azimuth time t is: The two-way slant distance history with double square roots is as follows: In the formula: r represents the coordinate of any target in the range direction, in meters; C represents the speed of sound in water, in meters per second; b i R represents the distance between the i-th receiving element and the transmitting element, in meters; i (y;r) represents the two-way slant range history of the i-th receiving array element, in meters.

3. The underwater acoustic image reconstruction method based on Legendre polynomial approximation as described in claim 1, characterized in that: The two-way slant distance history after the fourth-order Legendre polynomial approximation in step S2 is as follows: R i (y;r)=c0p0(y)+c1p1(y)+c2p2(y)+c3p3(y)+c4p4(y) In the formula: p0(y), p1(y), p2(y), p3(y), and p4(y) are the zeroth, first, second, third, and fourth order basis functions derived from the Rodrigues expression of the Legendre polynomial, respectively; and c is the coefficient of each order after approximation by the fourth order Legendre polynomial. n The formula for calculating (n=0,1,2,3,4) is: in, <R i (y;r),p n (y)> represents a two-way slant distance history with double square roots and basis functions of various orders p. n The inner product of (y)(n=0,1,2,3,4) over the interval [-1,1]; n denotes the order of the Legendre polynomial; The two-way slant range history approximated by the fourth-order Legendre polynomial is presented as a power series expression in terms of the azimuth slow time t. for: The expressions for the coefficients are as follows: Where: c0, c1, c2, c3, and c4 are the coefficients of each order after the fourth-order Legendre polynomial approximation of the accurate two-way slant distance history.

4. The underwater acoustic image reconstruction method based on Legendre polynomial approximation as described in claim 1, characterized in that: The baseband echo signal ss of the point target in step S3 i The expression for (τ,t;r) is: In the formula: p is the broadband signal emitted by the multi-subarray synthetic aperture sonar transmitting element; τ represents the range-direction fast time, in seconds; is the echo delay time of the i-th receiving element, in seconds; r represents the coordinates of any target in the range direction, in meters; C represents the speed of underwater acoustics, in meters per second; f c The center frequency of the transmitted signal is represented by Hertz; j is the imaginary unit.

5. The underwater acoustic image reconstruction method based on Legendre polynomial approximation as described in claim 1, characterized in that: The phase expression of the two-dimensional frequency domain system function in step S4 is: In the formula: f τ f t These represent the instantaneous range frequency and the Doppler frequency corresponding to the fast range time and the slow azimuth time, respectively, both in Hertz; f c The signal's center frequency is represented by Hertz; C represents the speed of sound in water by meters per second; R... cen_i k 1_i k 2_i k 3_i k 4_i These are the coefficients of the expression in power series form.

6. The underwater acoustic image reconstruction method based on Legendre polynomial approximation as described in claim 1, characterized in that: The phase expression of the two-dimensional frequency domain system function in step S5, in power series form with respect to the instantaneous frequency in the range direction. Φ i (f τ ,f t )≈φ a_i (f t )+φ rcm_i (f t ;r)f τ +φ src_i (f τ ,f t )f τ 2 ; Here, the coefficients v in the above formula a_i (f t ), φ rcm_i (f t ;r) The expressions are as follows: In the formula: f τ f t These represent the instantaneous range frequency and the Doppler frequency corresponding to the fast range time and the slow azimuth time, respectively, both in Hertz; f c The signal's center frequency is represented by Hertz; C represents the speed of sound in water by meters per second; R... cen_i k 1_i k 2_i k 3_i k 4_i These are the coefficients of the expression in power series form.

7. The underwater acoustic image reconstruction method based on Legendre polynomial approximation as described in claim 1, characterized in that: In step S6, the phase that undergoes secondary range compression in the two-dimensional frequency domain is extracted, and its expression is φ. src_i (f τ ,f t The filtering function H performs range pulse compression and secondary range compression processing on the echo signal of each receiving array element in the two-dimensional frequency domain. r_i (f τ )for: Here, P(f) τ f represents the spectrum of the broadband signal emitted by the transmitting elements of a multi-subarray synthetic aperture sonar, which is dimensionless; τ This represents the instantaneous frequency in the range direction corresponding to the fast time in the range direction, in Hertz; * indicates the conjugate operation; r s φ represents the distance to the reference target in the survey area, in meters; r represents the coordinate of any target in the range direction, in meters; src_i (f τ ,f t ) represents the square of the instantaneous frequency f in the range direction in step S5. τ 2 coefficient; f t This represents the Doppler frequency corresponding to the azimuth slow time, in Hertz; Represents expression P * (f τ )exp{-jφ src_i (f τ ,f t At the reference target distance r s The phase at that point.

8. The underwater acoustic image reconstruction method based on Legendre polynomial approximation as described in claim 1, characterized in that: The expression for the distance curvature curve extracted in the distance-Doppler domain in step S7 is as follows: Furthermore, the range migration ΔR within the range-Doppler domain can be obtained. i (f t The expression for r is: Here, r represents the coordinates of any target in the range direction, in meters; f t This represents the Doppler frequency corresponding to the azimuth slow time, in Hertz; φ rcmi (f t The expression for r is: In the formula: f c The signal's center frequency is represented by Hertz; C represents the speed of sound in water by meters per second; R... cen_i k 1_i k 2_i k 3_i k 4_i These are the coefficients of the expression in power series form.

9. The underwater acoustic image reconstruction method based on Legendre polynomial approximation as described in claim 1, characterized in that: In step S8, the phase φ that undergoes azimuth pulse compression in the range-Doppler domain is extracted. a_i (f t ); Filter function H used for azimuth pulse compression a_i (f t ;r) is: H a_i (f t ;r)=exp{-jφ a_i (f t )} Here, r represents the coordinates of any target in the range direction, in meters; f t This represents the Doppler frequency corresponding to the azimuth slow time, in Hertz; φ a_i (f t The expression for ) is: In the formula: f c The signal's center frequency is represented by Hertz; C represents the speed of sound in water by meters per second; R... cen_i k 1_i k 2_i k 3_i k 4_i These are the coefficients of the expression in power series form.

10. The underwater acoustic image reconstruction method based on Legendre polynomial approximation as described in claim 1, characterized in that: The result sS(f) after coherent superposition in step S9 t ;r) is: In the formula: sS i (f t ;r) represents the result of azimuth pulse compression performed by the i-th receiving array element in the range-Doppler domain; N is the total number of receiving array elements; r represents the coordinates of any target in the range direction, in meters; f t This represents the Doppler frequency corresponding to the azimuth slow time, and the unit is Hertz.

Images