Improvement of Stripmap SAR Autofocus Based on Mini- mum-Entropy Criterion

Autofocus is an essential part of the SAR imaging process. Multi-subaperture autofocus algorithm is a commonly used autofocus algorithm for processing SAR stripmap mode data. The multi-subaperture autofocus algorithm has two main steps, the first is to estimate the phase error gradient within the subaperture, the second is to splice the phase error gradient, that is, to remove the shift amount between the estimated adjacent subapertures’ error gradients. Previous gradient-splicing algorithms assume that the estimation of subaperture error is accurate, but when the estimation of subaperture phase error gradients is not accurate enough, these algorithm performance will be degraded. A new phase error gradient splicing algorithm is proposed in this paper. It roughly estimates the shift amount first, and then finely estimates the shift amount based on the minimum-entropy criterion, which can improve the robustness of splicing especially when the estimation of the phase error gradients of the subaperture is not accurate enough. To speed up the algorithm, a variable-step-size search method is used. Simulation and experimental results show that the algorithm has enough accuracy and still has good performance when other splicing algorithms doesn’t perform well.


Introduction
Synthetic aperture radar (SAR) is a kind of high-resolution imaging radar which is not limited by diurnal variation, clouds and other factors [1]. However, airborne SAR imaging has high requirements for the flight path of aircraft. If a SAR is going to record data, the aircraft carrying the radar needs to fly relatively smoothly. Due to the disturbance of air flow, it is impossible for the aircraft to fly in a straight line at a constant speed. So after acquiring SAR data, compensating the data is necessary. The original compensation only considers the motion error recorded by the inertial navigation system on the platform, but the compensation effect depends on the accuracy of the inertial navigation system.
Autofocus is an essential step in synthetic aperture radar data processing, which can compensate for uncorrected range cell migration and other errors. Multi-subaperture autofocus is widely used in strip mode SAR autofocus. Firstly, the SAR data is divided into some subapertures, and the phase error (including its gradients) in each subaperture is estimated. Then, gradients of the subaperture phase error (SPE) should be spliced together to obtain the phase error gradient of the whole aperture.
When using multi-subaperture autofocus methods to compensate the error of the strip mode SAR, some autofocus algorithms for spotlight mode SAR can be used to estimated the SPE. Map Drift (MD) [2] and Phase Gradient Autofocus (PGA) [3] are two classical autofocus algorithms which can be used in SAR spotlight mode and also can be used to estimated the subapertures' phase error. Based on these two classical algorithms, many other autofocus algorithms are derived, such as Multiple-Aperture Map Drift (MAM) [4], weighted PGA (WPGA) [5], Quality Phase-Gradient (QPGA) [6] and so on.
At this point, the image needs to be divided into several subapertures. There may be overlaps between adjacent subapertures. 0 i t is the center of azimuth time of the ith subaperture. In the azimuth dimension, multiply the ith subaperture by a chirp signal: where sq  is the squint angel, V is the speed of the platform carrying the radar,  is the wavelength of the radar carrier wave, 0 R is the distance from the range cells to the radar, 0 ( ) a K R is the Doppler frequency rate. By IFFT the subaperture's data which are in 2-D time domain, a sub image that only contains the irradiated targets within the ith subaperture will be obtained [11]. By using PGA or other algorithms for this sub image, the SPE and their gradients in the azimuth time domain of ith subaperture will be obtained. Fig.1 is a schematic diagram of the error gradients of each subaperture obtained by PGA in one simulation, and the error gradients of the adjacent subapertures on the overlaps have a significant shift.

Combine the SPE gradient
i L is the number of azimuthal pulses of the ith subaperture, l is the number of azimuthal pulses of the overlapping part of the adjacent subapertures. Generally speaking, the length of each subaperture is equal, i.e.
represent the actual shift amount of the ith subaperture's phase error gradients and (i-1)th subaperture's phase error gradients at this time. Therefore, before splicing the error gradients to obtain the whole aperture's errors, When using PCA algorithm, differences in {} i nL  is used to get the curvature. It assumes that in overlaps of adjacent subapertures, although the phase error gradients are different, the curvatures is the same. In this algorithm, can actually be estimated as follows: where 0 jl  , and j is the ordinal of the splice point.
When using the PGA-LS algorithm, in the overlap, the phase error gradients of the previous subaperture are subtracted from the phase error gradients of the latter subaperture, then the differences are averaged to estimate 1, ii −  , which can be expressed as: In the PGA-MD algorithm,

