A Correction Method to Systematic Phase Drift of a High Resolution Radar for Foreign Object Debris Detection

: Due to the small size and various types of foreign object debris (FOD), radar detection of FOD on airport runways is a great challenge, and there are often a large number of false alarms in the detection results. Arc-scanning synthetic aperture radar (AS-SAR) is an emerging method for detecting FOD targets, which achieves omnidirectional coverage with a very high azimuth resolution. However, this method faces a similar challenge. A direct way to reduce false alarms is to increase the detection threshold based on enhancing the target signal-to-noise ratio (SNR), and in this paper, the coherent accumulation of multiple images is used to improve the target SNR. The stable phase is also an important feature of the target distinguishing background. Therefore, it is important to maintain the stability of the target phase. Aiming at the systematic phase drift (SPD) caused by atmospheric disturbance and system hardware, a spatial and temporal model is established, a corresponding correction approach is proposed, and the performance of the correction approach is validated by ﬁeld experiments.


Introduction
Foreign object debris (FOD) may damage aircraft and equipment or threaten the safety of airport staff and passengers on the ground in active areas such as airport runways. Therefore, they have brought huge hidden risks to personnel and flight safety. FOD has been known to damage airplanes and puncture tires in recent years, which leads to huge economic losses [1][2][3]. At present, researchers have designed a lot of FOD detection equipment using photoelectric and radar technology [4][5][6][7]. Photoelectric technology can identify FOD targets but is susceptible to light intensity and extreme weather conditions such as rain and snow. Radar technology has the characteristics of a high detection rate, long detection distance, strong environmental adaptability, and all-day working, which is the main research field of current FOD detection. However, radar technology faces challenges such as high false alarm rates and difficult identification.
At present, typical radar FOD detection equipment mainly includes Tarsier radar [8], FOD Finder [9], and FODetect [10]. They all use a frequency-modulated continuous wave (FMCW) signal and have an antenna with a narrow beam. In addition, Feil et al. designed a radar FOD detection system using FMCW, which has a detection range of more than 110 meters (m) and can detect small targets such as nuts [11]. The Chengdu Saying company and the University of Electronic Science and Technology of China developed a millimeterwave FMCW radar with a bandwidth of 800 MHz to detect FOD in 2010 [12]. Shanghai Jiaotong University [13] and Beijing Jiaotong University [14,15] have also researched on SPD is estimated based on the establishment of the spatial and temporary model. Lastly, the correction of SPD is realized. The remainder of this paper is organized as follows. In Section 2, a description of the imaging model of AS-SAR and the corresponding FOD detection flow is provided. In Section 3, the correction approach of SPD is presented. The experimental results are given in Section 4. Then, the accumulation effect of the SPD correction is discussed in Section 5. Finally, the conclusion is given in Section 6.

AS-SAR Imaging Model
The linear FMCW signal adopted by AS-SAR is given as follows: where f c is the carrier frequency, T pu is the modulation period, k r is the chirp constant, and t is the fast time within a modulation period. Then, the radar received signal of a target is: where ∆t is the echo delay, ∆t = 2R/c. R is the distance from the target to the radar, and A is the amplitude attenuation coefficient. The AS-SAR receiver mixes the received signal with the reference to obtain the demodulated intermediate frequency (IF) signal, and the multiplier mixing expression is shown as follows: where t ∈ ∆t, T pu , s IF (t) is IF signal. The geometric model of AS-SAR is shown in Figure 1, where G is the surface plane; O is the center of the rotary axis, with the FOD target T f placed on the G plane; H is the radar height; Ap is the phase center of the radar antenna; L is the distance between Ap and the rotary axis; τ is slow time; ω is the angular velocity of the mechanical rotating arm, with the circular angle ϕ = ω(τ + t − t 0 ) corresponding to the Ap at time t + τ, where t 0 is the time corresponding to the closest distance between Ap and the target; and θ is the angle between the OT f line and the rotation plane of the radar rotating arm. Suppose the distance from O to T f is R and the distance from Ap to T f at time τ + t is R(τ + t).
Remote Sens. 2022, 14, x FOR PEER REVIEW 3 of 19 correction algorithm. Firstly, the phase control points are adaptively extracted through standardized features, and then the SPD is estimated based on the establishment of the spatial and temporary model. Lastly, the correction of SPD is realized. The remainder of this paper is organized as follows. In Section 2, a description of the imaging model of AS-SAR and the corresponding FOD detection flow is provided. In Section 3, the correction approach of SPD is presented. The experimental results are given in Section 4. Then, the accumulation effect of the SPD correction is discussed in Section 5. Finally, the conclusion is given in Section 6.

