Transient Scattering Echo Simulation and ISAR Imaging for a Composite Target-Ocean Scene Based on the TDSBR Method

: We propose an inverse synthetic aperture radar (ISAR) imaging algorithm for a composite target-ocean scene based on time-domain shooting and bouncing rays (TDSBR) method. To do this, we develop the TDSBR method to simulate the electromagnetic scattering echoes from a composite target-ocean scene. Then, we derive the ISAR imaging formulas based on electromagnetic scattering ﬁeld expressions. The demodulation, pulse compression, and azimuth inverse fast Fourier transform are performed on the scattering echo data to obtain ISAR images. Linear frequency modulation (LFM) pulses are usually used in traditional radar imaging. However, its long pulse width leads to more computation time for simulations. To increase efﬁciency, we replace LFM pulses with modulating Gaussian pulses. As a result, we modify the matched ﬁlter to achieve pulse compression under modulating Gaussian pulse excitations. Simulation results conﬁrm that this proposed method improves efﬁciency while ensuring image quality. This proposed technique fully incorporates the electromagnetic scattering mechanism and is accurate. It effectively combines the electromagnetic scattering algorithm and the radar imaging algorithm. Finally, we present ISAR images simulated by this method. These results conﬁrm the feasibility and efﬁciency of this technique.


Introduction
As active microwave imaging radar systems, synthetic aperture radar (SAR) and inverse synthetic aperture radar (ISAR) have been widely used in military and civil applications because they can be used continuously [1,2]. At the present, researchers pay more attention to the imaging algorithms, and few research works exist on the acquisition of the radar echo. It is undeniable, however, that the accurate and efficient acquisition of the radar echo data plays an essential role in designing radar imaging systems and the analytics of imaging algorithms. The echo data is needed for the system design to determine several parameters and evaluate its performance. During development, echo data is also required to verify the performance of each part of the system. Many researchers have studied imaging algorithms based on measured data [3][4][5][6][7][8][9][10]. Although the measured data is accurate, this approach is often very costly. For instance, the imaging of large targets requires a large anechoic chamber and a sizeable measured target model, which are both expensive to build and take a long time to construct. For the space-borne SARs, it is impossible to collect actual data before the system is developed successfully. Benefiting from the rapid development of computer technology, computational electromagnetic simulations are often used instead of experimental measurements to obtain the scattering characteristics of the targets. The electromagnetic simulation incorporates the scattering mechanism between the electromagnetic wave and the target. Such accurate simulation solutions address the issues attributed to the experimental measurement.
The research on microwave imaging based on computational electromagnetics has been developed over a long period of time. In terms of target imaging, R. Bhalla and H. Ling [11] presented an ISAR imaging technique using bistatic scattering field data based on the physical optics (PO) method. They also used a shooting and bouncing rays (SBR) method to calculate multiple scattering data. They showed that in cases where multiple scattering is the dominant scattering mechanism, ghost artifacts might appear in both the monostatic and bistatic images. In 2010, Miao Sui [12] developed a near-field electromagnetic scattering method to predict the high-frequency wideband EM scattering from the targets in the near-field region. They then proposed an approach for a twodimensional high-resolution ISAR image sequence generation and synthesis. For composite scattering imaging, Chen [13] proposed an efficient hybrid algorithm to simulate the scattering from a ship-like target on the sea surface. The scattering model is well suited to deal with the SAR imagery simulation of a large ship on the sea surface. Ye Zhao [14] developed a bistatic SAR imaging method for metallic targets in the ocean. The image intensity distribution for the ocean surface is analyzed, and the facet scattering model is used to give the individual returns from separate facets. In 2017, Feng Xu [15] presented a simplified small perturbation method (SPM) to evaluate the scattering of finite-length conducting rough surfaces. The simplified SPM greatly facilitated the derivation of the analytical expression of the SAR image. Currently, the radar echo is often obtained using the frequency-domain electromagnetic scattering method. High range resolution is then achieved by employing frequency sweep. However, the actual electromagnetic waves emitted by radar are mostly pulse signals with a definite frequency bandwidth and a higher carrier frequency. Compared with frequency-domain (FD) electromagnetic scattering methods, TD methods can avoid the time-consuming frequency-scanning process. The total wideband scattering results can be obtained by a single Fourier transform of TD scattering echoes, which can significantly improve the computational efficiency; thus, TD electromagnetic scattering methods are more suitable for simulating radar echoes.
There is extensive research on the TD electromagnetic scattering methods. There are many classical TD numerical methods, such as the finite-difference time-domain (FDTD) [16,17] and the TD integral equation [18,19]. These methods can be used to obtain the accurate scattering echo of the target. However, with the increase of the electrical size of a target, the computational resources and simulation time required by these numerical methods are sharply increased. Hence, they are not suitable for the scattering calculation of larger targets. In order to obtain the TD scattering echoes of large targets, researchers developed a series of TD high-frequency scattering methods based on FD high-frequency scattering methods. Among them, the TDSBR method effectively combines the advantages of TD geometrical optics (TDGO) and TD physical optics (TDPO) [20,21], and fully considers the existing reflection fields on the target. This method is widely used in the electromagnetic scattering of electrically large and complex targets. In 1993, Ling and Bhalla [22] proposed a TDSBR method based on the classical FDSBR algorithm. In [22], they then derived a closed-form TD ray tube integral expression to calculate the TD response of complex targets quickly. In 2015, Zhou [23] presented a TDSBR method based on the beam-tracing (BT) technique to analyze the transient scattering responses from large PEC objects illuminated by a pulsed plane wave. In 2019, Guo [24] extended the application scope of the TDSBR method from far-field scattering to near-field scattering and proposed a near-field TDSBR method for large and complex targets. In radar imaging, Guo [25,26] proposed ISAR imaging technology based on the TDPO and TDSBR methods, respectively. However, the previously mentioned TD methods focus on a single target. There are few studies on transient scattering echo simulations of composite scenes. In particular, to the best of our knowledge, there is no research on the radar imaging of composite scenes by using TD simulated scattering echoes.
Because of the high cost of experimental measurements and the difficulty in obtaining non-cooperative target image data, we propose an ISAR imaging method for a composite target-ocean scene by using simulated scattering echoes. As the electromagnetic wave emitted by radar is a TD pulse signal with a frequency bandwidth, we use the TD electromagnetic scattering method to simulate the scattering echo. The difficulty of a composite Remote Sens. 2022, 14, 1183 3 of 18 scene scattering echo simulation is that the coupling scattering between the target and the sea surface is very complex. In this paper, based on a ray-tracing model, we propose a TDSBR method for calculating the transient scattering echo from a composite target-ocean scene. In the TDSBR method, the coupling scattering field comes from the radiation field of the equivalent electromagnetic current generated by the electromagnetic wave repeatedly bouncing between the target and the sea surface. The scattering field expressions of the conducting target and the dielectric sea surface are derived. Then, we derive the ISAR imaging formulas based on the electromagnetic scattering echo expressions used in TDSBR. In addition, to improve the simulation efficiency, we replace the LFM pulse with a modulating Gaussian pulse. We then modify the matched filter to achieve pulse compression under a modulating Gaussian pulse excitation.

