# Composite Scattering from a PEC Object above Two-Layered Rough Soil Surfaces Based on an Efficient Bi-Iterative Model

Authors: Ke Li; Lixin Guo; Juan Li
Journal: International Journal of Antennas and Propagation
 (2022)

---

## Abstract

In the study, an efficient bi-iterative physical optics (PO) electromagnetic model is used to analyze the composite scattering from a PEC object above two-layered rough soil surfaces for the first time. The PO method is applied to calculate the direct scattered fields from the object and the upper layer rough surface and to deal with the mutual coupling effect between the object and the two-layered rough surfaces. During the mutual coupling process, a bi-iterative strategy is considered to solve the coupling interaction between the object and the upper rough surface and the coupling interaction between the upper and lower rough surfaces. We assume a 2-D scattering problem. The result is validated by comparing it with that obtained from the conventional method of moments (MoM). Comparison results from different object types show that the scattered curves calculated using the proposed method resemble the ones obtained from MoM. The comparison of computer resources also shows the efficiency of our method. Meanwhile, the normalized radar cross-section of the composite model is analyzed under different parameters.


---

## Body

## 1. Introduction

Actual natural backgrounds, such as the Earth’s surface (forest, snow, desert, and ocean scenarios), and all kinds of artificial material surfaces, can be regarded as random rough surface models. Scattering of electromagnetic (EM) waves from rough surfaces has been studied for over half a century since the 1950s [1, 2]. The natural surfaces are modeled as single layered rough surfaces in the early days. However, many scenarios mentioned above consist of two-layered or multilayer rough surfaces. Radar remote sensing from the Earth holds great realistic significance, such as retrieval of soil properties and finding camouflaged targets above the land. Thus, most related applications concerning surveillance to the soil, target tracking, and navigation communication in complex backgrounds are explored. For example, L-band soil moisture active and passive missions are identified as a priority for NASA’s Science Mission Directorate [3], L-band Soil Moisture and Ocean Salinity Missions are conducted by the European Space Agency [4], and an experiment of autonomous vessel tracking in simulated rough sea conditions is implemented by Villanova University [5].Apart from the experiments, theoretical simulation methods are also used to examine rough surfaces. Because of the flexibility and high precision, numerical methods and their improved methods are always used to calculate the composite scattering from the object and rough surfaces. For example, the method of moments (MoM) [6], the finite difference time domain (FDTD) method [7], and the finite element method (FEM) [8]. However, the numerical methods consume a lot of computer resources, which limit their scope of application. Notably, proper combinations, such as the combination of approximate method with numerical method can result in high computational efficiency while maintaining the desirable precision. Examples are the hybrid analytic-numerical algorithm combining the Kirchhoff approximation (KA) with the method of moments (MoM) [9], the method combining KA with the hybrid finite element-boundary integral [10], and the hybrid method combining physical optics (PO) with PO [11, 12]. However, these works are either ineffective or not applied further on most actual scenarios, such as forest covered with vegetation, land covered with snow, and even same visual soil with different moisture contents in layer depth. To date, some methods are used to calculate the properties of layered rough surfaces, such as iterative physical optics (IPO) model [13], physical optics approximation (POA) algorithm [14], or used to analyze the polarimetric scattering from two-layered random rough surfaces with and without buried objects, such as the KA model [15] and the steepest descent fast multipole method (SDFMM) [16]. However, research on composite scattering from the object above the layered rough surface is limited, except the work complemented by Wang and Li [17, 18]. This scarcity is due to the intricate coupling interaction between the object and the layered rough surface that complicates the problem further.The aforementioned problems motivate the investigation on scattering from the object above multilayered soil rough surfaces. In this study, the soil is modeled as a two-layered rough surface with a Gaussian spectrum [17] by using the Monte Carlo method, and each layer has different soil moisture contents. The objects are cylinder and square cylinder. A bi-iterative PO is adopted to solve the coupling interaction between the object and the upper rough surface and the coupling interaction between the upper and lower rough surfaces. The normalized radar cross-section (NRCS) of the composite model is validated by comparing the results with those obtained from the conventional MoM. In this study, all fields and currents are assumed to have a time-harmonic dependence in the form e−iωt, which is suppressed in this study.