1,
ii −  is estimated using the relative shift of the targets on images formed by the adjacent autofocused subapertures [9], which can be expressed as: In (5), i  is the relative shift of shown targets in the sub images formed by the (i-1)th subaperture and ith subaperture, the second term is the inherent shift of the target. (5) is a simplified expression of the algorithm, the original algorithm uses the least-square criterion to estimate the value.
Whether PCA or PGA-LS, the estimator of isn't accurate enough [9]. On the premise that the errors of subaperture has been accurately estimated, the subaperture can generate a clear sub image, so can be accurately estimated by MD algorithm. However, MD algorithm isn't robust enough. When the platform motion error is too large or the scene lack target with strong reflection, these sub images is hard to become clear. In this case, MD algorithm isn't the best gradient-splicing algorithm to estimate

Minimum-Entropy criterion
Minimum-Entropy criterion is derived from the theory of Image entropy. In a MN  gray image, ( , ) g m n represent the gray level of the point whose coordinate is (m,n). There are several definitions of image entropy [8], and the expression of image entropy in this paper is as follows: Assume that the total energy of the image (i.e. C) is constant, it can be inferred from (6) and (7) that the more average the energy distribution of the image, the greater the entropy of the image. For a SAR image, the worse the focusing effect, the more average the energy distribution of the image. Therefore, the entropy of compensated SAR image can be used as a measure of the accuracy of compensation for SAR data. The error estimation algorithms for SAR data based on minimum-entropy criterion overcomes the disadvantage of PGA relying on strong points and can be used in many scenarios. The minimum-entropy-based error estimation algorithms for SAR data can be divided into two kinds. The first kind estimate the error point by point [16] [21], which have high accuracy but large amount of calculation, and they are usually used to estimate the the whole aperture's error. The second kind approximate the real error function by polynomial [14] [22], and they can be used in multi-subaperture autofocus. Since the quadratic term in a polynomial has the greatest impact on focusing, only estimating the quadratic term of subaperture's phase error is an efficient method. To make the estimation faster, a Chebyshev polynomial fitting method is proposed to get the quadratic coefficients of the SPE [10].
For multi-subaperture autofocus algorithms, the estimation accuracy of the shift amount of the adjacent subaperture's error gradients will affect the splicing of the error gradients, thereby affecting the compensation for SAR data. So entropy also can be used as a measure of whether the shift amount is estimated accurately or not.

Proposed Gradient-Splicing Algorithm
This section may be divided by subheadings. It should provide a concise and precise description of the experimental results, their interpretation, as well as the experimental conclusions that can be drawn.

Rationale
Subaperture division needs to be considered before error estimation. It is affected by the algorithm used to estimate the error within the subaperture. It includes setting the length of the overlap between adjacent subapertures, i.e. l in this paper. When the polynomial coefficients of the SPE are estimated by the algorithm based on minimum-entropy criterion, the SPE gradient is smooth, which is different from that obtained by the PGA algorithm. Through some simulation and real data processing, it is found that in this case, the shift amount can still be well estimated without overlaps between subapertures. So the overlaps between the subaperture may be unnecessary in this case. In this paper, the cases of overlaps existing and not existing, i.e. l=0 and l>0, are both discussed.
After the phase error gradients of each subaperture are obtained, two or more adjacent subapertures are combined to form a "Combined Subaperture" (CS). The ith CS contains N subapertures ( iN  ), which is composed of the ( 1 iN − + )th to the ith subapertures. N is a function of the subaperture's ordinal number, which is expressed as:  If i>1, let the ith subaperture and ( ) 1 Ni− subapertures before it be combined to form the ith CS, whose length is  subaperture. Assume that the phase error gradients of the first i-1 subapertures has removed their shift amount, the next thing is to estimate should be subtracted H , which is expressed as: and IFFT the ith CS to get an image. Name this image "Img i". Img i is the result of coherent superposition of the two sub images, one of which is formed by the data of the ith subaperture and the other is formed by the data of the N(i)-1 subapertures before the ith. This algorithm is based on such a hypothesis, that is if one target's images in these two sub images will overlap to the maximum extent, and Img i will be clear. When  , these two images will be shifted gradually, so one target's energy disperse to a larger area, resulting in the increase of image entropy. If the shift amount is large enough, some targets will appear twice on Img i . Let However, it is difficult to determine the functional expression of () i Ex and whether Ex is a monotonic function when x< Ex by mathematical analysis. Direct one-dimensional search within a range is a viable method in this case. To narrow the search range,