AS-SAR Imaging Model
The linear FMCW signal adopted by AS-SAR is given as follows: where is the carrier frequency, is the modulation period, is the chirp constant, and is the fast time within a modulation period. Then, the radar received signal of a target is: where ∆ is the echo delay, ∆ = 2 / . is the distance from the target to the radar, and is the amplitude attenuation coefficient. The AS-SAR receiver mixes the received signal with the reference to obtain the demodulated intermediate frequency (IF) signal, and the multiplier mixing expression is shown as follows: where ∈ ∆ , , ( ) is IF signal.
The geometric model of AS-SAR is shown in Figure 1, where G is the surface plane; O is the center of the rotary axis, with the FOD target Tf placed on the G plane; is the radar height; Ap is the phase center of the radar antenna; is the distance between Ap and the rotary axis; is slow time; is the angular velocity of the mechanical rotating arm, with the circular angle = ( + − ) corresponding to the Ap at time + , where is the time corresponding to the closest distance between Ap and the target; and is the angle between the OTf line and the rotation plane of the radar rotating arm. Suppose the distance from O to Tf is and the distance from Ap to Tf at time + is ( + ). During imaging, because R L, R(τ + t) can be approximated through Taylor expansion as: In this paper, the L of the AS-SAR system is 1 m and the monitoring distance is more than 250 m, so Equation (4) is valid.
The imaging process coherently accumulates the target echo signals within the radar accumulation angle. To analyze the target echo, substituting ∆t = 2R(τ + t)/c into Equation (3), the s IF (t) yields: The corresponding S IF ( f ) is: A 1 is the amplitude attenuation coefficient, which is independent of phase. There are three exponential terms corresponding to the distance phase, residual video phase (RVP), and azimuthal Doppler phase. Therefore, imaging is the process of phase-matched filtering for the above items, in which RVP is ignored in imaging because its value change is too small [28,29,39].
The distance between the target and the radar axis is constant and there are different delays only for different azimuth angles. The frequency of the rotating arm circular angle ϕ is f ϕ . Then, S IF (f) is written as two-dimensional S(f, f ϕ ). Assume that the horizontal beam width of the antenna is ϕ ω . The bandwidth B ϕ of f ϕ is decided by the rotation speed of the turntable ω and the carrier frequency wavelength λ, which can be calculated using the following equation [39]: Additionally, S(f, f ϕ ) can be calculated as [28,29,39]: So, a phase-matched filter can be designed: Then, the filtering output is: where A 2 is the amplitude attenuation coefficient, which is independent of phase. The two-dimensional inverse fast Fourier transform (IFFT) is performed on S M (f, f ϕ ), and the temporal SAR image s I (R, ϕ) in polar coordinates can be obtained. The imaging flow is shown in Figure 2.