Description of TDSBR Method
In Figure 1, we consider the composite target-ocean scene illuminated by an incident plane wave. The scattering echo consists of three parts: target scattering, sea surface scattering, and coupled scattering. The scattering echoes from these three parts can be easily obtained using the TDSBR method. The TDSBR method consists of both TDGO and TDPO. The main goal of the TDSBR method is to imitate the propagation behavior of electromagnetic waves by using rays. As shown in Figure 1, the initial rays hit the surface of the target and the sea, respectively. When a ray intersects with the surface of an object, it reflects. The reflected rays continue to travel forward, intersecting objects until no more reflections occur. The TDGO method is used to calculate the intensity and phase value of the reflected rays. The object's surface that intersects with the rays is thought to be illuminated by electromagnetic waves. The scattering field of the illuminated surface is calculated using the TDPO method. wave emitted by radar is a TD pulse signal with a frequency bandwidth, we use the TD electromagnetic scattering method to simulate the scattering echo. The difficulty of a composite scene scattering echo simulation is that the coupling scattering between the target and the sea surface is very complex. In this paper, based on a ray-tracing model, we propose a TDSBR method for calculating the transient scattering echo from a composite target-ocean scene. In the TDSBR method, the coupling scattering field comes from the radiation field of the equivalent electromagnetic current generated by the electromagnetic wave repeatedly bouncing between the target and the sea surface. The scattering field expressions of the conducting target and the dielectric sea surface are derived. Then, we derive the ISAR imaging formulas based on the electromagnetic scattering echo expressions used in TDSBR. In addition, to improve the simulation efficiency, we replace the LFM pulse with a modulating Gaussian pulse. We then modify the matched filter to achieve pulse compression under a modulating Gaussian pulse excitation.

