Translational-Invariant Ultrahigh-Resolution Spaceborne Bistatic SAR Focusing and Image Coordinate Analysis Using Bistatic Polynomial Model

We have proposed a novel process flow for translational-invariant (TI) and ultrahigh-resolution (UHR) bistatic synthetic aperture radar (BSAR) operations, wherein a UHR BSAR image is obtained via coherent integration of sub-aperture images, and its two-dimensional (2D) coordinates are then interpreted. The main contributions of the proposed method are as follows. First, to address additional aliasing of the TI BSAR spectrum after sub-aperture division, linear range walk compensation and its restoration are subsequently performed. This efficiently extends the duration of the Doppler domain after sub-aperture division, accommodating the skewed shape of the TI BSAR spectrum. Second, spatial-variant focusing is achieved via the bistatic polynomial omega-K algorithm (BP-OKA), for which a key consideration is highlighted to perform BP-OKA in the sub-image integration process. In its absence, the same scatterers appeared at different positions in each sub-aperture image, causing failure of coherent accumulation. Finally, to retrieve the range–azimuth coordinates of targets within an illuminated scene, we established a high-order equation using bistatic polynomials, and the 2D image coordinates of targets were determined by solving the equation. This new concept is fully compatible with the proposed BP-OKA and can be used in various spaceborne TI BSAR scenarios. All these contributions were verified using numerical simulations. Accordingly, the proposed processing flow can provide focused SAR images with accurate image coordinate interpretations for TI and UHR BSAR operations.