FOD Detection Flow
There are many types of FOD on the runway, such as a piece of asphalt concrete or cement concrete, a wrench, a twisted metal strip or the fuel tank cap, and so on, and more of them are small in size and weak in electromagnetic scattering. Moreover, people, vehicles, and aircraft on runways may form strong scattering areas in the images, affecting the target detection performance. The classical detection process of SAR images is shown in Figure 3, which mainly includes constant false alarm rate (CFAR) detection and mathematical morphology (MM) filtering [40,41]. CFAR detector involves a kind of algorithm that uses the difference between the target grayscale and its surrounding neighborhood grayscale for separation, which is widely used for SAR target detection. As a local detection algorithm, the thresholds of the CFAR detector are determined adaptively. The process of MM filtering is used to remove clutter with large differences in size and element targets. In general, through the selection of the size and shape of structural elements, morphological filtering first corrodes the image to remove the small clutter whose size does not meet the element target and then expands the image to maintain the element target in the image. To remove large-scale clutter in the scene, morphological filtering also adopts a similar method to remove small clutter. First, the image with the target removed is obtained, and then the image with the large-scale clutter removed can be obtained by subtracting the image from the original image. After that, target classification can even be realized by feature extraction for detection results [42][43][44]. According to the actual situation of FOD detection, false alarms are the main challenge of classical detection algorithms owing to the lower SNR. AS-SAR can monitor airport runways and other areas for a long time and form sequence of images. Multi-image coherent accumulation is introduced to improve the target SNR and reduce false alarms while ensuring the detection rate. The practical detection flow is given in Figure 4. The effectiveness of this flow depends on the stability of the target phase, so SPD correction is required. At the same time, the stable phase is also an important feature of the target distinguishing background. Through the derivation of Equation (10), the phase of the target in the image is only related to the distance. Theoretically, the phase of the target reflects the information of its position, which should be independent of other factors.

FOD Detection Flow
There are many types of FOD on the runway, such as a piece of asphalt concrete or cement concrete, a wrench, a twisted metal strip or the fuel tank cap, and so on, and more of them are small in size and weak in electromagnetic scattering. Moreover, people, vehicles, and aircraft on runways may form strong scattering areas in the images, affecting the target detection performance. The classical detection process of SAR images is shown in Figure 3, which mainly includes constant false alarm rate (CFAR) detection and mathematical morphology (MM) filtering [40,41]. CFAR detector involves a kind of algorithm that uses the difference between the target grayscale and its surrounding neighborhood grayscale for separation, which is widely used for SAR target detection. As a local detection algorithm, the thresholds of the CFAR detector are determined adaptively. The process of MM filtering is used to remove clutter with large differences in size and element targets. In general, through the selection of the size and shape of structural elements, morphological filtering first corrodes the image to remove the small clutter whose size does not meet the element target and then expands the image to maintain the element target in the image. To remove large-scale clutter in the scene, morphological filtering also adopts a similar method to remove small clutter. First, the image with the target removed is obtained, and then the image with the large-scale clutter removed can be obtained by subtracting the image from the original image. After that, target classification can even be realized by feature extraction for detection results [42][43][44].

FOD Detection Flow
There are many types of FOD on the runway, such as a piece of asphalt concrete or cement concrete, a wrench, a twisted metal strip or the fuel tank cap, and so on, and more of them are small in size and weak in electromagnetic scattering. Moreover, people, vehicles, and aircraft on runways may form strong scattering areas in the images, affecting the target detection performance. The classical detection process of SAR images is shown in Figure 3, which mainly includes constant false alarm rate (CFAR) detection and mathematical morphology (MM) filtering [40,41]. CFAR detector involves a kind of algorithm that uses the difference between the target grayscale and its surrounding neighborhood grayscale for separation, which is widely used for SAR target detection. As a local detection algorithm, the thresholds of the CFAR detector are determined adaptively. The process of MM filtering is used to remove clutter with large differences in size and element targets. In general, through the selection of the size and shape of structural elements, morphological filtering first corrodes the image to remove the small clutter whose size does not meet the element target and then expands the image to maintain the element target in the image. To remove large-scale clutter in the scene, morphological filtering also adopts a similar method to remove small clutter. First, the image with the target removed is obtained, and then the image with the large-scale clutter removed can be obtained by subtracting the image from the original image. After that, target classification can even be realized by feature extraction for detection results [42][43][44]. According to the actual situation of FOD detection, false alarms are the main challenge of classical detection algorithms owing to the lower SNR. AS-SAR can monitor airport runways and other areas for a long time and form sequence of images. Multi-image coherent accumulation is introduced to improve the target SNR and reduce false alarms while ensuring the detection rate. The practical detection flow is given in Figure 4. The effectiveness of this flow depends on the stability of the target phase, so SPD correction is required. At the same time, the stable phase is also an important feature of the target dis- According to the actual situation of FOD detection, false alarms are the main challenge of classical detection algorithms owing to the lower SNR. AS-SAR can monitor airport runways and other areas for a long time and form sequence of images. Multi-image coherent accumulation is introduced to improve the target SNR and reduce false alarms while ensuring the detection rate. The practical detection flow is given in Figure 4. The effectiveness of this flow depends on the stability of the target phase, so SPD correction is required. At the same time, the stable phase is also an important feature of the target distinguishing background.