## 2. Theoretical Formulations of the Bi-Iterative Model

The bi-iterative model in this study is an extension of our work in [12]: a PEC object partially buried in a rough sea surface. Specifically, we developed it further on and expanded the method to a more complicated model than before with an object above layered rough surfaces. The complication is the single layer rough surface is extended to a two-layered rough surface. Thus, the coupling scattering is complicated further and the scattering mechanism should be explored between the two-layered rough surfaces. Furthermore, the numerical calculation of the rough surface often accounts for the most amount of the computation in the problem of composite scattering of the object and rough surfaces, so the attempt to develop a method to measure the scattering effectively, accelerate the calculation, and analyze the coupling scattering mechanisms of the object and rough surfaces (especially for the two-layer rough surfaces) is a valuable task. The detailed solution is expressed as follows:The geometry of the problem of interest is shown in Figure1. k^i is the incident wave unit vector, and k^s is the scattered wave unit vector. θi is the incident angle, θs is the scattered angle. Meanwhile, we assume that the simulated rough surfaces have large radii of the curvature to enable the use of the PO method.Figure 1
The geometry of a PEC object above a two-layered rough surface.The soil surface is modeled as a two-layered rough surface, and the profile of the two-layered rough surface divides the entire space into three spaces. They are labeled asΩ0, Ω1, and Ω2. So denotes the object in space Ω0,S1 denotes the upper layer surface with a profile ς1x, and S2 denotes the lower layer with a profile ς2x; So, S1, and S2 are used as subscripts in the subsequent current variables. S1 and S2 are generated by the Monte Carlo method. We assume that S1 and S2 are uncorrelated. h is the altitude of the object, and d is the average thickness of the S1 layer. h and d must be sufficiently large to avoid the profile spatial overlap. Ω0, Ω1, and Ω2 are characterized by relative permittivities εr0, εr1, and εr2. For the model in our study, εr0=1.


### 2.1. The Direct Scattering from SpaceΩ0

We use horizontal (HH) polarization as an example. The incident waveEincr=y^ϕincr excites equivalent induced currents on the surfaces located in space Ω0 with the assumption that So and S1 are independent, where ϕincr is Thorsos tapered plane wave [19]. The currents on So and S1 are represented as Jso1 and Js11,Ms11. The superscript number represents the order number of currents, and they are the same for subsequent variables. The first-order equivalent currents are obtained by the PO method and represented as(1)Jso1r=2n^o×Hincrr∈So,(2a)Js11r=1−R01Hn^1×Hincrr∈S1,(2b)Ms11r=1+R01HEincr×n^1r∈S1,where n^o is the surface normal vector to So, n^1 is the upper surface normal vector to S1, and the reflection coefficient R01H is expressed as follows:(3)R01H=εr0cosθl−εr1−εr0sin2θlεr0cosθl+εr1−εr0sin2θlcosθl=−k^i·n^1.We provide a detailed explanation to (3) because it will be used frequently in the following analysis with a regular pattern. R01H is the reflection coefficient on the boundaries of Ω0 and Ω1, and subscripts “0” and “1” are related to the parameters in Ω0 and Ω1, respectively. “01” indicates that the wave is incident in Ω0 and the refractive in Ω1. θl is the local incident angle on the boundary and is determined by the incident wave and normal vectors of boundary at the incident side.


### 2.2. The Mutual Coupling Strategy among the Boundaries