I. INTRODUCTION
Bistatic synthetic aperture radar (BSAR) is a type of synthetic aperture radar (SAR) system in which a transmitter and receiver are spatially separated on different platforms. The geometry of the BSAR system offers considerable reliability and flexibility in designing SAR missions and has many advantages over typical monostatic SAR systems: frequent monitoring, good electromagnetic concealment, acquisition of multiaspect scattering information, and capability of high-resolution and wide-swath imaging. These advantages The associate editor coordinating the review of this manuscript and approving it for publication was Brian Ng .
The translational-invariant (TI) BSAR configuration, where transmitting and receiving platforms follow the same track or parallel tracks with almost identical velocities, is a special case of the general BSAR system [18], [19], [20], [21], [22], [23], [24], [25], and [25]. The restricted behavior of TI BSAR facilitates the operation of two separate platforms. Using this advantage, spaceborne BSAR with a TI configuration has been operated for interferometry and tomography, e.g., TanDEM-X [7] and RADARSAT-2/3 [8]. Moreover in the future, more advanced TI BSAR systems will be developed, as multiple TI BSARs can be used to establish a multistatic SAR configuration, which is one of the future development directions of the SAR system for providing wide-swath and high-resolution images [41], [42], [43]. In particular, with the increasing demand for centimeterscale SAR images, ultrahigh-resolution (UHR) imaging is becoming imperative in TI spaceborne BSAR applications, for which dedicated SAR processing should be devised.
In the TI BSAR system, the baseline between the transmitter and receiver is almost constant over time, indicating that the spatial variance of the echo signal is negligible in the azimuthal direction. Using the azimuth-invariant signal property, several attempts at TI BSAR imaging have been made. In most cases, bistatic range histories were formulated with double hyperbolic range models (HRMs). In [18] and [19], after dip-moveout preprocessing, the BSAR data were converted into equivalent monostatic data. Subsequently, the data were focused on using monostatic algorithms. In [20], [21], [22], [23], [24], and [25], the two-dimensional (2D) frequency spectrum of the double HRMs was derived, and typical monostatic algorithms were modified to accommodate it. All the aforementioned TI BSAR approaches have been demonstrated to have sufficient accuracy in a rectilinear geometry, wherein both the transmitter and receiver move along straight paths. In traditional low-resolution spaceborne SAR systems, such a linear-orbit approximation is valid because the integration time is short. However, as the integration time increases with the achievement of higher resolutions, the curved paths of the two satellites exhibit large orbital arcs. Subsequently, the bistatic range histories deviate significantly from the double HRMs; thus, classical TI BSAR approaches fail to provide focused images. Accordingly, a focused approach for TI UHR spaceborne BSAR using general expressions of bistatic range histories should be developed.
In addition to the curved-path issue, Doppler spectrum aliasing is another problem that should be addressed for UHR spaceborne SAR focusing, as the total bandwidth of the radar echo signal significantly exceeds the pulse repetition frequency (PRF) in the UHR SAR system. Generally, SAR systems are operated with a PRF of several kHz, while the total Doppler bandwidth of the echo signal with UHR operations is tens of kHz. In previous studies [26], [27], [28], [29], [30], [31], [32], [33], [34], such PRF limitations were addressed via signal processing, which can be divided into two categories: sub-apertureand deramping-based techniques. Deramping-based techniques [26], [27], [28] involve convolutions between the entire raw data and a selected chirp signal along the azimuth direction, which is computationally expensive and requires a sufficient hardware memory capacity. In comparison, sub-aperture techniques are advantageous for spaceborne SAR applications [29], [30], [31], [32], [33], [34], because they are conceptually simple and facilitate parallel computing.
Sub-aperture techniques divide all received signals into azimuth sub-apertures, thereby avoiding PRF aliasing. The key aspect of this technique is the recovery of the full-aperture spectrum from sub-aperture data. In typical approaches [29], [30], [31], sub-aperture data are recombined in the signal domain before focusing, whereas in recent years, more advanced frameworks have been proposed, where sub-aperture images are coherently integrated in the image domain after sub-aperture focusing [32], [33], [34]. Subimage integrations can provide a full-aperture UHR image, as well as subsequent low-resolution images within a single received data. This allows consideration of additional applications, such as moving-target detection and clutter suppression, using subsequent sub-aperture images [35], [37]. In addition, sub-aperture images can be obtained by an onboard processor while receiving radar echoes, which facilitates the design of real-time processing [32], [34]. Accordingly, sub-image integration has attracted more attention than typical sub-aperture recombination in the signal domain. However, previous studies [32], [33], [34] have focused on monostatic SAR systems. Although accurate registration is the most important consideration for coherent accumulation of sub-images, it was not discussed in [32], [33], and [34].
To interpret BSAR images effectively, it is imperative to clarify how illuminated areas are mapped to 2D image coordinates. In fact, this mapping relationship depends on the range-Doppler property of focusing. The most well-known relationship is close-range and zero-Doppler time [1], [31], where targets are mapped to the range and azimuth time that are concurrently reached at zero Doppler frequency. However, this is only valid for HRM-based focusing, which limits its use in TI spaceborne BSAR scenarios, because the focusing quality is significantly degraded when the transmitter and receiver are separated by a large baseline and exhibit severe orbital arcs. Therefore, an image coordinate analysis method that is compatible with a focusing approach using general range formulas should be developed; however, related topics have not yet been studied.
Motivated by the aforementioned problems, in this study, we propose a full processing framework for TI UHR spaceborne BSAR imaging, in which three main issues are addressed: 1) Doppler spectrum aliasing, 2) curved-orbit geometry, and 3) interpretation of BSAR image coordinates. Our proposal includes the following contributions:

1) Doppler Extensions via Subsequent Linear Range Walk Compensation (LRWC) and Its Restorations in Sub-Aperture Divisions:
According to the signal properties in the BSAR system, the Doppler centroid varies with respect to the range frequencies, which skews the support region of the 2D spectrum.

2) Bistatic Polynomial Omega-k Algorithm (BP-OKA):
To address the curved-orbit issue, the polynomial omega-K algorithm (P-OKA) was proposed for monostatic SAR focusing [38], [39], [40]. The basic property of P-OKA, which provides high-quality images via range-variant focusing with an accurate range model, is well suited for azimuth-invariant TI BSAR signals. Nevertheless, the signal formulas for bistatic and monostatic SARs differ; thus, further modification is required for P-OKA to be compatible with the BSAR system. To this end, we have introduced a focusing scheme that is a modified form of monostatic P-OKA using bistatic polynomials, denoted as BP-OKA. In contrast to the double HRMs reported in [18], [19], [20], [21], [22], [23], [24], and [25], the bistatic polynomials can accurately express actual bistatic range histories even in large curved paths. Thus, compared with the conventional TI BSAR focusing method used in [18], [19], [20], [21], [22], [23], [24], and [25], BP-OKA can provide superior image quality by achieving complete range-variant processing. Moreover, we illustrate the key considerations required to perform BP-OKA within the framework of sub-image integration, without which subsequent subaperture images are incorrectly registered, resulting in the failure of coherent integration for UHR imaging.