The Correction Approach to SPD
SPD is caused by atmospheric disturbance and system noise. It has both two-dimensional space-varying and time-varying characteristics. There are always some points with fixed positions in the actual scene, with high SNR and stable scattering characteristics. According to the derivation of the imaging model, the phase change of these scattering points is mainly caused by SPD. Therefore, these scattering points can be used to estimate SPD and called stable phase control points (SPCP) in this paper. This section introduces an adaptively SPD correction approach based on the SPCP, which mainly includes three parts. First, a spatial and temporal model of SPD is established, and the corresponding SPD estimation method is given. Second, the screening method of the SPCP is shown, which is used to estimate the SPD. Third, the complete process flow of the correction approach is obtained.

Modeling and Estimation of SPD
Setting the values of the monitoring target in two sequential images as and , SPD ∆ can be calculated by the following equation: where (•) represents the calculated phase angle and * represents the complex conjugate. Considering that SPD is caused by the system itself and the atmosphere, ∆ can be where ∆ represents the phase caused by the system itself, ∆ represents the phase caused by the atmosphere, and ∆ is the noise phase. In general, ∆ is constant. The temperature and humidity of the atmosphere change in space over time, which causes the nonuniformity of electromagnetic characteristics, so that the propagation speed and direction of the signal change continuously when passing through the atmosphere. When the atmospheric changes along the range and azimuth are independent of each other, the first-order function model of SPD along the distance and azimuth is established [37,38]: where is the carrier frequency wavelength, is the ground range from target to radar, and is the angle between the target and the initial scanning direction of the radar. , , and are weighting coefficients. Equation (13) shows that as long as three stable strong scattering points with high SNR are found, the simultaneous equation can be established through their phase drift, and the weighting coefficient can be solved. Here, a stable strong scattering point was described as the SPCP. In fact, the number of SPCPs is much greater than 3. Therefore, the least squares (LS) method will be used to estimate the weighting coefficient. Let ⃑ =

The Correction Approach to SPD
SPD is caused by atmospheric disturbance and system noise. It has both two-dimensional space-varying and time-varying characteristics. There are always some points with fixed positions in the actual scene, with high SNR and stable scattering characteristics. According to the derivation of the imaging model, the phase change of these scattering points is mainly caused by SPD. Therefore, these scattering points can be used to estimate SPD and called stable phase control points (SPCP) in this paper. This section introduces an adaptively SPD correction approach based on the SPCP, which mainly includes three parts. First, a spatial and temporal model of SPD is established, and the corresponding SPD estimation method is given. Second, the screening method of the SPCP is shown, which is used to estimate the SPD. Third, the complete process flow of the correction approach is obtained.