1,
ii −  will be roughly estimated using the existing method first, and then the search will be performed within the neighborhood of the result of the rough estimation.

Preprints (www.preprints.org) | NOT PEER-REVIEWED | Posted: 27 April 2021
Suppose the phase error gradients of ith subaperture (before removing the shift amount) is as follows: The case of l>0 is considered first. The ordinal number of the pulses in the ith CS is i cn  is the corrected gradient of the ith subaperture's phase errors, which can be expressed as Now, start estimating can be expressed as follows: This step is similar to the error-gradient-splicing step in PGA-LS [9] [11], but the algorithm for estimating SPE is not limited to PGA. inaccurate now and will be further modified in the next steps. The gradients are spliced as follows: Suppose   is calculated as follows: This step is similar to the error gradient splicing step in PCA [10]. and {} i n  can be simplify as: After getting the rough estimator of ' i  is set, which is expressed as: One-dimensional search is used here to estimate ' i  . The center of the search range is generally set to 0. The radius of the search range will affect the efficiency and accuracy of the proposed algorithm. If the interval radius is too small, ' i  won't be included; if the interval radius is too large, the computational efficiency will decrease. So the search range should be as small as possible but contain the possible maximum value of  Considering the calculation accuracy and efficiency, the variable-step search method is used to search be the center of the kth search's range (k>0), the search step size be k d , and the search range is the following arithmetic progression: If the minimum value of  , so the (k+1)th search range is the arithmetic progression as follows: Let K denote the number of searches for ' i  . The value of K depends on the desired accuracy. The accuracy of the estimation of However, when only three elements are compared each time, due to the large search step in the first few searches,

()
i En x is easy to fall into the local minimum and cannot get the global minimum. To avoid this situation, the 0th search is added.
The 0th search also can be seen as the coarse estimation of    13:Accumulate the gradients to get the phase error of the whole aperture.