3) Imaging Coordinate Analysis Using Bistatic Polynomials:
To interpret image coordinates after focusing according to general expressions of bistatic ranges, we introduce a new method for defining rangeand-azimuth coordinates using a bistatic polynomial relationship. This new relationship is completely different from the conventional close-range and zero-Doppler time and is fully compatible with the bistatic polynomial-based focusing that can provide focused images in various TI BSAR scenarios even with large curved paths and long baseline geometries. The remainder of this paper is organized as follows. Section II presents the mathematical formulas of the TI spaceborne BSAR signal model with bistatic polynomials. Section III describes the proposed processing framework. In Section IV, numerical simulation results are presented to verify the efficacy of the proposed method. Finally, conclusions are presented in Section V.

II. SIGNAL FORMULATIONS
In this section, we illustrate the signal formulations of the TI UHR spaceborne BSAR geometry, as shown in Fig. 1; the origin of the coordinate system is located at the scene center O, and the x, y, and z axes represent the cross-track, along-track, and height directions, respectively. In the figure, the two satellites, including a transmitter S T and receiver S R , are separated by a baseline length b BI , and they move along parallel curved paths with almost identical velocities. θ T and θ R represent the squint angles of S T and S R , respectively, which establishes the bistatic angle θ BI . S T transmits a linear frequency-modulated signal, and S R receives an echo reflected by an illuminated scene during integration time T .
After range compression neglecting envelope terms, the received BSAR signal from a target P p = [x p , y p , z p ] T can be expressed in a range-frequency (f τ ) and azimuth-time (η) format, as follows: where c represents the speed of light and f c represents the radar operating frequency. In (1), R BI η; P p denotes the bistatic range histories of P p , which is the total distance from S T to P p and from P p to S R : Here, S T (η) and S R (η) are the position vectors of the transmitter and receiver with respect to η, respectively, and |·| denotes the norm operator. For obtaining accurate expressions of bistatic range histories in the curved-path geometry, polynomial range models higher than the fifth order were adopted in [13], [14], and [45]. Although the previous studies focused on translational-variant BSAR, the high-order polynomials also achieve sufficient accuracy in the TI BSAR configuration, even with large curved paths. Accordingly, 54464 VOLUME 11, 2023 Authorized licensed use limited to the terms of the applicable license agreement with IEEE. Restrictions apply. we expand (2) via the sixth-order polynomials using the Taylor series, as follows: where η c = T /2 represents the azimuth central time, r BI p = R BI η c ; P p represents the bistatic range of P p at η c , and k BI p,n represents the n th -order (n = 1, 2, . . . , 6) polynomial coefficients. Using (3), the bistatic Doppler frequencies of P p at f c can be obtained as follows: where represents the n th -order bistatic Doppler parameters of P p . The 2D spectrum of P p for the BSAR signal is obtained via the Fourier transform (FT) of (1) with respect to η, as follows: where f η represents the azimuth frequency, and BI f τ , f η represents the 2D BSAR phase spectrum of P p . According to principle of stationary points, BI f τ , f η ; P p in (6) can be computed as where η * f η represents the stationary phase point that satisfies d dη Using (3) and (8), we obtain the following relationship: In (9), η * f η can be obtained via the method of series reversion [10], and the result is where Y = c (f τ +f c ) − k BI p,1 , and .
(11) VOLUME 11, 2023 As indicated by (7) and (11), BI f τ , f η ; P p in (6) is defined by bistatic polynomials that change geometrically within the illuminated scene. Suppose that all the targets are located at the same height as the reference point of O; then, the geometric variance in (6) is only dependent on the crosstrack (x) and along-track (y) directions; the variations in x and y correspond to the range and azimuth variance, respectively. According to the TI spaceborne BSAR signal property, the azimuth variance of BI f τ , f η ; P p is sufficiently small to be neglected; subsequently, sole consideration in TI BSAR processing accommodates spatial variance of the phase spectrum in the range directions. Therefore, we can where P x = [x, 0, 0] T highlights targets distributed in crosstrack direction.

