Abstract

A hybrid model based on the geometrical optics and physical optics (GO-PO) method and facet-based asymptotical model (FBAM) for the composite EM scattering from sea surface with ship has been previously developed. But, the GO-PO method faces a difficult problem that the lower efficiency in testing the visibility of patches to the incident and reflected waves. Hence, to improve the computation efficiency of the GO-PO method, a kd-tree accelerated hybrid model is presented in this paper. The kd-tree method firstly carries out the intersection test between the ray and the bounding box where the patch is located before the operation between the ray and the patch. If the distance parameters of the bounding box do not meet the requirement, the ray is not necessary to enter the bounding box, which can reduce the needless intersection tests. Numerical results show that the proposed kd-tree accelerated hybrid model could improve the efficiency about three times and could be applied to the study of composite EM scattering characteristics of electrically large sea surface with a ship target.

1. Introduction

The study of composite electromagnetic (EM) scattering characteristics from sea surface with a complex target is always a difficult problem both in the EM experiment measurement and computer simulation. However, the coupling effect between rough sea surface and target, including multiple scattering, has a very important influence on the radar cross section (RCS) of the target and sea surface. The modeling of composite EM scattering from rough sea surface and target involves the scattering theories of rough surface and complex target, which makes the research of EM coupling mechanism and the modeling of coupling field more complicated. Therefore, the composite scattering modeling of target and rough sea surface possesses important significance both in theories and practical applications.

For the problems of the EM scattering modeling and the coupling mechanism of target and rough surface, many scholars have done a lot of related researches and proposed various relatively effective methods. The representative numerical algorithms are as follows: the generalized forward-backward method (GFBM) [1, 2], the finite element method (FEM) [3, 4], the finite-difference time-domain (FDTD) algorithm [5], the 3D multilevel UV method [6], and so on. Due to the problems of huge memory consumption and low computational efficiency, the numerical methods are difficult to handle the composite scattering modeling for the complex target over a larger sea surface, which have a big disparity to the engineering application of actual problem. Another kind of methods are based on the high frequency hybrid algorithm, such as the four path model (FPM) [7, 8], the iterative physical optics (IPO) method [9], the shooting and bouncing ray (SBR) algorithm [10], the hybrid Kirchhoff approximation and Method of Moment (KA-MoM) algorithm [11], the bidirectional analytic ray-tracing (BART) method [12], the extended propagation-inside-layer expansion (E-PILE) method [13], the geometrical optics and physical optics (GO-PO) method [14], and so on. Each of these approaches has its merit and demerit. Comparatively speaking, the GO-PO method is a reliable method to calculate the multiple EM scattering of the target, and it can also be applied to the calculation of coupling scattering of the target and rough surface. But, this method faces a difficulty that the lower efficiency in testing the visibility of patches to the incident and reflected waves. Aiming at the low efficiency problem of the EM scattering model, there are some acceleration schemes, among which the kd-tree method is a popular one. However, the kd-tree method is mainly used to deal with the scattering of the pure target. Relevant studies on the application of the kd-tree method to the composite model of a larger sea surface and complex target is still not sufficient and even very little.

Therefore, in order to improve the efficiency of the GO-PO method in calculating the scattering of the ship target and the coupling scattering between sea surface and ship target, this paper applies the kd-tree method [1517] to optimize and improve the GO-PO method [14] proposed in our earlier work. When applying the kd-tree method, the simulation scene should be divided into multilevel subbounding boxes to construct the kd-tree. Then, the intersection test between the ray and the bounding box, where the patch is located, is firstly implemented before the intersection test between the incident wave or reflected wave ray and the patch. If the distance parameters of the bounding box do not meet the requirement, the ray is not necessary to enter the bounding box, which can reduce the needless intersection tests. Therefore, the model combined with the kd-tree and GO-PO methods can greatly improve the computing speed, which can be applied to the study of composite EM scattering characteristics of electrically large sea surface with a complex ship target.

