1. A multilayer overlapping sub-aperture imaging method of a large squint time-varying parameter SAR is characterized by comprising the following steps:

step one, adopting a time-varying parameter system to transmit signals, receiving echo signals and performing pulse compression processing;

step two, only carrying out azimuth interpolation on the wave number domain signals obtained after pulse pressure processing;

step three, compressing the signals processed in the step two in the distance direction, and dividing the signals into N layers of overlapped sub-apertures along the azimuth direction, wherein the N layers of overlapped sub-apertures correspond to the N + 1-dimensional azimuth direction signals;

performing fast Fourier transform on the front N-dimensional azimuth signals one by one dimension to obtain an estimated value of a distortion coordinate of each dimension; solving an estimated value of a real coordinate according to the estimated value of the distortion coordinate, and compensating the space-variant secondary phase error by using the estimated value of the real coordinate; performing fast Fourier transform on the (N + 1) th dimensional azimuth signal to obtain an estimated value of an (N + 1) th dimensional distortion coordinate;

and step five, reconstructing the N + 1-dimensional azimuth direction signal processed in the step four to obtain an SAR image.

2. The imaging method according to claim 1, wherein the dividing into N layers of overlapping sub-apertures along the azimuth direction in step three is performed by:

according to the length M of a single sub-aperture₁By Δ₂For data extraction ratio, extract M₂A sub-aperture constituting a first layer overlapA sub-aperture; by Delta₂Δ₃For data extraction ratio, extract M₃A first layer of overlapping sub-apertures forming a second layer of overlapping sub-apertures; by Delta₂Δ₃Δ₄For data extraction ratio, extract M₄A second layer of overlapping sub-apertures, constituting a third layer of overlapping sub-apertures … …, and so on until N layers of overlapping sub-apertures are obtained;

wherein M is₁……M_NIs a positive integer, and M₁……M_NThe value of (a) ensures that the secondary phase of each layer of overlapping sub-apertures is less than pi/4.

3. The imaging method according to claim 1, wherein the dividing into N layers of overlapping sub-apertures along the azimuth direction in the third step, and then obtaining an N + 1-dimensional azimuth signal, the N + 1-dimensional azimuth signal being represented as:

where i' is a range sequence after the range fast Fourier transform, m₁represents-M₁/2,…,M₁Sequence of/2, m₂represents-M₂/2,…,M₂The sequence of/2, and so on; x is the number of_pAnd y_pIs the coordinate of the imaging target, U₁And U₂Geometrical distortions of the distance and azimuth directions caused by the wave front hypothesis, respectivelyVariable, U₃Are azimuth quadratic phase coefficients caused by wave front hypothesis, which are only determined by scene definite value parameters and scene target point coordinates x_p,y_pDetermination of,. DELTA.K_yRepresenting the distance number bandwidth, Δ, after azimuth interpolation₂Representing the first data extraction ratio, Δ₂Δ₃Indicating a second data extraction ratio, … … Δ₂Δ₃…Δ_N+1Expressing the extraction proportion of the Nth data; j is an imaginary symbol, u_xIs the azimuthal wavenumber interval.

4. The imaging method according to claim 1, wherein in step four, the solving of the estimated value of the true coordinate from the estimated value of the distorted coordinate is:

solving an estimated value of the real coordinate in the form of an analytic solution according to a higher-order signal estimation model, wherein the higher-order signal estimation model is expressed as:

where i and m represent fast and slow time sequences, respectively, j is an imaginary symbol, x_pAnd y_pIs the coordinate of the imaging target, U₁And U₂The geometrical distortion, U, of the distance and orientation directions respectively caused by the wave front hypothesis₃Is an azimuthal quadratic phase coefficient, K, caused by a wave front hypothesis_x(i, m) is the azimuth wave number, K_y(i) Distance-wise wave number;

the analytic solution is represented as:

wherein

Wherein N represents the nth dimension azimuth signal, and N belongs to [1, 2]；Yc represents the distance from the ground projection of the radar to the center of the scene at the moment of the center of the synthetic aperture, beta₀Showing the downward view angle of the synthetic aperture center moment, phi the downward-looking angle of the radar speed, delta the complementary angle of the ground oblique angle, H the flying height of the platform at the radar synthetic aperture center moment, y_QAs an estimate of the distance-wise distorted coordinates, x_QIs an estimate of the azimuthal distortion coordinate.

