CN104237885B - A kind of diameter radar image orientation secondary focusing method - Google Patents

Abstract

The invention belongs to diameter radar image azimuth focus technical field, particularly to a kind of diameter radar image orientation secondary focusing method.The present invention comprises the following steps: step 1, obtains the SAR image after imaging processing；Step 2, is divided into the image block of multiple same size by SAR image；As i=1, perform step 3；Step 3, setting iterative parameter, as j=1, performs step 4；Step 4, draws the elevation of the impact point that the pixel given is corresponding, and the pixel given is the jth pixel of the i-th image block of N number of image block；Calculate the radar equivalent speed that the pixel given is corresponding；Construct the quadratic phase error penalty function that the pixel given is corresponding；The pixel given of the i-th image block of N number of image block is carried out phase compensation；Judging whether the phase compensation procedure of all pixels of all image blocks completes, if completed, then all image blocks after phase compensation being merged into the SAR image after orientation is to secondary focusing.

Description

Synthetic aperture radar image azimuth secondary focusing method

Technical Field

The invention belongs to the technical field of synthetic aperture radar image azimuth focusing, and particularly relates to a synthetic aperture radar image azimuth secondary focusing method. The method comprises the steps of firstly establishing a polynomial model among radar equivalent speed, target azimuth, distance and elevation positions, then utilizing the polynomial model and a target positioning result to quickly calculate the radar equivalent speed, constructing a Doppler domain equivalent speed space-variant error compensation function according to the polynomial model and the target positioning result, and carrying out azimuth secondary focusing on a traditional SAR imaging result.

Background

Synthetic Aperture Radars (SAR) are widely used in military and civil fields, such as battlefield investigation, ocean monitoring, agricultural census, topographic mapping, etc., due to their all-weather high-resolution earth observation capability throughout the day.

The conventional SAR imaging processing algorithm mainly includes time-domain Back Projection (BP), Range-Doppler (RD), linear Scaling (CS), and modified algorithms of the above algorithms. The BP algorithm can be regarded as a direct application of a Digital Beam Forming (DBF) technology, but the operation efficiency is poor. The latter processing algorithms can be regarded as fast realization algorithms of the BP algorithm under certain assumed conditions (such as polynomial models of the equivalent speed of the radar along the distance direction). For high-resolution SAR imaging processing, CS algorithms are widely applied, for example, an ecs (extended CS) processing algorithm is adopted for TanDEM-XTMSP, and the registration processing of the existing interferometric synthetic Aperture Radar (InSAR) mainly uses one image as a reference image, and registers the rest images with the reference image, and applies the registration result to each image.

