Width estimation of a subsurface water-filled crack using internal multiple reflections from electromagnetic scattering model of a cylinder

Abstract

A width estimation method for a subsurface water-filled crack using internal multiple reflections is proposed. The energy transport and path length of internal multiple reflections are analyzed from electromagnetic scattering model of a dielectric cylinder using the extended ray theory. The energy transport shows that some multiple rays transport along the surface of the cylinder to the receiver, which form the creeping waves, and some multiple rays refract toward the receiver, focusing on the radial direction. The conclusion that the path lengths are multiples of a diameter distance is utilized to estimate the width of the cylinder. In order to distinguish the adjacent reflections with small time interval, wavefronts and resonances of internal multiple reflections are also analyzed, based on the fact that the creeping waves and refracted waves at the radial direction show different resonances in the time-frequency domain. In order to verify the proposed method, experiments both in free space and subsurface scene are given to estimate the width of water-filled cracks.

Keywords:

• Internal multiple reflections
• dielectric cylinder
• wavefronts
• width estimation

1. Introduction

The width estimation of subsurface water-filled cracks is important in the diagnosis of underground anomalies. If the cracks cannot be detected in time or estimated correctly, they would extend and may result in serious security accident [1–7]. In recent years, many techniques have been used for crack detection and estimation [2–4].The possibility of crack detection with GPR in asphalt pavement was proved by Ahmad [3]. In the study, the influencing parameters and limitations were analyzed for crack detection by GPR. However, the cracks are visible on the surface of the pavement. For a totally invisible subsurface crack, it has not been mentioned. The apertures of water-filled and air-filled fractures were determined by Koivisto [5]. In the study, the fractures were layered. The GPR signal polarity was utilized to distinguish fractures with different material. Moreover, the relationship of fracture aperture and vertical resolution of the antenna was investigated. These studies provide valuable methods for crack estimation. However, the internal scattering process has not been involved.

Usually, cracks are regarded as cylindrical shape to obtain analytic solutions in theoretic analysis [8–10]. The electromagnetic (EM) scattering of a dielectric cylinder with circular cross-section by means of transient time-domain technique has been made great progress [11–15]. However, for the case of a dielectric cylinder with electric size comparable to the working wavelength, the geometrical theory of diffraction interpretation and internal multiple reflections should be also considered in order to obtain a physical point of view of the scattering characteristics of such target [16–19]. An excellent research on such problem of EM scattering from a dielectric cylinder lying in a dielectric medium was studied by Ludwig and Leviatan [11]. In this research, the ray tracing principle for a dielectric cylinder illuminated by a plane wave was illustrated. The detailed process of the ray tracing was explained, which includes the incident ray, backscattered ray, creeping ray, and refracted ray, providing guidelines to intuitively understand the process of EM wave travelling in the dielectric cylinder.

With the deep exploration of the dielectric target, the internal multiple reflections have been attracted more interestingly [20–24]. An EM scattering process on a dielectric cylindrical object was analyzed based on an extended ray theory (ERT) by Ohmori et al. [20]. In the study, the transient scattering responses by a circular dielectric cylinder were constructed by the ERT, according to the specular reflections and internal multiple reflections. The rigorous solution of the EM scattering from an infinite dielectric cylinder was studied by Sasamori et al. [23] in terms of multiple inside refracting and creeping waves. In their study, the geometric optics expression, such as reflection, transmission, and diffraction, was derived by means of appropriate high-frequency approximation. Particularly, the physical meaning of the geometric optics was explained in detail.

The researches above gave high-frequency approximation models for the consideration of contribution of the internal multiple reflections of a dielectric cylinder. However, the propagation properties of the internal multiple reflections related to the size and dielectric feature of the cylinder have not been involved.

In this paper, a width estimation method for a subsurface water-filled crack using internal multiple reflections is proposed. The crack is modeled by a dielectric cylinder and the electromagnetic scattering model is analyzed using the ERT. The energy transport and path length of internal multiple reflections are deduced. The conjoint analysis of wavefronts and resonances is performed to show the characteristic of internal multiple reflections in time-frequency domain. Then, the time interval of adjacent wavefronts at different resonant frequencies is calculated to estimate width of the cylinder. Finally, experiments both in free space and subsurface scene are given to verify the proposed method.