5. The imaging method according to claim 1, wherein the time-varying parameter system is specifically: adjusting radar parameter carrier frequency and frequency in real time along the azimuth direction to counteract the distance wave number variation caused by the variation of downward viewing angle and slant ground declination;

the carrier frequency and the frequency modulation rate are expressed as:

where c is the speed of light, f_c0And gamma₀Respectively representing the carrier frequency of the initial transmission signal and the modulation frequency, beta, of the chirp signal₀The downward viewing angle at the moment of the center of the synthetic aperture is shown, beta (t) represents the instantaneous downward viewing angle from the platform to the center of the scene, and alpha (t) represents the instantaneous included angle formed by the projection of the instantaneous reference slope distance on the ground plane and OY (Y axis).

6. The imaging method of claim 1, wherein said pulse compression process includes a de-line tone and a remove remaining video phase operation.

Background

Synthetic Aperture Radar (SAR) is an active microwave imaging Radar, has all-weather and all-day earth surface information acquisition capability, plays an important role in aspects of military reconnaissance, resource exploration, topographic mapping and the like, and is an indispensable technical means in the field of earth observation at present. The large squint bunching SAR is an important working mode, the complementary angle (squint angle) of the included angle between the beam sight direction and the track direction can reach eighty degrees or even more than eighty degrees, scene information can be sensed in advance, and the reconnaissance capability of the radar can be obviously improved. However, the large squint SAR imaging has the problems of poor image resolution and low imaging accuracy. Under a large squint configuration, a traditional constant radar parameter system causes severe distortion of a multilayer Function (PSF) of a target Point in an SAR image, and spatial resolution deteriorates. The larger the squint angle is, the larger the range migration is, the more the time-space varying coupling of the echo is serious, and the imaging processing algorithm needs to have strong time-varying echo decoupling capacity because the two-dimensional coupling rule has time-space dependence. The traditional time domain imaging processing method has the capacity, but under a large squint configuration, the calculation amount is large, the efficiency is too low, and the real-time requirement of the SAR is not met; the processing precision of the traditional range-Doppler domain algorithm is reduced along with the increase of an oblique angle until the imaging cannot be performed, and the congenital defect exists; the traditional wave number domain Algorithm is represented by a Polar Format Algorithm (PFA), and realizes two-dimensional decoupling through two-dimensional interpolation of a wave number domain, but is limited by the assumption of a wave front plane, the farther the distance from the center of a scene is, the larger the approximation error is, the more limited the imaging width is, and the smaller the imaging width is. Therefore, there is a considerable need to develop SAR imaging systems and algorithm studies under large squints.

Disclosure of Invention

In view of the above, the invention provides a multilayer overlapping subaperture imaging method of a large squint time-varying parameter SAR, which can correct a distortion target point expansion function, improve imaging performance, omit an original distance interpolation operation, simplify an imaging processing flow, and improve operation efficiency; and the phase position can be compensated with higher precision, and the imaging width is improved.

The specific scheme of the invention is as follows:

a multilayer overlapping sub-aperture imaging method of a large squint time-varying parameter SAR comprises the following steps:

step one, adopting a time-varying parameter system to transmit signals, receiving echo signals and performing pulse compression processing;

step two, only carrying out azimuth interpolation on the wave number domain signals obtained after pulse pressure processing;

step three, compressing the signals processed in the step two in the distance direction, and dividing the signals into N layers of overlapped sub-apertures along the azimuth direction, wherein the N layers of overlapped sub-apertures correspond to the N + 1-dimensional azimuth direction signals;

performing fast Fourier transform on the front N-dimensional azimuth signals one by one dimension to obtain an estimated value of a distortion coordinate of each dimension; solving an estimated value of a real coordinate according to the estimated value of the distortion coordinate, and compensating the space-variant secondary phase error by using the estimated value of the real coordinate; performing fast Fourier transform on the (N + 1) th dimensional azimuth signal to obtain an estimated value of an (N + 1) th dimensional distortion coordinate;

and step five, reconstructing the N + 1-dimensional azimuth direction signal processed in the step four to obtain an SAR image.

Further, in the third step, the N layers of overlapping sub-apertures are divided along the azimuth direction, and the specific division process is as follows:

according to the length M of a single sub-aperture₁By Δ₂For data extraction ratio, extract M₂A sub-aperture forming a first layer of overlapping sub-apertures; by Delta₂Δ₃For data extraction ratio, extract M₃A first layer of overlapping sub-apertures forming a second layer of overlapping sub-apertures; by Delta₂Δ₃Δ₄For data extraction ratio, extract M₄A second layer of overlapping sub-apertures, constituting a third layer of overlapping sub-apertures … …, and so on until N layers of overlapping sub-apertures are obtained;

wherein M is₁……M_NIs a positive integer, and M₁……M_NThe value of (a) ensures that the secondary phase of each layer of overlapping sub-apertures is less than pi/4.

Further, in step three, the azimuth direction is divided into N layers of overlapping sub-apertures, and then an N + 1-dimensional azimuth direction signal is obtained, where the N + 1-dimensional azimuth direction signal is expressed as:

where i' is a range sequence after the range fast Fourier transform, m₁represents-M₁/2,…,M₁Sequence of/2, m₂represents-M₂/2,…,M₂The sequence of/2, and so on; x is the number of_pAnd y_pIs the coordinate of the imaging target, U₁And U₂The geometrical distortion, U, of the distance and orientation directions respectively caused by the wave front hypothesis₃Are azimuth quadratic phase coefficients caused by wave front hypothesis, which are only determined by scene definite value parameters and scene target point coordinates x_p,y_pDetermination of,. DELTA.K_yRepresenting the distance number bandwidth, Δ, after azimuth interpolation₂Representing the first data extraction ratio, Δ₂Δ₃Indicating a second data extraction ratio, … … Δ₂Δ₃…Δ_N+1Expressing the extraction proportion of the Nth data; j is an imaginary symbol, u_xIs the azimuthal wavenumber interval.

Further, the step four of solving the estimated value of the real coordinate according to the estimated value of the distorted coordinate is as follows:

solving the estimation value of the real coordinate in the form of analytic solution according to a higher-order signal estimation model, wherein the higher-order signal estimation model is expressed as:

where i and m represent fast and slow time sequences, respectively, j is an imaginary symbol, x_pAnd y_pIs the coordinate of the imaging target, U₁And U₂The geometrical distortion, U, of the distance and orientation directions respectively caused by the wave front hypothesis₃Is an azimuthal quadratic phase coefficient, K, caused by a wave front hypothesis_x(i, m) is the azimuth wave number, K_y(i) Distance is in wave number.

The analytical solution is expressed as:

wherein

Wherein N represents the nth dimension azimuth signal, and N belongs to [1, 2]；Y_cRepresents the distance, beta, of the projection of the radar on the ground to the center of the scene at the moment of the center of the synthetic aperture₀Showing the downward view angle of the synthetic aperture center moment, phi the downward-looking angle of the radar speed, delta the complementary angle of the ground oblique angle, H the flying height of the platform at the radar synthetic aperture center moment, y_QAs an estimate of the distance-wise distorted coordinates, x_QIs an estimate of the azimuthal distortion coordinate.

Further, the time-varying parameter system specifically includes: adjusting radar parameter carrier frequency and frequency in real time along the azimuth direction to counteract the distance wave number variation caused by the variation of downward viewing angle and slant ground declination;

the carrier frequency and the frequency modulation are represented as:

where c is the speed of light, f_c0And gamma₀Respectively representing the carrier frequency of the initial transmission signal and the modulation frequency, beta, of the linear frequency modulated LFM signal₀The downward viewing angle at the moment of the center of the synthetic aperture is shown, beta (t) represents the instantaneous downward viewing angle from the platform to the center of the scene, and alpha (t) represents the instantaneous included angle formed by the projection of the instantaneous reference slope distance on the ground plane and OY (Y axis).

Further, the pulse compression process includes the operations of de-line tone and removing the remaining video phase.

Has the advantages that:

(1) a multilayer overlapping subaperture imaging method of a large squint time-varying parameter SAR is characterized in that a multilayer overlapping subaperture imaging algorithm is provided by dividing the SAR into N layers of overlapping subapertures along the azimuth direction, so that the imaging width is improved; on the basis, the estimation value of the real coordinate is used for compensating the space-variant secondary phase error, so that the phase is compensated with higher precision, the problem that points far away from the center of the scene cannot be well focused is solved, and the imaging width is further improved.