III. PROPOSED METHOD A. OVERVIEW
The main objectives are image focusing and its coordinate analysis in a TI UHR spaceborne BSAR system. Conceptual illustrations and key points of the proposed method are highlighted in Fig. 2.
For focusing, the full-aperture data are divided into subapertures, and partial-resolution images are obtained from each sub-aperture dataset. These sub-aperture images are coherently integrated into the image domain to recover the full resolution. The attractive features of our focusing process are as follows. First, to resolve the additional aliasing introduced by the skewed shape of the sub-aperture BSAR spectrum, we propose extended sub-image integration that includes subsequent procedures, such as LRWC, Doppler extensions, and its restoration. This proposal efficiently extends Doppler durations from the original PRF and accommodates a skewed BSAR spectrum through a simple implementation. Second, by employing the azimuth-invariant signal property in TI BSAR, we implement BP-OKA, which completely performs range-variant focusing using the 2D phase spectrum of bistatic polynomials, as given by (12). In BP-OKA, bulk focusing is first conducted using multiplication (RFM). Subsequently, residual range-variant errors are eliminated by generalized Stolt interpolations (G-SI); the 2D frequency phase of bistatic polynomials is decomposed into range-variant and invariant terms, and then SI is performed. To coherently integrate the sub-images provided by BP-OKA, they should be accurately registered such that the same scatterers appear at the same pixel position. This holds when G-SI in sub-aperture focusing is conducted under the following conditions: 1) the range-invariant phase should be derived using the same bistatic polynomials, and 2) the resampling grid for G-SI should be the same for all Doppler frequencies.
In this study, we computed the parameters for G-SI using the bistatic polynomials of full apertures prior to performing subaperture BP-OKA. This preprocessing is a key consideration for performing BP-OKA within the framework of sub-image integration.
In addition to image focusing, we propose a method to clarify how illuminated areas are mapped to 2D imaging coordinates in general TI spaceborne BSAR scenarios. Using the range-variant properties of bistatic polynomials in the focusing process, we established a high-order equation for a certain target and concurrently retrieved its range and azimuth positions. The proposed method is fully compatible with the proposed focusing process, and accordingly, it can be used in various TI BSAR scenarios, including those with large curved paths and long baselines.
In Subsections III-B and III-C, we present details regarding the proposed image focusing method and its coordinate analysis, respectively.
represents the Doppler frequencies of full-aperture data, f BI Ref ,1 represents the reference Doppler centroid for a target positioned at the scene center (P Ref =[0, 0, 0] T ), and F FA represents the interval of f η . Note that F FA should be set to accommodate the skewed spectrum of the entire BSAR signal. Thus, F FA is larger than the original radar PRF of F 0 . Suppose that RFM is applied to (12); then, the residual 2D phase that introduces range-variant defocusing becomes To compensate for (14), G-SI is introduced, where (14) is first decomposed into rangevariant [V (P x )] and invariant terms [U f τ , f η ], and subsequently SI is performed. For phase decomposition, where r BI x = r BI x − r BI Ref is defined at the central time of the full-aperture data, i.e., η c , and r BI x = R BI (η c ; P x ) and r BI Ref = R BI (η c ; P Ref ) denote the bistatic ranges of P x and P Ref at η c , respectively. In (15), U f τ , f η is related to the G-SI operation and can be obtained via polynomial fitting with respect to the range frequency, i.e., f τ , using BI f τ , f η ; P x and r BI x , as follows: where χ w f η (w = 0, 1, 2, . . . , W ) represents the polynomial parameters of the range-invariant term U f τ , f η . Using U f τ , f η in (16), a new frequency for G-SI is defined as Because In (17), f G−SI τ,min and f G−SI τ,max are defined as respectively, where B r represents the radar operating bandwidth [39].
The key parameter in G-SI is χ w f η in (16), which defines f G−SI τ (f τ , f η ) according to phase decomposition of (14). As given by (16), χ w f η is obtained from BI f τ , f η ; P x and r BI x , and r BI x in the range-variant term is related to the range position of a target. This implies that if we compute χ w for each sub-aperture with different r BI x values, the target appears at different range positions for all the sub-aperture images. In addition, to have the same range bandwidth in the f G−SI τ domain, the same f G−SI τ,uniform values should be used for resampling in G-SI for all Doppler frequencies; otherwise, sub-aperture images have different range resolutions.
Here, prior to performing sub-aperture BP-OKA, we compute G-SI parameters using bistatic polynomials of fullaperture data in advance, where χ w f η is derived using the same r BI x value, and the resampling grid for G-SI, i.e., f G−SI τ,uniform in (18), is set to be same for all Doppler frequencies.
The predefined parameters, i.e., χ w f η and f G−SI τ,uniform , are then applied to each sub-aperture data for BP-OKA.

