Sunday, February 25, 2024

Fast Array Ground Penetrating Radar Localization by CNN-Based Optimization Method | IEEE Journals & Magazine | IEEE Xplore

A. Near-Field MUSIC

Fig. 1. - Configuration of a DOA problem with two variables.
Fig. 1. Configuration of a DOA problem with two variables.

Fast Array Ground Penetrating Radar Localization by CNN-Based Optimization Method | IEEE Journals & Magazine | IEEE Xplore


Publisher: IEEE
Abstract:
This article presents an optimization-based approach to overcome redundancy arising from the multivariables enumeration process in multiple signal classification (MUSIC). By incorporating Broyden–Fletcher–Goldfarb–Shanno (BFGS) optimization, the computational speed of the MUSIC algorithm is significantly improved while maintaining mathematical accuracy. The optimization techniques require reasonable initial values to start the iteration, while for single target imaging purposes, the initial values can be acquired by the boundary between the near field and the far field. To generate suitable initial values for the optimization, we employ a modified convolutional neural network (CNN) to approximate the boundaries between the near and far fields, which vary with array system properties. Besides, the proposed method introduces a method for the Hessian matrix and gradient initialization for the BFGS method. Using simulation results as training samples, the modified CNN successfully establishes boundary approximations. Simulation and experimentation confirm the feasibility of our proposed method, showing its advantages in both accuracy and computation speed compared to existing techniques.
 
Topic: Ground Penetrating Radar and Artificial Intelligence Technology: Innovation and Applications
Page(s): 4663 - 4673
Date of Publication: 24 January 2024
ISSN Information:
Publisher: IEEE

SECTION I. Introduction

The significance of Ground-Penetrating Radar (GPR) has grown markedly and become a valuable nondestructive and high-precision method. Nowadays, GPR is one of the most important tools for monitoring subsurface conditions and is widely used for various purposes [1], [2], [3], [4], [5], [6], [7], [8]. Recent notable work by Zhou et al. [2] designed a multisensor data acquisition (MDA) platform and proposed an associated scale probability-based pipeline mapping (SP-PM) algorithm for reducing the acquisition point requirement. In addition, Zhou et al. [3] proposed a Gaussian-process regression-based method for underground cable pursuit that shows better stability and accuracy compared to conventional methods.

Moreover, as the demand for objective recognition and localization in complicated subsurface conditions continues to surge, neural networks earn popularity for effectively processing the data. Recent work by Li et al. [4] introduced a new deep learning approach YOLOv4-hyperbola for high-precision tree root localization using simulated and real data as the training sets. In addition, Sun et al. [5] proposed a clutter removal network (CR-net). By training from real large-scale hybrid data, the method can remove the residual clutter and tackle the target echo deformation caused by the previous work.

For higher resolution imaging purposes, the application of MUSIC algorithms in recent array GPR signal processing is maintaining the trends. To improve the feasibility of the MUSIC algorithm, many improvements are made in two aspects—one is to improve the accuracy, and another is to reduce the computational complexity [9], [10]. There has been a growing trend in the application of numerical optimization techniques in direction of arrival (DOA) estimation or time domain-DOA (TDOA). These techniques are increasingly being used to enhance the accuracy and efficiency of DOA estimation methods [9], [10], [11], [12], [13], [14].

For instance, second-order source (SOS) localization [9], rank reduced (RARE) [10], and reduced dimension (RD) MUSIC [11] manipulate matrices operation to simplify the search procedure. Moreover, some works applied optimization techniques, for example, Qian [15] introduced Newton method to a manifold reconstruction-based ESPRIT algorithm for computational efficiency. Similarly, Cui et al. [16] proposed a DOA estimation algorithm based on PSO-Gauss-Newton. This hybrid approach leverages particle swarm optimization (PSO) to generate initial values.

In this article, optimization techniques and neural networks are combined to improve time efficiency. The proposed method has three parts. The first one is to use a numerical optimization Broyden–Fletcher–Goldfarb–Shanno (BFGS) quasi-Newton technique for removing the redundancy of the MUSIC algorithm.

The second part is to evaluate the initial guess. It is generated by the boundary area between the near field and the far field, which depends on the properties of the radar system itself. However, summarizing an explicit formulation of a function from a large number of array properties is mathematically challenging, instead, we use a modified convolutional neural network (CNN) technique to generate an approximation without loss of generality. This adapted CNN, distinct from conventional image-processing CNNs, features fewer convolution layers and a specialized structure for generating precise initial guesses.