(2) After the estimated value of the distorted coordinate is obtained, the estimated value of the real coordinate is solved in an analytic solution mode according to a higher-order signal estimation model, so that the compensation phase precision is higher, and the imaging width is further improved.

(3) The radar parameters are adjusted by adopting the time-varying parameter constitution, the spatial resolution performance is improved, the original distance interpolation operation is omitted, the imaging processing flow is simplified, and the operation efficiency is improved.

Drawings

Fig. 1 is a flow chart of a multilayer overlapping subaperture imaging method of a large squint time-varying parameter SAR according to the present invention.

Fig. 2 is a schematic diagram of the geometry of a large squint SAR.

Fig. 3 is a schematic diagram of the division of the azimuth N-layer overlapping sub-apertures.

Fig. 4 is a flow chart of sub-aperture compensation processing of the front N-dimensional azimuth signal.

FIG. 5 is a diagram illustrating simulation results for verifying superiority of the time-varying parameter system in an embodiment.

FIG. 6 is a schematic representation of the imaging widths of PFA, one layer of OSA and two layers of MOSA in the examples.

FIG. 7 is a diagram illustrating simulation results for verifying the superiority of the MOSA according to an embodiment.

Fig. 8 is a list of key parameters of a large squint SAR.

Detailed Description

The invention is described in detail below by way of example with reference to the accompanying drawings.

The invention discloses a multilayer overlapping subaperture imaging method of a large squint time-varying parameter SAR, which comprises the following specific steps as shown in figure 1:

step one, adopting a time-varying parameter system to transmit signals, receiving echo signals and performing pulse compression processing; performing pulse compression processing includes de-line tone Dechirp and removing remaining video phase RVP operations.

The geometry of the large squint SAR is shown in fig. 2. In the figure, a coordinate system XYZ is established by taking a scene central point as an origin O, XOY represents a ground plane, OZ, namely the direction of a Z axis deviates from the center of the earth and is vertical to the ground surface, and OY, namely the direction of a Y axis, is along the projection direction of the slant distance at the ground plane at the moment of the center of the synthetic aperture. The point P is any scattering point in the scene, the coordinate of the point P is (xp, yp,0), and the point A and the point S respectively represent the position of the platform at the central moment of the synthetic aperture and at the sampling moment of any azimuth direction. H denotes the platform flight height at the moment of the radar synthetic aperture center, t is the azimuth time (slow time), rc (t) denotes the instantaneous reference slope distance from the platform to the scene center, and rp (t) denotes the instantaneous slope distance from the platform to the scattering point P. α (t) represents the instantaneous angle formed by the projection of the instantaneous reference slope distance Rc (t) on the ground plane and OY. Beta (t) denotes the instantaneous downward view angle of the platform to the center of the scene, beta₀Representing the downward view at the moment of the center of the synthetic aperture. Phi denotes the angle of depression of the radar speed and delta denotes the complement of the ground squint angle.

Assuming that the chirp signal is a transmission signal, but unlike the traditional signal with fixed radar platform parameters, under the time-varying parameter system, the radar platform parameters (carrier frequency and modulation frequency) are related to azimuth sampling time and position, which are expressed as

Where t denotes slow time, f_c0And gamma₀Respectively representing the carrier frequency of the initial transmission signal and the modulation frequency of the LFM signal, f_cAnd gamma respectively represents the carrier frequency and the tuning frequency which change along with slow time, and alpha (t) represents the instantaneous included angle formed by the projection of the instantaneous reference slope distance Rc (t) on the ground plane and OY. Beta (t) denotes the instantaneous downward view angle of the platform to the center of the scene, beta₀Representing the downward view at the moment of the center of the synthetic aperture.

The transmitted signal may be represented as

Where j is an imaginary symbol, τ is the distance time (fast time), Tp is the pulse width, rect (-) is the window function, f_cAnd gamma denotes a carrier frequency and a modulation frequency, respectively, which vary with a slow time.

Demodulating the signal to baseband, the echo signal of any target point P can be expressed as

Where c is the speed of light and rp (t) represents the instantaneous slope distance of the platform to the scattering point P.