2) SUB-APERTURE FOCUSING
For sub-aperture focusing, we divide the full-aperture BSAR signal of (1) into L sub-apertures in the azimuth time domain as Ss BI,l f τ , η l ; P p , where represents the azimuth time of the l th sub-aperture, and η c,l and T l represent the central time and duration of η l , respectively. By taking the FT of Ss BI,l f τ , η l ; P p with respect to η l , we obtain SS BI,l f τ , f η,l ; P p , where represents the sub-aperture Doppler frequency, F 0 represents the original radar PRF, and f BI,l Ref,1 = −f c k BI,l Ref,1 /c represents the reference Doppler centroid at the l th sub-aperture.
In a Typical Monostatic SAR Geometry, The Sub-Aperture Spectrum is Fully Accommodated Within the Radar PRF, i.e., F 0 , Without Spectrum Aliasing. However, in the BSAR Case, the Doppler Centroids Change with Respect to f τ in the 2D Spectrum Domain (f τ -f η,l ): This property skews the spectral shape of the BSAR and introduces additional spectral aliasing, which is referred to as spectral folding, as shown in Fig. 3(a). Such spectral folding significantly degrades the image quality.
In most previous BSAR studies, LRWC was used to eliminate the effect of f τ -variant Doppler centroids, as follows: where By taking the azimuth FT of (23), the 2D spectrum SS 1 BI,l f τ , f η,l ; P p is obtained, where the spectrum folding phenomenon is corrected with the de-skewed shape of the spectrum, as shown in Fig. 3(b).
To maintain the signal property of TI BSAR, we extend the interval of sub-aperture Doppler frequencies from the radar PRF, i.e., F 0 , to F ext 0 with zero padding, and the result becomes SS 1 BI,l f τ , f ext η,l ; P p , as shown in Fig. 3(c), where In (25), F ext 0 should be set to a sufficiently large to accommodate the skewed shape of the sub-aperture Doppler spectrum. Subsequently, we restored the LRWC effect in the azimuth time domain using where Ss 1 BI,l f τ , η ext l ; P p represents the inverse FT (IFT) result of SS 1 BI,l f τ , f ext η,l ; P p , η ext l represents the azimuth time of the l th sub-aperture after the proposed Doppler extensions, and conj (·) denotes the conjugate operator. By taking the azimuth FT of (26), SS 2 BI,l f τ , f ext η,l ; P p is obtained; thus, the spectrum folding issue is addressed with the extended Doppler duration of F ext 0 while the azimuth-invariant signal property of the original TI-BSAR signal is maintained, as shown in Fig. 3(d). Consequently, by including (23)- (26) in the sub-image integration process, completely focused TI-BSAR images can be obtained via azimuth-invariant focusing of BP-OKA without spectrum aliasing.
For sub-aperture focusing with BP-OKA, RFM is first conducted, as follows: Then, G-SI frequencies for f ext η,l are retrieved using the predefined parameter χ w f η = f ext η,l : where Using (28), we can express (27) in the f G−SI τ domain as where r BI IM,p and η BI IM,p represent the range and azimuth positions of P p determined by the range-Doppler relationship To attain full-aperture resolutions after sub-image integration, the range of Doppler frequencies should be F FA , which includes the full-aperture Doppler spectrum of the TI BSAR signal. Here, we extend the f ext η,l of (31) to f FA η,l via zero padding as SS 2,G−SI By taking the 2D IFT of SS 2,G−SI BI,l f G−SI τ,uniform , f FA η,l ; P p , the l th sub-aperture image is obtained, as follows:

3) SUB-IMAGE INTEGRATIONS
For coherent integrations, Doppler centroid offsets among the sub-aperture images should be corrected. Let us consider subsequent 2D spectra of the l th 0 and l 1 = l 0 + 1 th sub-images 54468 VOLUME 11, 2023 Authorized licensed use limited to the terms of the applicable license agreement with IEEE. Restrictions apply. [i.e., I 0 BI,l 0 (τ, η) and I 0 BI,l 1 (τ, η)], as shown in Fig. 4(a). Although the Doppler centroids of these two images were different as f BI,l 0 Ref,1 and f BI,l 1 Ref, 1 , these values were defined in the same Doppler bin. This property overlapped with the 2D spectra after sub-image integration, yielding a failure of fullaperture spectrum recovery. Therefore, we corrected the Doppler frequency offsets for each sub-aperture image with respect to the reference image as follows: represents the l th Doppler centroid offset, and l Ref is the index for the reference subaperture image. After (34) was applied, the 2D spectra were shifted by f BI,l Ref , resulting in Fig. 4(b), and all the subaperture images could be coherently integrated as follows:

C. IMAGE COORDINATE ANALYSIS USING BISTATIC POLYNOMIALS
Because I 1 BI,FA (τ, η) in (35) is obtained via bistatic polynomial-based focusing, the conventional relationship of close-range and zero-Doppler time is invalid for its imaging coordinate analysis. Therefore, in this section, we clarify how  the TI BSAR images shown in (35) are interpreted geometrically after BP-OKA.
In BP-OKA, x-variant bistatic polynomials of P x are used. The relationship between the range-variant Doppler centroids, i.e., f BI DC r BI x , is then defined as where α and β represent the zeroth-and first-order polynomials of f BI DC r BI x , respectively. We aim to map P p to the range-azimuth image coordinates after BP-OKA using the VOLUME 11, 2023 54469 Authorized licensed use limited to the terms of the applicable license agreement with IEEE. Restrictions apply.  range-variant Doppler property in (36), as shown in Fig. 5, where the Doppler frequencies and bistatic ranges of P p reach f BI DC r BI x = r BI IM,p and r BI IM,p at η BI IM,p , simultaneously, as follows: Here, r BI IM,p = R BI η = η BI IM,p ; P p . As shown in Fig. 5, if we obtain η BI IM,p in (37), we can determine the range and azimuth positions of P p as r BI IM,p and η BI IM,p , respectively. Subsequently, the main focus of our imaging coordinate analysis is to obtain the η BI IM,p that satisfies (37). 54470 VOLUME 11, 2023 Authorized licensed use limited to the terms of the applicable license agreement with IEEE. Restrictions apply. Here, we substitute (3) and (4)    By numerically solving the high-order equation of (39), η BI IM,p = η BI IM,p + η c is obtained, and concurrently, r BI IM,p = R BI η BI IM,p ; P p is obtained. Fig. 6 shows a flowchart of the proposed process, wherein image focusing and coordinate analysis are conducted in the uhr and ti bsar geometries. As shown in Fig. 6, analysis of the computational complexity of image focusing facilitates estimation of the imaging time. In the proposed focusing process, a UHR image is obtained via simple summation of L sub-aperture images, wherein the parallel-computing process of sub-aperture focusing primarily determines the total computational cost. Thus, we focus on analyzing the computational complexity of sub-aperture focusing. Suppose that the L TH sub-aperture data are given by N η,l and N τ samples in the azimuth and range directions, respectively. After zero padding is applied twice in the azimuth direction, as given by (25) and (32), N η,l becomes N ext η,l and N FA η,l (I.E., N η,l < N ext η,l < N FA η,l ). In the proposed sub-aperture focusing, there are N η,l N τ and 2N ext η,l N τ phase multiplications, N τ N η,l -point fast fourier transforms (FFTS), 2N τ N ext η,lpoint FFTS, N τ N FA η,l -point FFTS, N ext η,l N τ -Point FFTS, and N ext η,l interpolations with 2 (M − 1), where M represents the interpolation kernel length. Therefore, the computational load Of the L TH sub-aperture focusing is