The rest of this paper is organized as follows. Section 2 describes the EM scattering model for the composite scene of sea surface and target. Section 3 firstly introduces the method of kd-tree construction for the composite model in this study, and then gives how to test the visibility factors based on kd-tree. Section 4 gives some numerical results, including the analysis of the scattering characteristics of the pure ship target, the scattering characteristics of the composite model, and the influence of the number of patches on the RCS. Moreover, the comparisons of the efficiency between the kd-tree accelerated hybrid method and the original one are shown. Section 5 forms the conclusion.

2. EM Scattering Model of Sea Surface and Target

As shown in Figure 1, the scattering contributions from the composite scene of sea surface and target mainly include sea surface scattering contribution, target scattering contribution (which includes primary and multiple scattering), and the coupling scattering contribution between the target and sea surface. The coupling effect includes the multiple scattering contribution of the target due to the sea surface reflection wave and the multiple scattering contribution of sea surface due to the target reflection wave. In Figure 1, is the unit vector of the incident wave; , , and are the unit vectors of the reflected wave of -th and -th patches on the target and -th patch on sea surface, respectively. , and are the visibility factors of -th and -th patches on the target and -th patch on sea surface to the incident wave. If one patch, such as -th patch on the target, can be irradiated by the incident wave, the visibility factor is valued by . On the contrary, is set as ; is the visibility factor of -th patch on the target to the reflected wave of -th patch on the target; is the visibility factor of -th patch on sea surface to the reflected wave of -th patch on the target; and is the visibility factor of -th patch on the target to the reflected wave of -th patch on sea surface.

Then, the total scattering field from the composite scene of sea surface and target is given by the superposition of the three kinds of scattering contribution. The scattering contribution from sea surface in this study is computed based on the facet-based asymptotical model (FBAM) proposed in our earlier work [18]. The scattering contribution from the ship target and the coupling scattering contribution between target and sea surface are computed by the GO-PO method; more details can be found in [14] proposed in our earlier work. Although the GO-PO method is an effective method to calculate the multiple scattering of targets and the coupling scattering between sea surface and target, it still faces the difficulty that the lower efficiency of testing the visibility factors , , , , and (which are shown in Figure 1) due to the enormous number of triangular patches on the composite scene and the complex structures of the ship target. In the following section, the kd-tree method is described to improve the efficiency of the GO-PO method.

3. Visibility Tests Based on Kd-Tree

3.1. Kd-Tree Construction

First of all, the ship target model is partitioned into small triangle patches with the help of CAD tools. Subsequently, the sea surface is generated based on the Monte Carlo method and decomposed by a mount of tilted square facets, which can be easily divided into triangle patches through computer programming. Then, the patches on the ship target below the waterline and the patches on the sea surface below the ship target are both discarded. So far, all triangle patches on the sea surface and ship target are integrated and uniformly numbered.

Then, the entire simulation scene will be divided into multilevel subbounding boxes to construct the kd-tree. The largest outside bounding box of the entire simulation scene is a cuboid described by the maximums and minimums of x, y, and z coordinates. This largest bounding box is marked by , which corresponds to the root node and contains all patches. Next, the cuboid bounding box is divided into two child cuboid bounding boxes and along an axial direction, which form a binary tree, namely, a root node is split into two internal nodes or leaf nodes. The splitting plane has the highest intersecting probability with the ray, and the splitting plane can be any plane perpendicular to x, y, and z axes. If the following two conditions, that is, the number of patches contained in a subbounding box is larger than the user-defined number and the depth of the kd-tree is less than the maximum depth, are both satisfied, the subbounding box needs to be further subdivided. The above subdivision process is carried out recursively until one of the abovementioned two conditions cannot be satisfied. In addition, if the patch is across the splitting plane, it should belong to the subbounding box where the center of the patch is located. Aiming at the composite model of the sea surface and ship target in this study, the kd-tree construction for the model is given, as shown in Figure 2. In Figure 2, is an internal node corresponding to the subbounding box above the sea surface. are the leaf nodes, in which there are no patches. are the internal nodes containing a certain amount of patches. It is important to note that contains all patches on the ship target, and contains all patches on the sea surface. Obviously, the two subbounding boxes need to be further subdivided, and the depth of the subdividing of and is set as twenty in this study. is then divided into two child cuboid bounding boxes and . Due to space limitations, further illustration is omitted. The detail of kd-tree construction is well described in Pharr and Humphreys’ book [16].

