Joint Estimation of Azimuth and Distance for Far-Field Multi Targets Based on Graph Signal Processing

: Target position estimation is one of the important research directions in array signal processing. In recent years, the research of target azimuth estimation based on graph signal processing (GSP) has sprung up, which provides new ideas for the Direction of Arrival (DoA) application. In this article, by extending GSP-based DOA to joint azimuth and distance estimation and constructing a fully connected graph signal model, a multi-target joint azimuth and distance estimation method based on GSP is proposed. Firstly, the fully connection graph model is established related to the phase information of a linear array. For the fully connection graph, the Fourier transform method is used to solve the estimated response function, and the one-dimensional estimation of azimuth and distance is completed, respectively. Finally, the azimuth and distance estimation information are combined, and the false points in the merging process are removed by using CLEAN algorithm to complete the two-dimensional estimation of targets. The simulation results show that the proposed method has a smaller mean square error than the Multiple Signal Classiﬁcation (MUSIC) algorithm in azimuth estimation under the condition of a low signal-to-noise ratio and more accurate response values than the MUSIC algorithm in distance estimation under any signal-to-noise ratio in multi-target estimation.


Introduction
Target position estimation is one of the important research directions in array signal processing, which had been a popular research object in the field of autonomous driving radar systems, military radars, and satellite navigation [1,2]. Recently, some new methods of Direction of Arrival (DoA) estimation based on graph signal processing have emerged. Graph signal processing uses a new data structure that studies the connection of things, which has shown excellent performance in many fields such as graph neural network and graph cuts [3]. Some related works can be found focusing on the use of graph signals to deal with DoA estimation problems in the radar array system [4][5][6], microphone and speakers [7,8], and the sonar array system [9]. Experiments show that the graph signal processing-based DoA methods have better performance than traditional algorithms such as Multiple Signal Classification (MUSIC) in a low signal-to-noise ratio environment [10,11].
However, these works use a kind of non-fully connected graph signal model. The nonfully connected graph signal model will produce a relatively large deviation in the estimation, because it does not use all the information between the array element and the target. What is more, these studies have not studied the distance estimation, since the non-fully connected graph signal method will produce a relatively large deviation in the estimation of the distance, which will lead to the distance being unable to be obtained [12,13].
To solve the problem of the non-fully connected graph signal model and the lack of distances estimation, a fully connected graph signal model is proposed in this article, based on which a multi-target azimuth and distance joint estimation method is proposed. By constructing fully connected graph structure data to introduce more element phase information into the adjacency matrix, the phase information between the array element and the target is fully utilized, which can achieve higher estimation accuracy. The expression of an adjacency matrix and the estimation algorithm of a fully connected graph signal are derived for azimuth estimation. At the same time, the fully connected graph signal algorithm is also used in distance estimation, which has not been done by existing GSP-based DoA estimation algorithms. What is more, combined with azimuth estimation and distance estimation, multi-target joint two-dimensional estimation can also be realized. However, false points appear in joint estimation, which is solved by the CLEAN algorithm [14] in this article. Finally, a series of experiments and Monte Carlo analysis are designed to verify the effectiveness of the proposed algorithm.
The rest of this article is arranged as follows: Section 2 reviews the radar array signal model and the subspace of signal and noise, then it introduces the graph signal modeling of the radar array model and derives two special estimation graph matrices from this graph matrix; Section 3 shows the simulation and Monte Carlo experiment results; Section 4 discusses and analyzes the estimation principle and simulation results of the GSP algorithm; Section 5 draws conclusions.

Array Signal Model
We assume that there are M radar array elements arranged equidistantly in space to form a linear uniform array. The transmitting array and receiving array are juxtaposed, with M − 1 equidistant transceiver antennas, and the array element spacing is d = λ/2, where λ is the wavelength. The radar array transmits a single frequency signal with a frequency interval of ∆f from the first array element to form a stepped frequency broadband signal, as shown in Figure 1. To solve the problem of the non-fully connected graph signal model and the lack of distances estimation, a fully connected graph signal model is proposed in this article, based on which a multi-target azimuth and distance joint estimation method is proposed. By constructing fully connected graph structure data to introduce more element phase information into the adjacency matrix, the phase information between the array element and the target is fully utilized, which can achieve higher estimation accuracy. The expression of an adjacency matrix and the estimation algorithm of a fully connected graph signal are derived for azimuth estimation. At the same time, the fully connected graph signal algorithm is also used in distance estimation, which has not been done by existing GSPbased DoA estimation algorithms. What is more, combined with azimuth estimation and distance estimation, multi-target joint two-dimensional estimation can also be realized. However, false points appear in joint estimation, which is solved by the CLEAN algorithm [14] in this article. Finally, a series of experiments and Monte Carlo analysis are designed to verify the effectiveness of the proposed algorithm.
The rest of this article is arranged as follows: Section 2 reviews the radar array signal model and the subspace of signal and noise, then it introduces the graph signal modeling of the radar array model and derives two special estimation graph matrices from this graph matrix; Section 3 shows the simulation and Monte Carlo experiment results; Section 4 discusses and analyzes the estimation principle and simulation results of the GSP algorithm; Section 5 draws conclusions.

Array Signal Model
We assume that there are radar array elements arranged equidistantly in space to form a linear uniform array. The transmitting array and receiving array are juxtaposed, with equidistant transceiver antennas, and the array element spacing is , where is the wavelength. The radar array transmits a single frequency signal with a frequency interval of from the first array element to form a stepped frequency broadband signal, as shown in Figure 1. For the estimation problem under far-field conditions, it can be considered that the received echo signal is a plane wave relative to the radar array, so we assume that the array receives the echo single-frequency continuous-time signal to estimate azimuth and distance information reference to [15] as below: For the estimation problem under far-field conditions, it can be considered that the received echo signal is a plane wave relative to the radar array, so we assume that the array receives the echo single-frequency continuous-time signal to estimate azimuth and distance information reference to [15] as below: phase delays for azimuth and distance estimation. The phase delay of the signal arriving at the m-th element is as follows: where θ i is the azimuth value of the i-th target of the multi-target, R i is the distance value of the i-th target in the multi-target, c is the speed of light propagating in vacuum, and d is the array element spacing. The phase delay of the echo signal includes the delay caused by the array structure itself and the delay caused by different frequency signals. These two phase delays are multiplied exponentially, so they can be simplified to the exponential addition of e, and the expression of phase delay is p reference to [15]: Simplify the exponential term in (3) to e −jω M τ θ i ,R i 1,M and write the array echo signal in vector form, and we can get the vector that uniquely describes the phase shift of the source, which is called the steering vector [15]: It can be seen from the formula that the steering vector is related to the target's azimuth θ, distance R, and the number of array elements M. We will discuss how to express this form in the next section. The actual received signal is a signal that has undergone Q = 128 times snapshots in the time domain, so it is a discrete uniformly sampled signal. It is known that combined with the steering vector in Equation (4), the final echo signal x(k) ∈ C M×M is the superposition form of M array element signal vectors: Because this article uses a MIMO radar array, so s m (k), m = 1, · · · , M represents the orthogonal array element signal of the m-th frequency; the specific form will be analyzed in the next section; and n(k) represents the zero-mean Gaussian white noise that is not coherent with time in the spatial domain. Its magnitude is ε n(k)n H (k) = σ 2 n Iδ(n), and it can be seen that the Gaussian white noise is only related to the time interval.