The third part contains the modification of the BFGS search. We initialize the Hessian matrix and first-order derivative estimation directly by leveraging the outcomes of multiple 1-D MUSIC searches. Moreover, BFGS may converge slowly near the optimal solution, which may result in potential on-grid/off-grid issues in the MUSIC algorithm, as discussed in [17]. To address this, we introduce a local 3×3 pixel search in the final optimization steps.

We conducted numerical simulations and a laboratory experiment to assess the viability of our proposed method and to benchmark it against the conventional fast MUSIC algorithm. The results unequivocally demonstrate that our method achieves nearly identical accuracy while significantly reducing computational time.

The rest of this article is organized as follows. Section II-A provides an introduction to the 2-D MUSIC technique designed for near-field applications, then delves into the mathematics of surfaces generated through the multiplication of the covariance matrix and the steering vector. Section II-B presents our modified CNN approach for generating initial guess. Section II-C introduces the modifications of the BFGS method to accelerate the algorithm. Section III-A is the ablation test. Section III-B encompasses a comparison between the boundary of near and far-field calculated by the conventional method and the proposed method. Sections III-C and III-D feature numerical simulation and experimentation based on GPR. Section III-E is the complexity comparison between classic and recent fast MUSIC algorithms. Sections IV comprises the discussion. Finally, Section V concludes this article.

SECTION II. Methodology

A. Near-Field MUSIC

According to the configuration of Fig. 1, a uniform linear array receiving signal can be donated as

S(t)=x(t)A(θ,r)+N(1)
View SourceRight-click on figure for MathML and additional features.where x(t)is the received signal of the reference element in the time domain, steering matrix A(θ,r) contains two unknowns, range r and angle θ, S(t) is the received signal matrix for all channels, andN is a noise matrix.

Fig. 1. - Configuration of a DOA problem with two variables.
Fig. 1. Configuration of a DOA problem with two variables.

After certain matrix operations, the spatial spectrum can be plotted by using the following formula:

F(r,θ)=A(θ,r)EHNENAH(θ,r)(2)
View SourceRight-click on figure for MathML and additional features.where F(θ,r) is the spatial spectrum function, ENis the noise subspace of X(t), and * donates the conjugate transpose [1], [12].

According to (2), where the target is located has the smallest amplitude, which is nonnegative. More explicitly, it can be written as

F(r,θ)=[1ejωτM]a1,1aM,1a1,MaM,M1ejωτM.(3)
View SourceRight-click on figure for MathML and additional features.

Or in a polynomial form

F(θ,r)=a1,1+a1,2conj(ejωτ2)++aM,MejωτMconj(ejωτM)(4)
View SourceRight-click on figure for MathML and additional features.where M is the number of antenna elements and τican be expressed as
τi=r2+(i1)2d22r(i1)dcosθrc.(5)
View SourceRight-click on figure for MathML and additional features.

Note that the steering matrix has no more unknowns except r and θ after the array system is set, the exponential terms in (5) can be precalculated by enumerating all possible combinations of r and θ. Based on (4), the imaginary part of conj(ejωτM) will be decimated by the summation and so for ai.j, the result can be considered as 2-D surfaces as shown in Fig. 2. Consequently, optimization techniques can be applied to search for the minimum value.

Fig. 2. - Example of a 2-D MUSIC result, the ridge appears in the marked area.
Fig. 2. Example of a 2-D MUSIC result, the ridge appears in the marked area.

A ridge appears to be gradually vertical to the angle axis, and the signal source is increasingly fitting the perfect far-field condition. It indicates a reasonable initial guess on the angle parameter. Mathematically, it can be explained by using the numerator of τi

limrr2+(i1)2d22r(i1)dcosθrr2r=0.(6)
View SourceRight-click on figure for MathML and additional features.

B. CNN-Based Initial Guesses Generation

In this article, the proposed method does not generate an entire spatial spectrum. First of all, the initial guess on the range parameter will be determined. Then a 1-D MUSIC will be applied to generate a two-variable initial guess to start the optimization process. An example of the mentioned process is shown in Fig. 3.

Fig. 3. - Process of initialization.
Fig. 3. Process of initialization.