3.2. Visibility Tests

Figure 3 gives the test flow of the visibility factors and . Accordingly, the tests of and are similar to those of , and the tests of are similar to those of . Once the kd-tree is constructed, the rays for all incident and reflected directions can be traced efficiently. Differently, the backward ray-tracing is used for all steps.

First of all, we should label the patches in the entire simulation scene that satisfy the condition ( is normal vector of the patch), for only these patches can possibly be illuminated by the incident wave. And, these patches should further be tested whether they are shaded by other patches, which also satisfy the condition . As Figure 4(a) shows, taking backward central ray-tracing of the -th patch, for example, the traversal begins with the leaf node that the -th patch belongs to. At the leaf node, we let the backward ray of the -th patch do the intersection test with the plane , where the -th patch is located by an intersection point ( and are the midpoints, is the normal vector, and represents the distance parameter between an intersection point and the starting point of the ray), whose distance parameter is denoted by . If is an interior node of the -th patch, the -th patch is shaded by the -th patch when the incident rays pass into the -th patch, and should be greater than zero. If the ray is intersected with one patch, the traversal is immediately stopped and the visibility factor of the -th patch to the incident wave is marked by . On the contrary, if there is no intersection at the leaf node, the traversal moves onto the next node. What needs to be emphasized is that the intersection test between the backward ray and the bounding box of the node is firstly done. If the distance parameter of the bounding box do not meet the requirement, namely, , as shown in Figure 4(b), the ray is not necessary to enter the bounding box. If , the ray enters the bounding box. Meanwhile, if the node is an interior node, the kd-tree is traversed down until a leaf node is encountered, and then the intersection test is processed consequently. For every node encountered in the process of traversing, we firstly let the ray do the intersection test with the bounding box of the node. Afterward, the traversal continues until the nil node is encountered or the ray is intersected with one patch. In the entire course, if the nil node is encountered, we consider that the -th patch is illuminated by the incident wave, and the visibility factor of the -th patch to the incident wave is marked by .

For the determination of the visibility factor , we let the -th patch do the test with the -th patch that is illuminated by the incident wave (i.e., ). During the test, if the direction of the reflected wave of the -th patch does not satisfy (i.e., the -th patch cannot be illuminated by the wave reflected by the -th patch), the test skips the -th patch to the next patch that is illuminated by the incident wave. Otherwise, the backward ray test with the -th patch is done. If the ray intersects the -th patch and the distance parameter of the intersection point is denoted by , it suggest that the -th patch is likely illuminated by the wave reflected by the -th patch unless it is shaded by the other patch between them. Accordingly, we implement the backward ray-tracing to find whether the ray is cut off by the other patch during the propagation. The traversal process is similar to the visibility test of . Except for the difference that if the -th patch is shaded by the -th patch (as shown in Figure 5(a), the distance parameter of the intersection point between the backward ray of the -th patch and the -th patch is denoted by ), should be less than . When there is no intersection at the leaf node that the -th patch belongs to, the traversal moves onto the next node. Then, the intersection test between the backward ray and the bounding box of the node is done. If the distance parameters of the bounding box do not meet the requirement, namely, , as shown in Figure 5(b), the ray is not necessary to enter the bounding box. If , the ray enters the bounding box and performs the further operations, as shown in Figure 5(c).

4. Results and Discussion