With the improvement of the resolution of the SAR image (especially the resolution in decimeter order), the assumptions of a walking-stopping model, an orbit hyperbolic model, atmospheric delay non-space-variation and the like which are usually adopted by the traditional imaging processing are not applicable any more. Liuyan et al, in the article "Echo models and imaging algorithm for high-resolution SAR on high-speed platform" (IEEE Transactions on Geoscience and Remote Sensing,2012,50(3): 933-. The influence of the high-order model of the orbit on the imaging is analyzed in the article "Processing of ultra high-Resolution spectroscopy sliding guided light SAR data on curved restriction" (IEEE Transactions on Aerospace and electronic Systems,2013,49(2):819 and 839), and an imaging Processing algorithm based on the high-order model is given. It should be noted that the equivalent speed of the satellite-borne radar is space-variant with azimuth time, radar slope and ground elevation. However, the frequency domain imaging algorithm cannot adjust the doppler modulation frequency (i.e. the radar equivalent velocity) according to the target position along the azimuth direction, which causes the generated SAR image to have a spatially varying defocus effect.

Disclosure of Invention

The invention provides a synthetic aperture radar image azimuth secondary focusing method, aiming at the problem that the influence of ground elevation on equivalent speed cannot be considered in the traditional SAR imaging processing algorithm, so that partial region defocusing of a high-resolution SAR image is caused. The method comprises the steps of firstly establishing a polynomial model among radar equivalent speed, target azimuth, distance and elevation positions, then utilizing the polynomial model and a target positioning result to quickly calculate the radar equivalent speed, constructing a Doppler domain equivalent speed space-variant error compensation function according to the polynomial model and the target positioning result, and carrying out azimuth secondary focusing on a traditional SAR imaging result.

In order to achieve the technical purpose, the invention is realized by adopting the following technical scheme.

A synthetic aperture radar image azimuth secondary focusing method comprises the following steps:

step 1, acquiring synthetic aperture radar echo data, carrying out SAR imaging processing on the synthetic aperture radar echo data to obtain an SAR image, and acquiring a prior digital elevation model capable of covering an imaging area;

step 2, dividing the SAR image into N image blocks with the same size, wherein the number of pixel points of each image block is M, N is a natural number greater than 1, and M is a natural number greater than 1; setting an iteration parameter i equal to 1, 2., and when i equal to 1, executing a step 3;

step 3, setting an iteration parameter j equal to 1,2, and executing step 4 when j equal to 1;

step 4, carrying out SAR image target positioning aiming at given pixel points to obtain the elevations of target points corresponding to the given pixel points, wherein the given pixel points are the jth pixel points of the ith image blocks of the N image blocks;

step 5, calculating the radar equivalent speed corresponding to the given pixel point by using the elevation of the target point corresponding to the given pixel point;

step 6, constructing a secondary phase error compensation function corresponding to the given pixel point by using the radar equivalent speed corresponding to the given pixel point;

step 7, performing azimuth Fourier transform on the ith image block of the N image blocks to obtain a distance Doppler domain signal corresponding to the ith image block of the N image blocks; performing phase compensation on the range-Doppler domain signal corresponding to the ith image block of the N image blocks by using the secondary phase error compensation function corresponding to the pixel point given in the step 6 to obtain a signal subjected to phase compensation; performing azimuth inverse Fourier transform on the signal subjected to the phase compensation to obtain an ith image block subjected to the phase compensation; setting the ith image block of the N image blocks as the ith image block subjected to phase compensation processing;

step 8, judging whether i is equal to N, if i is not equal to N, judging whether j is equal to M, if j is equal to M, enabling the value of i to be increased by 1, then returning to the step 3, and if j is not equal to M, enabling the value of j to be increased by 1, then returning to the step 4; and if j is equal to M, the value of j is increased by 1, and then the step 4 is returned, and if j is equal to M, the obtained 1 st image block to the Nth image block which are subjected to phase compensation processing are combined into the SAR image subjected to azimuth secondary focusing.

The invention is characterized by further improvement:

the specific substeps of the step 4 are as follows:

(4.1) setting the initial reference elevation value of a target point corresponding to the given pixel point as h₀The given pixel point is the jth pixel point of the ith image block of the N image blocks; let elevation variable h be h₀；

(4.2) according to the following SAR geometric positioning equation set, iteratively solving the coordinates of a target point corresponding to the given pixel point under the geocentric fixed coordinate system:

$\left\{\begin{array}{c}\frac{{2v}_m\left(t_m\right)\·[p_t-p_m\left(t_m\right)]}{{\mathrm{\λr}}_1}=0\\ |p_r-p_m\left(t_m\right)|=r_1\\ \frac{p_{t,x}^2+p_{t,y}^2}{{(R_e+h)}^2}+\frac{p_{t,z}^2}{{(1-f)}^2{(R_e+h)}^2}=1\end{array}\right.$

wherein, t_mIs the azimuth time, v_m(t_m) For the satellite velocity vector, λ is the wave of the synthetic aperture radar transmitted signalLength, r₁For synthetic aperture radar slant range, p_tRepresenting the coordinates of the target point corresponding to the given pixel point in the earth-centered-fixed coordinate system, p_t＝(p_t,x,p_t,y,p_t,z)^TThe superscript "T" being the transpose of a matrix or vector, p_t,x、p_t,yAnd p_t,zX-axis coordinate, Y-axis coordinate and Z-axis coordinate, p, of target point corresponding to given pixel point under earth center fixed coordinate system_m(t_m) For the phase center position, | p, of the receiving antenna of the synthetic aperture radar_t-p_m(t_m) I denotes the vector p_t-p_m(t_m) Of a mold of R_eEquator radius, f is the earth oblateness factor;

(4.3) converting coordinates of a target point corresponding to the given pixel point under the geocentric fixed connection coordinate system to a coordinate system where the prior digital elevation model is located;

(4.4) interpolating the elevation value h of the target point corresponding to the given pixel point according to the coordinate of the target point corresponding to the given pixel point in the coordinate system of the prior digital elevation model₁；

(4.5) setting a threshold h for setting the end of iteration_thres；

(4.6) calculating the elevation value h of the target point corresponding to the given pixel point obtained after interpolation₁Initial reference elevation value h of target point corresponding to given pixel point₀If h is different from the value of₁And h₀Satisfies | h₁-h₀|≤h_thresIf the elevation of the target point corresponding to the given pixel point is h, then executing the step 5; otherwise, let h be h₁Then, substep (4.2) is performed.

In step 5, the corresponding radar equivalent velocity v of the given pixel point is obtained according to the following formula_e：

$\begin{array}{c}{v}_e=v_e(t_{m,\mathrm{ref}},R_{\mathrm{ref}},h_{\mathrm{ref}})+k_{t_a}\·(t_m-t_{m,\mathrm{ref}})+k_R\·(R-R_{\mathrm{ref}})+\\ k_{R^2}\·{(R-R_{\mathrm{ref}})}^2+k_h\·(h-h_{\mathrm{ref}})+k_{h,R}\·(h-h_{\mathrm{ref}})\·(R-R_{\mathrm{ref}})\end{array}$

Wherein, t_m,refFor reference azimuth time, R_refFor radar reference of slope, h_refTo a reference elevation, v_e(t_m,ref,R_ref,h_ref) For a known radar equivalent velocity corresponding to a reference pixel point,k_R、k_hand k_h,RH is the elevation of a target point corresponding to the given pixel point for the set five coefficients; t is t_mR is the radar slant distance of a target point corresponding to the given pixel point at the azimuth moment, $R=\sqrt{R_{\mathrm{ref}}^2+{v_e}^2(t_{m,\mathrm{ref}},R_{\mathrm{ref}},h_{\mathrm{ref}})*{(t_m-t_{m,\mathrm{ref}})}^2}.$

in step 6, a given quadratic phase error compensation function corresponding to the pixel point is constructed by the following formula

Wherein R is the radar slant distance of a target point corresponding to a given pixel point, lambda is the wavelength of a synthetic aperture radar transmission signal, v_gRadar equivalent velocity, Δ v, used in SAR imaging processing for step 1_e＝v_e-v_g，f_dIs the doppler frequency.

The specific substeps of the step 7 are as follows:

(7.1) performing azimuth Fourier transform on the ith image block of the N image blocks to obtain a distance Doppler domain signal S (f) corresponding to the ith image block of the N image blocks_d)_i；

(7.2) using the quadratic phase error compensation function corresponding to the pixel point given in the step 6 to the range-Doppler domain signal S (f) corresponding to the ith image block of the N image blocks in the sub-step (7.1)_d)_iPerforming phase compensation to obtain a phase-compensated signal S_c(f_d)_i，S_c(f_d)_iComprises the following steps:

wherein,constructing a secondary phase error compensation function corresponding to the given pixel point constructed in the step 6;

(7.3) the phase compensated signal S obtained in substep (7.2)_c(f_d)_iPerforming azimuth inverse Fourier transform to obtain an ith image block subjected to phase compensation; and setting the ith image block of the N image blocks as the ith image block subjected to phase compensation processing.

The invention has the beneficial effects that: 1) the invention provides a method for azimuth secondary focusing processing, which effectively solves the problem that the influence of ground elevation on equivalent speed cannot be considered in the traditional SAR imaging processing algorithm, so that partial area defocusing of a high-resolution SAR image is caused; 2) the invention establishes a polynomial model of the equivalent speed of the satellite-borne SAR system, which is spatially-varied along with the azimuth, the distance and the target elevation, and can quickly and effectively calculate the equivalent speed spatially-varied error compensation quantity of the traditional PTA algorithm by utilizing the model, thereby reducing the computation quantity of the traditional PTA algorithm.

Drawings

FIG. 1 is a flow chart of a synthetic aperture radar image azimuth secondary focusing method of the present invention;

FIG. 2 is a flow chart of SAR image positioning in the present invention;

FIG. 3 is a schematic view of digital elevation of an observation scene in a simulation experiment;

FIG. 4a is a schematic diagram of an imaging processing result obtained by three methods for a point target 1 in a simulation experiment;