In the case of offering a reasonable initial guess on the range, a modified CNN is employed to avoid the complicated mathematical investigation. Compared to the conventional CNN for image processing, the proposed method has simplified convolution and max pooling layers and a more suitable structure for generating initial guesses.

Before the introduction of the proposed CNN, an essential explanation of the input data for training or prediction should be made. According to (4), a surface is generated by the summation of exponential terms, and each of them is multiplied by ai.j,in this case, we name these exponential terms ai.j terms and, associatively, their surface, ai.j surfaces. Theseai.j terms are fundamental elements for determining the range parameter and generating surfaces with a given domain of angle and range, and an example of a 4-elements uniform linear array is shown in Fig. 4.

Fig. 4. - Surfaces generated by the ${a}_{i,j}$ terms. All of them will only contain angle responses beyond a certain range.
Fig. 4. Surfaces generated by the ai,j terms. All of them will only contain angle responses beyond a certain range.

The boundary range between the near-field and the far-field of (4) is generated by interaction among these surfaces and as shown in Fig. 5. Thus, the input only contains the system properties and labeled ranges. The labels are marked manually for the training. The rule that labeling is followed is marking the transition area between the near field and the far field, as shown in Fig. 3

Fig. 5. - Schematic demonstration of the proposed network and the coefficient from the covariance matrix.
Fig. 5. Schematic demonstration of the proposed network and the coefficient from the covariance matrix.

The surfaces seen in this case are structurally different from the datasets used for conventional CNN for imaging processing purposes. Under this circumstance, the focus is on the range in which the angle response is weak. Therefore, the modification is made for the purpose, after surface generation, a summation kernel with Hanning window scanning along the range direction is added as a modification. The process can be written as

Lm=n=1N{0.50.5cos(2πnN))Am,ni,j}(7)
View SourceRight-click on figure for MathML and additional features.where Am,ni,jpresents element on the mth row and the nth column from the matrix generated by the ai,j surface. The result of the operation yields an array that will be convoluted by a modified kernel. The modified convolution kernel is presented to record the condition of the response from the angle and range. The process is shown in Fig. 6.

Fig. 6. - Convolutional kernel design for ${a}_{i,j}$surfaces by using an example dataset.
Fig. 6.

Convolutional kernel design for ai,jsurfaces by using an example dataset.

Compared to the 2-D convolutional method for imaging processing, the modified method can keep the information on the angle response with very little distortion and record the boundary between the near field and the far field.

The resulted arrays can be connected to a fully connecting layer (FCL) with a large number of inputs. To reduce the amount of the input before the FCL, we design a max pooling process to store the kth largest value of the convoluted array, which follows:

Θk+1=maxWs.t. W<Θk+1(8)
View SourceRight-click on figure for MathML and additional features.where Θk+1 is the kth largest value in the input arrayW and Θis the array after max pooling. Besides, most of the elements in Ware 0 as shown in Fig. 7, which contributes very little to the network training. Thus, the max pooling is essential.

Fig. 7. - Vector ${\rm W} $ in matrix representation.
Fig. 7. Vector W in matrix representation.

After the convolution and max pooling are done for all surfaces, an FCL will be linked for the final estimation of the range parameter.

The FCL structure is shown in Fig. 8, and the network outputs can be solved in an overdetermined system with the given the input datasets. The solution divided by the associate elements in covariance will be used for FCL training. The procedure can be considered as solving the following:

2ΓkRe{Z}=rktraining(9)
View SourceRight-click on figure for MathML and additional features.where Z contains all upper-triangle matrix elements in the covariance matrix of the received signal, rktraining is the kth label of the training data l, and 2Γk is given by
f(Θi,j)=Γkk(10)
View SourceRight-click on figure for MathML and additional features.where Θi,jrepresents the maxing pooling result from the ai.jsurface. It is the fundamental idea of the proposed CNN, to approximate fit the unknown function f(Θk).

Fig. 8. - FCL structure of the proposed method, two biases are added for better fitting precision.
Fig. 8. FCL structure of the proposed method, two biases are added for better fitting precision.

Above all, the overall structure of the modified CNN is schematically shown in Fig. 9. In this work, 66 simulations are carried out on a 4-element ULA array system with different element divisions and center frequencies. The input data for training are listed in Table I.

Fig. 9. - Structure of the proposed modified CNN.
Fig. 9. Structure of the proposed modified CNN.
TABLE I Input Parameters of the Training Data
Table I- Input Parameters of the Training Data

