# A New Shooting Bouncing Ray Method for Composite Scattering from a Target above the Electrically Large Scope Sea Surface

Authors: Peng Peng; Guo Lixin
Journal: Mathematical Problems in Engineering
 (2017)

---

## Abstract

A new shooting and bouncing ray (SBR) simulator based on the hybrid scheme of GO/PO/SDFM/EEC method is developed for the accurate prediction of composite scattering from a low altitude target above the electrically very-large-scale sea surface. It can adequately deal with the complex local electromagnetic interactions between the target and the large scope sea surface. The method is compared with the exact computational electromagnetic solver FEKO-MLFMM to validate its accuracy and efficiency. Then, it is applied to simulate the bistatic and monostatic scattering characteristics of an airplane above the electrically large sea surface at X-band, for different sea states. The results reveal the contributions from the target, sea surface, and interactions, which are of significance for radar target detection and remote sensing in real maritime environments.


---

## Body

## 1. Introduction

Composite electromagnetic scattering from a target above sea surface has attracted great attention in recent years, since it has a wide range of applications in target detection and tracking, remote sensing, radar imaging regime, and so forth [1–4]. In this issue, the EM scattering from the sea surface and the coupling interactions between the target and sea surface can greatly affect the overall radar target characteristics. In the present researches, many of the electromagnetic theories treat sea surface as a sort of rough surface and handle the EM simulations with various numerical methods, such as the most used method of moment (MoM) [5], finite difference time domain (FDTD) [6], finite element method (FEM) [7], and multilevel fast multipole algorithm (MLFMA) [8]. But they are generally limited by the computational requirements in terms of time and memory for handling the electrically small scale EM simulation issue, which is far from the practical needs. In comparison, high frequency approximation methods have the obvious advantage of electrically large EM simulations. Shooting and bouncing ray (SBR) method [9–12] is a popular and effective high frequency technique in consideration of the multipath coupling mechanisms, which can offer convincing physical insights and adequate accuracies. The traditional SBR method comprises a ray tracing process on the geometrical optics (GO) and physical optics (PO)/Kirchhoff approximation (KA) method to handle the coupling interactions among the meshed facets of the simulation model. In this process, a fine mesh size (smaller than 1/10 of the incident wavelength) for both the target and the sea surface is required. However, it can cause a heavy computation burden, which bottlenecks SBR by not being able to maintain both the accuracy and the efficiency in calculating the composite scattering from a target at the practical maritime scene with an electrically very large size and complex scattering mechanisms.The purpose of this paper is to fix this problem by an improved SBR simulator. In the simulator, a semideterministic facet model with two-scale formulas [13, 14] is employed to describe the scattering from the electrically very large sea surface. A fast bidirectional ray tracing method [15] combined with the KD-tree [16] acceleration technique is employed to describe the interactions among the target and sea facets, which can keep both the accuracy and the efficiency in the ray tracing process. An equivalent edge current (EEC) method is used as a modification to consider the diffraction field contributions from the edged structures. The accuracy and efficiency of the simulator are approved, and it is well applied in studying the composite scattering from an airplane target above the electrically large scope sea in the practical scene at different sea states.


## 2. Theoretical Model and Formulations

### 2.1. The Semideterministic Facet Model (SDFM)

The practical sea surface profile is a complex random system with composite surface structures. The physical structure of the large-scale sea surface can be generated by the Monte Carlo method [17], based on Elfouhaily’s sea spectrum and spreading function [18]. In the semideterministic facet model (SDFM), the sea surface is approximated by the large planar facets with short (capillary) wave superimpositions. The short waves can be expressed in the form of sinusoidal wave expansions, as shown in Figure 1.Figure 1
Sea surface illustration with SDFM.According to the Bragg theory, only the wave component traveling along the radar sight at the Bragg resonant frequency makes the major contribution to the radar scattering. The Bragg wave is modeled by sinusoidal waves traveling towards and away from the radar sight as in Figure1 with the following form: (1)ζcapr′=Bκccos⁡κc·r′-ωct,wherein the wave amplitude is B(κc)=S(κc)/ΔS. ΔS is the facet area. S(kc) is the capillary spectrum in the higher part of Elfouhaily’s sea spectrum. κc is the spatial capillary wave vector, whose direction is along the projection line ql of the scattering vector q=k(ks-ki) on the tilted plane of the rough sea facet.