Description of TDSBR Method
In Figure 1Error! Reference source not found., we consider the composite targetocean scene illuminated by an incident plane wave. The scattering echo consists of three parts: target scattering, sea surface scattering, and coupled scattering. The scattering echoes from these three parts can be easily obtained using the TDSBR method. The TDSBR method consists of both TDGO and TDPO. The main goal of the TDSBR method is to imitate the propagation behavior of electromagnetic waves by using rays. As shown in Figure 1Error! Reference source not found., the initial rays hit the surface of the target and the sea, respectively. When a ray intersects with the surface of an object, it reflects. The reflected rays continue to travel forward, intersecting objects until no more reflections occur. The TDGO method is used to calculate the intensity and phase value of the reflected rays. The object's surface that intersects with the rays is thought to be illuminated by electromagnetic waves. The scattering field of the illuminated surface is calculated using the TDPO method. Several efficient ray-tracing techniques have been developed to find the illuminated surface of an object. In this paper, the surface of the target and sea is subdivided into triangular patches. As shown in Figure 2, we introduce the forward and backward ray tracing technology, as was previously done in [27], into the TDSBR. Forward ray tracing is used to obtain the illuminated triangle; then, backward ray tracing is used to judge whether the triangle adjacent to the illuminated triangle is illuminated. Several efficient ray-tracing techniques have been developed to find the illuminated surface of an object. In this paper, the surface of the target and sea is subdivided into triangular patches. As shown in Figure 2, we introduce the forward and backward ray tracing technology, as was previously done in [27], into the TDSBR. Forward ray tracing is used to obtain the illuminated triangle; then, backward ray tracing is used to judge whether the triangle adjacent to the illuminated triangle is illuminated. Remote Sens. 2022, 14, x FOR PEER REVIEW 4 of 19 Next, we calculate the scattering fields of the illuminated triangles. Before calculating the scattering fields, the TDGO method is used to obtain the incident fields of the triangles. The electric and magnetic fields of the initial incident plane wave with the excitation pulse ( ) p t are as follows: and  Next, we calculate the scattering fields of the illuminated triangles. Before calculating the scattering fields, the TDGO method is used to obtain the incident fields of the triangles. The electric and magnetic fields of the initial incident plane wave with the excitation pulse p(t) are as follows: where e i 0 and h i 0 are amplitude vectors.k i 0 denotes the normalized incidence direction. δ(t) is the impulse function. c is the speed of light in free space. The symbol " * " represents the convolution operation.
For the mth reflection, the incident fields In the above equation, r c m is the center of triangular patch. ω =k i m −k s and β = ω − (ω ·n m )n m . ∆v i = v i+1 − v i . ν i is the position vector of the vertices of Tri m in the global coordinate system, and v 4 = v 1 . ε(t) is a step function.
To verify the accuracy and efficiency of the TDSBR method, results from TDSBR are shown and compared with the full-wave solutions that were simulated using the multilevel fast multipole algorithm (MLFMA) obtained using FEKO commercial software. In addition, the results of the frequency domain SBR (FDSBR) are given as a reference. In the TDSBR method, the excitation pulse is the modulating Gaussian pulse: where f 0 is the carrier frequency. T = 2τ. τ = 4/ f b . f b denotes the bandwidth. Figure 3 shows the definition of the incident and the observation direction in the global coordinate system. Figure 4 shows the monostatic (θ i = 45 • , ϕ i = 0 • ) and bistatic (θ i = 45 • , ϕ i = 0 • , θ s = 45 • , ϕ s = 180 • ) scattering results of a composite cube-ocean scene (see Figure 4a). The side length of the cube is 0.2 m. The height of the cube's center is 0.15 m. The size of the rough sea surface is 1 m × 1 m. Figure 5 shows the monostatic (θ i = 10 • , ϕ i = 90 • ) and bistatic (θ i = 45 • , ϕ i = 0 • , θ s = 45 • , ϕ s = 180 • ) scattering results of a composite missileocean scene (see Figure 5d). The size of the missile is about 0.33 m × 0.1 m × 0.026 m. The height of the missile's center is 0.1 m. The size of the rough sea surface is 1 m × 1 m. The power spectrum of the sea surface is PM spectrum. The wind speed is 3 m/s. The relative dielectric constant of the sea is 55.9 − 37j. The rough sea surface is generated by the Monte Carlo method. The carrier frequency f 0 is 6 GHz and the bandwidth is 4 GHz.
In Figures 4 and 5, the wideband scattering fields obtained by MLFAM and FDSBR are converted into TD by the inverse fast Fourier transform (IFFT) and compared with those obtained from TDSBR. The transient results by TDSBR are transformed to the FD using the fast Fourier transform (FFT). It can be seen from the figures that the results from TDSBR and FDSBR agree well, but are slightly different from those using MLFMA; TDSBR and FDSBR belong to high frequency algorithms and the error between these results and those from MLFMA is due to the use of high frequency approximation. For wideband radar cross section (RCS), we give the root mean square error (RMSE) between TDSBR-FFT and MLFMA. The RMSE is expressed as: where X and Y represent the two sets of data to be compared. The "n" is the number of elements in the data.  In Figure 4 and 5, the wideband scattering fields obtained by MLFAM and FDSBR are converted into TD by the inverse fast Fourier transform (IFFT) and compared with those obtained from TDSBR. The transient results by TDSBR are transformed to the FD using the fast Fourier transform (FFT). It can be seen from the figures that the results from TDSBR and FDSBR agree well, but are slightly different from those using MLFMA; TDSBR and FDSBR belong to high frequency algorithms and the error between these results and those from MLFMA is due to the use of high frequency approximation. For wideband radar cross section (RCS), we give the root mean square error (RMSE) between TDSBR-FFT and MLFMA. The RMSE is expressed as: where X and Y represent the two sets of data to be compared. The "n" is the number of elements in the data. Table 1 compares the computational costs between the three methods. The experiment platform is a desktop computer with Quad 3.2 GHz CPU and 16 GB memory. The simulation time is estimated based on the single thread using Windows 10. As we can see, the ratio of computation time of TDSBR and FDSBR is about 1:5, and MLFMA takes 10 6 times longer than TDSBR.     Table 1 compares the computational costs between the three methods. The experiment platform is a desktop computer with Quad 3.2 GHz CPU and 16 GB memory. The simulation time is estimated based on the single thread using Windows 10. As we can see, the ratio of computation time of TDSBR and FDSBR is about 1:5, and MLFMA takes 10 6 times longer than TDSBR.