We assume then−1 th-order equivalent induced currents Json−1 and Js1n−1,Ms1n−1 scatter fields in space Ω0. We also introduce the equivalent currents on the lower surface of S1 and denote them as Js1′n−1,Ms1′n−1. Js1′n−1,Ms1′n−1 and Js2n−1,Ms2n−1 scatter fields in space Ω1. From the boundary conditions, we can reason out that(4)Js1′n−1=−Js1n−1;Ms1′n−1=−Ms1n−1n≥2.An overview of the coupling iterative procedure is described as follows: (a) the scattered fields excited byJs1n−1andMs1n−1 travel to So, and the n th-order currents Json on So are obtained; (b) the scattered fields excited by Json−1 travel to S1, and the n th-order currents Js1−1nandMs1−1n on the upper surface of S1 are obtained; (c) the scattered fields excited by Js1′n−1andMs1′n−1 travel to S2, and the n th-order currents Js2nandMs2n on S2 are obtained; (d) the scattered fields excited by Js2n−1andMs2n−1 travel to S1, and the n th-order currents Js1−2nandMs1−2n on the lower surface of S1 are obtained. The n th equivalent currents are obtained by the PO method. The detailed coupling procedure is shown in Figure 2.Figure 2
The mutual coupling strategy among the boundaries.In Figure2, the primary and secondary currents are marked in black and blue, respectively. The solid line with an arrow represents the coupling relation described as (a)–(d). The dashed line with an arrow represents the direction of the current updating, which indicates that the secondary currents will replace the last currents and can be used to calculate the next order currents. The newly induced currents on the upper surface of S1 are expressed as follows:(5)Js1n=Js1−1n−Js1−2n;Ms1n=Ms1−1n−Ms1−2nn≥2.Then, equivalent currents on the lower surface ofS1 can be obtained from (4). The currents on the upper and lower surfaces will iterate in Ω0 and Ω1, respectively. This model is the so-called bi-iterative model. The initial currents are the directly induced currents excited by the tapered plane wave in (1) and (2). Notably, Js21=0andMs21=0.The expression of the induced equivalent currents on the surfaces can be written as follows:(6)Jsonr=2n^o×L1Js1n−1,Ms1n−1r∈So,(7a)Js1−1nr=1−R01Hn^1×L1Json−1r∈S1,(7b)Ms1−1nr=1+R01HK1Json−1×n^1r∈S1,(8a)Js2nr=1−R12Hn^2×L2Js1′n−1,Ms1′n−1r∈S2,(8b)Ms2nr=1+R12HK2Js1′n−1,Ms1′n−1×n^2r∈S2,(9a)Js1−2nr=1−R10Hn^1′×L2Js2n−1,Ms2n−1r∈S1,(9b)Ms1−2nr=1+R10HK2Js2n−1,Ms2n−1×n^1′r∈S1.Equations (6)–(9) correspond to (a)–(d), respectively. The definition of reflection coefficients R01H, R12H, and R10H in (7)–(9) is similar to that in (3), all the corresponding parameters can be taken into the expression regularly. However, the local incident angle is determined by the local incident wave and normal vectors on the boundaries. The local incident wave vectors are the directions from the source point to the field point, and they replace k^i in (3). n^2 is the surface normal vector to S2, n^1′ is the lower surface normal vector to S1. L1, K1, L2, and K2 are linear operators, and they can be found in [12].


### 2.3. Far Zone Scattering of the Composite Model

The total equivalent currents on the surfaces inΩ0 can be expressed as follows:(10)Jsor=Jso1r+Jso2r+⋯+Jsonrr∈So,(11a)Js1r=Js11r+Js12r+⋯+Js1nrr∈S1,(11b)Ms1r=Ms11r+Ms12r+⋯+Ms1nrr∈S1.In our simulations, an error functionυ of the object is chosen as a rule to control the convergence in this model. The iteration begins with the first-order currents Jso1r, and terminates with Jsor reaching stable. The error function [20] is defined as follows:(12)υ=Jsonr∑m=1nJsomr,where the notation of · indicates the Euclidean 2-norm of a vector. When υ<0.001, the interaction process is assumed to be stable and will stop. Then, the far zone fields Eso and Es1 in Ω0 can be calculated through the stable currents on the object and the upper layer rough surface by using Huygens’s principle. The NRCS [19] can be calculated by(13)σθs=limr⟶∞2πrEs1θs+Eoθs2Pinc.Pinc is the power impinging upon the composite model, and it can be found in [21]. The process of vertical (VV) polarization is analogous to that of HH polarization, and the incident wave Hincr=y^ϕincr. It will be omitted in this paper.