FIG. 4b is a schematic diagram of an imaging result obtained by three methods for the point target 2 in the simulation experiment;

fig. 4c is a schematic diagram of an imaging processing result obtained by three methods for the point target 3 in the simulation experiment.

Detailed Description

The invention will be further described with reference to the accompanying drawings in which:

referring to fig. 1, a flowchart of a synthetic aperture radar image azimuth secondary focusing method of the present invention is shown. The synthetic aperture radar image azimuth secondary focusing method comprises the following steps:

step 1, acquiring synthetic aperture radar echo data, carrying out SAR imaging processing on the synthetic aperture radar echo data to obtain an SAR image, and acquiring a prior digital elevation model capable of covering an imaging area, specifically, in step 1, acquiring auxiliary parameters (including pulse repetition frequency PRF, radar carrier frequency fc, azimuth signal bandwidth Ba, signal bandwidth Br, distance sampling rate Fr) of the SAR imaging processing while acquiring the SAR image

Step 2, dividing the SAR image into N image blocks with the same size, wherein N is a natural number greater than 1, in the embodiment of the invention, the more pixel points of each divided image block are, the higher the compensation precision is, but the operation complexity is increased, so that the pixel points of each divided image block need to be selected in a compromise mode between precision and complexity, and in the embodiment of the invention, the number of the pixel points of each divided image block is 32 or 64.

An iteration parameter i is set to 1, 2.

And 3, setting an iteration parameter j equal to 1, 2.

And 4, carrying out SAR image target positioning on the given pixel point to obtain the elevation h of a target point corresponding to the given pixel point, wherein the given pixel point is the jth pixel point of the ith image block of the N image blocks.

Fig. 2 is a flow chart of positioning the SAR image in the present invention. In the embodiment of the invention, the SAR image target positioning is carried out on the given pixel points by adopting a positioning method assisted by a prior DEM (ground elevation information (DEM) of the positioning method is obtained through an SRTM DEM, the accuracy of the prior DEM is not high, so that the ground actual elevation closest to the prior DEM is solved by adopting an iterative method, namely the actual position and the actual elevation of the corresponding pixel points are solved). The method comprises the following specific substeps:

(4.1) setting the initial reference elevation value of a target point corresponding to the given pixel point as h₀The given pixel point is the jth pixel point of the ith image block of the N image blocks; let elevation variable h be h₀。

(4.2) according to the following SAR geometric positioning equation set, iteratively solving the coordinates of a target point corresponding to the given pixel point under the geocentric fixed coordinate system:

$\left\{\begin{array}{c}\frac{{2v}_m\left(t_m\right)\·[p_t-p_m\left(t_m\right)]}{{\mathrm{\λr}}_1}=0---\left(a\right)\\ |p_r-p_m\left(t_m\right)|=r_1---\left(b\right)\\ \frac{p_{t,x}^2+p_{t,y}^2}{{(R_e+h)}^2}+\frac{p_{t,z}^2}{{(1-f)}^2{(R_e+h)}^2}=1---\left(c\right)\end{array}\right.$

wherein, t_mIs the azimuth time, v_m(t_m) For satellite velocity vectors, λ is the synthetic aperture radarUp to the wavelength of the transmitted signal, r₁For synthetic aperture radar slant range, p_tRepresenting the coordinates of the target point corresponding to the given pixel point in the earth-centered-fixed coordinate system, p_t＝(p_t,x,p_t,y,p_t,z)^TThe superscript "T" being the transpose of a matrix or vector, p_t,x、p_t,yAnd p_t,zX-axis coordinate, Y-axis coordinate and Z-axis coordinate, p, of target point corresponding to given pixel point under earth center fixed coordinate system_m(t_m) For the phase center position, | p, of the receiving antenna of the synthetic aperture radar_t-p_m(t_m) I denotes the vector p_t-p_m(t_m) Of a mold of R_eThe radius of the equator is, f is an earth oblateness factor, the equation (a) in the equation set is an SAR image imaging Doppler geometric equation, the equation (b) in the equation set is an SAR image imaging slant distance equation, and the equation (c) in the equation set is an SAR image imaging earth ellipsoid equation.

And (4.3) converting coordinates of a target point corresponding to the given pixel point under the geocentric fixed coordinate system to a coordinate system (geodetic coordinate system) of the prior digital elevation model.

(4.4) interpolating the elevation value h of the target point corresponding to the given pixel point according to the coordinate of the target point corresponding to the given pixel point in the coordinate system of the prior digital elevation model₁。

(4.5) setting a threshold h for setting the end of iteration_thres。

(4.6) calculating the elevation value h of the target point corresponding to the given pixel point obtained after interpolation₁Initial reference elevation value h of target point corresponding to given pixel point₀If h is different from the value of₁And h₀Satisfies | h₁-h₀|≤h_thresIf the elevation of the target point corresponding to the given pixel point is h, then executing the step 5; otherwise, let h be h₁Then, substep (4.2) is performed.

Step 5, calculating the equivalent speed v of the radar corresponding to the given pixel point_e。

The method comprises the following specific substeps:

in the embodiment of the invention, the radar equivalent speed v corresponding to the given pixel point is obtained according to the following formula_e：

$\begin{array}{c}{v}_e=v_e(t_{m,\mathrm{ref}},R_{\mathrm{ref}},h_{\mathrm{ref}})+k_{t_a}\·(t_m-t_{m,\mathrm{ref}})+k_R\·(R-R_{\mathrm{ref}})+\\ k_{R^2}\·{(R-R_{\mathrm{ref}})}^2+k_h\·(h-h_{\mathrm{ref}})+k_{h,R}\·(h-h_{\mathrm{ref}})\·(R-R_{\mathrm{ref}})\end{array}$

Wherein, t_m,refFor reference azimuth time, R_refFor radar reference of slope, h_refTo a reference elevation, v_e(t_m,ref,R_ref,h_ref) For a known radar equivalent velocity corresponding to a reference pixel point,k_R、k_hand k_h,RH is the elevation of a target point corresponding to the given pixel point for the set five coefficients; t is t_mR is the radar slant distance of a target point corresponding to the given pixel point at the azimuth moment, $R=\sqrt{R_{\mathrm{ref}}^2+{v_e}^2(t_{m,\mathrm{ref}},R_{\mathrm{ref}},h_{\mathrm{ref}})*{(t_m-t_{m,\mathrm{ref}})}^2}.$

step 6, constructing a secondary phase error compensation function corresponding to the given pixel point

Specifically, construction of the given pixel point correspondence quadraticPhase error compensation function

Where R is the radar slant distance of the target point corresponding to the given pixel point, λ is the wavelength of the synthetic aperture radar transmitted signal, v_gRadar equivalent velocity, Δ v, used in SAR imaging processing for step 1_e＝v_e-v_g，f_dIs the doppler frequency.

And 7, performing secondary focusing processing (phase compensation processing) on the given pixel points.

The method comprises the following specific substeps:

(7.1) performing azimuth Fourier transform (FFT) on the ith image block of the N image blocks to a range Doppler domain to obtain a range Doppler domain signal S (f) corresponding to the ith image block of the N image blocks_d)_i；

(7.2) utilizing the quadratic phase error compensation function corresponding to the pixel point given in step 6The distance Doppler domain signal S (f) corresponding to the ith image block of the N image blocks in the sub-step (7.1)_d)_iPerforming phase compensation to the phase compensated signal S_c(f_d)_iComprises the following steps:

(7.3) the phase compensated signal S obtained in substep (7.2)_c(f_d)_iPerforming inverse Fourier transform (IFFT) to complete the height of the given pixel pointAnd (5) compensating the space-variant error to obtain the ith image block subjected to phase compensation. And setting the ith image block of the N image blocks as the ith image block subjected to phase compensation processing.