### 2.2. The Hybrid GO/PO/SDFM Scheme and Formulations

The improved SBR method is based on a hybrid GO/PO/SDFM scheme. The calculation process is illustrated in Figure2.Figure 2
The geometry of the facet-based hybrid scheme.Both the target and the physically modeled sea surface are meshed by triangular facets for ray tracing. The incident plane wave is modeled as a grid of ray tubes. The ray tubes are tracing among facets according to the GO principles. Each illuminated facet is working at a re-radiation pattern like a small antenna. The scattering far field upon each facet is calculated according to the corresponding high frequency scattering mechanisms. Each ray is traced until it does not have any intersection with any facets. The total scattering field is contributed by the scattering from all the facets.The target surface is assumed to have a (perfect electric conductor, PEC) structure, where the scattered fields are solved by the Stratton-Chu equation in PO integral approximation as given by [19](2)EPO=jke-ikR4πR∬η·k^s×k^s×Jr′eikk^s·r′ds′.k^s is the unit vector of the scattering vector ks. R is the far-field distance. J(r′) is the current density upon surface. The integral in (2) can be solved by Gordon’s method [20]. The target facets’ size is requested to be smaller than the 1/10 incident wavelength.The scattering from sea facets is handled by SDFM as given by [14](3)ESDFM=k21-εeikRi4πRFpq∬ζr′e-iq·r′dr′.ζ(r′) upon sea facets has a sinusoidal form as in (1). The integral in (3) can be analytically solved by using Euler’s formula and the expansion of a plane wave in terms of Bessel functions as follows [21]:(4)I·=∬ζr′e-iq·r′dr′=ΔS2nze-iq·rcBκc+∑n=-∞∞-inJnqzBκc+I0κc++Bκc-∑n=-∞∞-inJnqzBκc-I0κc-.In the Bessel series, only the dominant terms n=0,±1 are considered. Jn(·) is the n-order Bessel function. At the facet local coordinate in Figure 1, I0(κc) is given by(5)I0κc=e-i1+nωctsinc⁡Δx21+nκcx-qx-qzzx·sinc⁡Δy21+nκcy-qy-qzzy+e-i1-nωctsinc⁡Δx21-nκcx+qx+qzzx·sinc⁡Δy21-nκcy+qy+qzzy.


### 2.3. The Improved Ray Tracing Techniques

The SDFM allows using large sea facets, which can substantially reduce the number of facets. However, when target facets have different mesh size from sea facets, some of the ray tubes may split. Discarding these rays may lead to calculation errors. To solve this problem and improve the ray tracing efficiency, bidirectional ray tracing combined with the KD-tree acceleration technique is employed.As illustrated in Figure3, all meshed facets are numbered and stored in the computer memory beforehand. The grey area is assumed to be the region shot by the ray from facet 0. Several steps are employed to find the IDs of the illuminated facets. Firstly, a forward ray is used to carry out an intersection test to find one of the facets in the illuminated region (i.e., facet 22). Then, a binary tree is built as in Figure 4, where facet 22 is the root and three neighboring facets are its children. Next, a backward ray is used to check the visibility of the three children facets (i.e., facets 12, 21, and 23); in this situation, they are also illuminated. Thereafter, the backward ray tracing process continues for these three children facets until facets outside the illuminated area appear (such as facet 10). Then, the backward ray tracing process is terminated. In each forward tracing process, intersection tests are executed with the KD tree of sea surface built to determine whether the facet is illuminated by the ray. Finally, the total scattering field can be obtained by vector summing of the scattering contributions from all the illuminated facets.Figure 3
The geometry of the bidirectional ray tracing process.Figure 4
A binary tree built according to the ray tracing.


### 2.4. Edge Diffraction Consideration with EEC