## 2.1. The Direct Scattering from SpaceΩ0

We use horizontal (HH) polarization as an example. The incident waveEincr=y^ϕincr excites equivalent induced currents on the surfaces located in space Ω0 with the assumption that So and S1 are independent, where ϕincr is Thorsos tapered plane wave [19]. The currents on So and S1 are represented as Jso1 and Js11,Ms11. The superscript number represents the order number of currents, and they are the same for subsequent variables. The first-order equivalent currents are obtained by the PO method and represented as(1)Jso1r=2n^o×Hincrr∈So,(2a)Js11r=1−R01Hn^1×Hincrr∈S1,(2b)Ms11r=1+R01HEincr×n^1r∈S1,where n^o is the surface normal vector to So, n^1 is the upper surface normal vector to S1, and the reflection coefficient R01H is expressed as follows:(3)R01H=εr0cosθl−εr1−εr0sin2θlεr0cosθl+εr1−εr0sin2θlcosθl=−k^i·n^1.We provide a detailed explanation to (3) because it will be used frequently in the following analysis with a regular pattern. R01H is the reflection coefficient on the boundaries of Ω0 and Ω1, and subscripts “0” and “1” are related to the parameters in Ω0 and Ω1, respectively. “01” indicates that the wave is incident in Ω0 and the refractive in Ω1. θl is the local incident angle on the boundary and is determined by the incident wave and normal vectors of boundary at the incident side.


## 2.2. The Mutual Coupling Strategy among the Boundaries

We assume then−1 th-order equivalent induced currents Json−1 and Js1n−1,Ms1n−1 scatter fields in space Ω0. We also introduce the equivalent currents on the lower surface of S1 and denote them as Js1′n−1,Ms1′n−1. Js1′n−1,Ms1′n−1 and Js2n−1,Ms2n−1 scatter fields in space Ω1. From the boundary conditions, we can reason out that(4)Js1′n−1=−Js1n−1;Ms1′n−1=−Ms1n−1n≥2.An overview of the coupling iterative procedure is described as follows: (a) the scattered fields excited byJs1n−1andMs1n−1 travel to So, and the n th-order currents Json on So are obtained; (b) the scattered fields excited by Json−1 travel to S1, and the n th-order currents Js1−1nandMs1−1n on the upper surface of S1 are obtained; (c) the scattered fields excited by Js1′n−1andMs1′n−1 travel to S2, and the n th-order currents Js2nandMs2n on S2 are obtained; (d) the scattered fields excited by Js2n−1andMs2n−1 travel to S1, and the n th-order currents Js1−2nandMs1−2n on the lower surface of S1 are obtained. The n th equivalent currents are obtained by the PO method. The detailed coupling procedure is shown in Figure 2.Figure 2
The mutual coupling strategy among the boundaries.In Figure2, the primary and secondary currents are marked in black and blue, respectively. The solid line with an arrow represents the coupling relation described as (a)–(d). The dashed line with an arrow represents the direction of the current updating, which indicates that the secondary currents will replace the last currents and can be used to calculate the next order currents. The newly induced currents on the upper surface of S1 are expressed as follows:(5)Js1n=Js1−1n−Js1−2n;Ms1n=Ms1−1n−Ms1−2nn≥2.Then, equivalent currents on the lower surface ofS1 can be obtained from (4). The currents on the upper and lower surfaces will iterate in Ω0 and Ω1, respectively. This model is the so-called bi-iterative model. The initial currents are the directly induced currents excited by the tapered plane wave in (1) and (2). Notably, Js21=0andMs21=0.The expression of the induced equivalent currents on the surfaces can be written as follows:(6)Jsonr=2n^o×L1Js1n−1,Ms1n−1r∈So,(7a)Js1−1nr=1−R01Hn^1×L1Json−1r∈S1,(7b)Ms1−1nr=1+R01HK1Json−1×n^1r∈S1,(8a)Js2nr=1−R12Hn^2×L2Js1′n−1,Ms1′n−1r∈S2,(8b)Ms2nr=1+R12HK2Js1′n−1,Ms1′n−1×n^2r∈S2,(9a)Js1−2nr=1−R10Hn^1′×L2Js2n−1,Ms2n−1r∈S1,(9b)Ms1−2nr=1+R10HK2Js2n−1,Ms2n−1×n^1′r∈S1.Equations (6)–(9) correspond to (a)–(d), respectively. The definition of reflection coefficients R01H, R12H, and R10H in (7)–(9) is similar to that in (3), all the corresponding parameters can be taken into the expression regularly. However, the local incident angle is determined by the local incident wave and normal vectors on the boundaries. The local incident wave vectors are the directions from the source point to the field point, and they replace k^i in (3). n^2 is the surface normal vector to S2, n^1′ is the lower surface normal vector to S1. L1, K1, L2, and K2 are linear operators, and they can be found in [12].