IV. RESULTS
The effectiveness of the proposed method was verified using numerical simulations. To obtain the simulated dataset, we considered the TI BSAR geometry with curved paths, as shown in Figs. 1 and 7. The transmitter S T and receiver S R were separated in the cross-track (x) and along-track (y) directions by a baseline length of b BI = 381.9 km, and they moved along parallel curved paths at the same height of 550 km with identical velocities of 7583 m/s. The squint angle of the transmitter was θ T = 0.01 • , and that of the receiver was θ R = 21.88 • ; this resulted in a large bistatic angle of θ BI = 22.73 • . A total of 441 point targets were deployed within a 3.2 km × 2.8 km scene, as shown in Fig. 8. For UHR data collection, the targets were illuminated by a staring spotlight, where the transmitted radar bandwidth was 1 GHz and the Doppler bandwidth was 27 kHz. The simulation parameters are presented in Table 1.
To evaluate the key contributions of our proposals, in Section IV-A, we examine the focusing results of BP-OKA combined with extended sub-image integrations, and in Section IV-B, we verify the proposed analysis for 2D image coordinates using bistatic polynomials.

A. FOCUSING ANALYSIS USING BISTATIC POLYNOMIALS 1) COUMPATION OF G-SI PARAMETERS FOR FULL APERTURES
For the 441 point targets shown in Fig. 8, echo data were collected from the UHR and TI BSAR operations with the parameters presented in Table 1. After range compression, we divided the data into 10 sub-apertures (l = 1, 2, . . . , 10) in the azimuth directions, and concurrently, G-SI parameters [χ w f η and f G−SI τ,uniform ] were computed for full-aperture Doppler frequencies. Using predefined parameters, subaperture images were obtained via BP-OKA and coherently integrated to recover the full azimuth resolutions. Fig. 9 presents the subsequent coherent integration results of the sub-images for the central target. As shown, the azimuth resolutions of the target improved as the amount of subaperture synthesis increased. Fig. 10(a) shows the final focusing result of the proposed focusing process; its 2D frequency spectrum is shown in Fig. 10(b). For more specific verification, we enlarged the impulse response functions (IRFs) of five targets [P A , P B , P C , P D , P E ], as marked with yellow circles in Figs. 8 and 10(a). Enlarged IRFs of the five targets are shown in Fig. 11. Likewise, for quantitative analysis, image quality parameters of 441 targets were measured, including the resolution, peak sidelobe ratio (PSLR), and integrated sidelobe ratio (ISLR). The means and standard deviations of the measurements are presented in Table 2. As indicated by the results, all the targets were completely focused, verifying that the proposed method can achieve spatially variant focusing in UHR and TI BSAR geometries.

2) EFFECT OF PRIOR COMPUTATION OF G-SI PARAMETERS ON SUB-IMAGE INTEGRATIONS
In this subsection, we examine the effect of predefined parameters, i.e., χ w f η and f G−SI τ,uniform , on the final focusing performance of sub-image integrations.
In contrast to the proposed processing chain shown in Fig. 6, in each sub-aperture processing, we computed χ w f ext η,l in (29) using different r BI x,l = r BI x,l − r BI Ref,l (l = 1, 2, . . . , 10) values, where r BI x,l = R BI η c,l ; P x and r BI Ref,l = R BI η c,l ; P Ref denote bistatic ranges of P x and P Ref at the central time of the l th sub-aperture data, i.e., η c,l . Next, f G−SI τ (f τ , f ext η,l ) was defined by (28), and different f G−SI τ,uniform (l) (l = 1, 2, . . . , 10) grids were used for the l th subaperture G-SI, where and Fig. 12 shows the enlarged IRFs of the five targets [P A ,P B ,P C ,P D ,P E ] after sub-image integration without consideration of prior computations of the G-SI parameters. As shown in Fig. 12(c), the central target P C was always located at the center of the image after RFM, attaining a full-aperture azimuth resolution with coherent integration. In contrast, because of different G-SI conditions, other targets [P A ,P B ,P D ,P E ] appeared at different range pixels for every sub-aperture image, resulting in the failure of the coherent accumulation of sub-aperture images, as shown in Figs. 12(a), (b), (d), and (e). These results verified that as shown in Fig. 6, G-SI parameters, i.e., χ w f η and f G−SI τ,uniform , should be defined a priori for full apertures before sub-aperture focusing. Otherwise, the sub-aperture images cannot be properly integrated, because of incorrect registration results.