This example is given here to validate the improved kd-tree accelerated method. The RCS of a cube on a PEC square plane is shown in Figure 6. The result is simulated in the backscattering case and for VV polarization. The wavelength of the incident wave is set by . The side length of the square plane is and the cube is . From Figure 6, a good agreement is observed between the improved method and the original one. Compared with the original method, the improved method is the coalition of the original method and the kd-tree. And, the patches that have influence on the scattering are the same for the two methods. The kd-tree added in the improved method does not affect the accuracy of the result.

Next, we analyse the performance and efficiency of the improved kd-tree accelerated method in calculation of the composite scattering from the sea surface and ship target. The geometric structure and dimensional parameters of the ship target and the coordinate system for the composite scattering of the sea surface and ship are shown in Figure 7, where is the incident angle, incident azimuth angle, scattering angle, and scattering azimuth angle. If there are no special instructions, the following results are computed at 5 GHz and VV polarization. For the bistatic scattering, the incident direction is , , and the scattering azimuth angle is . The 2D sea surface is generated using the Elfouhaily et al. [19] spectrum; the sampling number is Nx = Ny = 256, the interval is dx = dy = 0.76 m, the size is Lx = Ly = dx (Nx − 1) = 193.8 m, the number of triangle patches on the sea surface is 2  (Nx − 1)  (Ny − 1) = 130050, the wind is upwind at a speed of 5 m/s, and the relative dielectric constant of the sea water is calculated by the Klein and Swift model [20] at 20°C and 35% salinity.

Figure 8 gives the bistatic RCS of the pure ship target: (a) is with different methods and (b) with different number of patches. From Figure 8(a), one can observe that the result evaluated by the improved method agrees well with the result given by the original method. From Figure 8(b), it is worth noting that, with the increasing number of patches, the RCS of the ship increases on the whole. In addition, the RCS has a bigger change when the number of patches increases from 4340 to 9576. And, it has a little bit less change when the number of patches increases from 9576 to 27527. It is because that, in the GO-PO method, as long as the meshing patch can give a relatively precise description of the illumination states, the patches need not to be partitioned smaller.

To further contrast the efficiency between the improved kd-tree accelerated method and the original one, Table 1 gives the time consumed by the two methods in the simulation of bistatic scattering from the pure ship target with different number of patches, which correspond to the three curves in Figure 8(b). For the original method, all the patches of the model need to be tested the visibility. With the improved one, only the patches in the specific bounding box (i.e., its distance parameter meets the corresponding requirement) need to be tested. As shown in Table 1, the original method consumes much time and the improve method increases efficiency about three times.

Figure 9 gives the comparisons of the bistatic RCS of total field from the ship and sea surface: (a) is with different number of patches on the ship and (b) with different number of patches on sea surface. In Figure 9(a), the number of patches on the ship target for the three curves is 4340, 9576, and 27527, respectively. The number of patches on sea surface is 130050. In Figure 9(b), the number of patches on sea surface for the three curves is 32258, 79202, and 130050, respectively. It should be noted that the red dash line is with Nx = 128, dx = 0.76 m, Lx = 96.52 m, and the number of patches is 32258, and the green dotted line is with Nx = 200, dx = 0.76 m, Lx = 151.24 m, and the number of triangle patches is 79202. The number of patches on the ship target is 9576.

From Figure 9(a), one can observe that, with the increasing number of patches on the ship, the backscattering result increases on the whole, but the forward scattering result basically remains unchanged. It is because that, for the composite scattering of the sea surface and ship, the scattering in the forward direction is dominated by the contribution of sea surface and the scattering in the backward direction is dominated by the contribution of the ship target and coupling field. In Figure 10 that displayed in the following text, one can clearly observe such a phenomenon. From Figure 9(b), with the increasing number of patches on the sea surface, the result increases on the whole except that at the range of θs = −55∼−45 and near θs = −90, the result is almost unchanged. This is because that the simulations of the three lines are with the same ship target but with different size of sea surfaces. And, the scattering is dominated by the contribution of the ship target at the range of θs = −55∼−45 and near θs = −90.