Simulation about l = 0
In this simulation, the proposed algorithm will be compared with other existing gradient splicing algorithm. In order to make the simulation scene closer to the real scene, hundreds points are randomly placed in the simulation scene, and hundreds points are arranged in letters. The point's reflection coefficient is set as 1 or 0.8 or 0.6, which results in different echo power of each point. The echo is produced in this simulation scene, and is added with Gaussian white noise with power of -15dB.
The simulation parameters are listed on Table 1: In the simulation, trajectory of radar carrier on Y-axis and Z-axis is shown in Figure  6(b). In real signal processing, SAR data need to be compensated by IMU data first, but this simulation doesn't involve the compensation about IMU, so the motion error is set little, which can be consider as the residual error, but it is enough to cause the decline of image quality. The Y position function of the trajectory is:  In the scene, the phase and amplitude of each point in the letters are seriously affected by the surrounding points, so it is relatively difficult to extract the error directly from the data using the traditional autofocus algorithm. Considering the minimum-entropy error estimation algorithm which estimate the quadratic term of error is robust and widely used, it is used to estimate the subaperture's phase error here. minimum-entropy-based algorithms for estimating higher-order errors have also been tried. Although these estimation method seems estimate the SPE and the gradients more accurately, after getting the final imaging result it can be observed that their effect is worse, so they aren't used here. The quadratic coefficients can be directly used as the curvature, so the data can be compensated by "PCA". Because there is no overlapping between the subapertures, when using PCA, the subaperture's phase error gradient line needs to extend an azimuthal pulse.
Though better than PGA, the subaperture's error estimated by the minimum-entropy-based algorithms is still not very accurate. For example, Figure 7 shows two images obtained by respectively using the minimum-entropy-based algorithm and the PGA algorithm to compensate the sixth subaperture's data, and Figure 7(a) is clearer than Figure 7(b); however, in Figure 7(a), the letter "N" is still a little fuzzy. Next, the algorithm proposed in this paper is used to calculate the shift amount.  The number of sampling points in azimuth direction is 6144. The SAR image's entropy after compensation in various cases is shown in Table 2. From the table, the entropy of the compensated image is related to N 0 . When 0 =2, although compensation has effect, but the effect is not good enough. When 0 =5, ˆi  is more accurate than that of 0 =2, and the image entropy after compensation is close to that without motion error. The gradient functions of subapertures before and after splicing and the error function curve are shown in Figure 8. From the Figure 9, it can be seen that when N 0 =5, the compensated image is almost the same as the error-free image. The image whose phase error second derivative has been estimated by the polynomial fitting method.
The fitting methods are tried. When the polynomial fitting or interpolation method is used to estimate the second derivative of the of the whole aperture phase error, because there are many kinds of polynomial fitting or interpolation methods, and their results are different, it is difficult to ensure the accuracy of the results. In this simulation, When the fitting method is used, every subaperture's length is 512 and 0 l = . After estimating the phase error of each subaperture, 12 quadratic coefficients are obtained. The number of polynomial terms is set to 12. In [10], overlapping should be set between adjacent subapertures, so the case of 256 l = is also simulated. Because of overfitting, although the phase error gradient appears smooth, the performance of using polynomial fitting method is far worse than that of the proposed algorithm. In addition, polynomials with more terms are also used for fitting, but their final compensation performances are all not well. These may due to that polynomial fitting brings uncontrollable error. After changing polynomial fitting into interpolation, the effect is not improved. So, the phase error estimating algorithm which only estimates the quadratic term of the SPE, the method of fitting / interpolation is not suitable for obtaining the whole aperture phase error in this scene.
Next, MD will be used to splice gradients. Let i S be the ith subaperture data after autofocus, i C be the cross correlation matrix of the compensated sub images which are formed by (i-1)th and ith subaperture. When MD algorithm is used for estimator, i  should be extract from i C . i C can also be expressed as follows: where means Hadamard product, and IFFT works in azimuth. However, due to the slight deformation of the compensated image (stretched or compressed in the azimuth direction) and the inaccurate focusing of the image corresponding to the subaperture, it is difficult to obtain the accurate relative shift between images. Figure 11(c) is the image representation of 6 C (the brighter the pixel, the higher the corresponding value), and it has been upsampled 20 times, so it has 10240×2048 pixels. Under the ideal condition, the peak's position of each column of i C should be almost the same, that is, in Figure 11(c) the brightest pixel in each column are concentrated on the same row. However, in fact the brightest pixels in each column are not concentrated on one row, which is also shown on Figure 11(d). It affect the estimator of i  . Though it seems that i  can still be extracted from 6 C , MD needs to be further modified to fulfill it. In conclusion, the proposed splicing algorithm has the best performance in this simulation, .

Simulation about l > 0
When l>0, i.e. there are overlaps between adjacent subapertures, the algorithm is still accurate. This simulation is designed to proves it. The simulation parameters are the same as those in Table 1, except that the azimuth beam angle is changed to 1°, and there are only point targets in the scene. The two changes greatly reduce the difficulty of SPE estimation.
The Y position function of the radar carrier's flight path is:   PGA+the proposed algorithm ( 0 =4) 8.4220179 After using various algorithms to compensate SAR data, the entropy of the resulting image is listed in Table 3. The entropy of the image compensated by "PGA+the proposed algorithm ( 0 =4)" differs by about 0.0039 from that of the image compensated by the PGA-MD algorithm, but the entropy of the image compensated by PGA-LS is 1.6716499 more than that of PGA-MD, and the entropy of the image compensated by PCA is 1.8044046 more than that of PGA-MD. The images of the data which are compensated by PGA-MD, "PGA + the proposed algorithm" respectively and the image formed by the data without error are shown in Figure 13. To compare on the azimuth dimension and the range dimension, one target is chosen in the scene. This point are surrounded by a red dotted line. The sections of the point obtained by the "PGA+the proposed algorithm" and PGA-MD almost coincide at this point, and are very close to the section of the SAR image without motion error and noise, which are shown in Figure 13(e)(f)(g)(h).
From the results of splicing, image entropy and the point-target analysis, the proposed algorithm has similar performance to MD and can splice accurately in this simulation.