2. Estimated method based on EM scattering model of a cylinder

2.1. Analysis of EM scattering model for a dielectric cylinder in free space

The EM scattering model can be obtained by analyzing the ray tracing paths based on ERT for a dielectric cylinder with circular cross-section of radius a illuminated by a plane wave. The cylinder is homogeneous and isotropic with relative permittivity ε_r, and relative permeability μ_r. The relative permeability is equal to that of in free space so that μ_r = μ₀. The cylinder is long enough compared to the radius such that the boundary of circular ends of the cylinder can be ignored. In this paper, the plane wave illuminates normally to the axis of the cylinder with the direction of the electric field perpendicular to the axis of the cylinder. It should be noted that the wavelength of the incident wave is comparable to the electrical size of the cylinder such that the diffraction and transmission cannot be neglected. Backscattered field with co-polarization is received to obtain strong echoes. The incident wave is represented by a system of rays, striking the surface at the points of impact [25]. Although the plane wave illuminates the cylinder normally, each ray strikes the point of impact at a local incident angle θ_i from the incident ray to the normal. It is obvious that the local incident angle does not exceed 90^°, as shown in Figure 1(a).

Figure 1. (a) A dielectric cylinder impinged by a plane wave. (b) Illustration of ray path for a dielectric cylinder impinged by a plane wave. (c) Propagation processes of internal multiple reflections.

Figure 1(b) shows a ray tracing schematic diagram. The specularly reflected ray (SR) identifies the scattering event that the EM wave is reflected totally on the surface of the cylinder. Moreover, the incident wave traverses the interior of the cylinder and forms multiple internal reflected rays (IR_m, m = 1, 2, …). is the number of wave traveling from one side to the other inside the cylinder [20].

According to the geometrical optics ray path schema and the ERT [26], the received EM waves are given in terms of contributions from the specularly reflected rays and the multiple internal reflected rays. We focus on the latter to investigate the internal multiple propagation properties. Firstly, the global coordinate system XYZ is set up, which satisfies the right-hand screw rule and let Z-axis perpendicular to the surface of the paper. A ray striking the point ϕ_0u is given as an example, as shown in Figure 1(c). The subscript u denotes the point on the surface and . The point of impact behaves locally like an infinitesimal plane interface, wherein the tangent line of this point is set up as x-axis along the clockwise direction and the outward normal vector of this point is regarded as y-axis and the direction of z coincides with that of Z. The ray at a local incident angle θ_i will partly reflect out as SR ray and partly refract inside the cylinder subject to the Snell’s law, i.e. . Here, θ_i and θ_r are local incident angle and refracted angle in the local coordinate system xyz.

The first time refracted ray is denoted as . When the ray reaches the point ϕ_1u inside the surface of cylinder, it partly reflects internally and the rest refracts outside the cylinder, which are denoted as and , respectively. Similarly, the reflected ray of again propagates as two parts at the reflected point ϕ_2u, denoted as and , respectively. As seen from the process, the ray propagating m times inside the cylinder will result in m-time output refractions. The propagation process will continue until the energy vanishes.

Firstly, we attempt to identify the energy transport of internal multiple reflections. For each local coordinate system xyz determined by the incident or reflected point ϕ_vu, v = 0, 1, 2, …, u = 0, 1, 2, …, the relationship between the local coordinate system xyz and the global coordinate system XYZ is shown as follows:(1)

where α is the angle from Y to y, relating to the local incident angle θ_i and its corresponding refracted angle θ_r. According to the geometric relationship, α is derived as α = mπ + θ_i – 2 mθ_r.

The incident electric field and magnetic field in the global XYZ coordinate system are given as(2)

where E₀ is the magnitude of incident field and η₀ is the intrinsic impedance in free space. Θ_i is from the phases determined from the path lengths of each scattering event. Assume the EM wave reflects (m−1) times and then the m-time refracted electric field and magnetic field are given by(3a) (3b)