Modeling and Estimation of SPD
Setting the values of the monitoring target in two sequential images as α 1 e jϕ 1 and α 2 e jϕ 2 , SPD ∆ϕ can be calculated by the following equation: where angle(·) represents the calculated phase angle and * represents the complex conjugate.
Considering that SPD is caused by the system itself and the atmosphere, ∆ϕ can be where ∆ϕ s represents the phase caused by the system itself, ∆ϕ a represents the phase caused by the atmosphere, and ∆ϕ n is the noise phase. In general, ∆ϕ s is constant. The temperature and humidity of the atmosphere change in space over time, which causes the nonuniformity of electromagnetic characteristics, so that the propagation speed and direction of the signal change continuously when passing through the atmosphere. When the atmospheric changes along the range and azimuth are independent of each other, the first-order function model of SPD along the distance and azimuth is established [37,38]: where λ is the carrier frequency wavelength, r g is the ground range from target to radar, and ϕ g is the angle between the target and the initial scanning direction of the radar. β 0 , β 1 , and β 2 are weighting coefficients. Equation (13) shows that as long as three stable strong scattering points with high SNR are found, the simultaneous equation can be established through their phase drift, and the weighting coefficient can be solved. Here, a stable strong scattering point was described as the SPCP. In fact, the number of SPCPs is much greater than 3. Therefore, the least squares (LS) method will be used to estimate the weighting coefficient. Let β = [β 0 , β 1 , β 2 ]. Equation (13) can be expressed as a matrix as follows: where e is the noise matrix, and h = 4π For multiple SPCPs, the matrix is expressed as an equation The LS algorithm can be used, and Then,∆ϕ is estimated as follows:∆ Due to the time continuity of temperature and humidity changes causing SPD, the estimated SPD can be filtered by a Kalman filter. The state transition matrix can be established, such as the equation: where q k−1 is Gaussian noise and F 1 is the state transfer matrix. The observation equation is established as follows: According to actual, the standard deviation of SPD change caused by noise shall not exceed 2 degrees, so the variance matrix Q corresponding to q is set to 4. Considering that the error of SPCPs' SPD measured by AS-SAR is better than 3 degrees when processing two adjacent frames of images, so the variance matrix V corresponding to v is set to 9. In addition, it is assumed that SPD changes slowly with temperature and humidity, and the observed value is SPD, so F 1 = 1 and F 2 = 1. Then, time-domain filtering of the SPCP can be realized.