C. Modified BFGS Method for 2-D MUSIC

Knowing the initialization can be generated from the range parameter and a 1-D MUSIC search, we can start the optimization by using various techniques. In this work, the BGFS quasi-Newton method is applied to the optimization to avoid calculating the full Hessian matrix [13].

The entire process has two advantages: one is the computational speed. Both optimization-based MUSIC and the modified CNN vastly reduce the calculation, which offers the opportunity to use single target time domain MUSIC for localizing coherent targets. Another is that the initial values are no longer required as they are generated by the proposed modified CNN. Though the BFGS algorithm uses an approximation of the Hessian matrix, the initial estimation of the Hessian matrix can be directly calculated by using the results of multiple 1-D MUSIC searches. It means the derivatives can be approximately calculated by using the adjacent pixels. The procedure is shown in Fig. 10. Thus, knowing that

B(0)p(0)=F(x(0))(11)
View SourceRight-click on figure for MathML and additional features.and approximated Hessian matrix B0 should be numerically close to the Hessian matrix H. In 2-D MUSIC case, the Hhas two variables, which can be expressed as
B(0)=H(θ(0),r(0))=2F(θ(0),r(0))2F(θ(0),r(0))θ22F(θ(0),r(0))θr2F(θ(0),r(0))θr2F(θ(0),r(0))r2(12)
View SourceRight-click on figure for MathML and additional features.where x(0)1 is the initial value acquired by other 1-D MUSIC search and p(0)=x(0)1x(0). And F(x(0)) is the first-order derivative approximated from x(0)1and x(0)as follows:

Fig. 10. - Calculation of the Hessian matrix by using multiple local 1-D MUSIC search.
Fig. 10.

Calculation of the Hessian matrix by using multiple local 1-D MUSIC search.

F(x(0))==F(θ(0),r(0))=F(θ(0),r(0))θF(θ(0),r(0))rF(θ(0),r(0))F(θ(0)1,r(0)1)θ(0)θ(0)1F(θ(0),r(0))F(θ(0)1,r(0)1)r(0)r(0)1(13)
View SourceRight-click on figure for MathML and additional features.we can solve B(0) if any of the diagonal elements in (11) is known. For example, to solve the second-order derivate of F(θ,r)versus θ, we can use an approximation of (11)
2F(θ(0),r(0))θ2=F(θ(0),r(0))θx02F(θ(0),r(0))θx01θ(0)2θ(0)1(14)
View SourceRight-click on figure for MathML and additional features.where F(θ(0),r(0))θ|x02 stands for the approximated first-order derivate resulting from x(0)2and x(0). Besides, for simplicity, the 1-D searches can be done in a small scale as shown in Fig. 9.

The only to-be-estimated parameter left is α(0). We can set α(0) to the division of the range. However, in this article, we do not go into detail work on α(0) in the iterations, because the 1-D search α(k+1)=argminF(x(k)+a(k)p(k)) yields accurate results. Nonetheless, α(k+1)=argminF(x(k)+a(k)p(k)) search may cause redundancy around the solution as shown in Fig. 11. The BFGS method converges to the optimal solution at a quadratic rate, which causes small steps that may raise the on-grid off-grid issue in the MUSIC algorithm [14]. Thus, in the last few steps of the optimization, a threshold will be added to examine the distance between x(k) andx(k1). Both distances along the angle parameter and range parameter should be larger than the size of a pixel, in this case r(k)r(k1)>0.05 and θ(k)θ(k1)>1o. When neither of the adjustments is fitted, the BGFS will be stopped and a 2-D MUSIC will be used for a 3×3 pixel search to determine the optimal solution. The 2-D MUSIC follows the equation given in (2). Above all, the proposed BFGS method is shown in a pseudoalgorithm way in Algorithm 1.

Algorithm 1: Modified BFGS Method for 2-D MUSIC.

Input:

Covariance matrix RXX, L

1:Generate initial range parameter r(0) by modified CNN
2:Find x(0)=(r(0),θ(0)) by 1-D MUSIC search
3:Calculate B(0) and P(0) by using adjacent pixels
4:Set convergency condition: x(k-1) and x(k) in same pixels after L steps
5:For k=0,1,2… until converges
6:Solve B(k)p(k)=F(x(k))
7:Update sk=αp(k)
8:Update: x(k+1)=skx(k)
9:Update: yk=F(x(k+1))F(x(k))
10:Update: Bk+1=Bk+ykyTkyTkykBkykyTkBTksTkBksk
11:End for
12:Start 3×3 MUSIC search for the target
13:End
Fig. 11. - Proposed 3×3 pixels search demonstration.
Fig. 11.