where is from the phase and will be discussed later. R₁₁ is the reflection coefficient at each reflected point ϕ_vu, v = 1, 2, …, u = 1, 2, …. τ₀₁ and τ₁₀ are the transmission coefficients from free space to the cylinder at the incident point ϕ_0u and the transmission coefficients from the cylinder to free space at the output point ϕ_vm, v = 1, 2, …, respectively. Here, the point of impact behaves locally like an infinitesimal plane interface. With this respect, the reflection coefficient and transmission coefficients are considered as the same with those in plane interface, subject to the local coordinate system xyz. Thus(4)

Define , and then (4) can be written as .

Then the transmission coefficient τ₀₁ is given as(5)

Similarly, the reflection coefficient R₁₁ and transmission coefficient τ₁₀ are derived as follows:(6)

The geometrical divergence coefficient D is deduced from the reference [25] and given directly as .

Thus (3) can be rewritten as(7a) (7b)

Then the energy transport expressed by Poynting vector in the global coordinate system XYZ via (7), combining (1), can be written as(8)

In (8), it can be found that due to θ_i < 90°. Hence, when m is odd, the energy transport direction is between and or and . This means that for the backscattering (i.e. ) received responses, the energy from the direction cannot be received at all but part of the energy from or direction may be received due to the creeping waves along the tangent of the cylinder. However, when m is even, the energy propagation direction is between and or and . This indicates that the energy from the direction can be totally received and part of the energy from or direction may also be received. Hence, for the odd m-time internal multiple reflections, the received energy from reflected point ϕ_mu may be smaller than that for the next point ϕ_(m+1)u, although the received energy from point ϕ_mu should be larger than that for the point ϕ_(m+1)u due to the path loss.

Then, we consider the path length of waves propagating inside the cylinder. As illustrated in Figure 1(c) the path length l_n from the point ϕ_vu to the point ϕ_(v+1)u can be written as(9)

where d is the diameter of the cylinder and d = 2a. It is obvious that the path length is affected by the local incident angle, the permittivity and the diameter. When the EM wave is normal incident, i.e. θ_i = 0^°, which result in the refracted angle is also zero, the path length equals to the diameter d. For other incident angles, the refracted angle is non-zero and the path length is the chord length, which is shorter than the diameter. When the local incident angle is unchanged, the refracted angle will vary as the permittivity changes. In this paper, the diameter d is to be estimated and thus it is unchanged. Then the path length can be calculated by analyzing the coefficient of d, namely, the ratio of path length l_n to d.

In Figure 2(a), the ratio of path length to the diameter of the cylinder with respect to incident angles decreases as the incident angle increases. For the case of the relative permittivity larger than 20, the ratio can approach 1. This indicates that the path length is almost the distance of a diameter. As the permittivity decreases, the ratio deviates from 1 to a certain extent. However, the deviation extent is small. For instance, the ration is still close to 0.87 at a large incident angle for the case of . Similar conclusion can be also obtained for the ratio with respect to permittivities, as shown in Figure 2(b). For the case of , the ratio is beyond 0.9 at each incident angle. From an overall view, it can be concluded that for most cylindrical material of , the internal propagating length can be regarded as the multiples of the diameter.

Figure 2. (a) Ratio of path length to the diameter (L/d) with respect to incident angles. (b) (L/d) with respect to relative permittivities. (c) (t/τ) with respect to incident angles. (d) (t/τ) with respect to relative permittivities.

Now we consider the propagation time computed by from the phase determined from the scattering process. The following equation is used as the time delay from the propagation of the EM wave inside the cylinder [24]:(10)

where indicates the real part. The time delay in (10) is normalized by the timeτ. Let for the EM wave propagating the distance of diameter d = 2a of the cylinder. Here, c is the speed of light in free space.

In (10), the path length of the ray propagating m times inside the cylinder can be given by(11)

where Φ_i(ϕ_0u) is determined from the phase of the incident ray at the point ϕ_0u while Φ_s(ϕ_mu) from the phase at the point ϕ_mu to the observation point. s(ϕ_mu) is the length of the ray propagating along the surface from the refracted point ϕ_mu to the observation point. l_n is the distance from the incident point ϕ_0u to the output point ϕ_mu.