Up to now, the scattering from all the target and sea facets is considered. But in fact, the scattering far field also includes diffractions from the edged structures. They are considered by EEC method as follows:(6)Ed=ik04πe-jk0rr∫ls^×s^×t^Ie+s^×t^Imeik0s^·r′dt′.Ie and Im have the following form:(7)Ie=i2Einc·t^DeEECksin2⁡βi+i2ηHinc·t^DemEECksin⁡βi,Im=i2ηHinc·t^ksin⁡βisin⁡βsDmEEC.Einc and Hinc are the incident electric and magnetic field. t^ is the unit vector along edges. The integral is made along edges. DeEEC, DemEEC, and DmEEC are the EEC diffraction coefficients [22]. Their values are mainly determined by the angles between the neighboring facets.


## 2.1. The Semideterministic Facet Model (SDFM)

The practical sea surface profile is a complex random system with composite surface structures. The physical structure of the large-scale sea surface can be generated by the Monte Carlo method [17], based on Elfouhaily’s sea spectrum and spreading function [18]. In the semideterministic facet model (SDFM), the sea surface is approximated by the large planar facets with short (capillary) wave superimpositions. The short waves can be expressed in the form of sinusoidal wave expansions, as shown in Figure 1.Figure 1
Sea surface illustration with SDFM.According to the Bragg theory, only the wave component traveling along the radar sight at the Bragg resonant frequency makes the major contribution to the radar scattering. The Bragg wave is modeled by sinusoidal waves traveling towards and away from the radar sight as in Figure1 with the following form: (1)ζcapr′=Bκccos⁡κc·r′-ωct,wherein the wave amplitude is B(κc)=S(κc)/ΔS. ΔS is the facet area. S(kc) is the capillary spectrum in the higher part of Elfouhaily’s sea spectrum. κc is the spatial capillary wave vector, whose direction is along the projection line ql of the scattering vector q=k(ks-ki) on the tilted plane of the rough sea facet.


## 2.2. The Hybrid GO/PO/SDFM Scheme and Formulations

The improved SBR method is based on a hybrid GO/PO/SDFM scheme. The calculation process is illustrated in Figure2.Figure 2
The geometry of the facet-based hybrid scheme.Both the target and the physically modeled sea surface are meshed by triangular facets for ray tracing. The incident plane wave is modeled as a grid of ray tubes. The ray tubes are tracing among facets according to the GO principles. Each illuminated facet is working at a re-radiation pattern like a small antenna. The scattering far field upon each facet is calculated according to the corresponding high frequency scattering mechanisms. Each ray is traced until it does not have any intersection with any facets. The total scattering field is contributed by the scattering from all the facets.The target surface is assumed to have a (perfect electric conductor, PEC) structure, where the scattered fields are solved by the Stratton-Chu equation in PO integral approximation as given by [19](2)EPO=jke-ikR4πR∬η·k^s×k^s×Jr′eikk^s·r′ds′.k^s is the unit vector of the scattering vector ks. R is the far-field distance. J(r′) is the current density upon surface. The integral in (2) can be solved by Gordon’s method [20]. The target facets’ size is requested to be smaller than the 1/10 incident wavelength.The scattering from sea facets is handled by SDFM as given by [14](3)ESDFM=k21-εeikRi4πRFpq∬ζr′e-iq·r′dr′.ζ(r′) upon sea facets has a sinusoidal form as in (1). The integral in (3) can be analytically solved by using Euler’s formula and the expansion of a plane wave in terms of Bessel functions as follows [21]:(4)I·=∬ζr′e-iq·r′dr′=ΔS2nze-iq·rcBκc+∑n=-∞∞-inJnqzBκc+I0κc++Bκc-∑n=-∞∞-inJnqzBκc-I0κc-.In the Bessel series, only the dominant terms n=0,±1 are considered. Jn(·) is the n-order Bessel function. At the facet local coordinate in Figure 1, I0(κc) is given by(5)I0κc=e-i1+nωctsinc⁡Δx21+nκcx-qx-qzzx·sinc⁡Δy21+nκcy-qy-qzzy+e-i1-nωctsinc⁡Δx21-nκcx+qx+qzzx·sinc⁡Δy21-nκcy+qy+qzzy.