3) COMPARISION WITH CONVENTIONAL TI BSAR FOCUSING
We compared the image quality performance of the proposed and conventional methods in [23] under the simulation conditions presented in Table 1. Sub-aperture recombination was first performed in the signal domain by exploiting [31], and then the conventional focusing approach of [23] derived from double HRMs was applied to the recombined data. Fig. 13 shows the five IRFs, i.e., [P A , P B , P C , P D , P E ], after conventional focusing, and the means and standard deviations of image quality parameters for 441 targets are presented in Table 3. In contrast to the result of the proposed focusing method, as shown in Fig. 11 and Table 2, targets are significantly blurred-particularly in the azimuthal direction-except for the central target. This result indicates 54472 VOLUME 11, 2023 Authorized licensed use limited to the terms of the applicable license agreement with IEEE. Restrictions apply.    that double HRMs are not sufficiently accurate to represent the actual range histories of the TI BSAR geometry with large curved paths, preventing complete spatially variant focusing for the entire illuminated area. Consequently, the conventional method is unsuitable for TI and UHR spaceborne BSAR applications, whereas the proposed method provides superior image quality with complete focusing results.

B. IMAGE COORDINATE ANALYSIS USING BISTATIC POLYNOMIALS
In this subsection, the proposed image coordinate analysis method based on bistatic polynomial focusing is investigated. The range and azimuth locations that reached the peak values for the IRFs of the 441 point targets in Fig. 10(a) were measured; all the IRFs were oversampled 20 times, and their estimates were obtained by numerically solving the proposed high-order equation of (39). Figs. 14(a) and (b) show the absolute errors between the measurements and estimates in the range and azimuthal directions, respectively, and their means and standard deviations are presented in Table 4. For both the range and azimuth directions, the location errors were close to zero, with a small bias and deviation. This implies that after bistatic polynomial focusing, targets were correctly mapped to the 2D image coordinates using the relationship in (37). Consequently, the proposed mapping rule allows accurate geolocation analysis after bistatic polynomial focusing, even in a TI BSAR geometry with large baseline and orbital arcs.

V. CONCLUSION
A novel processing flow of image focusing, and its coordinate analysis are proposed for TI and UHR BSAR operations, where a UHR image is obtained via the coherent accumulation of sub-aperture images, and its 2D image coordinates are interpreted. The proposed processing flow makes the following contributions. First, to address the additional aliasing of the TI BSAR spectrum after sub-aperture division, we introduced a Doppler extension, in which LRWC and its restorations were subsequently performed. This extended concept efficiently addresses the aliasing issue while maintaining the azimuth-invariant signal properties of the TI BSAR. Second, we introduced BP-OKA, which completely achieved spatially variant focusing in the TI BSAR geometry, even with large curved paths and a long baseline. For exploiting BP-OKA in the framework of sub-image integration, a key consideration was highlighted, where G-SI parameters are defined for full-aperture Doppler frequencies prior to performing sub-aperture imaging. In its absence, the same scatterers appear at different range positions in each subaperture image, causing failure of coherent accumulation. Finally, to retrieve the range-azimuth coordinates of the targets within an illuminated scene, we established a highorder equation using bistatic polynomials, and their 2D image coordinates were interpreted. This new concept is compatible with the proposed focusing process and can be used in various spaceborne TI BSAR scenarios. Numerical simulations demonstrated that the proposed processing flow can achieve complete spatially variant focusing with accurate image coordinate interpretations in UHR and TI BSAR geometries, and all the aforementioned contributions were verified.