## 2.3. Far Zone Scattering of the Composite Model

The total equivalent currents on the surfaces inΩ0 can be expressed as follows:(10)Jsor=Jso1r+Jso2r+⋯+Jsonrr∈So,(11a)Js1r=Js11r+Js12r+⋯+Js1nrr∈S1,(11b)Ms1r=Ms11r+Ms12r+⋯+Ms1nrr∈S1.In our simulations, an error functionυ of the object is chosen as a rule to control the convergence in this model. The iteration begins with the first-order currents Jso1r, and terminates with Jsor reaching stable. The error function [20] is defined as follows:(12)υ=Jsonr∑m=1nJsomr,where the notation of · indicates the Euclidean 2-norm of a vector. When υ<0.001, the interaction process is assumed to be stable and will stop. Then, the far zone fields Eso and Es1 in Ω0 can be calculated through the stable currents on the object and the upper layer rough surface by using Huygens’s principle. The NRCS [19] can be calculated by(13)σθs=limr⟶∞2πrEs1θs+Eoθs2Pinc.Pinc is the power impinging upon the composite model, and it can be found in [21]. The process of vertical (VV) polarization is analogous to that of HH polarization, and the incident wave Hincr=y^ϕincr. It will be omitted in this paper.


## 3. Results and Discussion