Experimental data processing
The object of this processing are the SAR data of an airport in China and its surrounding environment. Through the accumulating the data of IMU, the diagram of the change of aircraft speed with time on each axis is obtained as Figure 14. The SAR data can be compensated roughly by these data. However, due to some reasons, the IMU data isn't accurate, which affects the compensation based on motion parameters.

Figure 15. The SAR image without autofocus
In this experiment, the RD algorithm is used for imaging. The result of coarse imaging is shown in the Figure 15. Number of azimuthal sampling points is 32768. The data that have been corrected for range cell migration is divided into 32 subapertures, i.e.

1024( 1,2,
,32) . Considering that the buildings and some targets with strong reflection are concentrated between the 1601st range cell and the 3648th range cell in the image, we use the part between these two range cells to estimate the phase error to reduce the calculation of estimation.
From the view of the image without autofocus, there is no isolated strong point in some azimuth coordinates, and there are some mountain slopes in the SAR image, which reflect radar waves so strongly that they are brighter than some man-made targets. Hence it is difficult to estimate the SPE by using the PGA algorithm, and the minimum-entropy-based algorithm is chosen to estimate the SPE.  is estimated accurately.
Firstly, the quadratic coefficients of the errors of the subapertures are estimated, and the gradients of the error are obtained. Next, the algorithm proposed in this paper is used for compensation. On the 0th search, 0 0.003 d = and 5 n = . The number of searches for ' i  is set as 25, and 0 is set as 2, 3, 4 and 5. Figure 16 shows the results of some intermediate steps of the autofocus algorithm. Figure 16(a) is the sub image formed by compensated first 3 subapertures' data, Figure 16(b) is the sub image formed by the compensated fourth subaperture, Figure 16(c) is the result of superposition of Figure 16(a) and Figure 16(b) when 3,4  is estimated accurately. Figure 17(b) shows the obtained phase error gradient function curve of the whole aperture. The SAR image compensated by the proposed algorithm ( 0 =5) is shown in Figure 18 PCA's splicing method also can be used, but the error in the subaperture is still estimated based on the minimum-entropy criterion, and only the quadratic term is estimated. Figure 17(a) shows the phase error gradient function curve of the whole aperture which is obtained by PCA's splicing method. There are slight differences between Figure  17(a) and Figure 17(b). In this experiment, comparing PCA with the proposed algorithm, we can see that the results are close, because the gradient change is relatively flat. But the effect of the new algorithm is slightly better. In azimuth direction, the energy of the target's first side lobe is smaller, and the main lobe and side lobe of the targets are separated better when using the proposed algorithm, which are reflected in Figure 19. Because there is a lot of irregular reflectors in this SAR image and their reflection is relatively strong, it is difficult to calculate 1, i i −  by MD algorithm in this experiment. For example, Figure 20(a) is the image presentation of the correlation matrix which is obtained by correlating the first and the second subapertures in azimuth direction and then upsampling 10 times. It has 10240×4096 pixels. The peaks' position of each column (i.e. each range cell) in Figure 20(a) is not on the same row, which is shown in Figure 20 After processing SAR data with each algorithm, the entropy of the resulting image is listed in Table 4. It is worth noting that when using the proposed algorithm and 0 = 3, the entropy of the image is smaller than that of 0 =4. This may be due to the fact that entropy does not fully reflect the compensation effect.

Table 4. Entropy in Various Cases
For the autofocus of SAR strip mode imaging, it is necessary to select not only the appropriate SPE estimation algorithm, but also the appropriate error gradient splicing algorithm. The proposed algorithm is based on the minimum-entropy criterion and its principle is simple. In the two simulations and a real-data processing, compared with other algorithms, the algorithm is proved to be accurate and robust. However, the algorithm also has some limitations that will slightly affect the use of the algorithm. It is suitable for cases where the accurate phase errors in SAR data are difficult to obtain. Considering the advantages of the algorithm, the algorithm deserves further study.