Here, Φ_i(ϕ_0u) is given as the average distance from the top point of the circular cylinder to the end of the lit region, i.e. Φ_i(ϕ_0u) = a/2. Φ_s(ϕ_mu)and s(ϕ_mu) are determined according to the path propagation and energy transport in Figure 1(c). When m is odd, Φ_s(ϕ_mu) ≈ 0 and ; When m is even, while s(ϕ_mu) ≈ 0. Then the time delay between the odd m-time reflection and its next even m-time reflection was analyzed, i.e.(12)

According to (12), the time delay normalized to the time of propagation of a diameter distance is plotted as the incident angles or permittivities change, as shown in Figure 2(c) and (d). Figure 2(c) shows that the normalized time delay is close to 1 at all incident angles for most permittivities. Even for the case of , the normalized time delay can reach around 1.2. Further on, for much smaller permittivities, the normalized time delay will beyond 1 to a greater extent at small incident angles (θ_i = 0^°, 24^°), as shown in Figure 2(d). In practical application, most geological scatters possess the relative permittivity much larger than 4. Particularly, for the water-filled crack, the relative permittivity is approximate 81 and the normalized time delays at all incident angles are quite close to 1, which means that the time delay of the EM wave traveling inside the cylinder between the odd m-time reflection and its next even m-time reflection is the time the EM wave propagating a diameter distance of the cylinder.

In order to verify the analysis, an example is given by the commercial software XFDTD based on the Finite Difference Time Domain method. A cylinder with diameter , relative permittivity ε_r = 81 is illuminated by a Gaussian plane wave with center frequency 2.5 GHz. The backscattering field is shown in Figure 3.

Figure 3. Backscattering field of a cylinder with and .

In Figure 3, it is found that the scattering events in Figure 1 are well reflected in the time-domain responses. The SR shows the maximum magnitude and then the internal multiple reflections are reflected by subsequent peaks, as labeled as IR_m. Based on the previous analysis, the normalized time delay is given as(13)

where is the time of the m-time output refraction to the observer. Here the time indicates the time of EM wave propagating a diameter distance of the cylinder. The theoretic time delays and the simulated time delays are calculated according to (10) and (13), respectively. The results are listed in Table 1.

Table 1. Theoretic and simulated values of normalized time delay.

m	2	3	4	5	6
Theoretic	1.12	2.11	3.11	4.11	5.11
Simulated	1.00	2.15	3.15	4.27	5.29

The simulated values show good agreement with the theoretic values. In fact, this also corresponds to the previous analysis that the time delay of the EM wave between the odd m–time reflection and its next even m–time reflection is the time the wave propagating a diameter distance of the cylinder. Besides, it is found that the magnitude of IR4 is much smaller than that of IR5. This can be explained by (8) that some energy transports along the opposite direction of the observer such that only some creeping energy can be received for the case of IR4. However, for the case of IR5, part of the energy transports along the received direction and also creeps along the surface to be received. From the time-domain responses in Figure 3, it can be observed that for the even m-time internal reflection, such as IR1, IR3, and IR5, there are gradually decreasing oscillations behind the primary peak. This may be caused by the effect of dielectric material on the EM wave propagation. However, for the odd m-time internal reflection, the refracted outside wave may be filtered by the creeping propagation along the surface of the cylinder.

2.2. Electromagnetic scattering model for a subsurface dielectric cylinder

The dielectric cylinder discussed above is located horizontally at depth h in a lossless ground with a relative permittivity ε_s. A plane wave illuminates the cylinder normally from free space. The analysis of waves acting on the cylinder is similar to that in free space. However, the initial waves impinging on the cylinder are not the incident waves in free space any more, but the refracted waves from free space to the ground. Meanwhile, the received waves also undergo the propagation from the ground to free space. With this respect, the received field can be revised as(14)