Proposed 3×3 pixels search demonstration.

SECTION III. Simulation and Experimentations

We embarked on an analysis of the accuracy of the modified CNN. Our evaluation begins with an ablation analysis of biases and the max pooling process. Then we delved into an assessment of far-field and near-field assumptions. We contrasted multiple widely employed criteria for defining the far field, such as the intuitive “two times the wavelength” rule and the more nuanced Fresnel distance. This comparative examination helps us determine which assumption aligns better with our specific research context. To address the complexities of real-world applications, we introduced a two-target GPR simulation and a verification of MIMO GPR calibration. These enable us to assess their robustness and suitability for demanding imaging conditions.

In the end, a complexity analysis of the target searching stage of MUSIC and its variants was conducted to prove the efficiency of the proposed work.

A. Ablation Analysis and Forwarding Prediction

This section aims to show the necessity of using the proposed components in the network for efficiency and precision. Because the network has a simple structure, the analysis is only carried out on biases and max pooling. Totally, the following five structures of the network are in the analysis:

  1. proposed : Max pooling+ FCL(6×6) +2 biases;

  2. structure 1: FCL(20×20)+2 biases;

  3. structure 2: Max pooling+ FCL(6×1) +1 biases;

  4. structure 3: Max pooling+ FCL(6×6×6) +3 biases;

  5. structure 4: Max pooling+ FCL(6×6) +1 biases.

The results shown in Figs. 12 and 13 reveal that the loss of any component lowers the network accuracy and the extra components merely improve the performance. It indicates the proposed method effectively reduces input data dimensionality without sacrificing performance. Moreover, when dealing with a large training dataset or fine-grained range division, this technique proves highly advantageous.

Fig. 12. - Comparison between the larger and smaller size FCL.
Fig. 12.

Comparison between the larger and smaller size FCL.

Fig. 13. - Absolute error of the structures.
Fig. 13.

Absolute error of the structures.

In the aspect of prediction capability, we use a new set of labeled data for the analysis. The input data parameters are shown in Table II. The results shown in Fig. 14 prove the feasibility of the proposed network for forwarding prediction. In this case, due to the potential errors introduced during manual labeling, we are unable to confidently assert that the proposed max pooling technique exhibits superior performance. This situation prompts us to consider that neither of the networks can perfectly align with the labeled data. Nevertheless, the margin of error remains under 0.1m, except for the 21st input. Thus, the overall performance is suitable for the boundary estimation.

TABLE II Input Parameters of the Training Data for Prediction
Table II- Input Parameters of the Training Data for Prediction
Fig. 14. - Prediction on 1.2 GHz system.
Fig. 14.

Prediction on 1.2 GHz system.

Fig. 15. - Comparison among the commonly used near-field and far-field assumptions.
Fig. 15.

Comparison among the commonly used near-field and far-field assumptions.

B. Comparison Between Conventional Near-Field Assumption

This section addressed the significance of adopting the proposed initialization method for optimization due to differences in the boundary or transition zone of the MUSIC algorithm compared to conventional approaches.

The boundary between the near field and the far field is often approximated based on the wavelength of the electromagnetic waves involved. Typically, the near field is considered to extend up to a distance of approximately two wavelengths from the source, while the far field starts beyond this distance. It is worth noting that the largest distance in the first Fresnel zone is practically determined by

fFresnel=λD2(15)
View SourceRight-click on figure for MathML and additional features.where D is the distance between two antenna elements. On the other hand, the Fraunhofer distance demarcates the far field, and its minimum distance is characterized by
fFresnel=2ζ2λ(16)
View SourceRight-click on figure for MathML and additional features.where ζis the largest electric length of the antenna elements. For simplification purposes, half-wavelength and quarter-wavelength antennas are employed, and their electric lengths are assumed to be equivalent to their geometrical sizes. The comparison result is shown in Fig. 15. The above results indicate the completely different behavior that MUSIC has versus the frequency and element division.

C. GPR Simulation of Two Buried Objects