## 2.3. The Improved Ray Tracing Techniques

The SDFM allows using large sea facets, which can substantially reduce the number of facets. However, when target facets have different mesh size from sea facets, some of the ray tubes may split. Discarding these rays may lead to calculation errors. To solve this problem and improve the ray tracing efficiency, bidirectional ray tracing combined with the KD-tree acceleration technique is employed.As illustrated in Figure3, all meshed facets are numbered and stored in the computer memory beforehand. The grey area is assumed to be the region shot by the ray from facet 0. Several steps are employed to find the IDs of the illuminated facets. Firstly, a forward ray is used to carry out an intersection test to find one of the facets in the illuminated region (i.e., facet 22). Then, a binary tree is built as in Figure 4, where facet 22 is the root and three neighboring facets are its children. Next, a backward ray is used to check the visibility of the three children facets (i.e., facets 12, 21, and 23); in this situation, they are also illuminated. Thereafter, the backward ray tracing process continues for these three children facets until facets outside the illuminated area appear (such as facet 10). Then, the backward ray tracing process is terminated. In each forward tracing process, intersection tests are executed with the KD tree of sea surface built to determine whether the facet is illuminated by the ray. Finally, the total scattering field can be obtained by vector summing of the scattering contributions from all the illuminated facets.Figure 3
The geometry of the bidirectional ray tracing process.Figure 4
A binary tree built according to the ray tracing.


## 2.4. Edge Diffraction Consideration with EEC

Up to now, the scattering from all the target and sea facets is considered. But in fact, the scattering far field also includes diffractions from the edged structures. They are considered by EEC method as follows:(6)Ed=ik04πe-jk0rr∫ls^×s^×t^Ie+s^×t^Imeik0s^·r′dt′.Ie and Im have the following form:(7)Ie=i2Einc·t^DeEECksin2⁡βi+i2ηHinc·t^DemEECksin⁡βi,Im=i2ηHinc·t^ksin⁡βisin⁡βsDmEEC.Einc and Hinc are the incident electric and magnetic field. t^ is the unit vector along edges. The integral is made along edges. DeEEC, DemEEC, and DmEEC are the EEC diffraction coefficients [22]. Their values are mainly determined by the angles between the neighboring facets.


## 3. Numerical Simulations and Discussions

### 3.1. Validations

The total scattering fieldEsca can be obtained by the vector sum of the scattered field from all the illuminated facets and their edges. The composite scattering characteristics can be evaluated by the scattering coefficient, as defined by (8)σ=limr→∞⁡4πr2Esca2Einc2.As an initial investigation, the scattering from an individual6λ×6λ sea facet with a sinusoidal wave expression of ξ=0.015cos⁡(2πx/Λ) as shown in Figure 1 is calculated. Λ=λ/(2sin⁡θi) is the spatial wavelength satisfying the Bragg resonance condition. The dielectric constant of the sea water is frequency-dependent, which is calculated according to the Klein and Swift model [23] at 20°C and 32.5% of salinity. The incident frequency is 1 GHz. The incident angle is θi = 30°, φi = 0°. The scattering upon the sea facet is analytically solved by (4) and (5). The simulation is compared with the exact numerical solver, MLFMM, in the commercial software FEKO. The results are shown in Figure 5.Figure 5
Bistatic scattering from an individual sea facet.
(a)
HH polarization
(b)
VV polarizationIt can be seen that the results from the SDFM agree well with the FEKO-MLFMM results for both types of polarization. In Table1, SDFM shows much lower computational complexity and higher efficiency than the MLFMM solver in FEKO, so it has an obvious priority to be applied in the evaluation of scattering from the electrically very large sea surface.Table 1
Comparison of the calculation time and memory.
Solution
SDFM
FEKO-MLFMM
Time
10.423 s
2460.248 s
Memory
23.8 MB
5814.8 MBNext, we investigate the bistatic scattering characteristics of an electrically small composite model of a 0.1 m × 0.1 m × 0.1 m cubic target above a 1 m × 1 m flat plane. The cube is set at 0.5 m high above the plane facet.Figure6 shows the HH (a) and VV (b) polarization bistatic scattering coefficients at a frequency of 10 GHz. The incident angle is 45∘. The simulation results are also compared with the results from FEKO-MLFMM, which show a good agreement. But the FEKO-MLFMM solver will consume much more computation time and memory as shown in Table 2, so it is very hard to accomplish electrically very large simulations in the real scene.Table 2
Comparison of the calculation time and memory.
Solution
Hybrid scheme
FEKO-MLFMM
Time
128.343 s
17411.34 s
Memory
64.6 MB
12697.6 MBFigure 6
Bistatic scattering from a cube above a planar facet.
(a)
HH polarization
(b)
VV polarization