ISAR Imaging Algorithm
Section 2 shows that the total echo of the composite scene is the sum of the echoes from all of the illuminated triangular patches. In this section, the ISAR imaging algorithm is derived based on the scattering field expression (12) of Section 2.
As it is shown in Figure 6, we consider a plane wave scanning a small angle ∆ϕ alonĝ ϕ near ϕ 0 . In this example, the angular velocity is ω A . Thus, the azimuth angle ϕ is wherek and Collecting radar returns at different look angles to form the cross-range profile.
Assuming the polarization direction of the receiving antenna is p , the received echo is Assuming the polarization direction of the receiving antenna isp, the received echo is where w a (η) is a window function.
In general, |β| and ω · ∆v i in Equation (16) are not equal to zero for most observation directions. Thus, I(t, η) can be written as the superposition of the step function: In general, the transmitted signal p(t) with carrier frequency f 0 can be expressed as: where g(t) is the baseband signal and T is the pulse width.
In order to avoid the integral operation in the convolution operation, we use the first derivative of p(t) as the excitation impulse. The echo expression can then be rewritten as: We then derive the imaging formula under small-bandwidth small-angle approximation based on (25). First, the echo needs to be demodulated to the baseband domain. This procedure can be completed by multiplying the echo signal by the phase factor exp(−j2π f 0 t). The echo after demodulation can be written as: Then, by taking the range Fourier transform of s 0 (t, η), we can get The next step is matched filtering. For the LFM pulse, the matched filter is usually set to: The spectrum G( f ) of LFM is a window function of complex-type. However, since the spectrum of the modulating Gaussian pulse is a Gaussian-type function, the matching filter is unable to achieve the pulse compression effect. Therefore, we redesign the matched filter as follows: where w r ( f ) is a window function. It can be seen from the filter's expression that the filter is suitable for modulating the Gaussian pulse and for the LFM pulse.
The output of the matched filter is For small values of ω A η, cos(ω A η) can be approximated by 1 and sin(ω A η) can be approximated by ω A η. Therefore,k i 0 (η) is reduced tô and A(η) can be approximated as Therefore, S 1 ( f , η) can be simplified as: By further applying the small-bandwidth approximation: we can get where λ 0 = c/ f 0 . Finally, we take the range and azimuth inverse Fourier transform (IFT) of S 1 ( f , η). We can obtain the final ISAR image expression: where p r and p a are sinc-type functions. The image intensity of each triangular patch can be obtained by (38). The imaging process for one of the triangles is the same for all the triangles. The overall image is the sum of the image of all the illuminated triangles.

Numerical Examples
We present the ISAR images generated using simulated echo data from TDSBR. A modulating Gaussian pulse is used as the excitation pulse. For comparison, we also provide images generated using the LFM pulse and the FDSBR method to validate the accuracy and efficiency of the proposed method. The expression of the LFM pulse is: where K is the linear frequency and K = 5 × 10 15 Hz/s in this paper.