The part was conducted to simulate a simplified circumstance of a stationary GPR array for lunar drilling navigation. A real case is the Chang-E 5 launcher, which used array GPR for imaging the Moon's subsurface condition as a part of the predrilling process [18]. In this operation, the batteries of the launcher, drilling system, and array system are limited; thus, fast imaging is heavily required.

The simulation has two buried idealistic point reflectors located in (0.25 m, 0.43 m) and (0.5 m, 0.87 m), associatively, (0.5 m, 60°) and (1 m, 60°). The simulation assumes the array system is idealistic, omnidirectional, and closely located on the ground. The GPR stationary ULA has a center frequency of 1 GHz and contains four elements with a division of 0.12 m. The subsurface contains one homogenous layer and the relative electric permittivity is 3, following the Moon regolith data [18]. Since GPR are commonly wideband systems that violate the basic assumption of the MUSIC algorithm, we recommend a low-complexity algorithm given by a mature GPR layer detection work [1]. By following the method, two targets can be differed and localized by the MUSIC algorithm. Different from the spatial smooth process in the given work, we directly use TDOA approaches with the MUSIC estimator in (2) as the phase shift is addressed. To avoid the potential coherency caused by targets, the MUSIC estimator conducts one target evaluation that can be done by enumerating the possible combinations for the collected delays.

After evaluating the covariance matrices, the proposed network generates 1.4 m as an initial guess for the range parameter. The results are compared with conventional back-projection imaging and the original quasi-Newton method demonstrates the advantages in Figs. 16 and 17.

Fig. 16. - Simulation and the results of the two-target localization. (a) Simulation condition. (b) Result of the back projection imaging. (c) First target search process. The red dots are generated by the proposed method and the background is the spatial spectrum generated by the conventional 2-D scan. (d) Second target search process.
Fig. 16.

Simulation and the results of the two-target localization. (a) Simulation condition. (b) Result of the back projection imaging. (c) First target search process. The red dots are generated by the proposed method and the background is the spatial spectrum generated by the conventional 2-D scan. (d) Second target search process.

Fig. 17. - Traces of the final search steps of the proposed BFGS (blue dots) and original BFGS (red dots).
Fig. 17.

Traces of the final search steps of the proposed BFGS (blue dots) and original BFGS (red dots).

D. Laboratory Experiment

This part examined the proposed method for a MIMO GPR case. MIMO GPR system Yakumo is an SFCW system with eight antenna pairs, working from 65 MHz to 1.5 GHz [19], [20]. It has been used and achieved many successful experiments since 2007. However, due to the aging of the cables, many of them are required to be replaced. Thus, the replaced cables may suffer from the uncertain delay, which should be compensated. We conducted the metal sphere imaging to evaluate if the cable calibration was suitable for use. A focused target in the correct position means the correctness of the calibration. However, back projection (BP) imaging has a weak capability of focusing, therefore, MUSIC methods were used. However, imaging from all Tx-Rx combinations by using original MUSIC is time-consuming. Consequently, we introduced the proposed method for time efficiency. The experiment setup is shown in Fig. 18. The antenna configuration is shown in Fig. 19. The results are shown in Fig. 20.

Fig. 18. - Experiment setup.
Fig. 18.

Experiment setup.

Fig. 19. - Yakumo antenna configuration.
Fig. 19.

Yakumo antenna configuration.

Fig. 20. - Metal sphere imaging by (a) BP, (b) MUSIC, (c) RD-MUSIC, and (d) proposed method.
Fig. 20.

Metal sphere imaging by (a) BP, (b) MUSIC, (c) RD-MUSIC, and (d) proposed method.

In this case, the Yakumo is intentionally put upside down, and the mental sphere is 0.8 m above the geometry midpoint.

For visualization purposes, we only used the combination of Tx1 and all Rx to conduct the BP imaging, RD-MUISC, and the proposed fast MUISC in this section. We assume the imaging center at Rx5 and the sphere is located at (0.1, 0.8, 0) in meter. The results in figure show the accuracy of the proposed method, which matches the result of MUSIC and its variants with an error of 0.02 m. Note that neither of the given algorithm gives a result that exactly matches (0.1, 0.8, 0), as the distance in the Z direction exist, which affects localization.

E. Computational Complexity Analysis

This part compared classical and recent fast MUSIC algorithm with the proposed method. For one target localization, the complexity of the subspace method has a O(M3) base in general, where M is the snapshot of the signal. The complexity of the final stage from target localization generated by (2) are all O(M2J), where J is the number of array elements.