### 3.2. The Composite Scattering from an Airplane above the Sea Surface

In the following simulations, the SBR simulator is used for the prediction of composite scattering from an airplane above an electrically large scope sea. The geometry of the airplane model is shown in Figure7. It has a full length of 16 m, wing span of 12 m, and a height of 2.5 m. The sea surface is physically modeled by the Monte Carlo method. The calculation region is 128 m × 128 m. The target is located 5 m high above the sea. The incident frequency is set to 10.0 GHz (X-band) HH polarization. The composite target and sea model are extremely electrically large which is almost impossible for the exact numerical solver.Figure 7
The geometry of an airplane model.Figure8 shows the bistatic scattering characteristics of the composite model, target, and sea surface. The sea state is evaluated by the wind speed at 10 m high of the sea surface U10 in Elfouhaily’s spectrum. In this simulation, U10=5 m/s. The incident angle is θi = 45°, φi = 90°. The scattering angle is θs = −90°~90°. The sea scattering and composite scattering results are obtained by averaging 50 surface realizations.Figure 8
The bistatic scattering characteristics of the airplane, sea surface, and composite model.It can be seen that the scattering from the composite model is influenced by the scattering from the target and sea surface as well as the target-surface interactions. In the specular directions, the scattering from the sea surface dominates the total scattering. In the backward directions, the composite scattering is mainly contributed by the target scattering and the coupling scattering components. Further simulations are made on the bistatic and monostatic scattering at different sea states as shown in Figure9.Figure 9
The composite scattering of the airplane above the sea surface at different sea states. (a) Bistatic scattering. (b) Monostatic scattering.
(a)
(b)From the comparisons in Figure9, one can see that, at the low level sea state (U10 = 3 m/s), the scattering in the backward directions is apparently strengthened by the coupling interactions between the target and the sea surface. At the high sea state (U10 = 10 m/s), the fluctuation of sea surface is stronger with the wind speed increasing. The scattering intensities in the specular directions are attenuated, but the diffuse scattering intensities are enhanced. The phenomenon is also observed in another related study about the composite scattering from the target-sea model [24].


## 3.1. Validations