ISAR Images of Targets
At first, we consider the ISAR images of a dihedral corner reflector with a side length of 1 m. The plane wave is polarized in HH. The carrier frequency f 0 = 10 GHz and the bandwidth is set to 1 GHz. The center azimuth angle ϕ 0 = 0 • . Based on the azimuth resolution, the look angle of radar varies from −0.05 rad to 0.05 rad around ϕ 0 . Therefore, the resolution is 0.15 m in both range and azimuth. Figure 7 shows the geometric model of the dihedral corner reflector.

Numerical Examples
We present the ISAR images generated using simulated echo data from TDSBR. A modulating Gaussian pulse is used as the excitation pulse. For comparison, we also provide images generated using the LFM pulse and the FDSBR method to validate the accuracy and efficiency of the proposed method. The expression of the LFM pulse is: where K is the linear frequency and K = 5 × 10 15 Hz/s in this paper.

ISAR Images of Targets
At first, we consider the ISAR images of a dihedral corner reflector with a side length    Figure 9shows the ISAR images and profile maps in range and azimuth under different incident angles 0 θ . As shown in Figure 9, the images at different incident angles are consistent. For a regular dihedral angle reflector, no matter what the incident angle is, the electromagnetic wave propagation distance caused by double bounces is the same. This is equivalent to the backscattering of an imaginary scatterer at the corner of the reflector. Figure 10shows the ISAR images at 0 45 θ =° by using the LFM pulse and the FDSBR method. The two images shown in Figure 10are almost identical to those shown in Figure 9a.   Figure 9 shows the ISAR images and profile maps in range and azimuth under different incident angles θ 0 . As shown in Figure 9, the images at different incident angles are consistent. For a regular dihedral angle reflector, no matter what the incident angle is, the electromagnetic wave propagation distance caused by double bounces is the same. This is equivalent to the backscattering of an imaginary scatterer at the corner of the reflector. Figure 10 shows the ISAR images at θ 0 = 45 • by using the LFM pulse and the FDSBR method. The two images shown in Figure 10   Next, we simulate the ISAR images of a complex missile target as shown in Figure  11. The size of the missile is 6.4 m × 2.12 m × 0.5 m. The carrier frequency is 0 f = 10 GHz.   Next, we simulate the ISAR images of a complex missile target as shown in Figure  11. The size of the missile is 6.4 m × 2.12 m × 0.5 m. The carrier frequency is 0 f = 10 GHz.   Next, we simulate the ISAR images of a complex missile target as shown in Figure  11. The size of the missile is 6.4 m × 2.12 m × 0.5 m. The carrier frequency is 0 f = 10 GHz. Next, we simulate the ISAR images of a complex missile target as shown in Figure 11. The size of the missile is 6.4 m × 2.12 m × 0.5 m. The carrier frequency is f 0 = 10 GHz. The incident angle is θ 0 = 45 • and the center azimuth angle is ϕ 0 = 0 • . The polarization is VV. The resolution is 0.15 m in both range and azimuth. The incident angle is 0 45 θ =° and the center azimuth angle is 0 0 φ =°. The polarization is VV. The resolution is 0.15 m in both range and azimuth.   Figure 13shows the normalized ISAR images and profile maps of the missile target using a modulating Gaussian pulse, a LFM pulse, and the FDSBR method. As shown, these images are in perfect agreement. By comparing Figure  12a and 13a, it can be found that the distribution of bright spots on the image corresponds to the distribution of the pulse peaks of the echo in the range profile. From the images of Figure 13, we can see the main bright spots located around the nose, wings, and tail of the missile. Since the warhead has a hemispherical shape and a substantial backscatter contribution, it forms a bright spot in the image. In the image of wings, the appearance of one of the strips is also enhanced. This brightness occurs because the wings have a specific thickness, as shown in Figure 12. The plate of the wing parallel to the yoz plane is exposed to the electromagnetic wave and contributes to the backscattering. Therefore, the main reason for the formation of the brighter strip is that it contains the scattering contribution from the wing's plate parallel to the yoz plane.   Figure 12 shows the normalized transient scattering echoes and wideband scattering fields of the missile at ϕ 0 = 0 • . Figure 13 shows the normalized ISAR images and profile maps of the missile target using a modulating Gaussian pulse, a LFM pulse, and the FDSBR method. As shown, these images are in perfect agreement. By comparing Figures 12a  and 13a, it can be found that the distribution of bright spots on the image corresponds to the distribution of the pulse peaks of the echo in the range profile. From the images of Figure 13, we can see the main bright spots located around the nose, wings, and tail of the missile. Since the warhead has a hemispherical shape and a substantial backscatter contribution, it forms a bright spot in the image. In the image of wings, the appearance of one of the strips is also enhanced. This brightness occurs because the wings have a specific thickness, as shown in Figure 12. The plate of the wing parallel to the yoz plane is exposed to the electromagnetic wave and contributes to the backscattering. Therefore, the main reason for the formation of the brighter strip is that it contains the scattering contribution from the wing's plate parallel to the yoz plane.   Figure 13shows the normalized ISAR images and profile maps of the missile target using a modulating Gaussian pulse, a LFM pulse, and the FDSBR method. As shown, these images are in perfect agreement. By comparing Figure  12a and 13a, it can be found that the distribution of bright spots on the image corresponds to the distribution of the pulse peaks of the echo in the range profile. From the images of Figure 13, we can see the main bright spots located around the nose, wings, and tail of the missile. Since the warhead has a hemispherical shape and a substantial backscatter contribution, it forms a bright spot in the image. In the image of wings, the appearance of one of the strips is also enhanced. This brightness occurs because the wings have a specific thickness, as shown in Figure 12. The plate of the wing parallel to the yoz plane is exposed to the electromagnetic wave and contributes to the backscattering. Therefore, the main reason for the formation of the brighter strip is that it contains the scattering contribution from the wing's plate parallel to the yoz plane.    Figure 13shows the normalized ISAR images and profile maps of the missile target using a modulating Gaussian pulse, a LFM pulse, and the FDSBR method. As shown, these images are in perfect agreement. By comparing Figure  12a and 13a, it can be found that the distribution of bright spots on the image corresponds to the distribution of the pulse peaks of the echo in the range profile. From the images of Figure 13, we can see the main bright spots located around the nose, wings, and tail of the missile. Since the warhead has a hemispherical shape and a substantial backscatter contribution, it forms a bright spot in the image. In the image of wings, the appearance of one of the strips is also enhanced. This brightness occurs because the wings have a specific thickness, as shown in Figure 12. The plate of the wing parallel to the yoz plane is exposed to the electromagnetic wave and contributes to the backscattering. Therefore, the main reason for the formation of the brighter strip is that it contains the scattering contribution from the wing's plate parallel to the yoz plane.