The signal is then Dechirp and RVP removed, and the processed signal is mapped from the analog domain to the digital domain, which may be represented as

Where i and m denote the fast time sequence and the slow time sequence(s), respectively₁(τ, t) represents the analog signal, s₂(i, m) represents a digital signal), Nr and Na are the number of sampling points in the distance direction and the azimuth direction, respectively, and R_c(m) and R_p(m) each represents R_c(t) and R_pNumber of (t)Field representation, F_sIs the sampling rate.

Based on the wavefront plane assumption and the geometry shown in FIG. 2, equation (4) is projected into the wavenumber domain and the wavenumber K is oriented in azimuth_xSum distance wavenumber K_yPerforming a two-dimensional taylor expansion to obtain a higher order signal model, expressed as:

wherein x is_pAnd y_pIs the coordinate of the imaging target, U₁And U₂The geometrical distortion, U, of the distance and orientation directions respectively caused by the wave front hypothesis₃Are azimuth quadratic phase coefficients caused by wave front hypothesis, which are only determined by scene definite value parameters and scene target point coordinates x_p,y_pDetermining; i and m represent fast time sequences and slow time sequences, respectively. Azimuth wave number K_x(i, m) and distance wavenumber K_y(i) As the following formula (6), the distance direction wave number only changes with i due to the adoption of a time-varying parameter system

Where α (m) is a digital domain representation of α (t).

And step two, only carrying out azimuth interpolation on the wave number domain signals obtained after the pulse pressure processing.

Because a time-varying parameter system is adopted, the echo signals expressed by the formula (5) are distributed in a trapezoidal shape in a wave number domain, distance interpolation correction is not needed, azimuth interpolation processing can be directly carried out, and the azimuth interpolation processing is realized through the following mapping

The signal after azimuth interpolation is expressed as

Wherein, K_x' (m) represents the azimuth wave number after azimuth interpolation, as shown in the following formula:

wherein u is_xIs the azimuthal wavenumber interval.

And step three, performing range compression on the signals processed in the step two, and dividing the signals into N layers of overlapped sub-apertures along the azimuth direction, wherein the N layers of overlapped sub-apertures correspond to N + 1-dimensional azimuth direction signals, and specifically the N layers of overlapped sub-apertures correspond to the N + 1-dimensional azimuth direction signals.

Divide into N layers of overlapping subaperture along the azimuth, the concrete division process is:

according to the length M of a single sub-aperture₁By Δ₂For data extraction ratio, extract M₂A sub-aperture forming a first layer of overlapping sub-apertures; by Delta₂Δ₃For data extraction ratio, extract M₃A first layer of overlapping sub-apertures forming a second layer of overlapping sub-apertures; by Delta₂Δ₃Δ₄For data extraction ratio, extract M₄A second layer of overlapping sub-apertures, constituting a third layer of overlapping sub-apertures … …, and so on until N layers of overlapping sub-apertures are obtained.

Wherein M is₁……M_NIs a positive integer, and M₁……M_NShould ensure that the secondary phase of each layer of overlapping sub-apertures is less than pi/4. As shown in equation (25), the number of sub-apertures M of the first layer₂Second order phase in third order phaseIs limited to within pi/4.

The interpolated signal is compressed in the range direction, i.e. FFT in the range direction, as shown in equation (8), to obtain

Where i' is the range sequence after the range FFT, Δ K_yRepresenting the distance number bandwidth after azimuth interpolation and obtaining the estimated value of the distance distortion coordinate

Then, the substrate is divided into N layers of overlapping sub-apertures along the azimuth direction, as shown in fig. 3. M₁Denotes the length of the first layer sub-aperture, Δ₂Representing the first data extraction ratio, M₂Denotes the number of subapertures of the first layer, Δ₂Δ₃Indicating the second data extraction ratio, and so on. m is₁represents-M₁/2,…,M₁Sequence of/2, m₂represents-M₂/2,…,M₂The sequence of/2, and so on. The azimuth sequence m can be written as

m＝m₁+Δ₂m₂+Δ₂Δ₃m₃+…+Δ₂…Δ_N+1m_N+1 (12)

By substituting equation (12) for equation (10), a signal divided into N overlapping subapertures in the azimuth direction can be obtained, which is expressed as:

step four, performing fast Fourier transform on the front N-dimensional azimuth signals one by one dimension to obtain an estimated value of each dimension of distortion coordinates; solving an estimated value of a real coordinate in an analytic solution mode according to a higher-order signal estimation model, and compensating a space-variant quadratic phase error; and performing fast Fourier transform on the (N + 1) th dimensional azimuth signal to obtain an estimated value of the (N + 1) th dimensional distortion coordinate.

The front N-dimensional azimuth signal is sub-aperture compensated, i.e. compensated for the space-variant Quadratic Phase Error (QPE) by using the estimated value.

The flowchart of step four is shown in fig. 4. First, to ensure m₁The dimension signal can be well focused, and the length M of the first layer sub-aperture is selected to be proper₁In the second phase of formula (13)Is limited to within pi/4. For the first dimension direction (m)₁Direction), the first dimension direction FFT-ed signal is represented as:

wherein m'₁Represents m₁The first estimate, i.e. the first rough estimate, of the first layer of the directional distorted coordinates from the directional FFT-ed sequence is:

where ρ is_x1Is the azimuthal first coarse resolution. According to distorted coordinatesTo true coordinatesThe mapping of (2) is derived from equation (17) to equation (24), and an analytic solution form is obtained, from equations (11) and (15), an estimated value of distance to real coordinates and a first rough estimate of azimuth to real coordinates are obtained, as follows:

wherein f is₁And f₂To solve the analytic solution expression of the real coordinates, the specific derivation is as follows:

firstly, the methodAndis as follows

Wherein the content of the first and second substances,Y_crepresenting the distance of the radar projection on the ground to the center of the scene at the moment of the center of the synthetic aperture,from the formulae (11) and (15), the following formula is obtained

In the general formula (17)Is substituted for the first expression in formula (18) to obtain

R_p0＝R_c0-y_Qsinβ₀ (19)

Substitution of the above formula (19) intoIn the expression of (1), the expression (18) becomes

Therefore, the temperature of the molten metal is controlled,can be expressed as

According to R_p0And (19) to obtain the following equation

Then, formula (21) is substituted for formula (22), resulting in the following equation

Then, the value can be obtained by the formula (23)As shown in the following formula (24), after that,can be obtained by equation (21). To this end, from distorted coordinatesTo true coordinatesThe mapping of (2) is deduced.

Thus, an expression of an analytical solution is obtained:

wherein

Wherein the content of the first and second substances,Y_crepresents the distance, beta, of the projection of the radar on the ground to the center of the scene at the moment of the center of the synthetic aperture₀Showing the downward view angle of the synthetic aperture center moment, phi the downward-looking angle of the radar speed, delta the complementary angle of the ground oblique angle, H the flying height of the platform at the radar synthetic aperture center moment, y_QAs an estimate of the distance-wise distorted coordinates, x_QIs an estimate of the azimuthal distortion coordinate.

The analytical solution expressions for different dimensions are:

wherein N represents an nth dimension azimuth signal, and N is less than or equal to N +1 (the form of an expression of each dimension azimuth signal analytical solution is the same, and an estimated value of an azimuth distortion coordinate changes along with the change of dimensions);

is obtained from equation (16)Performing phase compensation, writing the compensated signal as

To ensure m₂The dimension signal can be well focused, and the proper number M of the first layer sub-apertures is selected₂So that the second phase in the third phase of equation (25)Is limited to within pi/4. For the second dimension azimuth direction (m)₂Direction), the second dimension direction FFT-ed signal is represented as:

wherein m'₂Represents m₂The second rough estimate of the direction-distorted coordinates from the sequence after direction FFT is:

where ρ is_x2For the second coarse resolution of the azimuth direction, the estimated value of the distance-to-real coordinates and the second coarse estimation of the azimuth-to-real coordinates are obtained from equations (11) and (27), as shown below

Wherein the content of the first and second substances,is an estimate of the distance to the real coordinates (which may also be referred to as a second estimate of the distance to the real coordinates),is a second rough estimate of true coordinates of the bearing direction, f₃And f₄Solving the analytic solution expression of the real coordinate to obtainPerforming phase compensation, writing the compensated signal as

In the same way, carry out m₃…m_NDimension is processed by sub-aperture, and the signal after N times of compensation is written as