Table 2 gives the time consumed by the improved kd-tree accelerated method and the original one in the simulation of composite scattering from the sea surface and target with different number of patches, which correspond to the curves in Figure 9. As shown in Table 2, the improve method increases efficiency about three times.

The following simulations give the scattering characteristics of the ship and sea surface. The number of the triangle patches on the ship target is 9576 and that on the sea surface is 130050. Figure 11 shows the change of RCS of the pure ship target with incident angle and incident azimuth angle. The result is simulated in the backscattering case. From Figure 11, one can obviously see the RCS distribution of the ship in different incident directions. Figure 12 gives the change of RCS with incident angle at different incident planes. The curves plotted at planes of , , and , respectively, represent the incident EM wave slanting at the bow, stern, and side of the ship. It is observed that, when the incident wave slant at different directions, the peak values of the RCS will be distributed at different incident angles due to the different visibility of the patches to the incident wave and reflected wave and the different shade relationship between the patches.

Figure 13 displays the change of bistatic RCS of the ship with scattering angle and scattering azimuth angle. In more details, Figure 14 gives the change of RCS with scattering angle at different scattering planes. And, the larger RCS can be observed in the scattering planes of and , and the RCS is relatively small in the scattering plane of .

Figure 10 gives each contribution of bistatic scattering from the ship target on a sea surface, including RCS of the target, sea surface, coupling field, and total field. From Figure 10, it is worth noting that the presence of ship significantly enhances the bistatic scattering in the backward direction, especially at the range of θs = −55∼−45 and near θs = −90. In the forward direction, the scattering is dominated by the sea surface contribution.

Figure 15 gives the simulations for the change of RCS of coupling field and total field from the composite scene with different wind speeds (the wind is at a speed of 5 m/s and 10 m/s). From Figure 15, it is worth noting that the effect of wind speed on coupling field is not obvious. The apparent reason for such a situation is that with the increasing wind speed, the roughness of the sea surface increases, but the slopes of facets on sea surface along the x direction have a little change when the wind direction, incident direction, and observation direction are equal to 0°, so the change of the visibility of patches on ship and sea surface to the incident wave and reflection wave is not obvious. For the total field, the level of RCS increases when the wind speed increases except in the mirror direction, the total field decreases with the increase of the wind speed.

Figure 16 displays the change of bistatic RCS of coupling field and total field from the ship and sea surface with different polarization. From Figure 16, one can observe that the coupling field at VV polarization is larger than that at HH polarization. For the total field, the forward scattering result at HH polarization is larger than that at VV polarization, and the backscattering result has opposite situation.

5. Conclusion

A kd-tree accelerated hybrid model combined with the GO-PO method and FBAM is presented to calculate the scattering of the ship target and coupling scattering between sea surface and ship target. Compared with other calculation models, this model not only ensures the calculation accuracy, but also improves the efficiency by accelerating the program of testing the visibility of patches to the incident and reflected waves. Therefore, the hybrid model is relatively effective and can be applied to the study of composite EM scattering characteristics of electrically large sea surface with a ship target. However, the EM coupling mechanism between sea surface and target is very complex and difficult to be fully expounded. So, some insight into the coupling scattering mechanisms still remains to be continued.

Data Availability

The simulation data used to support the findings of this study are included within the article.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

This work was supported in part by the National Natural Science Foundation of China (grant nos. 61801416, 61701428, and 61861043), in part by Shaanxi Provincial Science and Technology Department (grant nos. 2019JQ-237 and 2019JQ-120), in part by Shaanxi Provincial Education Department (grant no. 18JK0872), in part by the Scientific Research Foundation of Yanan University (grant no. YDBK2016-16), and in part by the Open Foundation of Fudan University Key Laboratory for Information Science of Electromagnetic Waves (MOE) (grant no. EMW201910).