Discrete Analysis of Array Signal Model
Equation (5) describes the received signal form of multiple targets but requires a more detailed mathematical description form; on the basis of the array signal analysis, let us analyze the following parameters of the MIMO radar array in this article to derive the specific form of the received signal x(k) [ In this article, we consider the joint estimation of azimuth and distance. The corresponding situation is the single-transmit and single-receive situation of the radar array, which is in line with the assumptions of traditional radar systems and MIMO radar arrays with co-located antennas, and the number of transmitting antennas M equal to the number of receiving antennas N, so according to the above parameters and Equation (5), it can be known that the baseband signal expression of the receiving antenna array reference to [16] is: where X ∈ C M×Q represent x(k) in (5), which is the snapshot signal samples collected by M receiving antennas; Σ = diag {σ i } I i=1 is the diagonal matrix of the RCS scattering coefficient of the far-field targets; W ∈ C M×Q is the sampling noise term of n(k) in Equation (5), usually additive white gaussian noise (AWGN). The receiving steering vector B includes the two-way phase delay, so the expression is the same as the steering vector st (θ i ,R i ,M) in Equation (4): The transmitting steering vector C is the phase delay during transmission; at this time, the signals of all the array elements are in the same phase, and the phase delay is 0, so the expression of the receiving steering vector is: Therefore, it can be known that the result of summing for the steering vector in Equation (5) is equal to CΣB T ∈ C M×M in Equation (6). The advantage of MIMO radar array over traditional phased array radar is that it can transmit signals orthogonal to each other antennas; so, for the transmitted signal S ∈ C M×Q , its expression is: Among them, Orth Q is the Hadamard matrix form of M array elements, and its size is M * Q. F Q is the Q snapshots time sampling of transmitting array elements; f s is the sampling frequency.
Due to the orthogonality of S, we can get a diagonal matrix I M = (1/Q)SS H . Although there are non-zero elements on the off-diagonal line, the following experiments prove that it has no effect on the accuracy of the result. By right multiplication (1/Q)S H to Equation (6), the matched filter output result reference to [16] as below: where X rec = (1/Q)XS H ∈ C M×M represents the mathematical model of the reception received by the radar array, and Z = (1/Q)WS H is the vectorization of the white noise W. However, because the MIMO array radar form of single-transmit and single-receive is used in this article, which transmits orthogonal signals, so the non-diagonal elements in X rec are all 0. In order to facilitate the theoretical derivation, we will only take the diagonal values from X rec and form a vector: The final received signal X rec ∈ C M×1 will be used in the estimation algorithm in the subsequent Sections.

Echo Signal Covariance Matrix
For the second-order statistics of the echo signal x m (k) in (5), we can find the covariance matrix R xx reference to [17]: Covariance is one of the most commonly used second-order statistics in array signal processing. σ 2 i refers to the intensity coefficient of the i-th target echo signal. The significance of finding the covariance matrix for the echo signal is because the expression R xx contains information about multiple targets [9]. Multi-angle information and eigenvalue decomposition can be carried out. Perform EVD eigenvalue decomposition reference to [17] on the covariance matrix R xx to get: The eigenvalue matrix of the eigenvector Λ = {λ 1 , λ 2 , · · · , λ m } and the covariance matrix Q = [q 1 , q 2 , · · · , q m ] can be decomposed, but there is no direct description of the azimuth and distance information in the eigenvector Q, but the orthogonality of the eigenvector Q can provide information about the signal subspace and the noise subspace [17]. By traversing the parameter θ and arranging the pairs in descending order of eigenvalues R xx , the signal and noise subspace decomposition form can be obtained [18]: When the traversed azimuth is aligned with the target, there will be a maximum eigenvalue containing I targets in the eigenvalue block matrix Λ S ∈ C I×I , and its corresponding eigenvector Q S ∈ C M×I includes all the return guide vectors, so the noise subspace Q S ∈ C M×(M−I) is theoretically the same as all return waveguide vectors are orthogonal.

Graph-Based Joint Estimation Method and Algorithm Derivation
In this Section, we will introduce how to use the graph signal model to modeling the stepped frequency radar array, construct the framework of the full connection graph signal adjacency matrix through the expression of the steering vector of the radar receiving signal, and give the expression of the elements of the graph signal adjacency matrix in the case of azimuth, distance, and joint estimation. Secondly, through the data model in Section 2.2, the method of solving adjacency matrix elements and the response function of graph signal algorithm are deduced and analyzed, and a general algorithm for solving graph signal is proposed. Finally, through the construction of search domain, the graph signal solution algorithm of azimuth, distance, and joint estimation is obtained, which solves the problem of target estimation.

Fully Connected Graph of Array Signal Model
The graph data structure can be used for network modeling [19], graph processing [20], and biological computing [21]; in recent years, graph signals have been proposed as tools for parallelization and vectorization techniques [22].
In the previous Section, we analyzed the structure of radar frequency array and its data model. To study the graph adjacency matrix of radar array, we need to analyze the relationship between the spatial arrangement order and phase of graph signals. The radar signal array is arranged at equal intervals in the spatial domain, and the echo signal received by each unit has different phase due to different position and frequency. The actual radar array structure is shown in the figure below: The radar step frequency array joint estimation model is shown in Figure 2, and the corresponding graph signal structure is constructed through the radar signal array structure of Figure 2, as shown in Figure 3. The m-th element of the radar array corresponds to the m-th node in the graph signal, which constitute the point set V = {v 1 , · · · , v m }; inspired from the previous work of nofull connected graph signal algorithm [18], the expression of the phase difference of the m-th element relative to the n-th element in the array shown below according to Equation (4): Remote Sens. 2022, 14, x FOR PEER REVIEW 7 of 30 structure of Figure 2, as shown in Figure 3. The -th element of the radar array corresponds to the -th node in the graph signal, which constitute the point set ; inspired from the previous work of nofull connected graph signal algorithm [18], the expression of the phase difference of the -th element relative to the th element in the array shown below according to Formula (4): is different from the previous work [23], we use a fully connected way to construct the graph signal. The proposed graph signal matrix structure diagram is shown in Figure 3: · · · , α θ i ,R i m,n . This article is different from the previous work [23], we use a fully connected way to construct the graph signal. The proposed graph signal matrix structure diagram is shown in Figure 3.
It can be seen from the graph signal model that the size of the adjacency matrix depends on the number of radar elements M and step frequency points f M . Now we know that the number of points of the graph signal is M; the number of edges is M * M; and the expression of the edges is Equation (15). Then, the arrangement structure of adjacency matrix A θR is as follows: It can be seen from the graph signal model that the size of the adjacency matrix depends on the number of radar elements and step frequency points . Now we know that the number of points of the graph signal is ; the number of edges is ; and the expression of the edges is Formula (15). Then, the arrangement structure of adjacency matrix is as follows: It can be seen that the diagonal elements of are all 0, and the off-diagonal elements of each row are the delay values of the corresponding radar element relative to the diagonal element. However, these delay values cannot be obtained directly and need to be solved by specific matrix expressions, which will be explained in subsequent Sections.
Therefore, the joint estimation of the adjacency matrix is .
The above is the process of constructing the stepped frequency radar array graph adjacency matrix. The edges of the adjacency matrix contain the product of the delay caused by azimuth and distance estimation. In other words, the graph matrix contains estimates for two special cases, when the radar array only transmits single-frequency signals, in the weights of the edges degenerate to . Then, the expression of the edge becomes as follow: It can be seen that the diagonal elements of A θR are all 0, and the off-diagonal elements of each row are the delay values of the corresponding radar element relative to the diagonal element. However, these delay values cannot be obtained directly and need to be solved by specific matrix expressions, which will be explained in subsequent Sections. Therefore, the joint estimation of the adjacency matrix is A = A θR .
The above is the process of constructing the stepped frequency radar array graph adjacency matrix. The edges of the adjacency matrix contain the product of the delay caused by azimuth and distance estimation. In other words, the graph matrix contains estimates for two special cases, when the radar array only transmits single-frequency signals, ω i in the weights of the edges degenerate to ω. Then, the expression of the edge becomes as follow: At this time, the distance delay information caused by the frequency difference is missing, and then it degenerates into the graph signal matrix A = A θ for the azimuth estimation.
When the radar transmits the step frequency signal, but the number of array elements changes from M to 1, τ i in the weight of the edge degenerates to τ. Note that it does not mean that the delay of τ becomes a constant. Then, the expression of the edge becomes as follow: At this time, the delay caused by the distance is expressed in ω i,j , but the phase delay information brought by the array is missing; then, it degenerates into a graph signal matrix A = A R for distance estimation.

Graph Signal Estimation Algorithm
In this Section, we will derive the principle of solving the graph signal based on the definition of the adjacency matrix A θR (14), the definition of the graph signal edge α θ i ,R i m,n (13), the echo signal expression X rec (10), and the known conditions of the signal subspace Q S ∈ C M×(M−I) (8), to derive the principle of solving graph signal and give a general solution algorithm. Finally, for the false point problem of joint estimation, the CLEAN algorithm is used to eliminate false points.

1.
Graph signal estimation principle According to Equation (10) and the definition of the blue line signal in Figure 3, it can be known that X rec is a superposition form of orthogonal signals transmitted by M elements, so the received signal of the m-th element is: Therefore, the relationship of the phase field between two different array elements can be expressed as: This means that the received signal of the m-th element can be represented by the linear combination of the echo signal of the array element of M = m, so we assume that the signal vector Y ∈ C M×1 is the linear combination transformation of X rec : We take the first element in the Y vector and get the following formula: According to the Equations (6)-(8), the expression of X rec (m) are expands to: Substituting Equation (25) into Equation (24), we can get: The remaining terms of the vector Y satisfy the above expression, which means that we get the equation for solving the elements of the adjacency matrix A θ,R as below: This is a very important formula, and, in this way, it means that we can find each element in A θ,R with the linear algebraic expression of X rec . * X rec −1 only by writing the steering vector st(θ i , R i , m) rather than complex initialization settings (see the next Section for the specific algorithm). At the same time, X rec becomes a special eigenvector Q A of the adjacency matrix A θ,R corresponding to the unit eigenvalue, which inspired us to use the method similar to the covariance matrix R xx to extract the position and distance information in X rec , and it will be introduced in the next Section.

2.
General solving algorithm of graph signal From the theoretical derivation in the previous Section, we have noticed that the adjacency matrix A θ,R can also have special properties like the covariance matrix R xx after eigen-decomposition. The graph Fourier transform (GFT) can be applied to X rec to extract azimuth and distance information. Specifically, by performing subspace eigendecomposition on the adjacency matrix A θ,R , we can get the eigenvalue decomposition of the adjacency matrix A θ,R and the graph Fourier transform of the steering vector st θ i ,R i reference to [24]: where Y is the GFT of the X rec , and Q H A is the graph Fourier transform operator, which contains the steering vector st θ i ,R i . It can know that the steering vector st θ i ,R i is one or more eigenvectors of the adjacency matrix A θ,R , so when θ and R is aligned to the target, there will be non-zero positions corresponding to the unit eigenvalues of the steering vector st θ i ,R i .
Then, the adjacency matrix A θ,R is further subspace decomposed to obtain [25]: The eigenvalue matrix is Λ = diag M , · · · , 1 M , 0, · · · , 0 , the noise subspace V N ∈ C M×(M−I) satisfies that the covariance matrix is time-independent and completely orthogonal to the steering vector st θ i ,R i [26], so it can think that for the correct estimation result, an expression that measures the degree of orthogonality between the steering vector st θ i ,R i and the noise subspace can be derived as below: Generally, the most correct DOA result is found by searching the peak value, and the multiple root principle of eigenvalue matrix limits that the result of F GSP (θ i , R i ) is generally not too large. Therefore, the purpose of searching the peak value is realized by removing I maximum eigenvalues and taking the reciprocal: final expression for measuring the degree of orthogonality between steering vector st θ i ,R i and noise subspace can be obtained as follows: As shown in the equation, remove the eigenvectors corresponding to the I largest eigenvalues, then obtain the joint estimation result F GSP (θ i− , R i− ).

Joint azimuth and distance Estimation algorithm
Through the derivation of the above Sections, we now have obtained the graph matrix structure A θR (16) of the radar stepped frequency array and clarified the key expression for solving the elements α θ i ,R i m,n of graph matrix (27). Through the graph Fourier transform, the optimization function (32) of the azimuth and distance estimation algorithm is derived. Before implementing the specific algorithm, there are some prerequisites that need to be clarified:  (10) and (11) to construct a read received signal. In this article, for the convenience of calculation, X rec is compressed into a X rec ∈ C M×1 vector form. • α θ i ,R i m,n in adjacency matrix A θR is obtained by inversing calculation using Equation (27), by performing different inversion operations by distinguishing the three cases of the row number m in adjacency matrix A θR : where X i rec is the search domain corresponding to the i-th target, and the expression is exactly the same as X rec .
is the reciprocal of the two-norm Fourier transform Q H A of X rec with the k largest eigenvalues removed, where Q H A = EVD(A θR ) according to Equation (32). Now through the above preliminary conditions, the algorithm flow calculation can be performed, as described in Algorithm 1. By constructing the transmitting and receiving steering vectors C and B by using Equation (4), adding Gaussian white noise, the real received signal X rec is obtained. Through k targets, construct k search domains, calculate the corresponding graph matrix A θ i ,R i for each search domain grid, calculate the response function F GSP (θ i− , R i− ) of the current grid through Equation (32), and finally get the estimated result.
Receive steering vector B according to (7). Determine real received signal expression X rec according to (10). (29). (3) Response function solving step: Sort the eigenvalues by Q H A and delete I maximum eigenvalue.
The above is the general solution form of the graph signal solution estimation algorithm, but for distance, azimuth, and joint estimation, it is necessary to clarify the search domain R I or R I+1 on the basis of Algorithm 1; different response functions F GSP are distinguished according to different search domain R I or R I+1 , as shown in Algorithm 2. , θ i ∈ R I+1 , apply the Algorithm 1 to determine F GSP (θ i− , R).

Algorithm 2: Various search domains of graph signal algorithm
(2) Establish search domain θ ← θ max −θ min 2 , R i ∈ R I+1 apply the Algorithm 1 to determine According to Algorithm 2, we can know that when the joint estimation algorithm traverses the search space R I+1 , it needs to fix one type of parameter first, and after the other type of parameter is traversed to complete the estimation, then the other type of parameter is fixed for traversal. In this way, the phase coupling in the α θ i ,R i m,n can be eliminated, and the joint estimation response function F GSP (θ i− , R i− ) can be obtained.

CLEAN Algorithm for False Point Elimination
It can be seen from the Case c from the Algorithm 2, the solution method of the graph signal for joint estimation is to construct the received signal in the azimuth and distance, respectively, and then estimate them separately; so, for the multi-target case (I > 1), the final joint estimation is: There will be I * 2 peaks in the response graph of the estimation result, which is also one of the shortcomings of the joint estimation of the graph signal.
The CLEAN algorithm was first applied in the field of radio astronomy to reconstruct the sky brightness distribution closer to the real one [14]. It can effectively reduce the imaging sidelobe in the field of radar imaging. With the use of the CLEAN algorithm one can achieve the effect of eliminating false point targets [26]. Taking inspiration from [27], we can find the maximum peak value P 0 ; save its azimuth and distance {θ 0 , R 0 }; and, according to this information, find the corresponding F GSP (θ i− , R) and F GSP (θ, R i− ) response result functions and obtain the remaining signal by subtracting the corresponding function from the F GSP (θ i− , R i− ). Repeat the above steps until the remaining I maximum values of F GSP (θ i− , R i− ) [26]. In the end, only M peak target point will appear in the response map, as described in Algorithm 3.

Simulations and Results
In this Section, a series of numerical simulations are carried out to verify the effectiveness and performance of our proposed fully connected graph signal structure compared to traditional methods and non-fully connected graph structures. This article conducts verification experiments on three estimation problems and verifies and analyzes the azimuth, distance, and joint estimation, respectively.

Azimuth Estimation Results and Monte Carlo Analysis
For the azimuth estimation situation, consider a radar array composed of M = 8 collocated transceiver antennas with a uniform interval of d = λ/2. The carrier frequency is f 0 = 10 GHz, and the snapshot time Q = 128, under the 30 dB SNR scenario. In the dualtarget azimuth estimation experiment, it is assumed that there are two near-field targets located at {−30 • , 45 • }, respectively; in the three-target experimental azimuth estimation, it was assumed that there are three far-field targets located at {−45 • , −20 • , 30 • }. In addition, we assumed that the medium through which electromagnetic waves propagate is uniform, low-loss, and the relative permittivity ε 0 = 3, and the noise model in the numerical simulation is built with complex additive Gaussian white zero mean noise. Figure 4 shows the result response results of the graph signal azimuth to the dual and multi-target estimation. Figure 4b,d are the dimensionality reduction plan views of Figure 4a,c, respectively. The blue curve is the result of the azimuth estimation of the full GSP algorithm; the red curve is the result of the estimation of the nofull GSP algorithm according the previous work [23]; and the yellow curve is the original azimuth parameter. From Figure 4b of the dual target estimation, it can be seen that the blue line represented by the fully GSP can achieve the −3 db sidelobe effect at the −40 db response strength and can accurately estimate the target point {−30 • , 45 • }; compared with the red line of the nofull GSP algorithm, the full GSP has lower side lobes and more accurate results. From Figure 4d of the multi-target estimation, it can be seen that the full GSP can also accurately estimate the azimuth parameters of the target point {−45 • , −20 • , 30 • }. At this time, the F GSP (θ i− ) of nofull GSP has been incorrectly estimated, which proves that the full GSP algorithm is better in azimuth than the nofull GSP algorithm.
can accurately estimate the target point ; compared with the red line of the nofull GSP algorithm, the full GSP has lower side lobes and more accurate results. From Figure 4d of the multi-target estimation, it can be seen that the full GSP can also accurately estimate the azimuth parameters of the target point . At this time, the of nofull GSP has been incorrectly estimated, which proves that the full GSP algorithm is better in azimuth than the nofull GSP algorithm.  The experiment also carried out Monte Carlo analysis 500 times to analyze the influence of the signal-to-noise ratio of the azimuth estimation performance. For further comparison and analysis, the experiment quoted the classic MUSIC estimation algorithm [27]. It can be seen in Figure 5a that when the signal-to-noise ratio is lower than 0 dB, the signals of both full GSP and nofull GSP algorithm are better than MUSIC performance of the algorithm; in Figure 5b, it can be seen that because of the complexity of multi-target estimation algorithm, we set the search boundary under high SNR. Experiments show that the nofull GSP algorithm cannot accurately estimate the multi-target, and the search results appear on the search boundary or even exceed the boundary; the fully GSP algorithm still shows better performance than the MUSIC algorithm when the SNR is lower than 0 dB.

Distance Estimation Results and Monte Carlo Analysis
For the distance estimation situation, considering a single radar array element emitting N f = 8 frequency point signals to form a bandwidth signal, the experimental scene is that near-field dual-targets are located at {9750 m, 10, 040 m}, and near-field multi-targets are located at {9650 m, 10, 050 m, 10, 200 m}. Figure 6 shows the results of the estimated signal distance to the dual-target and multitarget estimation. From Figure 6b of dual-target estimation, it can be seen that the full GSP can still accurately estimate the distance parameter to the target point {9750 m, 10, 040 m}, However, the nofull GSP algorithm has been unable to estimate the distance of double targets, and the response function F GSP (R i− ) has been unable to find the accurate maximum value. From Figure 6d of multi-target estimation, it can be seen that the full GSP algorithm can also accurately estimate the distance parameters {9650 m, 10, 050 m, 10, 200 m}; but for the nofull GSP algorithm, its F GSP (R i− ) has a corresponding mean value similar to Gaussian white noise. Experiments show that the full GSP algorithm can also achieve the distance estimation. Next, the full GSP algorithm is analyzed by Monte Carlo experiment and compared with the performance of MUSIC algorithm. The experiment also carried out Monte Carlo analysis 500 times to analyze the influence of the signal-to-noise ratio of the azimuth estimation performance. For further comparison and analysis, the experiment quoted the classic MUSIC estimation algorithm [27]. It can be seen in Figure 5a that when the signal-to-noise ratio is lower than 0 dB, the signals of both full GSP and nofull GSP algorithm are better than MUSIC performance of the algorithm; in Figure 5b, it can be seen that because of the complexity of multi-target estimation algorithm, we set the search boundary under high SNR. Experiments show that the nofull GSP algorithm cannot accurately estimate the multi-target, and the search results appear on the search boundary or even exceed the boundary; the fully GSP algorithm still shows better performance than the MUSIC algorithm when the SNR is lower than 0 dB.

Distance Estimation Results and Monte Carlo Analysis
For the distance estimation situation, considering a single radar array element emitting frequency point signals to form a bandwidth signal, the experimental scene is that near-field dual-targets are located at , and near-field multi-targets are located at . Figure 6 shows the results of the estimated signal distance to the dual-target and multi-target estimation. From Figure 6b of dual-target estimation, it can be seen that the full GSP can still accurately estimate the distance parameter to the target point , However, the nofull GSP algorithm has been unable to estimate the distance of double targets, and the response function has been unable to find the accurate maximum value. From Figure 6d of multi-target estimation, it can be seen that the full GSP algorithm can also accurately estimate the distance parameters ; but for the nofull GSP algorithm, its has a corresponding mean value similar to Gaussian white noise. Experiments show that the full GSP algorithm can also achieve the distance estimation. Next, the full GSP algorithm is analyzed by Monte Carlo experiment and compared with the performance of MUSIC algorithm.
Angle RMSE(degree) Angle RMSE(degree) The distance estimation has also undergone 500 Monte Carlo analysis. Considering that the nofull GSP algorithm cannot perform the distance estimation, only the result of the full GSP and the MUSIC algorithm is compared here. As can be seen from Figure 7a,b, full GSP algorithm is better than the MUSIC algorithm under any SNR condition, and the estimation result is more accurate. The specific reasons will be analyzed in Section 4.

Joint Estimation Results
Finally, for the case of joint estimation, consider a radar array composed of M = 8 collocated transceiver antennas. Each array element transmits ∆f frequency interval signals to form a broadband signal. Under the 30 dB SNR condition, experiment scene 1 is single-target estimation, and the target parameter is {42 • , 9900 m}; experiment scene 2 is multi-target estimation, and the target parameter is {30 • , 10, 200 m},{−45 • , 9750 m}. The rest of the conditions are the same as the above experiment. Figure 8 shows the result of the joint estimation for single target. It can be seen from Figure 8a that the full GSP algorithm can accurately estimate the azimuth and distance parameters to the target point. The target point's parameters {42 • , 9900 m} is less than or equal to −35 dB, achieving good two-dimensional estimation results. Experiments have also proved that the graph signal method has good response results in azimuth, distance, and joint estimation. Figure 8b shows the slice diagram of the full GSP estimation algorithm in azimuth and distance dimension. Since F GSP (θ i− , R). * F GSP (θ, R i− ) will only produce one target point under single target joint estimation and will not produce false points, the next experiment shows multi-target joint estimation and CLEAN algorithm to eliminate false points.  The distance estimation has also undergone 500 Monte Carlo analysis. Considering that the nofull GSP algorithm cannot perform the distance estimation, only the result of the full GSP and the MUSIC algorithm is compared here. As can be seen from Figure 7a,b, full GSP algorithm is better than the MUSIC algorithm under any SNR condition, and the estimation result is more accurate. The specific reasons will be analyzed in Section 4. The distance estimation has also undergone 500 Monte Carlo analysis. Considering that the nofull GSP algorithm cannot perform the distance estimation, only the result of the full GSP and the MUSIC algorithm is compared here. As can be seen from Figure 7a,b, full GSP algorithm is better than the MUSIC algorithm under any SNR condition, and the estimation result is more accurate. The specific reasons will be analyzed in Section 4. or equal to −35 dB, achieving good two-dimensional estimation results. Experiments have also proved that the graph signal method has good response results in azimuth, distance, and joint estimation. Figure 8b shows the slice diagram of the full GSP estimation algorithm in azimuth and distance dimension. Since will only produce one target point under single target joint estimation and will not produce false points, the next experiment shows multi-target joint estimation and CLEAN algorithm to eliminate false points.
(a) (b)   Figure 9 shows the result of joint estimation of azimuth and distance for multi-point target. It can be seen from Figure 9a that the fully connected graph signal can perform accurate joint estimation of the azimuth and distance for multi-target points, but there will be the false target point phenomenon, and as shown in the figure, the target points {30 • , 9750 m}, {−45 • , 10, 020 m} are false, which is related to the graph signal solving algorithm itself. The sidelobe value near the target point parameter {30 • , 10, 020 m}, {−45 • , 9750m} is less than or equal to −35 dB, achieving a good two-dimensional estimation result. It can be seen that in the case of multiple targets, false points will appear, which need to be processed through the CLEAN algorithm, and the false points are eliminated through the construction of the Point Spread Function to obtain the correct result.  Figure 10 shows the results of eliminating false points in the complete GSP joint estimation results using the CLEAN algorithm according to Algorithm 3. In the process of clean algorithm, Figure 10a is the normalized response diagram of joint estimation of graphic signals after using CLEAN algorithm. It can be seen that CLEAN algorithm retains the distance and phase information of real point targets and eliminates false points accurately. Figure 10(b) is a sectional view of the azimuth and distance of the GSP algorithm. It can be seen that the fully connected map signal can accurately estimate the azimuth and distance information of multiple targets.  of graphic signals after using CLEAN algorithm. It can be seen that CLEAN algorithm retains the distance and phase information of real point targets and eliminates false points accurately. Figure 10b is a sectional view of the azimuth and distance of the GSP algorithm. It can be seen that the fully connected map signal can accurately estimate the azimuth and distance information of multiple targets.
(a) (b) Figure 9. GSP Joint estimation of azimuth and distance results for multi-targets (a) azimuth and distance dimensional response map (b) slice diagram response map. Figure 10 shows the results of eliminating false points in the complete GSP joint estimation results using the CLEAN algorithm according to Algorithm 3. In the process of clean algorithm, Figure 10a is the normalized response diagram of joint estimation of graphic signals after using CLEAN algorithm. It can be seen that CLEAN algorithm retains the distance and phase information of real point targets and eliminates false points accurately. Figure 10(b) is a sectional view of the azimuth and distance of the GSP algorithm. It can be seen that the fully connected map signal can accurately estimate the azimuth and distance information of multiple targets.

Discussion
In this Section, we will discuss and analyze the estimation principle and simulation results of GSP, which include the performance analysis and the principle interpretation of various algorithms; discussion on the nofull GSP, GSP and MUSIC algorithm complexity piquancy func/(dB) piquancy func/(dB)

Discussion
In this Section, we will discuss and analyze the estimation principle and simulation results of GSP, which include the performance analysis and the principle interpretation of various algorithms; discussion on the nofull GSP, GSP and MUSIC algorithm complexity and computational load comparison; and discussion on the solution of unknown estimated target number k.

Results and Performance Analysis of GSP Estimation Method
In Figures 5 and 7, we have done 500 Monte Carlo experiments on the azimuth estimation and distance estimation respectively, and compared the accuracy of GSP, Nofull GSP and music algorithm.
In the azimuth estimation experiment, we can observe that GSP algorithm is better than Nofull GSP algorithm in any SNR condition, because GSP algorithm build the adjacency matrix A θR by using the phase relationship between all array elements, so when the search domain parameters R I are aligned with the target I, the graph Fourier operator Q A forms the same linear phase relationship with the echo signal X rec , which satisfies Equation (27). At this time, the I maximum eigenvalues of the graph Fourier operator Q A correspond to the I targets of the echo signal, Therefore, GSP algorithm show better performance than Nofull GSP algorithm.
We can also observe that the performance of GSP algorithm is not as good as MUSIC algorithm when the SNR is greater than 0, because GSP algorithm uses the form of 2πf when constructing the steering vector, and in the distance search domain θ i ∈ R I under far-field conditions, the ratio of the numerator denominator of the exponential term of the guidance vector Equation (17) will approach a variable that changes greatly θ i with, as shown below: which makes the elements of the adjacency matrix A θR prone to periodic repetition, and the linear relationship between Q A and X rec of echo signal after graph Fourier transform is weakened, resulting in the poor performance of GSP algorithm than MUSIC algorithm under high SNR condition.
In the distance estimation experiment, because the nofull GSP algorithm has been unable to estimate the information about distance, we only compare the performance of GSP and MUSIC algorithm. It can be observed that GSP algorithm is better than MUSIC algorithm in the any SNR condition, because in the distance estimation, GSP algorithm uses the form of ∆f when constructing the steering vector, while in the distance search domain R i ∈ R I under far-field conditions, the proportion of the numerator denominator of the exponential term of the guidance vector Equation (19) will approach a small constant C, as shown below: This small constant C avoids the periodic repetition of the elements of the adjacency matrix A θR , making it easier for the graph Fourier operator Q A to form the same linear phase relationship with the echo signal X rec ; Moreover, GSP algorithm constructs an adjacency matrix on R I search dimensions to estimate the target. Compared with one search estimation of MUSIC algorithm, GSP algorithm has better orthogonality to the signal subspace Q S and noise subspace Q N . At the same time, graph Fourier operator Q A can better fit the accurate information of multi targets. Therefore, in distance estimation, GSP algorithm performs better than MUSIC algorithm under the condition of any SNR; At the same time, the construction of multiple search dimensions of GSP algorithm can ensure that the signal subspace Q S and noise subspace Q N are still orthogonal under the low SNR condition, which also explains that the performance of GSP algorithm is better than MUSIC algorithm in the case of low SNR in azimuth estimation.
The Monte Carlo experimental data results of azimuth and range algorithms are shown in Tables 1 and 2. From Table 1, it can be seen that in azimuth estimation, the performance of GSP algorithm is better than that of nofull GSP algorithm under the condition of any SNR; When the SNR is lower than 0 dB, the performance of GSP algorithm is better than MUSIC algorithm. It can be seen from Table 2 that in distance estimation, GSP algorithm has better performance than MUSIC algorithm under the condition of any SNR. As we can see from Table 3, in azimuth estimation, GSP algorithm has 45.6% performance improvement compared with MUSIC algorithm under the condition of low SNR; In distance estimation, GSP algorithm has 47.4% performance improvement compared with music algorithm under the condition of any SNR.

Analysis of Time Complexity and Computational Load
According to Section 4.1, the characteristic of GSP algorithm is that when the search domain R I is aligned with the target parameters, it will get a huge response value; It can be seen from Algorithm 1, the solution process of GSP algorithm includes solving the element value of the adjacency matrix A ∈ C N * N ; solving the eigenvector of the adjacency matrix according to Equation (29); removing I maximum values and solving the response function F GSP (32), and the above processes are completed iteratively in the search domain grid G = 100; The solution process of nofull GSP and GSP algorithm is exactly the same, except for a little difference in the process of graph modeling. Therefore, the computational complexity of nofull GSP algorithm is the same as GSP, and their computational complexity is deduced as follows: (1) Computing the element of the adjacency matrix A θR by the received signal X rec is O(N * (N − 1)). (2) Computing the EVD of adjacency matrix A θR is O N 3 .
(3) Sort Q A X rec , and delete the maximum value of the first I items, then counting down to get the response value result F GSP is O(N * I + 1). (4) The above process is about the iteration of the I-th power of the search domain G, and its algorithm complexity is O G I . Therefore, the algorithm complexity of nofull GSP and GSP algorithm is From the above analysis, it can be seen that the full GSP algorithm constructs the search domain R I and uses the increased computing cost to improve the performance of the estimation algorithm, which is indeed a defect in the design of the full GSP algorithm.
And the solution process of MUSIC algorithm includes: solving the covariance matrix R xx ; solving the eigenvector of the covariance matrix R xx ; removing I maximum values, then solving the noise subspace Q S and calculating the response value F MUSIC . The above processes are completed iteratively in the search domain grid G = 100. Therefore, the computational complexity of MUSIC algorithm is deduced as follows: (1) Computing the covariance matrix R xx by the received signal X rec is O N 2 .
(2) Computing the EVD of the covariance matrix R xx is O N 3 .
(3) Sort R xx , and delete the maximum value of the first I items, solving the noise subspace Q S and calculating the response value F MUSIC is O(N * I + 2). (4) The above process is about the iteration of search domain G, and its algorithm complexity is O(G).
By comparing Equations (36) and (37), it can be seen that the algorithm complexity of nofull GSP and GSP algorithm includes the high-order power O G I of the search grid, so their computational complexity mainly depends on the establishment and traversal of the search element domain G (related to the number of targets I). When the number of targets I increases, the increase of computation will be exponential, Therefore, the complexity of nofull GSP and GSP algorithm without parallel optimization is worse than MUSIC algorithm, the complexity of these algorithms is shown in Table 4 as below: Table 4. Computational Complexity result about MUSIC, nofull GSP and GSP algorithm.

Algorithm Computational Complexity
Therefore, we need to speed up the GSP and nofull GSP algorithm through the scheme of parallel computing to reduce the computing load and eliminate the influence of O G I . It can be seen that the search domain R I will produce a large number of core algorithm processes, which will show a great computational load in the process of serial operation. Fortunately, this large-scale small matrix is very suitable for parallel acceleration. We can greatly improve the computing speed of GSP and nofull GSP algorithm and reduce the computing time of GSP and nofull GSP algorithm by constructing the highdimensional tensor of adjacency matrix and making it into a data form suitable for MATLAB parallel accelerated computing. The schematic diagram of MATLAB parallel computing acceleration is shown in Figure 11 as below: Remote Sens. 2022, 14, x FOR PEER REVIEW 24 of 30 targets increases, the increase of computation will be exponential, Therefore, the complexity of nofull GSP and GSP algorithm without parallel optimization is worse than MU-SIC algorithm, the complexity of these algorithms is shown in Table 4 as below: Therefore, we need to speed up the GSP and nofull GSP algorithm through the scheme of parallel computing to reduce the computing load and eliminate the influence of . It can be seen that the search domain will produce a large number of core algorithm processes, which will show a great computational load in the process of serial operation. Fortunately, this large-scale small matrix is very suitable for parallel acceleration. We can greatly improve the computing speed of GSP and nofull GSP algorithm and reduce the computing time of GSP and nofull GSP algorithm by constructing the highdimensional tensor of adjacency matrix and making it into a data form suitable for MATLAB parallel accelerated computing. The schematic diagram of MATLAB parallel computing acceleration is shown in Figure 11 as below: In this experiment, AMD Ryzen 3600× processor and 32 g memory are used for accelerated experimental simulation under Windows10 system. The unoptimized GSP algorithm, optimized GSP algorithm and MUSIC algorithm are based on the condition that the search domain gird is 100. The results are as follows: As can be seen from Figure 12, under the search grid setting of , the operation time of the optimized GSP and nofull GSP algorithm is 465 times faster than that of the unoptimized GSP and nofull GSP algorithm, which is close to the MUSIC algorithm.
Thanks to the support of MATLAB for tensor parallel computing, in the complexity of nofull GSP and GSP algorithm can be approximately regarded as the lower order In this experiment, AMD Ryzen 3600× processor and 32 g memory are used for accelerated experimental simulation under Windows10 system. The unoptimized GSP algorithm, optimized GSP algorithm and MUSIC algorithm are based on the condition that the search domain gird is 100. The results are as follows: As can be seen from Figure 12, under the search grid setting of G = 100, the operation time of the optimized GSP and nofull GSP algorithm is 465 times faster than that of the unoptimized GSP and nofull GSP algorithm, which is close to the MUSIC algorithm.
Thanks to the support of MATLAB for tensor parallel computing, O G I in the complexity of nofull GSP and GSP algorithm can be approximately regarded as the lower order term of O(G), which reduces the computational load of nofull GSP and GSP algorithm to an acceptable range.

Estimation Method of Unknown Target Number K
In the experiment in Section 3, we all carry out verification simulation and a performance test on the GSP algorithm assuming that the number of targets is known. However, in the actual scene, we could not know the number of specific targets , so we need to conduct experimental analysis when the number of targets is unknown.
Firstly, the Unknown Number of Source Estimation and Blind Signal Estimation is a key problem in the field of Array Signal Processing. Without the number of targets, it is impossible to further measure and estimate the targets. For the typical target parameter estimation problem, when the number of targets should be known or the number of targets should be obtained through the Blind Signal Estimation algorithm in advance, the target parameters should be further estimated and measured. In the field of Blind Signal Estimation, the basic idea of source number estimation is to calculate the covariance matrix through the echo signal . By solving the SVD eigenvalues of the covariance matrix, the larger eigenvalues in the eigenvalues are accumulated to estimate the number of sources [28][29][30].

Estimation Method of Unknown Target Number K
In the experiment in Section 3, we all carry out verification simulation and a performance test on the GSP algorithm assuming that the number of targets I is known. However, in the actual scene, we could not know the number of specific targets I, so we need to conduct experimental analysis when the number of targets I is unknown.
Firstly, the Unknown Number of Source Estimation and Blind Signal Estimation is a key problem in the field of Array Signal Processing. Without the number of targets, it is impossible to further measure and estimate the targets. For the typical target parameter estimation problem, when the number of targets should be known or the number of targets K should be obtained through the Blind Signal Estimation algorithm in advance, the target parameters should be further estimated and measured. In the field of Blind Signal Estimation, the basic idea of source number estimation is to calculate the covariance matrix R xx (n) through the echo signal x(n). By solving the SVD eigenvalues of the covariance matrix, the larger K eigenvalues in the eigenvalues are accumulated to estimate the number of sources K [28][29][30].
Secondly, we can know from the general solving algorithm of graph signal in Section 2.4.2 that when we calculate the eigen decomposition of adjacency matrix A, the target number I maximum eigenvalues will be generated. It should be noted that the number of maximum eigenvalues does not depend on the preset K, but depends on the number of actual targets I, which means that when the number of search variables K in the search area increases to the number of actual targets I, there will be I maximum eigenvalues generated, what we need to do is to increase the number of search variables K. To sum up, the full GSP algorithm can indeed estimate the size of the real K by increasing the number of K while the number of targets K is unknown, but this process is the inverse process of the Unknown Number of Source Estimation algorithm. The experimental results can only verify that the full GSP algorithm has the characteristics of estimating the number of targets K but cannot verify that the GSP algorithm is a necessary condition of the Unknown Number of Source Estimation algorithm.
For the case of unknown target number I, the model Equation (6) needs to be modified as follows: The received signal X rec form of unknown target number K is generated through the above formula. Assuming that the current estimated target number is K, GSP estimation experiment is carried out. All conditions of this experiment are the same as those in Section 3, Figure 13 shows the response results, judgment threshold and algorithm principle of GSP algorithm under different target number K as below: The received signal form of unknown target number is generated through the above formula. Assuming that the current estimated target number is , GSP estimation experiment is carried out. All conditions of this experiment are the same as those in Section 3, Figure 13 shows the response results, judgment threshold and algorithm principle of GSP algorithm under different target number as below: We discuss different situations for different estimated number k as below: 1. When the estimates number , the search domain of GSP algorithm is , and function will not have sharp peaks and will cause confusion in the estimation results. As shown in (a) and (b) of Figure 13, the response of the search result is very small. This is because even if the parameter of the search domain is aligned with

1.
When the estimates number K < I, the search domain of GSP algorithm is R K , and F GSP function will not have sharp peaks and will cause confusion in the estimation results. As shown in Figure 13a,b, the response of the search result is very small. This is because even if the parameter R K of the search domain is aligned with any K in the target I, I maximum eigenvalues will appear as long as K < I. and deleting the eigenvalues of the first K large numbers in Q A X rec will not make F GSP have a huge response value, so it is impossible to obtain accurate estimation when K < I.

2.
When the estimated number K = I, the search domain of GSP algorithm R K is equal to R I ,F GSP function will produce K sharp peaks, and the size of each spike is greater than 10e4. At this time, it can be preliminarily determined that the estimated number of targets is K. As shown in Figure 13c, the searched parameter R K is equal to the actual number of targets I, resulting in I large eigenvalues. Therefore, when K = I, deleting K eigenvalues will make Q A X rec produce a huge response. At this time, it can be preliminarily determined that K is the target number.

3.
When the estimated number K > I, the search domain of GSP algorithm R K has one more dimension than R I , the F GSP function will also only produce I spikes, and no matter how the K value is increased, the F GSP response result remains unchanged. At this time, it can be determined that the estimated target is I. As shown in Figure 13d, at this time, the parameter R K of the search domain is aligned with the target I, whether the remaining K − I dimension is aligned with the target or not, it will only make Q A X rec produce I maximum values. At this time, even deleting the first K maximum values will not affect the estimation results F GSP , because the remaining eigenvalues have nothing to do with the target parameters, even if deleted, it can still produce a huge corresponding effect on. At this time, it can be determined that the previously estimated K is the target number I.

Conclusions
In this article, a joint estimation method based on graph signal for far-field targets was proposed, and a general solution to the problem of graph signal estimation was given. The proposed method constructs a fully connected graph signal matrix through the radar array model and archives the azimuth and distance joint estimation with sufficient estimation accuracy. In this article, the derivation and analysis showed that the estimation method based on fully connected graph signal makes full use of the phase information between array elements and has better performance than nofull GSP algorithm and MUSIC algorithm. The experimental results show that, in azimuth estimation, the GSP algorithm has a 45.6% performance improvement compared with the MUSIC algorithm under the condition of low SNR. In distance estimation, the GSP algorithm has 47.4% performance improvement compared with the MUSIC algorithm under the condition of any SNR under the far-filed conditions. Numerical experiments show that the GSP algorithm is superior to the MUSIC algorithm, because it exchanges the computational load for the improvement of algorithm performance.
For the problem of computational load, the unoptimized GSP algorithm shows a huge computational load in computational complexity. The main reason is that the construction of the search domain of GSP algorithm leads to exponential terms in its computational complexity. Through MATLAB parallel computing acceleration experiment, it was proved that the speed of the optimized GSP algorithm was increased by 465 times, which greatly improves the engineering applicability of the GSP algorithm. In the future, our work will focus on the accelerated application of GSP algorithm in Spark distributed computing framework.
For the assumption that the target number k of GSP algorithm is unknown, we drew a conclusion through theoretical derivation and experimental analysis: when the number of targets K is unknown, we can estimate the K dimension by gradually increasing the search domain R K . If the current response value is less than 10e4, continue to increase the K value to estimate the next search domain. If the current response value is greater than 10e4, and further increasing the K value will not affect the size of the response value, it can be considered that the K with the first response value reaching 10e4 is the number of real targets I.
Author Contributions: K.L. and J.J. conducted the algorithm design; Z.Y. and N.X. made a MATLAB implementation of the proposed algorithm and formulated the proposed algorithm. Z.Y. contributed to prepare and analyze the experimental data and the results. All authors were involved in modifying the article, the literature review, and the discussion of the results. All authors have read and agreed to the published version of the manuscript.