The Screening Method of the SPCP
In the process of SPD calculation, the selection of the SPCP is very important. For amplitude images, three features to screen SPCPs are defined: local image contrast F l , amplitude dispersion F d , and correlation coefficient F c . The screening process flow of the SPCP is shown in Figure 5. , (22 The larger is, the more stable the amplitude information is. Given that the calcu lation of only uses the amplitude information, it has the advantages of a small amoun of calculation and easy extraction.
For , the correlation coefficient of each pixel of the sequence of images needs to b calculated according to the following equation: where and are sliding neighborhood window sizes. Then, the mean ( , ) of eac ( , ) is calculated, and is acquired by The larger is, the stronger the correlation. involves the neighborhood of pixels which has a large amount of calculation, but the feature is stable and less disturbed b noise.
Each extracted feature is formed into a vector, ⃑ = [ , , ]. Considering that th values of each component of the feature vector are not uniform and have different effect on the linear classifier, it is necessary to normalize it. According to the principle of Ma halanobis distance measurement [45], this paper solves the problem of different dimen Let the sequence of images be X n , 1 ≤ n ≤ N, where N is the number of images. x n (i, j) represents the image X n amplitude value at pixel point (i, j). Then, F l , F d , and F c are defined as follows: For F l , this paper first acquires the neighborhood image slice I i of each pixel x n (i, j) in X n and calculate the mean value m n b (i, j) and standard deviation σ n b (i, j), then The larger F l is, the larger the SNR and the more stable the corresponding pixel. For F d , the amplitude mean value m(i, j) and standard deviation σ(i, j) of each pixel x n (i, j) in the sequence of images are calculated, and then it is calculated by The larger F d is, the more stable the amplitude information is. Given that the calculation of F d only uses the amplitude information, it has the advantages of a small amount of calculation and easy extraction.
For F c , the correlation coefficient of each pixel of the sequence of images needs to be calculated according to the following equation: where I and J are sliding neighborhood window sizes. Then, the mean m c (i, j) of each c n (i, j) is calculated, and F c is acquired by The larger F c is, the stronger the correlation. F c involves the neighborhood of pixels, which has a large amount of calculation, but the feature is stable and less disturbed by noise.
Each extracted feature is formed into a vector, F = [F l , F d , F c ]. Considering that the values of each component of the feature vector are not uniform and have different effects on the linear classifier, it is necessary to normalize it. According to the principle of Mahalanobis distance measurement [45], this paper solves the problem of different dimension distribution through covariance normalization. Covariance matrix Σ F can be obtained by where U is the unitary matrix, µ F is the mean value of feature vectors, and T is the transpose operation. Then, the normalized feature vector F u can be represented by the following equation: In the process of classification, a linear classification decision criterion is established: where W T and ω 0 are weighting coefficients, which can be obtained through known calibration target training. In the training process, to ensure the computational efficiency and the expansibility of the approach, the number of SPCPs in this paper is limited to 0.5-1% of the total number of image points. Then, on this basis, the training obtains the weighting coefficients when the number of SPCPs falling in fixed areas such as guardrails and buildings accounts for the highest number of total SPCPs. At the same time, the robustness and consistency of radar performance are good and the SNR meets the monitoring conditions; the classification process has expansibility. Figure 6 shows the complete process flow of the correction approach. It consists of three stages:

•
The screening of the SPCP. This stage can be completed in real time or offline. To prevent the FOD target from being used as an SPCP, it needs empty scene data. It processes the sequence of amplitude images through feature extraction and classification and sends the obtained SPCPs to the second stage.

•
The estimation of the SPD. In the actual processing process, using the SPCPs provided in the first stage and based on the SPD model, the SPD between two complex images is estimated in real time by the LS algorithm. Finally, the Kalman filter is carried out based on the multiframe SPD (greater than 2) to obtain the SPD for correction.

•
The correction of the SPD. Expected correction is achieved by subtracting the SPD, and the sequence of complex images after phase correction are produced and output.

•
The estimation of the SPD. In the actual processing process, using the SPCPs provided in the first stage and based on the SPD model, the SPD between two complex images is estimated in real time by the LS algorithm. Finally, the Kalman filter is carried out based on the multiframe SPD (greater than 2) to obtain the SPD for correction.

•
The correction of the SPD. Expected correction is achieved by subtracting the SPD, and the sequence of complex images after phase correction are produced and output.

Results
To demonstrate the effectiveness of the proposed correction approach to SPD, the experimental results based on the measured AS-SAR data are presented. Figure 7 shows the AS-SAR and its working field scene, and the corresponding image is shown in Figure 8. The 16-bit sampling data of radar echo is processed directly using the imaging process shown in Figure 2 to obtain the image s I . Considering the large value of s I amplitude, which is not conducted to display. Therefore, Figure 8 and subsequent AS-SAR images are obtained using the function 20 × log 10 (s I ). Thus, AS-SAR is able to effectively image the monitored area on the runway. As described above, the screening of SPCPs is very important. These points are extracted according to Section 3.2, as shown in Figure 9. Compared with the actual scene, the selected SPCPs are mainly fixed targets such as fences and houses. It can be assumed that their phase change is independent of the position change and is mainly caused by SPD. Thus, they can be used to estimate the whole scene SPD.

Results
To demonstrate the effectiveness of the proposed correction approach to SPD, the experimental results based on the measured AS-SAR data are presented. Figure 7 shows the AS-SAR and its working field scene, and the corresponding image is shown in Figure  8. The 16-bit sampling data of radar echo is processed directly using the imaging process shown in Figure 2 to obtain the image . Considering the large value of amplitude, which is not conducted to display. Therefore, Figure 8 and subsequent AS-SAR images are obtained using the function 20 ( ). Thus, AS-SAR is able to effectively image the monitored area on the runway. As described above, the screening of SPCPs is very important. These points are extracted according to Section 3.2, as shown in Figure 9. Compared with the actual scene, the selected SPCPs are mainly fixed targets such as fences and houses. It can be assumed that their phase change is independent of the position change and is mainly caused by SPD. Thus, they can be used to estimate the whole scene SPD.     To estimate the SPD more accurately, it is necessary to analyze the phase change history of SPCPs, as shown in Figure 10. All initial phases are aligned to 0, and the duration of continuous monitoring is 12 h. Figure 10a shows the phase change history of SPCPs with the same range and different azimuths. The phase difference between them is small. Figure 10b shows the phase change history of SPCPs with different ranges and the same azimuth. The farther they are from the radar, the greater the phase change. To estimate the SPD more accurately, it is necessary to analyze the phase change history of SPCPs, as shown in Figure 10. All initial phases are aligned to 0, and the duration of continuous monitoring is 12 h. Figure 10a shows the phase change history of SPCPs with the same range and different azimuths. The phase difference between them is small. Figure 10b shows the phase change history of SPCPs with different ranges and the same azimuth. The farther they are from the radar, the greater the phase change.
In this paper, the cylindrical targets shown in Figure 11 are analyzed, and their height and diameter are equal. The target is arranged on the runway, as shown in Figure 12, and Table 1 gives correspondence between cylinder size and target number.
To estimate the SPD more accurately, it is necessary to analyze the phase change history of SPCPs, as shown in Figure 10. All initial phases are aligned to 0, and the duration of continuous monitoring is 12 h. Figure 10a shows the phase change history of SPCPs with the same range and different azimuths. The phase difference between them is small. Figure 10b shows the phase change history of SPCPs with different ranges and the same azimuth. The farther they are from the radar, the greater the phase change. In this paper, the cylindrical targets shown in Figure 11 are analyzed, and their height and diameter are equal. The target is arranged on the runway, as shown in Figure 12, and Table 1 gives correspondence between cylinder size and target number.  In this paper, the cylindrical targets shown in Figure 11 are analyzed, and their height and diameter are equal. The target is arranged on the runway, as shown in Figure 12, and Table 1 gives correspondence between cylinder size and target number.   In the actual experimental process, the SPD is estimated using the LS algorithm and Kalman filter. Figure 13 gives an estimated SPD image of two AS-SAR images. Consistent with the actual measurement results, SPD is more affected by range. Then, a correction is achieved by subtracting the SPD. Figure 14 shows the phase change history of target 1 and clutter. Through comparison, it is found that the target phase changes more than 1.2 rad after 16 frames before the correction, while the target phase changes less than 0.5 when the SPD is corrected. Moreover, the phase change history of clutter is random. The correction approach proposed in this paper can effectively maintain the stability of the target  In the actual experimental process, the SPD is estimated using the LS algorithm and Kalman filter. Figure 13 gives an estimated SPD image of two AS-SAR images. Consistent with the actual measurement results, SPD is more affected by range. Then, a correction is achieved by subtracting the SPD. Figure 14 shows the phase change history of target 1 and clutter. Through comparison, it is found that the target phase changes more than 1.2 rad after 16 frames before the correction, while the target phase changes less than 0.5 when the SPD is corrected. Moreover, the phase change history of clutter is random. The correction approach proposed in this paper can effectively maintain the stability of the target phase. In addition, the stable phase can also be used to distinguish the target from the background.

Discussion
According to Section 2.2, false alarms are the main challenge of classical detection algorithms owing to the lower SNR. Coherent accumulation is needed in the FOD detection process to improve the SNR. The stable phase guarantees improved coherent accumulation performance. Figure 15 shows a comparison of coherent accumulation images. Compared with the original single image, the background noise in the scene is reduced after coherent accumulation. When more images are accumulated, the better the SNRs are. At the same time, the accumulation effect of the SPD correction is better than that of no correction.
The quantitative analysis results are shown in Figure 16 and Table 2. Figure 16 gives a comparison of target SNRs between uncorrected and corrected samples. Figure 16a shows the comparison of the target SNR after coherent accumulation using 4 images. The actual value curve of the target SNR, which is corrected, is mostly above its uncorrected value. Their corresponding mean values of the target SNRs are 20.68 dB and 21.58 dB, and the effect of coherent accumulation is better after correction. When the number of accumulated images becomes 8 and 16, it shows similar characteristics, as shown in Figure  16b

Discussion
According to Section 2.2, false alarms are the main challenge of classical detection algorithms owing to the lower SNR. Coherent accumulation is needed in the FOD detec tion process to improve the SNR. The stable phase guarantees improved coherent accu mulation performance. Figure 15 shows a comparison of coherent accumulation images Compared with the original single image, the background noise in the scene is reduced after coherent accumulation. When more images are accumulated, the better the SNRs are At the same time, the accumulation effect of the SPD correction is better than that of n correction.
The quantitative analysis results are shown in Figure 16 and Table 2. Figure 16 give a comparison of target SNRs between uncorrected and corrected samples. Figure 16 shows the comparison of the target SNR after coherent accumulation using 4 images. Th actual value curve of the target SNR, which is corrected, is mostly above its uncorrected value. Their corresponding mean values of the target SNRs are 20.68 dB and 21.58 dB, and the effect of coherent accumulation is better after correction. When the number of accu mulated images becomes 8 and 16, it shows similar characteristics, as shown in Figur 16b,c. For coherent accumulation using 8 images, the corresponding mean values are 26.5

Discussion
According to Section 2.2, false alarms are the main challenge of classical detection algorithms owing to the lower SNR. Coherent accumulation is needed in the FOD detection process to improve the SNR. The stable phase guarantees improved coherent accumulation performance. Figure 15 shows a comparison of coherent accumulation images. Compared with the original single image, the background noise in the scene is reduced after coherent accumulation. When more images are accumulated, the better the SNRs are. At the same time, the accumulation effect of the SPD correction is better than that of no correction. dB, 4.27 dB, and 3.86 dB, respectively. Table 2 gives the SNR results of 24 targets under different conditions, while the values with better SNRs are bold. Coherent accumulation can effectively improve the target SNR. Moreover, after SPD correction, the SNR of most targets improved significantly. Therefore, this approach makes the phase stability of the same target in multi-frame images better and proves the effectiveness of SPD correction from the side. The quantitative analysis results are shown in Figure 16 and Table 2. Figure 16 gives a comparison of target SNRs between uncorrected and corrected samples. Figure 16a shows the comparison of the target SNR after coherent accumulation using 4 images. The actual value curve of the target SNR, which is corrected, is mostly above its uncorrected value. Their corresponding mean values of the target SNRs are 20.68 dB and 21.58 dB, and the effect of coherent accumulation is better after correction. When the number of accumulated images becomes 8 and 16, it shows similar characteristics, as shown in Figure 16b,c. For coherent accumulation using 8 Table 2 gives the SNR results of 24 targets under different conditions, while the values with better SNRs are bold. Coherent accumulation can effectively improve the target SNR. Moreover, after SPD correction, the SNR of most targets improved significantly. Therefore, this approach makes the phase stability of the same target in multi-frame images better and proves the effectiveness of SPD correction from the side. Table 2. The SNR of targets (dB).

Conclusions
In this paper, maintaining the stability of the target phase is the key to coherent accumulation in FOD detection flow. Through the analysis of the AS-SAR imaging model, it is concluded that the theoretical phase of the target is only related to its position, and its change is the SPD, which is mainly caused by atmospheric disturbances and system hardware. Based on SPD modeling, a complete process flow of the correction approach is proposed. It adaptively screens SPCPs through the normalized features of local image contrast, amplitude dispersion, and core translation coefficient, and then makes full use of spatial filtering and temporal filtering to ensure the accuracy of correction. The real data results show that the SNR of most targets after the correction has improved significantly. Therefore, this approach makes the phase stability of the same target in multi-frame images better, which indirectly proves the effectiveness of SPD correction. However, the approach in this paper also has some shortcomings. First, the time model is relatively simple and does not introduce the relevant information of temperature and humidity changes. Second, the number of SPCPs screened is large, which affects the calculation efficiency. At the same time, the number of SPCPs at the edge of the runway is relatively small, which will also affect the spatial filtering effect of SPD. Therefore, it is necessary to further study and improve the spatial and temporal model of SPD.