The total scattering fieldEsca can be obtained by the vector sum of the scattered field from all the illuminated facets and their edges. The composite scattering characteristics can be evaluated by the scattering coefficient, as defined by (8)σ=limr→∞⁡4πr2Esca2Einc2.As an initial investigation, the scattering from an individual6λ×6λ sea facet with a sinusoidal wave expression of ξ=0.015cos⁡(2πx/Λ) as shown in Figure 1 is calculated. Λ=λ/(2sin⁡θi) is the spatial wavelength satisfying the Bragg resonance condition. The dielectric constant of the sea water is frequency-dependent, which is calculated according to the Klein and Swift model [23] at 20°C and 32.5% of salinity. The incident frequency is 1 GHz. The incident angle is θi = 30°, φi = 0°. The scattering upon the sea facet is analytically solved by (4) and (5). The simulation is compared with the exact numerical solver, MLFMM, in the commercial software FEKO. The results are shown in Figure 5.Figure 5
Bistatic scattering from an individual sea facet.
(a)
HH polarization
(b)
VV polarizationIt can be seen that the results from the SDFM agree well with the FEKO-MLFMM results for both types of polarization. In Table1, SDFM shows much lower computational complexity and higher efficiency than the MLFMM solver in FEKO, so it has an obvious priority to be applied in the evaluation of scattering from the electrically very large sea surface.Table 1
Comparison of the calculation time and memory.
Solution
SDFM
FEKO-MLFMM
Time
10.423 s
2460.248 s
Memory
23.8 MB
5814.8 MBNext, we investigate the bistatic scattering characteristics of an electrically small composite model of a 0.1 m × 0.1 m × 0.1 m cubic target above a 1 m × 1 m flat plane. The cube is set at 0.5 m high above the plane facet.Figure6 shows the HH (a) and VV (b) polarization bistatic scattering coefficients at a frequency of 10 GHz. The incident angle is 45∘. The simulation results are also compared with the results from FEKO-MLFMM, which show a good agreement. But the FEKO-MLFMM solver will consume much more computation time and memory as shown in Table 2, so it is very hard to accomplish electrically very large simulations in the real scene.Table 2
Comparison of the calculation time and memory.
Solution
Hybrid scheme
FEKO-MLFMM
Time
128.343 s
17411.34 s
Memory
64.6 MB
12697.6 MBFigure 6
Bistatic scattering from a cube above a planar facet.
(a)
HH polarization
(b)
VV polarization


## 3.2. The Composite Scattering from an Airplane above the Sea Surface

In the following simulations, the SBR simulator is used for the prediction of composite scattering from an airplane above an electrically large scope sea. The geometry of the airplane model is shown in Figure7. It has a full length of 16 m, wing span of 12 m, and a height of 2.5 m. The sea surface is physically modeled by the Monte Carlo method. The calculation region is 128 m × 128 m. The target is located 5 m high above the sea. The incident frequency is set to 10.0 GHz (X-band) HH polarization. The composite target and sea model are extremely electrically large which is almost impossible for the exact numerical solver.Figure 7
The geometry of an airplane model.Figure8 shows the bistatic scattering characteristics of the composite model, target, and sea surface. The sea state is evaluated by the wind speed at 10 m high of the sea surface U10 in Elfouhaily’s spectrum. In this simulation, U10=5 m/s. The incident angle is θi = 45°, φi = 90°. The scattering angle is θs = −90°~90°. The sea scattering and composite scattering results are obtained by averaging 50 surface realizations.Figure 8
The bistatic scattering characteristics of the airplane, sea surface, and composite model.It can be seen that the scattering from the composite model is influenced by the scattering from the target and sea surface as well as the target-surface interactions. In the specular directions, the scattering from the sea surface dominates the total scattering. In the backward directions, the composite scattering is mainly contributed by the target scattering and the coupling scattering components. Further simulations are made on the bistatic and monostatic scattering at different sea states as shown in Figure9.Figure 9
The composite scattering of the airplane above the sea surface at different sea states. (a) Bistatic scattering. (b) Monostatic scattering.
(a)
(b)From the comparisons in Figure9, one can see that, at the low level sea state (U10 = 3 m/s), the scattering in the backward directions is apparently strengthened by the coupling interactions between the target and the sea surface. At the high sea state (U10 = 10 m/s), the fluctuation of sea surface is stronger with the wind speed increasing. The scattering intensities in the specular directions are attenuated, but the diffuse scattering intensities are enhanced. The phenomenon is also observed in another related study about the composite scattering from the target-sea model [24].


## 4. Conclusions

In this article, a new SBR simulator is developed to calculate the composite scattering characteristics of a low altitude target above an electrically large scope sea surface. The accuracy and efficiency of the method are validated by the exact numerical solution. The simulator is employed for the scattering from an airplane model above the large scope sea at both low and high level sea states. The results indicate contributions from the target, sea surface, and interactions at different observation angles. In fact, some more complex sea structures such as the crest wave, whitecaps, and foams existing on the real sea surface will also influence the composite scattering. These need to be further considered in a future study.



---

*Source: 6201745-2017-11-29.xml*