Note that the computational complexity of the BFGS method is H(7K2+7K) [21], where K is the number of variables and H is the number of iterations. Modified BFGS contains a 1-D search, which causes O(M)complexity, and the initial operations for (11), (12), and (13) and the final stage 3×3 search cause the iteration number to be (HLU),where L is the local steps of 2-D search for gradient initialization and U is the total steps of 2-D search for final localization in Fig. 11. In this case, the local 1-D search is ignored because it contributes very little to the complexity. In addition, RARE, RD-MUSIC, and MUSIC have different computational complexity; however, in this work, we only evaluated the higher order term of these previous works, because the proposed method has an obvious order difference in the localization stage of (2), where the lower order complexity terms are negligible.

According to Table III, the complexity order is much reduced since the M-based multiplication is replaced by K-based ones in the BFGS method, where K = 2 for angle and range estimation.

TABLE III Comparison of Computational Complexity
Table III- Comparison of Computational Complexity

SECTION IV.

Discussion

The proposed method assumes localization involves one target since multiobject optimization in DOA estimation is problematic because of the potential coherency. Another issue to be addressed is the range initial value generated from the modified CNN. The input data of the modified CNN are manually marked, which may cause (7) unsolvable or has large error. As a result, the training stage will fail to suit accuracy. Thus, future study is needed to improve the following aspects.

  1. The shape of the object should be taken into consideration, where some improved MUSIC has already been done.

  2. The underground condition needs to be evaluated in a more realistic way with multilayer structures and a more complicated clutter source.

  3. Multiobject optimization should be developed for the DOA estimation cases. The technique should be robust in terms of coherency and advantageous in terms of the capability of global searching.

Though potential improvements are required for the proposed method to be fully employed, the novelty of this work is essential for future use of the optimization techniques in the subspace method. A self-generated initial value based on the system properties is much more convenient compared to evaluating a possible target location without any preknowledge. The modification employed for BFGS initial derivates and Hessian matrices is efficient. In the given simulations, the localization does not involve many nested loops. However, in practice, multiple variables are unknowns, for example, relative permittivity is usually undetermined because of the site condition, which requires an enumeration process to evaluate relative permittivity for a focused image. Moreover, the MUSIC for a wideband system requires calculations at every frequency point, where the proposed method can vastly reduce the calculation. Under such circumstances, the step reduction in the proposed method is advantageous.

SECTION V. Conclusion

This article introduces an optimization-driven solution to address redundancy issues that arise during the multivariable enumeration process in the MUSIC method. By integrating quasi-Newton optimization techniques, we achieve a significant enhancement in the computational efficiency of the MUSIC algorithm. Nevertheless, the optimization process necessitates appropriate initial values for successful iteration. In the context of single-target imaging, these initial values can be obtained by delineating the boundary between the near field and the far field. To establish suitable starting points for the optimization, we employ a modified CNN that approximates the boundary between the near and far fields. This boundary varies based on the properties of the array system. Our approach leverages simulated results as training data for the modified CNN, which effectively approximate these boundaries. Through both simulation and experimentation, we validate the effectiveness of our proposed technique. Our method demonstrates notable advantages in both accuracy and computational speed compared to existing approaches.




References