In this section, the NRCS of the composite model from a 2D PEC object and 1D underlying two-layered soil surfaces is calculated and analyzed for HH and VV polarizations. The numerical results are averaged by 30 Monte Carlo realizations unless noted otherwise.Scattering from a PEC cylinder and two-layered soil surfaces calculated by the bi-iterative PO method for one composite realization is shown in Figure3. The numerical results by the MoM [22] are presented to validate our code. The incident frequency is f=0.435GHz, the incident angle is θi=60∘, the length of soil surface is L=102.4λ, and the tapering parameter is g=L/4. The root mean square (RMS) height and correlation length for the upper and lower surfaces are δ1=0.1λ, l1=1.0λ, δ2=0.2λ, and l2=1.5λ. The cylinder radius is r=2λ. Other parameters are d=4λ and h=4λ. The soil moisture values are mv1=0.08cm3/cm3 and mv2=0.12cm3/cm3 for the upper and lower layers, and the corresponding relative permittivities are εr1=8.2,1.2 and εr2=11.1,1.4. The soil mixture model is the one by Peplinski et al. [15, 23]. As shown in Figure 3, the results by the two different methods match very well for most observation angles under different polarizations.Figure 3
Comparison of NRCS from a cylinder above two-layered soil surfaces using the bi-iterative PO and MOM with (a) HH and (b) VV polarizations.
(a)(b)The computational time and memory requirement of the two methods for one composite realization of Figure3 are shown in Table 1. The memory by our method is only 0.028% of that by MOM. The computational times are approximately 9.81%–19.6% of those by MoM for HH and VV polarizations. This finding illustrates that our method remarkably saves computer resources. All simulations are obtained on a computer with a 2.33 GHz processor (Intel (R) Core(TM)2 Quad CPU Q8200), 3.5 GB memory, and Visual Fortran 6.5 compiler.Table 1
CPU time and required memory in Figure3.
PolarizationMethodMemory (MB)Time (s)HHMOM2693.93214Bi-iterative PO0.7521VVMOM2693.93107Bi-iterative PO0.7521To further show the validity of our method, the D-RCS [24] from a cylinder above two-layered soil surfaces that is the same model as in Figure 3 is calculated by our method and is compared with MoM for one realization in Figure 4. As shown in Figure 4, the results by the two different methods match very well for most observation angles under different polarizations.Figure 4
The comparison of D-RCS from a cylinder above two-layered soil surfaces using the bi-iterative PO and MOM with (a) HH and (b) VV polarizations.
(a)(b)The dependence of the roughness of Gaussian rough surface on NRCS analyzed with HH polarizations is shown in Figure5. When the soil moistures are set as mv1=0.12cm3/cm3, and mv2=0.25cm3/cm3 for the upper and lower layer surfaces, and the corresponding relative permittivities are εr1=11.1,1.4 and εr2=21.6,2.1 at f=0.435GHz. θi=20∘, the other parameters are the same as those in Figure 3. Figure 5(a) shows the comparison for a fixed correlation length of l1=l2=1.0λ and a varying RMS height δ1=δ2 for both surfaces. An obvious peak appears in the specular direction for δ1=δ2=0.01λ, and the specular peak decreases with the increase in RMS height. The reason is that the rough surfaces become rough and the incoherent scattering becomes strong. Figure 5(b) shows the comparison for a fixed RMS height of δ1=δ2=0.2λ and a varying correlation length l1=l2 for both surfaces. The rough surface becomes smooth with the increase in the correlation length. As a result, a large amount of energy will be concentrated in the specular directions.Figure 5
Angular distribution of NRCS for different roughness, (a) differentδ,and (b) different l.
(a)(b)The backscattering for different soil moisture contents of the upper layer from the rough surface and composite models with VV polarization is shown in Figure6. f=0.435GHz and mv1 are set as 0.08, 0.12, and 0.2cm3/cm3 for the upper layer, and the corresponding relative permittivities are εr1=8.2,1.2, 11.1,1.4, and 17.4,2.1. The lower layer soil moisture mv2 is fixed to 0.25cm3/cm3 with the relative permittivity εr2=21.6,2.1. The backscattered field decreases approximately linear with the increase in the incident angle for the two-layered soil surface. However, the backscattering field oscillates with the increase in the incident angle for the composite model. This result is attributed to the multipath effect between the object and the underlying rough surface. In addition, the scattering field is increased in nearly all the observation angles with the increase in soil moisture because of the decrease in the loss tangent of the upper layer soil.Figure 6
Backscattering for different soil moisture contents of the upper layer for (a) two-layered soil surfaces and (b) the composite model.
(a)(b)


## 4. Conclusions

The composite scattering from a PEC object above a two-layered rough soil surface is investigated in this study. A bi-iterative PO method is used to deal with the problem. The NRCS of the composite model is validated by comparing with that from the conventional MoM. The comparison results show the validity of the adopted method. The comparison of computer resources also shows the efficiency of our method. The simulated results show that the coupling effect is remarkable in the composite model, and the soil moisture of the upper surface affects the NRCS considerably. This study can provide a fast and efficient theoretical reference for target detection above the soil and radar remote sensing retrieval of the soil moisture.In our future study, the complex 3D model of an object and layered rough surfaces will be considered for practical applications.



---

*Source: 6157022-2022-12-31.xml*