And step five, reconstructing the processed N + 1-dimensional azimuth signal to obtain an SAR image, namely, the reconstructed one-dimensional azimuth data of each range gate together to form an SAR two-dimensional image.

Fast Fourier Transform (FFT) is performed on the N + 1-dimensional azimuth signal:

the signal of the formula (30) is applied in the N + 1-dimensional azimuth direction (m)_N+1Direction) is subjected to FFT to obtain

At this time, the required one-dimensional azimuth data is stored in the N + 1-dimensional matrix as shown in equation (31). Therefore, vectorizing the FFT data, reconstructing one-dimensional azimuth data, and defining a new output azimuth variable as m':

where ρ is₁…ρ_xN+1The azimuth direction first to N +1 coarse resolutions.

And obtaining an SAR imaging image according to the reconstructed one-dimensional azimuth data azimuth signal.

The invention provides a multilayer overlapping subaperture imaging method of a time-varying parameter SAR. In the radar system level, aiming at the problem of large distortion of a target PSF of a constant parameter system, a target PSF is corrected by adopting time-varying adjustment radar parameters (time-varying parameters for short), a two-dimensional orthogonal PSF is realized under a large squint configuration, and the spatial resolution performance is improved. In the aspect of an imaging processing Algorithm, on the basis of a PFA (pulse-frequency domain Algorithm), in order to solve the problem that points far away from the center of a scene cannot be well focused, a higher-order signal estimation model is adopted, an overlapping Sub-Aperture Algorithm (OSA) is combined to compensate a space-variant phase error, and a Multi-layer overlapping Sub-Aperture imaging Algorithm (MOSA) is provided by changing the traditional orientation to multiple layers of one layer of Sub-Aperture OSA to N layers, so that the imaging width is further improved.

Examples

In order to verify the superiority of the time-varying parameter SAR system, the parameters in fig. 8 are used, and compared with the point target imaging performance under the conventional Back Projection Algorithm (BPA) and the variable parameter system BPA, the simulation result is shown in fig. 5. Fig. 5 (a) and (b) are plots of carrier frequency and modulation frequency, respectively, as a function of time; (c) and (d) are result diagrams of the simulation of the origin respectively, and it can be seen that under constant parameters, the azimuth direction and the range direction do not satisfy the orthogonality, and under a time-varying parameter system, the azimuth direction and the range direction satisfy the orthogonality, and the maximum resolution of the variable parameters is improved by 28%; (e) and (f) respectively representing the imaging results of two adjacent point targets, so that the imaging result of the conventional BP can hardly distinguish the two point targets, but the imaging result of the variable parameter BP can distinguish the two point targets, thereby proving the superiority of a time-varying parameter system and having better imaging resolution performance (improvement).

To verify the superiority of MOSA in terms of imaging width, taking MOSA of two overlapping sub-apertures in azimuth as an example, the simulation results are shown in fig. 6 and 7 compared with the one-layer OSA algorithm. In fig. 6, the areas within the black solid lines of (a), (b), and (c) represent the imaging widths of PFA, one OSA, and two MOSA, respectively, and it can be seen that the imaging widths of two MOSA layers are significantly improved compared to one OSA layer and are all significantly better than PFA. A point A (300, -100) within the imaging width of the first layer of OSA and the second layer of MOSA and a point B (290, -400) outside the imaging width of the first layer of OSA and within the depth of the second layer of MOSA are selected for point target simulation analysis, and the simulation result is shown in figure 7. In fig. 7, (a), (b) and (c) are the imaging results of the two algorithms at point a, because the imaging results are good within the imaging width of the two algorithms, and the azimuth peak sidelobe ratio can reach-13 dB; (e) and (f) and (g) are imaging results of two algorithms at the point B, the azimuth peak-to-side lobe ratio of the two-layer MOSA is-13.2 dB, the azimuth peak-to-side lobe ratio of the one-layer OSA is-10.7 dB, the imaging performance of the two-layer MOSA in the azimuth direction is obviously improved, and the superiority of the MOSA algorithm is proved.

The above embodiments only describe the design principle of the present invention, and the shapes and names of the components in the description may be different without limitation. Therefore, a person skilled in the art of the present invention can modify or substitute the technical solutions described in the foregoing embodiments; such modifications and substitutions do not depart from the spirit and scope of the present invention.