where Γ₀₁ and Γ₁₀ are transmit coefficients from free space to the ground and from ground to free space, respectively. k_s is the wavenumber in the ground. h_s is the distance of the refracted waves traveling in the ground and the range of h_s is from h to (h + a). Here, ξ^′ is used to replace ξ in the case of free space and . The replacement is also performed in D in (7) and (8), which is then denoted as D^′ in (14). The path length is given as . In fact, substitute this path length into (10) and it is found that the result will not change. The actual magnitude of time-domain responses from the subsurface cylinder will be much smaller than that in free space due to a great deal of loss from waves reflected directly from the ground surface.

In order to highlight the scattering responses from the cylinder, the impulse responses of the cylinder will be deconvoluted from the system responses. The detailed process is provided in reference [27].

3. Wavefronts and resonances analysis for dielectric cylinders

If the time interval between the m-time internal reflection and (m+1)-time internal reflection is not large enough to be distinguished in time-domain, other information, such as frequency, may be needed. Hence, time-frequency analysis based on the smooth-pseudo-Wigner-Ville (SPWV) distribution is used to reveal the scattering mechanisms [28]. It is known that wavefronts in time-frequency plane will be displayed as vertical lines corresponding to the occurrence of the scattering events (i.e. the peaks in time-domain). The resonances will be displayed as horizontal lines indicating that the scattering events interacted with each other.

Examples are given as follows to illustrate wavefronts for dielectric cylinders. The responses from cylinders with diameters of 4 and 8 cm, relative permittivities of 10, 33, 55, and 81, are computed by XFDTD, respectively. The incident Gaussian plane wave with center frequency 2.5 GHz illuminates a cylinder normally. The SPWV distribution is used for further analysis. In this subsection, the SPWV is performed with time smoothing window of 500, frequency smoothing window of 1125 and threshold of 0.01% for each time-domain response.

As shown in Figures 4(a) and 5(a), the wavefronts are coupled seriously due to small permittivity. The waves travel faster inside the cylinder with smaller permittivity comparing to those with greater permittivity, and then the received wave is considered to be a superposition of several waves. However, the wavefronts become gradually clearer as the permittivity increases, as shown in (b), (c) and (d) in Figures 4 and 5, respectively. In the time-frequency plane, the first wavefront with high energy corresponds to the specular reflection. Then the subsequent wavefronts show the internal multiple reflections. Hence, the number of wavefronts indicates the times of internal multiple reflections and the time interval between wavefronts can be used to estimate the diameter of the cylinder. The estimated results are listed in Tables 2 and 3, where Δt_m is the time interval between m-time and (m + 1)-time propagation inside the cylinder, which is reflected in terms of wavefronts in SPWV distribution. d_e is the estimated diameter (or width) of the cylinder according to the time interval.

Figure 4. (a) The SPWV distribution of the electric field from cylinder with d = 4 cm and . (b) From the cylinder with d = 4 cm, . (c) From the cylinder with d = 4 cm, . (d) From the cylinder with d = 4 cm, .

Figure 5. (a) The SPWV distribution of the electric field from cylinder with d = 8 cm and . (b) From the cylinder with d = 8 cm, . (c) From the cylinder with d = 8 cm, . (d) From the cylinder with d = 8 cm, .

Table 2. Width estimation of the cylinder with diameter 4 cm.

εr	m	1	2	3	4	5	6
10	△tm (ns)	–	1.00	–	–	–	–
	de (cm)	–	9.49	–	–	–	–
33	△tm (ns)	–	–	0.78	1.85	1.65	1.90
	de (cm)	–	–	4.07	9.66	8.62	9.03
55	△tm (ns)	1.36	1.21	1.07	1.21	0.61	1.14
	de (cm)	5.50	4.89	4.33	4.89	2.48	4.60
81	△tm (ns)	1.36	1.15	1.21	1.43	1.08	1.68
	de (cm)	4.56	3.80	4.06	4.76	3.60	5.60

Table 3. Width estimation of the cylinder with diameter 8 cm.

εr	m	1	2	3	4	5
10	△tm (ns)	–	0.94	0.76	–	–
	de (cm)	–	8.92	7.21	–	–
33	△tm (ns)	1.60	1.59	1.44	1.53	1.44
	de (cm)	8.36	8.30	7.52	7.99	7.52