1.
S. Zhao and I. L. Al-Qadi, "Super-resolution of 3-D GPR signals to estimate thin asphalt overlay thickness using the XCMP method", IEEE Trans. Geosci. Remote Sens., vol. 57, no. 2, pp. 893-901, Feb. 2019.
2.
X. Zhou et al., "Underground pipeline mapping from multipositional data: Data acquisition platform and pipeline mapping model", IEEE Trans. Geosci. Remote Sens., vol. 61, 2023.
3.
X. Zhou, Q. Chen, S. Lyu and H. Chen, "Mapping the buried cable by ground penetrating radar and Gaussian-process regression", IEEE Trans. Geosci. Remote Sens., vol. 60, 2022.
4.
S. Li, X. Cui, L. Guo, L. Zhang, X. Chen and X. Cao, "Enhanced automatic root recognition and localization in GPR images through a YOLOv4-based deep learning approach", IEEE Trans. Geosci. Remote Sens., vol. 60, 2022.
5.
H.-H. Sun, W. Cheng and Z. Fan, "Learning to remove clutter in real-world GPR images using hybrid data", IEEE Trans. Geosci. Remote Sens., vol. 60, 2022.
6.
W. Luo, Y. H. Lee, L. F. Ow, M. L. M. Yusof and A. C. Yucel, "Accurate tree roots positioning and sizing over undulated ground surfaces by common offset GPR measurements", IEEE Trans Instrum. Meas., vol. 71, 2022.
7.
H. Liu et al., "Evaluation of the antenna parameters for inspection of hidden defects behind a reinforced shield tunnel using GPR", Tunnelling Underground Space Technol., vol. 140, Jun. 2023.
8.
H. Liu, F. Ding, J. Li, X. Meng, C. Liu and H. Fang, "Improved detection of buried elongated targets by dual-polarization GPR", IEEE Geosci. Remote. Sens., vol. 20, 2023.
9.
K. Abed-Meraim, Y. Hua and A. Belouchrani, "Second-order near- field source localization: Algorithm and performance analysis", Proc. 30th Asilomar Conf. Signals Syst. Comput., pp. 723-727, 1996.
10.
J. He, M. N. S. Swamy and M. O. Ahmad, "Efficient application of music algorithm under the coexistence of far-field and near-field sources", IEEE Trans. Signal Process., vol. 60, no. 4, pp. 2066-2070, Apr. 2012.
11.
X. Zhang, W. Chen, W. Zheng, Z. Xia and Y. Wang, "Localization of near-field sources: A reduced-dimension MUSIC algorithm", IEEE Commun. Lett., vol. 22, no. 7, pp. 1422-1425, Jul. 2018.
12.
Y. Ma, Y. Zeng and S. Sun, "A deep learning based super resolution DoA estimator with single snapshot MIMO radar data", IEEE Trans. Veh. Technol., vol. 71, no. 4, pp. 4142-4155, Apr. 2022.
13.
N. Borijindargoon, B. P. Ng and S. Rahardja, "MUSIC-like algorithm for source localization in electrical impedance tomography", IEEE Trans. Ind. Electron., vol. 66, no. 6, pp. 4661-4671, Jun. 2019.
14.
Y. S. Yoon and M. G. Amin, "High-resolution through-the-wall radar imaging using beamspace MUSIC", IEEE Trans. Antennas Propag., vol. 56, no. 6, pp. 1763-1774, Jun. 2008.
15.
C. Qian, "A simple modification of ESPRIT", IEEE Signal Process. Lett., vol. 25, no. 8, pp. 1256-1260, Aug. 2018.
16.
X. Cui, R. Zhou, H. Chen, Y. Zhang, S. Li and J. Zhang, "A new DOA estimation algorithm based on PSO-Gauss-Newton", Proc. IEEE 9th Int. Conf. Inf. Commun. Netw., pp. 81-85, 2021.
17.
Z. Yang, L. Xie and C. Zhang, "Off-grid direction of arrival estimation using sparse Bayesian inference", IEEE Trans. Signal Process., vol. 61, no. 1, pp. 38-43, Jan. 2013.
18.
J. Feng, M. A. Siegler and M. N. White, "Shallow regolith structure and obstructions detected by lunar regolith penetrating radar at chang'E-5 drilling site", Remote Sens., vol. 14, no. 14, Jul. 2022.
19.
L. Yi, L. Zou, K. Takahashi and M. Sato, "High-resolution velocity analysis method using the ℓ-1 norm regularized least-squares method for pavement inspection", IEEE J. Sel. Topics Appl. Earth Observ. Remote Sens., vol. 11, no. 3, pp. 1005-1015, Mar. 2018.
20.
L. Zou, L. Yi and M. Sato, "On the use of lateral wave for the interlayer debonding detecting in an asphalt airport pavement using a multistatic GPR system", IEEE Trans. Geosci. Remote Sens., vol. 58, no. 6, pp. 4215-4224, Jun. 2020.
21.
A. Rodomanov and Y. Nesterov, "Greedy quasi-newton methods with explicit superlinear convergence", Soc. Ind. Appl. Math. J. Optim., vol. 31, pp. 785-811, 2021.

 

No comments:

Post a Comment

Breakthrough in Satellite Error Correction Improves Space Communications

Typical LEO Architecture and Segments Spectra of some LEO Link Losses Breakthrough in Satellite Error Correction Improves Space Communicatio...