Step 8, judging whether i is equal to N, wherein N is the number of image blocks divided by the SAR image in the step 2, if i is not equal to N, judging whether j is equal to M, wherein M is the number of pixel points of each divided image block, if j is equal to M, the value of i is increased by 1, then returning to the step 3, and if j is not equal to M, the value of j is increased by 1, and then returning to the step 4; and if j is equal to M, ending the process of azimuth secondary focusing of the synthetic aperture radar image, and combining the obtained 1 st image block to the Nth image block which are subjected to phase compensation processing into the SAR image subjected to azimuth secondary focusing.

In the embodiment of the invention, after the SAR image subjected to azimuth secondary focusing is obtained, geocoding is carried out on the SAR image subjected to azimuth secondary focusing, and the geocoded SAR image is obtained. In the geocoding process of the SAR image, processing steps such as geometric positioning, map projection, resampling correction and the like are completed on the basis of imaging processing, and finally a map image for people is formed.

After the SAR image subjected to geocoding is obtained, geometric correction is carried out on the SAR image subjected to geocoding, and the SAR image subjected to geometric correction is obtained. After the geometric correction is completed, the further corrected SAR image can be obtained by geometrically distorting the image (such as shading, overlapping and the like) through a GCP point (group control Points) and a priori Digital Elevation Model (DEM).

The effect of the present invention will be further described with reference to simulation experiments.

1) Simulation conditions are as follows:

the effectiveness of the verification algorithm is analyzed through simulation experiments, and the processing performance of the imaging processing method is compared with that of the existing imaging processing method. The simulation experiment input conditions are as follows: the satellite-borne SAR system with the resolution of 0.4m is adopted for simulation, and the specific parameters are as follows:

it should be noted that, in the table, the PRF is an equivalent pulse repetition frequency of the azimuth multi-channel SAR system or the beaming mode SAR system after doppler deblurring processing, and the simulated antenna pattern is rectangular.

The satellite orbit number is as follows:

referring to fig. 3, a schematic diagram of digital elevation of an observation scene in a simulation experiment is shown. In the simulation experiment, the observation scene is selected from the southwest mountain area of Xian, and the ground elevation information of the observation scene is obtained through an SRTM DEM (space shuttle radar topographic mapping task data elevation model). In fig. 3, the cross identifier represents a point target for simulation. The horizontal axis represents longitude in degrees, the vertical axis represents latitude in degrees, the color depth of the bar in fig. 3 represents the height of the location point, and the darker the color of the bar, the lower the height of the corresponding location point.

2) And (3) simulation result analysis:

after the simulated echo data are obtained, the imaging processing is respectively carried out by adopting a traditional CS algorithm (compressed sensing algorithm), an original PTA (accurate terrain and Aperture change) compensation algorithm and the invention (improved PTA compensation algorithm based on an equivalent velocity model). Three point targets located in the same range bin are selected and analyzed, and the ground elevations of the three point targets are 917.5m (corresponding to the point target 1 in fig. 3), 647.1m (corresponding to the point target 2 in fig. 3) and 776.3m (corresponding to the point target 3 in fig. 3). In the conventional CS algorithm, the radar equivalent velocity of the point target 2 is employed. Referring to fig. 4a, a schematic diagram of an imaging processing result obtained by three methods for the point target 1 in the simulation experiment is shown. Referring to fig. 4b, a schematic diagram of an imaging processing result obtained by three methods for the point target 2 in the simulation experiment is shown. Referring to fig. 4c, a schematic diagram of an imaging processing result obtained by three methods for the point target 3 in the simulation experiment is shown. In fig. 4a to 4c, the horizontal axis represents the azimuth distance of the target point in meters (m), and the vertical axis represents the intensity in decibels (dB); the circles represent the conventional CS algorithm (corresponding to "uncompensated" in the legend), the lines represent the raw PTA compensation algorithm (corresponding to "raw PTA" in the legend), and the asterisks represent the present invention (corresponding to "modified PTA" in the legend). As can be seen from FIGS. 4a to 4c, the elevation space-variant will cause the problems of defocusing and side lobe elevation of the image in the conventional CS imaging algorithm, wherein the azimuth resolution of each point target is 0.41m, 0.38m and 0.39m, and the peak-to-side lobe ratio is-10.06 dB, -13.27dB and-11.58 dB. The invention can effectively correct elevation space-variant errors, the processing performance of the invention is equivalent to that of the original PTA algorithm, the azimuth resolution of each point target is 0.38m, and the peak side lobe ratio is-13.27 dB.