ISAR Images of Composite Target-Ocean Scene
In the following example, we generate ISAR images of a composite target-ocean scene. First, we present the ISAR images of a PEC missile target above a rough sea surface. As shown in Figure 14, the missile's dimensions are the same as in Figure 11. For this case, the size of the sea surface is 20 m × 20 m, and h is the distance from the missile to the sea. The power spectrum of the sea surface is PM spectrum. The wind speed is 5 m/s. The relative dielectric constant of the sea is 55.9 − 37j. The incident angle is θ 0 = 60 • and the center azimuth angle is ϕ 0 = 0 • . The polarization is HH. The carrier frequency is f 0 = 10 GHz. The resolution is 0.15 m in both range and azimuth.

ISAR Images of Composite Target-Ocean Scene.
In the following example, we generate ISAR images of a composite target-ocean scene. First, we present the ISAR images of a PEC missile target above a rough sea surface. As shown in Figure 14, the missile's dimensions are the same as in Figure 11. For this case, the size of the sea surface is 20 m × 20 m, and h is the distance from the missile to the sea. The power spectrum of the sea surface is PM spectrum. The wind speed is 5 m/s. The relative dielectric constant of the sea is 55.9 37j − . The incident angle is 0 60 θ =° and the center azimuth angle is 0 0 φ =°. The polarization is HH. The carrier frequency is 0 f = 10GHz. The resolution is 0.15 m in both range and azimuth. According to the ray paths shown in Figure 1, the scattered echoes from different parts of the composite scene can be easily obtained. It is helpful to analyze the influence of the rough background on the images. Figure 15Error! Reference source not found. shows the ISAR images generated by total echo, differential echo, target echo, and coupling echo, respectively, at h = 3 m. The differential echo is obtained by subtracting the sea surface echo from the total echo. Equivalently, the differential echo is also the sum of the target echo and the coupling echo. The target echo is derived from the scattering contribution of the path ○ 1 in Figure 1, and the coupling echo corresponds to the path ○ 3 in Figure 1Error! Reference source not found.. These images are normalized and plotted in logarithmic scale with a dynamic range. As can be seen from Figure 15Error! Reference source not found.b, the coupling image is located behind the target image in the range axis. This occurs because the rays exhibit multiple bounces between the target and the sea surface. These bounces delay the echo in time and cause the coupling image to be delayed According to the ray paths shown in Figure 1, the scattered echoes from different parts of the composite scene can be easily obtained. It is helpful to analyze the influence of the rough background on the images. Figure 15 shows the ISAR images generated by total echo, differential echo, target echo, and coupling echo, respectively, at h = 3 m. The differential echo is obtained by subtracting the sea surface echo from the total echo. Equivalently, the differential echo is also the sum of the target echo and the coupling echo. The target echo is derived from the scattering contribution of the path 1 in Figure 1, and the coupling echo corresponds to the path 3 in Figure 1. These images are normalized and plotted in logarithmic scale with a dynamic range. As can be seen from Figure 15b, the coupling image is located behind the target image in the range axis. This occurs because the rays exhibit multiple bounces between the target and the sea surface. These bounces delay the echo in time and cause the coupling image to be delayed in the range axis. The coupling image is the mirror image of the target. Because the sea surface is rough, this modulation effect causes the coupling image to be divergent. in the range axis. The coupling image is the mirror image of the target. Because the sea surface is rough, this modulation effect causes the coupling image to be divergent.    Figure 16Error! Reference source not found. shows the ISAR images at h = 5 m. With the increase of the missile height, the coupling path between the target and the sea surface becomes longer. By comparing Figure 15and 16, it can be found that, with the increase of h, the coupling image becomes more dispersed on the range axis and the image intensity becomes weaker. In addition, to verify the accuracy and efficiency of the proposed method, we also present the ISAR images shown in Figure 17Error! Reference source not found. by using the LFM pulse and the FDSBR method. According to this comparison, the images obtained by using the three different excitation sources (modulating Gaussian pulse, LFM pulse, and FDSBR) are in perfect agreement. Section 4.3 gives a detailed analysis of the image similarity and the simulation efficiency.  In addition, to verify the accuracy and efficiency of the proposed method, we also present the ISAR images shown in Figure 17 by using the LFM pulse and the FDSBR method. According to this comparison, the images obtained by using the three different excitation sources (modulating Gaussian pulse, LFM pulse, and FDSBR) are in perfect agreement. Section 4.3 gives a detailed analysis of the image similarity and the simulation efficiency. Figure 16Error! Reference source not found. shows the ISAR images at h = 5 m. With the increase of the missile height, the coupling path between the target and the sea surface becomes longer. By comparing Figure 15and 16, it can be found that, with the increase of h, the coupling image becomes more dispersed on the range axis and the image intensity becomes weaker. In addition, to verify the accuracy and efficiency of the proposed method, we also present the ISAR images shown in Figure 17Error! Reference source not found. by using the LFM pulse and the FDSBR method. According to this comparison, the images obtained by using the three different excitation sources (modulating Gaussian pulse, LFM pulse, and FDSBR) are in perfect agreement. Section 4.3 gives a detailed analysis of the image similarity and the simulation efficiency.  Next, we simulate the ISAR images of a ship-ocean scene, as shown in Figure 18. The ship's size is 160 m × 20 m × 30 m. The size of the sea surface is 300 m × 300 m. The power spectrum of the sea surface is PM spectrum. The wind speed is 5 m/s. The relative dielectric constant of the sea is 55.9 − 37j. The incident angle is θ 0 = 45 • , and the center azimuth angle is ϕ 0 = 150 • . The polarization is HH. The carrier frequency is f 0 = 10 GHz. The resolution is 0.3 m in both the range and azimuth. Next, we simulate the ISAR images of a ship-ocean scene, as shown in Figure 18. The ship's size is 160 m × 20 m × 30 m. The size of the sea surface is 300 m × 300 m. The power spectrum of the sea surface is PM spectrum. The wind speed is 5 m/s. The relative dielectric constant of the sea is 55.9 37j − . The incident angle is 0 45 θ =°, and the center azimuth angle is 0 150 φ =°. The polarization is HH. The carrier frequency is 0 f = 10GHz.
The resolution is 0.3 m in both the range and azimuth. Error! Reference source not found. Figure 19 shows the ISAR images generated using the total echo, the differential echo, the target echo, and the coupling echo, respectively. These images are normalized and plotted with a logarithmic scale. By comparing the differential image (Figure 19b) with the target image (Figure 19c), it can be found that the strong scattering points in the image are caused by the dihedral corner structure on the target. The size of the coupling image is larger than that of the target image due to the rough modulation of the sea surface and multiple bounces of the ray. In Figure 19d, due to the substantial coupling scattering contributions of the dihedral corner structure formed by the hull and the sea surface, the ship's outline can be clearly seen from the coupling image. Due to the ship's shielding to the electromagnetic wave, only one side of the ship can be shown in this image. Figure 20shows the normalized ISAR images by using the LFM pulse and the FDSBR method. According to the comparison, the image in Figure  19a and the two images in Figure 20are in perfect agreement.   Figure 19 shows the ISAR images generated using the total echo, the differential echo, the target echo, and the coupling echo, respectively. These images are normalized and plotted with a logarithmic scale. By comparing the differential image (Figure 19b) with the target image (Figure 19c), it can be found that the strong scattering points in the image are caused by the dihedral corner structure on the target. The size of the coupling image is larger than that of the target image due to the rough modulation of the sea surface and multiple bounces of the ray. In Figure 19d, due to the substantial coupling scattering contributions of the dihedral corner structure formed by the hull and the sea surface, the ship's outline can be clearly seen from the coupling image. Due to the ship's shielding to the electromagnetic wave, only one side of the ship can be shown in this image. Figure 20 shows the normalized ISAR images by using the LFM pulse and the FDSBR method. According to the comparison, the image in Figure 19a and the two images in Figure 20 are in perfect agreement. rough modulation of the sea surface and multiple bounces of the ray. In Figure 19d, due to the substantial coupling scattering contributions of the dihedral corner structure formed by the hull and the sea surface, the ship's outline can be clearly seen from the coupling image. Due to the ship's shielding to the electromagnetic wave, only one side of the ship can be shown in this image. Figure 20shows the normalized ISAR images by using the LFM pulse and the FDSBR method. According to the comparison, the image in Figure  19a and the two images in Figure 20are in perfect agreement. target. The size of the coupling image is larger than that of the target image due to the rough modulation of the sea surface and multiple bounces of the ray. In Figure 19d, due to the substantial coupling scattering contributions of the dihedral corner structure formed by the hull and the sea surface, the ship's outline can be clearly seen from the coupling image. Due to the ship's shielding to the electromagnetic wave, only one side of the ship can be shown in this image. Figure 20shows the normalized ISAR images by using the LFM pulse and the FDSBR method. According to the comparison, the image in Figure  19a and the two images in Figure 20are in perfect agreement.  Table 2 presents the computation time of the images by using the modulating Gaussian pulse, the LFM pulse, and the FDSBR method. As it is shown, the computation time required for the imaging simulation using the modulating Gaussian pulse as the incident pulse is significantly less than that for the LFM pulse and the FDSBR method. The pulse width of the modulating Gaussian pulse is smaller than that of the LFM pulse under the same frequency bandwidth. Hence, the simulation time will be less. The imaging research based on the modulating Gaussian pulse, as presented in this paper, is a promising test for the imaging algorithm. Additionally, the TD method can also avoid the time-consuming frequency-sweep compared with the FD method, thus saving significant simulation time. As the size of the scene increases, the advantage of the modulating Gaussian pulse in computational efficiency becomes more significant.  Table 3 shows the similarity between the images under the modulating Gaussian pulse excitation and the images obtained by the other two methods. This similarity is quantified by the correlation coefficient shown in the following equation: The closer the correlation coefficient is to one, the more similar the two images become. As can be seen from Table 3, all similarities exceed 90%. The similarity between the modulating Gaussian pulse and the LFM pulse is close to one. This is because the matching filter used in this paper does not depend on the pulse. The images are nearly independent for both of the given pulse forms. The similarity between the results of modulating Gaussian pulse and the FDSBR method is slightly lower than that between the modulating Gaussian pulse and the LFM pulse. This is mainly caused by the difference between the imaging algorithm based on the FD echo and the imaging algorithm based on the TD echo, including the window function, FFT, and IFFT.

Conclusions
This paper proposes an accurate and efficient ISAR imaging method for the composite target-ocean scene based on the TDSBR method. The TDSBR method is used to obtain the TD electromagnetic scattering echoes of the composite target-ocean scene. Then the ISAR images are generated by performing demodulation, pulse compression, and the azimuth inverse fast Fourier transform on the scattering echo data. It is shown by a comparison with the precise full-wave method MLFMA that the performance of TDSBR can achieve good accuracy in simulating the scattering echo. Based on the echo expression, an image generation method is then proposed. This method includes echo demodulation, pulse compression, and the azimuth IFFT. To enhance efficiency, we replace LFM with a modulating Gaussian pulse and improve the matched filter to achieve pulse compression under a modulating Gaussian pulse excitation. Numerical examples verify the imaging accuracy and the efficiency of the proposed method.