55	△tm (ns)	1.90	2.05	1.98	2.05	1.93
	de (cm)	7.69	8.31	7.99	8.29	7.81
81	△tm (ns)	2.35	2.44	2.51	2.50	2.44
	de (cm)	7.83	8.13	8.37	8.33	8.13

It can be found that for the case of in Figure 4, most wavefronts cannot be observed and the estimated results are failed. This can be explained by the fact that the range resolution in such dielectric material is quite close to the diameter 4 cm of the cylinder [5]. Thus, the wavefronts are difficult to be distinguished. However, for the diameter 8 cm of the cylinder with the same permittivity, some wavefronts can be distinguished, as shown in Figure 5(a). Hence, the width can be estimated and the results are close to the true value, as listed in Table 3.

Comparing the results in Tables 2 and 3, it can be also seen that for the larger size and permittivity of the cylinder, the estimated result is better. In general, the results in Table 3 are better than those in Table 2. This is because the EM wave propagates longer path length in the cylinder with larger size and thus the time interval between wavefronts is large enough to be distinguished. In Table 2, some estimated results are approximately two times the true width (e.g. for the case of ). From the time-domain responses in Figure 4(b), it is found that the odd m-time refracted wave is too weak to be detected due to the small size of the cylinder. In this case, the wavefronts presented in the time-frequency plane show two adjacent even m-time refracted waves, which indicates two times the diameter distance. In fact, the energy from the odd m-time refracted wave is not strong for other cases. However, the wavefronts from odd m-time propagation inside the cylinder and those from even m-time propagation present alternately at different resonant frequencies. As an example given in Figure 5(c), the wavefronts from m = 3, 5, 7 in time-domain responses distribute at resonant frequencies between 0.75 and 1.22 GHz. The wavefronts from m = 4, 6, 8 distribute at resonant frequencies between 1.22 and 1.75 GHz. With this respect, the time interval between a wavefront at lower resonant frequencies and its next wavefront at higher resonant frequencies can be utilized to estimate the width of the cylinder.

4. Experiment results

4.1. Experiment in free space

To further verify the proposed method, real data are acquired by the experimental GSSI radar according to the scenario shown in Figure 6(a). The incident wave is the derived Gaussian impulse with center frequency 1 GHz. Water with , filled by a cylindrical plastic bottle, made of polycarbonate, a medium nearly transparent to radar waves [29]. In practical application, city pipes such as water pipelines are usually produced by PVC material. The water-filled PVC is also tested. The geometrical parameters are detailed displayed in Figure 6(a). In order to imitate the free space environment, the target is sited on a foam platform, which has negligible effect on our measurement. To avoid the interference of the background, the background measurement without a target is obtained firstly and then the measurement with a target is tested. The singular value decomposition (SVD) method is used in the deconvolution process to calculate the impulse response of the target, as described in detail in [27]. Here, we focus on the wavefront and resonance of the impulse response of the target, instead of the measurement signal.

Figure 6. (a) Experiment scenario for water-filled crack in free space. (b) The SPWV distribution of the impulse responses from water-filled plastic bottle. (c) From water-filled PVC pipe.

In this scenario, the background (shown in Figure 6(a)) is tested as background 1, where the direct wave appears at 3.24 ns in the time-frequency image, which is not shown in the paper. In our study, the plastic bottle is needed to fill with water. In order to extract the impulse response of water, an empty plastic bottle is sited on the foam to measure the response and we call it as background 2. The mean square error (MSE) between background 2 and background 1 is calculated as 4.7 × 10⁻³. The order of magnitude of the error is acceptable in our measurement system. As expected, the effect of the plastic bottle on the water can be ignored. Then water filled in the plastic bottle and in PVC pipe is tested, respectively. Note that the bottle or pipe is laid with its axis perpendicular to the direction of the electric field.

In this subsection, the SPWV distribution is performed with time smoothing window of 51, frequency smoothing window of 127 and threshold of 0.05% for each time-domain response.