In conclusion, the method can quickly and effectively calculate the equivalent speed space-variant error compensation quantity of the traditional PTA algorithm and reduce the calculation quantity of the traditional PTA algorithm. The effectiveness of the invention is verified by simulation experiment results.

It will be apparent to those skilled in the art that various changes and modifications may be made in the present invention without departing from the spirit and scope of the invention. Thus, if such modifications and variations of the present invention fall within the scope of the claims of the present invention and their equivalents, the present invention is also intended to include such modifications and variations.

Claims (4)

1. A synthetic aperture radar image azimuth secondary focusing method is characterized by comprising the following steps:

step 1, acquiring synthetic aperture radar echo data, carrying out SAR imaging processing on the synthetic aperture radar echo data to obtain an SAR image, and acquiring a prior digital elevation model capable of covering an imaging area;

step 2, dividing the SAR image into N image blocks with the same size, wherein the number of pixel points of each image block is M, N is a natural number greater than 1, and M is a natural number greater than 1; setting an iteration parameter i equal to 1, 2., and when i equal to 1, executing a step 3;

step 3, setting an iteration parameter j equal to 1,2, and executing step 4 when j equal to 1;

step 4, carrying out SAR image target positioning aiming at given pixel points to obtain the elevations of target points corresponding to the given pixel points, wherein the given pixel points are the jth pixel points of the ith image blocks of the N image blocks;

the specific substeps of the step 4 are as follows:

(4.1) setting the initial reference elevation value of a target point corresponding to the given pixel point as h₀The given pixel point is the jth pixel point of the ith image block of the N image blocks; let elevation variable h be h₀；

(4.2) according to the following SAR geometric positioning equation set, iteratively solving the coordinates of a target point corresponding to the given pixel point under the geocentric fixed coordinate system:

$\left\{\begin{array}{c}\frac{2v_m\left(t_m\right)\·\[p_t-p_m\left(t_m\right)\]}{{\mathrm{\λr}}_1}=0\\ |p_t-p_m\left(t_m\right)|=r_1\\ \frac{p_{t,x}^2+p_{t,y}^2}{{(R_e+h)}^2}+\frac{p_{t,z}^2}{{(1-f)}^2{(R_e+h)}^2}=1\end{array}\right.$

wherein, t_mIs the azimuth time, v_m(t_m) Is the satellite velocity vector, lambda is the wavelength of the synthetic aperture radar transmitted signal, r₁For synthetic aperture radar slant range, p_tRepresenting the coordinates of the target point corresponding to the given pixel point in the earth-centered-fixed coordinate system, p_t＝(p_t，x，p_t，y，p_t，z)¹The superscript "T" being the transpose of a matrix or vector, p_t，x、p_t，yAnd p_t，zX-axis coordinate, Y-axis coordinate and Z-axis coordinate, p, of target point corresponding to given pixel point under earth center fixed coordinate system_m(t_m) For the phase center position, | p, of the receiving antenna of the synthetic aperture radar_t-p_m(t_m) I denotes the vector p_t-p_m(t_m) Of a mold of R_eEquator radius, f is the earth oblateness factor;

(4.3) converting coordinates of a target point corresponding to the given pixel point under the geocentric fixed connection coordinate system to a coordinate system where the prior digital elevation model is located;

(4.4) target points corresponding to the given pixel pointsInterpolating the coordinates in the coordinate system of the prior digital elevation model to obtain the elevation value h of the target point corresponding to the given pixel point₁；

(4.5) setting a threshold h for ending the iteration_thres；

(4.6) calculating the elevation value h of the target point corresponding to the given pixel point obtained after interpolation₁Initial reference elevation value h of target point corresponding to given pixel point₀If h is different from the value of₁And h₀Satisfies | h₁-h₀|≤h_thresIf the elevation of the target point corresponding to the given pixel point is h, then executing the step 5; otherwise, let h be h₁Then performing substep (4.2);

step 5, calculating the radar equivalent speed corresponding to the given pixel point by using the elevation of the target point corresponding to the given pixel point;

step 6, constructing a secondary phase error compensation function corresponding to the given pixel point by using the radar equivalent speed corresponding to the given pixel point;

step 7, performing azimuth Fourier transform on the ith image block of the N image blocks to obtain a distance Doppler domain signal corresponding to the ith image block of the N image blocks; performing phase compensation on the range-Doppler domain signal corresponding to the ith image block of the N image blocks by using the secondary phase error compensation function corresponding to the pixel point given in the step 6 to obtain a signal subjected to phase compensation; performing azimuth inverse Fourier transform on the signal subjected to the phase compensation to obtain an ith image block subjected to the phase compensation; setting the ith image block of the N image blocks as the ith image block subjected to phase compensation processing;

step 8, judging whether i is equal to N, if i is not equal to N, judging whether j is equal to M, if j is equal to M, enabling the value of i to be increased by 1, then returning to the step 3, and if j is not equal to M, enabling the value of j to be increased by 1, then returning to the step 4; and if j is equal to M, the value of j is increased by 1, and then the step 4 is returned, and if j is equal to M, the obtained 1 st image block to the Nth image block which are subjected to phase compensation processing are combined into the SAR image subjected to azimuth secondary focusing.

2. The synthetic aperture radar image azimuth secondary focusing method according to claim 1, wherein in step 5, the radar equivalent velocity v corresponding to the given pixel point is obtained according to the following formula_e：

$\begin{array}{c}{v}_e=v_e(t_{m,\mathrm{ref}},R_{\mathrm{ref}},h_{\mathrm{ref}})+k_{t_a}\·(t_m-t_{m,\mathrm{ref}})+k_R\·(R-R_{\mathrm{ref}})+\\ k_{R^2}\·{(R-R_{\mathrm{ref}})}^2+k_h\·(h-h_{\mathrm{ref}})+k_{h,R}\·(h-h_{\mathrm{ref}})\·(R-R_{\mathrm{ref}})\end{array}$

Wherein, t_m，refFor reference azimuth time, R_refFor radar reference of slope, h_refTo a reference elevation, v_e(t_m，ref，R_ref，h_ref) For a known radar equivalent velocity corresponding to a reference pixel point,k_R、k_hand k_h，RH is the elevation of a target point corresponding to the given pixel point for the set five coefficients; t is t_mR is the radar slant distance of a target point corresponding to the given pixel point at the azimuth moment,

3. the synthetic aperture radar image azimuth secondary focusing method according to claim 1, wherein in step 6, the secondary phase error compensation function corresponding to the given pixel point is constructed by the following formula

Wherein R is the radar slant distance of a target point corresponding to a given pixel point, lambda is the wavelength of a synthetic aperture radar transmission signal, v_gRadar equivalent velocity, v, used in SAR imaging processing for step 1_eFor radar equivalent velocity, Δ v, corresponding to given pixel points_e＝v_e-v_g，f_dIs the doppler frequency.

4. The synthetic aperture radar image azimuth secondary focusing method according to claim 1, wherein the specific sub-steps of the step 7 are as follows:

(7.1) performing azimuth Fourier transform on the ith image block of the N image blocks to obtain a distance Doppler domain signal S (f) corresponding to the ith image block of the N image blocks_d)_i；

(7.2) using the quadratic phase error compensation function corresponding to the pixel point given in the step 6 to the range-Doppler domain signal S (f) corresponding to the ith image block of the N image blocks in the sub-step (7.1)_d)_iPerforming phase compensation to obtain a phase-compensated signal S_c(f_d)_i，S_c(f_d)_iComprises the following steps:

wherein,constructing a secondary phase error compensation function corresponding to the given pixel point constructed in the step 6;

(7.3) the phase compensated signal S obtained in substep (7.2)_c(f_d)_iPerforming azimuth inverse Fourier transform to obtain an ith image block subjected to phase compensation; and setting the ith image block of the N image blocks as the ith image block subjected to phase compensation processing.