Figure 6(b) shows the SPWV distribution of water filled by the plastic bottle, where the impulse response of water is calculated by SVD in terms of background 2 (including the empty plastic bottle). The SPWV distribution for water-filled PVC pipe is shown in Figure 6(c). It is found that the SPWV distributions for water filled by both plastic bottle and PVC pipe are similar. As demonstrated in the above subsection, the m-time internal reflection and (m+1)-time reflection can be distinguished. According to the time interval between a wavefront at lower resonant frequencies and its next wavefront at higher resonant frequencies, the width of the water-filled cylindrical plastic bottle and PVC pipe can be estimated and the results are listed in Table 4, which are close to the true values. It should be noted that there exist interferences of energy at two ends of the timeframe. This is caused by the truncation of the impulse response and the influence of the truncation is not considered in this paper.

Table 4. Estimated width for cracks in free space.

d (cm)	m	1	2	3
3.0	△tm (ns)	0.78	1.17	0.82
	de (cm)	2.59	3.91	2.73
3.5	△tm (ns)	1.03	0.97	1.07
	de (cm)	3.42	3.23	3.57

4.2. Experiment in subsurface

In order to simulate the subsurface water-filled crack, a water-filled plastic bottle with width d = 5.7 cm is located horizontally in the dry sand with relative permittivity 4.25 at depth h = 14 cm. The GSSI system is the same with that shown in free space at height H₀ = 25 cm beyond the sand surface. The experiment scenario is set as shown in Figure 7(a). Furthermore, the plastic bottle is replaced by a water-filled PVC pipe with width d = 3.5 cm. As analyzed in Section 2, the responses of water are deconvoluted from the system responses and the subsequent wavefront analysis is based on the deconvoluted impulse responses of water.

Figure 7. (a) Experiment scenario for subsurface water-filled crack. (b) The SPWV distribution of the impulse responses from water-filled PVC pipe. (c) From water-filled plastic bottle with d = 5.7 cm.

It is found in Figure 7(b) and (c) that the SPWV distribution for subsurface water-filled crack is similar to that in free space. In general, there exists less high frequency in Figure 7 than that in Figure 6. This is because the ground can be considered as a filter to smooth some high-frequency oscillation. However, comparing the SPWV distribution of water-filled PVC and water-filled plastic bottle, more low frequencies appear for the latter. It can be explained by the fact that the EM wave travels longer distance in water due to the larger width of the cylindrical bottle. In this case, the time interval between the m-time internal reflection and (m+1)-time reflection increases, which makes the resonant period increases as well and decreases the resonant frequency. As the same method proposed above, the widths are estimated and the results are listed in Table 5. The results show good agreement with the true values.

Table 5. Estimated width for subsurface cracks.

d (cm)	m	1	2	3
3.5	△tm (ns)	0.93	1.21	0.92
	de (cm)	3.10	4.03	3.07
5.7	△tm (ns)	1.84	1.66	1.75
	de (cm)	6.13	5.53	5.83

5. Conclusion

In order to estimate a subsurface water-filled crack, a method based on internal multiple reflections is proposed. The crack is modeled by a dielectric cylinder and the electromagnetic scattering model is analyzed using the extend ray theory. The energy transport and path length of internal multiple reflections are analyzed. The energy transport shows that odd m-time multiple rays transport along the surface of the cylinder to the receiver, which form the creeping waves, and even m-time multiple rays refract toward the receiver, focusing on the radial direction. The conclusion that the path lengths inside the cylinder are multiples of a diameter distance is utilized to estimate the width of the cylinder. Usually, the odd m-time reflection and its next even m-time reflection are difficult to be distinguished due to the small time interval. With this respect, wavefronts and resonances of internal multiple reflections are also analyzed, based on the fact that the odd m-time reflection and its next even m-time reflection show different resonances in the time-frequency domain. Hence, the time interval can be extracted by multiple reflections distributed at different resonant frequencies. In order to verify the proposed method, experiments both in free space and subsurface scene are given to estimate water-filled cylindrical cracks and the results show good agreement with the true values.

Disclosure statement

No potential conflict of interest was reported by the authors.

Funding

This work was supported by the National Natural Science Foundation of China [grant number 61371186]，[grant number 61162007]; Natural Science Foundation of Guangxi Province, China [grant number 2013GXNSFFA019004].