A direct positioning method based on orthogonality of cross-spectral subspace

Abstract

The present invention discloses a direct positioning method based on the orthogonality of cross-spectrum subspace. First, a reference node is selected, its received signal is set as the reference signal, and all received signals are segmented; secondly, the cross-spectrum between the segmented received signal and the reference signal is calculated, and the covariance matrix of the effective part of the obtained cross-spectrum is subjected to eigenvalue decomposition; then, the cost function under the cross-spectrum is constructed by using the orthogonality of the subspace; finally, the position of the radiation source is determined by searching for the superposition of these cost functions on the grid points of the target area. Compared with other direct positioning methods, this method has higher accuracy and lower complexity. In a multipath environment, this method can also maintain a high estimation performance. Since each monitoring node is only equipped with a single positioning sensor, it has a lower hardware cost than the direct positioning method based on the antenna array.

Description

Technical Field

[0001] The present invention relates to the field of passive positioning, and in particular to a direct positioning method based on cross-spectral subspace orthogonality. Background Art

[0002] Passive positioning technology has become an increasingly important positioning method due to its strong concealment and high positioning accuracy. It has been widely studied and applied in both military and civilian fields. According to the topological structure of the positioning station, radiation source positioning can be divided into single-station and multi-station radiation source positioning. The former is easy to implement but has low accuracy; the latter is complex to implement but has better performance. Positioning stations can be divided into array-based positioning stations and single-sensor-based positioning stations, among which single-sensor-based positioning stations have lower costs. In actual application scenarios, due to the complexity and diversity of channel propagation and the environment around the receiver, there is often multipath propagation between the source signal and the observation station. As an important part of passive positioning, the two-step positioning method first needs to measure intermediate parameters, such as direction of arrival (DOA), time of arrival (TOA), time difference of arrival (TDOA), and received signal strength (RSS), and then use these parameters to establish the positioning equation and calculate the target position. In multi-station positioning, the time difference based method is widely used due to its low complexity and real-time availability. Its key technology is time delay estimation (TDE). In two-step positioning, the measurement process of the intermediate parameters is independent of the calculation of the target position, which leads to the lack of spatial geometric constraints between the target and the observation station.

[0003] Unlike two-step positioning, the direct position determination (DPD) method does not require the estimation of intermediate parameters, and its accuracy is better than the traditional two-step method at low signal-to-noise ratio. However, the DPD algorithm often faces the problem of high algorithm complexity caused by high-dimensional matrix operations. In addition, in order to avoid the expensive hardware cost of the array-based DPD method, a TDOA-based DPD method has emerged, in which each observation station is equipped with a single sensor. With the increasing complexity of the electromagnetic environment, the transmission of signals in a multipath environment will cause many passive positioning algorithms to produce large errors. Therefore, there is still a lot of room for improvement in the passive positioning methods for multipath environments. Summary of the invention

[0004] The technical problem to be solved by the present invention is to provide a direct positioning method based on the orthogonality of the cross-spectrum subspace in view of the defects involved in the background technology. Based on the TDOA positioning scenario, by acquiring the cross-spectrum between multipath signals, the factors containing the time delay difference information from different nodes to the radiation source are separated, a cost function is established, and the estimated position of the radiation source that meets the subspace orthogonality is found by using the spectrum peak search method. The cross-spectrum of possible related signals is processed by spatial smoothing technology, so that the positioning result can be made more accurate. In addition, a low-cost single-sensor distributed node is used to receive the signal from the radiation source, and by calculating the cross-spectrum of other signals and the reference signal, direct processing of the received signal is avoided, thereby reducing the amount of calculation; and the complexity is reduced by data segmentation and screening of the effective part of the cross-spectrum.

[0005] The present invention adopts the following technical solutions to solve the above technical problems:

[0006] A direct positioning method based on cross-spectral subspace orthogonality comprises the following steps:

[0007] Step 1), receive signals from a radiation source p at an unknown location at L monitoring nodes with known locations, transmit the L received signals to a central site, L ≥ 3, and select one of the monitoring nodes as a reference node;

[0008] Step 2), evenly segment the L received signals, and divide each received signal into K segments; use the received signal of the reference node as the reference signal; for each of the other L-1 received signals, perform cross-correlation between the k-th segment received data and the k-th segment received data of the reference signal, k = 1, ..., K, and obtain (L-1) × K cross-correlation functions; perform discrete Fourier transform on the obtained cross-correlation functions in turn to obtain their corresponding cross-spectra;

[0009] Step 3), eigendecomposing the covariance matrix of the effective part of the cross spectrum to obtain its noise subspace and signal subspace, and constructing the cost function under the cross spectrum according to the orthogonality of the noise subspace and the signal subspace;

[0010] Step 4), determine the positioning area and the search step, grid the positioning area, perform peak search on the cost function values ​​on all grid points, and the coordinates corresponding to the peak are the estimation results for the radiation source.

[0011] As a further optimization scheme of the direct positioning method based on the orthogonality of the cross-spectral subspace of the present invention, the specific content of step 1) is as follows:

[0012] Let the coordinates of the reference node be q1, and the coordinates of the other L-1 monitoring nodes be q l, l = 2, ..., L, the discrete form of the received signal of each monitoring node in the multipath condition is:

[0013]

[0014] Where x1(n) represents the discrete form of the received signal of the reference node in the case of multipath, x l (n) represents the discrete form of the received signal of the non-reference node in the case of multipath, M represents the number of multipaths, and N l represents the length of the received signal, λ 1m , lm They represent x1(n), x l (n) amplitude coefficient; τ 1m is the time delay between the mth multipath component in the reference signal and the transmitted signal; τ lm is the time delay between the mth multipath component in the lth received signal and the transmitted signal; s(n) represents the transmitted signal of the unknown radiation source, s(n-τ 1m ) represents the time delay τ 1m The emission signal of the unknown radiation source, s(n-τ lm ) represents the time delay τ lm The emission signal of the unknown radiation source; σ1(n), σ l (n) are the reference node q1 and the monitoring node q l The noise generated when receiving the signal is considered to be zero-mean additive Gaussian noise.

[0015] As a further optimization scheme of the direct positioning method based on cross-spectral subspace orthogonality of the present invention, the step 2) comprises the following specific steps:

[0016] Step 2.1), assume that there is no multipath component in the reference signal, and replace x1(n) and x l (n) can be rewritten as

[0017]

[0018] Where s1(n) = λ 11 s(n-Δτ l ), μ lm =λ lm / λ 11 , Δτ lm =τ lm -τ1 is the time delay difference between other multipath signals and the reference signal, τ1=||q1-p|| / c, c represents the propagation speed of electromagnetic waves;

[0019] Since data segmentation of each received signal does not change the correlation between the received signals, and in order to reduce the amount of subsequent matrix operations, the length of N l The received data is evenly divided into K parts, and the length of each part is N0=N l / K; x1(n), x l The kth segment x of (n) 1k (n k ), x lk (n k ) is expressed as

[0020]

[0021] Among them, s1(n k ),σ1(n k )、s1(n k -Δτ lm ),σ l (n k ) represent s1(n), σ1(n), s1(n-Δτ lm ),σ l paragraph k of (n);

[0022] Step 2.2), let τ be the time interval variable, then receive data x lk (n k ) and reference data x 1k (n k )'s cross-correlation function The definition of

[0023]

[0024] Among them, E(·) represents expectation, (·) H represents the conjugate transpose; let the signal and noise be uncorrelated, Expand, Re-write as:

[0025]

[0026] in, is s1(n k ) and Also s1(n k ) and s1(n k -Δτ lm )’s cross-correlation function;

[0027] Step 2.3), since x 1k (n k )=s1(n k )+σ1(n k),have in and They are s1(n k ), x lk (n k ) and σ 1k (n k )'s power spectrum; Perform discrete Fourier transform and get x lk (n k ) and x 1k (n k )

[0028]

[0029] Among them, is considered as a noise term that follows a Gaussian distribution.

[0030] As a further optimization scheme of the direct positioning method based on cross-spectral subspace orthogonality of the present invention, the step 3) comprises the following steps:

[0031] Step 3.1), select any cross-spectrum from the cross-spectrum obtained in step 2), record the first position corresponding to the cross-spectrum greater than the preset minimum amplitude threshold TH1 as the starting point, and record the last position greater than TH1 as the end point, determine the number of selected points and the interval according to the length between the starting point and the end point, so that the distances between the points contained in the starting point and the end point are equal, and ensure that the distances between the selected points do not exceed the preset maximum distance threshold TH2, thereby obtaining a uniformly selected valid part of the cross-spectrum, and apply the determined starting point, end point, number of selected points and interval under the cross-spectrum to all cross-spectra obtained in step 2), and obtain the valid part y of all cross-spectra. lk (n d ),Right now

[0032]

[0033] Among them, y lk is the observation vector of the kth segment of the lth received data, Indicates the angular frequency corresponding to the subscript, k d is the index of the lower frequency limit, D is the interval of data selection, D>M, N d is the length of the observation vector, 1≤n d ≤N d ;

[0034] y lk The vector form is in:

[0035]

[0036] (·) T represents transpose;

[0037] Step 3.2), solve the covariance matrix of the effective part of each cross spectrum to obtain the covariance matrix R of the effective part of the cross spectrum of the kth segment of the lth received signal and the kth segment of the reference signal lk :

[0038]

[0039] in,

[0040] In order to balance the effect of noise, The average G 11 (ω) Replacement

[0041]

[0042] Because A lk and a lkm Regardless of the number of segments k, they can be rewritten as

[0043]

[0044] Add the covariance matrices of all segments to get their average covariance matrix R l :

[0045]

[0046] When there is no coherent signal, P l is non-singular, R l Perform eigendecomposition and obtain N d eigenvalues ​​and corresponding eigenvectors, and obtain

[0047]

[0048] Among them, R l N d The eigenvalues ​​are arranged from large to small. By R l The signal subspace spanned by the eigenvectors corresponding to the first M eigenvalues ​​of By R l After N d -The noise subspace spanned by the eigenvectors corresponding to the M eigenvalues, They are R l The first M eigenvalues ​​and the last N d -A diagonal matrix of M eigenvalues;

[0049] Step 3.3), in the ideal case without noise, according to the orthogonal relationship between the signal and noise subspaces, the following relationship is derived:

[0050]

[0051] Then there is

[0052]

[0053] Step 3.4), choose to estimate the location of the radiation source through peak search and establish the cost function of any point q0 in the selected area in:

[0054]

[0055] As a further optimization scheme of the direct positioning method based on the orthogonality of the cross-spectrum subspace of the present invention, in the step 3), in order to avoid the coherence of the signal, the cross-spectrum is processed by a spatial smoothing method to obtain a smoothed covariance matrix and then the feature decomposition is performed on it; that is, the cross-spectrum in the frequency domain is spatially smoothed to obtain a smoothed covariance matrix, which is used to replace R l :

[0056] definition for The d-th observation vector of is expressed as:

[0057]

[0058] Define Y lk is the combination of all observation vectors, and its expression is:

[0059]

[0060] Where D is the number of selected observation vectors;

[0061] The smoothed covariance matrix is Performing feature decomposition, we get

[0062]

[0063] in, yes The signal subspace of yes The noise subspace of N d The eigenvalues ​​are arranged from large to small. They are The first M, the last N d -A diagonal matrix consisting of M eigenvalues;

[0064] Use the obtained Replace the The cost function for the case where the signal has coherence is obtained:

[0065] Compared with the prior art, the present invention adopts the above technical solution and has the following technical effects:

[0066] The present invention is based on the TDOA positioning scenario and utilizes the distributed positioning nodes of a single sensor to receive signals, which is beneficial to reducing hardware costs; by utilizing the acquisition of the cross-spectrum between signals, data segmentation, and extraction of the effective part of the cross-spectrum, the dimension of matrix operations can be reduced, thereby achieving the effect of reducing the complexity of calculation and improving the real-time performance of positioning; the present invention is a DPD method that does not require parameter estimation and is applicable to both multipath environments and non-multipath environments. Compared with other DPD methods and TDOA methods, the present invention has higher positioning accuracy and more stable positioning performance. BRIEF DESCRIPTION OF THE DRAWINGS

[0067] Figure 1 A flow chart of a direct positioning method based on cross-spectral subspace orthogonality provided by the present invention;

[0068] Figure 2 A schematic diagram of a positioning scenario of the method of the present invention;

[0069] Figure 3 Schematic diagram of the spatial smoothing process in step 3;

[0070] Figure 4 CDF diagram of positioning error of the method of the present invention and other positioning methods for simulation data under different multipath numbers;

[0071] Figure 5 The measured data scene graph and the thermal map of the positioning area provided by the present invention;

[0072] Figure 6 This is a CDF diagram of positioning errors of the method of the present invention and other positioning methods for measured data. DETAILED DESCRIPTION

[0073] The technical solution of the present invention is further described in detail below in conjunction with the accompanying drawings:

[0074] The present invention can be implemented in many different forms and should not be considered to be limited to the embodiments described herein. On the contrary, these embodiments are provided to make this disclosure thorough and complete, and will fully express the scope of the present invention to those skilled in the art. In the accompanying drawings, components are enlarged for clarity.

[0075] The detailed process of the direct positioning method based on cross-spectral subspace orthogonality provided by the present invention is as follows: Figure 1 As shown in the figure: First, select the reference node, set its received signal as the reference signal, and segment all received signals; second, calculate the cross-spectrum between the segmented received signal and the reference signal, and perform eigenvalue decomposition on the covariance matrix of the effective part of the obtained cross-spectrum; then, use the orthogonality of the subspace to construct the cost function under the cross-spectrum; finally, determine the location of the radiation source by searching for the superposition of these cost functions on the grid points in the target area. The specific implementation is as follows:

[0076] Step 1: The positioning node receives the signal from the radiation source and transmits it to the central site:

[0077] like Figure 2 As shown in the figure, assume that there is a radiation source p with an unknown position, the number of distributed monitoring nodes is L (L ≥ 3), let the coordinates of the reference node be q1, and the coordinates of the non-reference monitoring node be q l , l = 2, ..., L, the discrete form of the received signal of each monitoring node in the multipath condition is:

[0078]

[0079] Where x1(n) represents the discrete form of the received signal of the reference node in the case of multipath, x l (n) represents the discrete form of the received signal of the non-reference node in the case of multipath, M represents the number of multipaths, and N l represents the length of the received signal, λ 1m and λ lm Indicates x1(n) and x in a multipath environment l (n) amplitude coefficient; τ 1m is the time delay τ between the mth multipath component in the reference signal and the transmitted signal lm is the time delay between the mth multipath component in the lth received signal and the transmitted signal; s(n) represents the transmitted signal of the unknown radiation source, s(n-τ 1m ) represents the time delay τ 1m The emission signal of the unknown radiation source, s(n-τ lm ) represents the time delay τ lm The emission signal of the unknown radiation source; σ1(n) and σ l (n) is the noise generated when the corresponding node receives the signal, which is considered to be zero-mean additive Gaussian noise.

[0080] Step 2: Evenly segment the L received signals, and divide each received signal into K segments; use the received signal of the reference node as the reference signal; for each of the other L-1 received signals, perform cross-correlation between its k-th segment received data and the k-th segment received data of the reference signal, k = 1, ..., K, and obtain (L-1) × K cross-correlation functions; perform discrete Fourier transform on the obtained cross-correlation functions in turn to obtain their corresponding cross-spectra;

[0081] Assume that there is no multipath component in the reference signal, and replace x1(n) and x l (n) can be rewritten as

[0082]

[0083] Where s1(n) = λ 11 s(n-Δτ l ), μ lm =λ lm / λ 11 , Δτ lm =τ lm -τ1,l=2,...,L is the time delay difference between other multipath signals and the reference signal, τ1=||q1-p|| / c, c represents the propagation speed of electromagnetic waves, s1(n-Δτ lm ) represents the time delay Δτ lm The emission signal of the unknown radiation source;

[0084] Since data segmentation of each received signal does not change the correlation between the received signals, and in order to reduce the amount of subsequent matrix operations, the length of N l The received data is evenly divided into K parts, and the length of each part is N0=N l / K; kth segment x of reference signal and non-reference signal 1k (n k ) and x lk (n k ),l=2,...,L;k=1,...,K is expressed as

[0085]

[0086] Among them, s1(n k ),σ1(n k ), s1(n k -Δτ lm ) and σ l (n k ) represent s1(n), σ1(n), s1(n-Δτ lm ) and σ l paragraph k of (n);

[0087] Let τ be the time interval variable, then the received data x lk (n k ) and reference data x 1k (n k )'s cross-correlation function The definition of

[0088]

[0089] Among them, E(·) represents expectation, (·) H represents the conjugate transpose. Assuming that the signal and noise are uncorrelated, Expand, Re-write as:

[0090]

[0091] in, is s1(n k ) and Also s1(n k ) and s1(n k -Δτ lm ) is the cross-correlation function. Since x 1k (n k )=s1(n k )+σ1(n k ),have in and They are s1(n k ), x lk (n k ) and σ 1k (n k ) is the power spectrum of Perform discrete Fourier transform and get x lk (n k ) and x 1k (n k )

[0092]

[0093] in, is considered as a noise term that follows a Gaussian distribution.

[0094] Step 3: Perform eigendecomposition on the covariance matrix of the effective part of the cross-spectrum to obtain its noise subspace. In order to avoid the coherence of the signal, the spatial smoothing method is used to process the data. According to the orthogonality of the noise subspace and the signal subspace, the cost function under the cross-spectrum is constructed:

[0095] When there is no coherence between each signal segment, select any cross-spectrum from the cross-spectrum obtained in step 2, record the first position corresponding to the cross-spectrum greater than the preset minimum amplitude threshold TH1 as the starting point, and record the last position greater than TH1 as the end point. Determine the number of selected points and intervals according to the length between the starting point and the end point, so that the distances between the points contained in the starting point and the end point are equal, and ensure that the distances between the selected points do not exceed the preset maximum distance threshold TH2, thereby obtaining a uniformly selected valid part of the cross-spectrum, and apply the determined starting point, end point, number of selected points and interval under the cross-spectrum to all cross-spectra obtained in step 2 to obtain the valid part y of all cross-spectra. lk (n d ):

[0096]

[0097] where y lk is the observation vector of the kth segment of the lth received data, Indicates the angular frequency corresponding to the subscript, k d is the index of the lower frequency limit, D (D>M) is the interval of data selection, N d is the length of the observation vector. lk The vector form is

[0098]

[0099] in,

[0100]

[0101] in(·) T represents the transpose. Therefore, we get y lk The covariance matrix R lk :

[0102]

[0103] in, It should be noted that data segmentation does not change the data correlation, which means that the cross-spectrum of each data segment and the reference data should be consistent in theory. 1k (n k ) should also be consistent across all data segments. However, additive noise can lead to differences in the cross-spectrum of each segment. To balance the effect of noise, The average G 11 (ω) Replacement

[0104]

[0105] At this time A lkand a lkm Regardless of the number of segments k, they can be rewritten as

[0106]

[0107] Add the covariance matrices of all segments to get their average covariance matrix R l :

[0108]

[0109] When there is no coherent signal, P l is non-singular, obtain N d eigenvalues ​​and corresponding eigenvectors, and obtain

[0110]

[0111] Among them, N d The eigenvalues ​​are arranged from large to small. is the signal subspace spanned by the eigenvectors corresponding to the first M eigenvalues, It is after N d -The noise subspace spanned by the eigenvectors corresponding to the M eigenvalues, and The first M eigenvalues ​​and the last N d -A diagonal matrix composed of M eigenvalues. According to the orthogonal relationship between the signal and noise subspaces, the following relationship is derived:

[0112]

[0113] Then there is

[0114]

[0115] In fact, due to the possible signal coherence, P l It is not necessarily non-singular. This problem can be solved by obtaining a smooth covariance matrix, that is, spatially smoothing the cross spectrum in the frequency domain to obtain a smooth covariance matrix to replace R l The specific forward space smoothing process is as follows: Figure 3 As shown, yes The d-th observation vector of is expressed as:

[0116]

[0117] Define Y lk For all combinations of observation vectors:

[0118]

[0119] The smoothed covariance matrix is Performing feature decomposition, we get

[0120]

[0121] in, yes The signal subspace of yes The noise subspace of N d The eigenvalues ​​are arranged from large to small. They are The first M, the last N d -A diagonal matrix consisting of M eigenvalues.

[0122] Considering the existence of noise and interference, and There is an approximate orthogonal relationship between them, namely And the point corresponding to the radiation source The result is closest to 0. Therefore, the cost function of any point q0 in the selected area is established

[0123] in,

[0124]

[0125] Due to the inversion process, the cost function of the point corresponding to the position of the radiation source is the largest among all cost function values.

[0126] Step 4: Select the positioning area and grid it, perform a spectrum peak search on the cost function value at each grid point, and obtain the estimated results of the radiation source coordinates:

[0127] The positioning area and search step are determined, and the area is gridded. The peak value of the cost function value on all grid points is searched, and the coordinates corresponding to the peak value are the estimated results of the radiation source.

[0128] Figure 4 This is a CDF diagram of the positioning error of the method of the present invention and other positioning methods for simulation data under different multipath numbers, where the number of nodes is L=4, the simulation signal is a 16QAM modulated signal, and its sampling rate is f s =125MHz, bandwidth is B=40MHz, data length is N l =12288, the length of each segment after segmentation is N0=256, the number of segments K=48, and the length of the effective part of the cross spectrum taken in each segment is N d=20, the interval is D=6, the number of one-dimensional searches is T=100, the signal-to-noise ratio SNR=6dB, and P=100 experiments are performed. The specific comparison algorithms are the TDOA algorithm based on generalized cross-correlation (TDOA-GCC), the TDOA algorithm based on normalized cross-spectrum (TDOA-CS), and the direct positioning method based on determinant (DPD-DET). Figure 4 (a) is the CDF curve comparison without multipath (M=1), Figure 4 (b) is a CDF comparison diagram in the presence of multipath (M=3). It can be seen from the figure that in multipath and non-multipath environments, the error CDF curve of the method of the present invention is closest to the vertical axis. In addition, the other three algorithms will produce large errors, which often means positioning failure, while the errors of the method of the present invention are relatively small. By comparison Figure 4 (a) and Figure 4 (b), it can be found that multipath has almost no effect on the proposed method, while other methods are affected to varying degrees. Therefore, in both multipath and non-multipath environments, the positioning accuracy of the proposed method is better than that of other methods.

[0129] Figure 5 The measured data scene diagram and the thermal map of the positioning area provided by the present invention, wherein Figure 5 (a) Figure 5 (b) are the physical pictures of the radiation source and four monitoring nodes. Figure 5 (c) shows the distribution of nodes and radiation sources and the heat map using node positioning. The signal received by the monitoring node from the radiation source is the source of the measured data of the present invention. The node receives the 16QAM modulated signal from the radiation source, which has a length of 30720, a bandwidth of 20MHz, a center frequency of 720MHz, and a sampling rate of f s =125MHz, take N0=2048, K=15, T=400. Figure 5 In (c), the irregular shaded part represents the estimated radiation source position and the surrounding situation. It can be seen that the estimated position of the radiation source with the darkest color is very close to the actual radiation source position, which verifies the effectiveness of the method of the present invention.

[0130] Figure 6 This is a CDF diagram of the positioning error of the method described in the present invention and other positioning methods for measured data. Figure 5 The same method was used and P = 20 experiments were performed. Figure 6 As shown in Figure 3, the overall trend of the CDF curve is consistent with the simulation test results, which proves the superiority of the proposed method in terms of estimation accuracy and practicality.

[0131] It will be understood by those skilled in the art that, unless otherwise defined, all terms (including technical and scientific terms) used herein have the same meaning as generally understood by those skilled in the art in the art to which the present invention belongs. It should also be understood that terms such as those defined in common dictionaries should be understood to have meanings consistent with the meanings in the context of the prior art, and will not be interpreted with idealized or overly formal meanings unless defined as herein.

[0132] The specific implementation methods described above further illustrate the objectives, technical solutions and beneficial effects of the present invention in detail. It should be understood that the above description is only a specific implementation method of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the scope of protection of the present invention.

Claims

1. A direct positioning method based on the orthogonality of cross-spectral subspaces, It is characterized in that The following steps are involved: Step 1), receive signals from a radiation source p at an unknown location at L monitoring nodes with known locations, transmit the L received signals to a central site, L ≥ 3, and select one of the monitoring nodes as a reference node; Let the coordinates of the reference node be q 1 , and the coordinates of the other L-1 monitoring nodes are q l , l = 2, ..., L, the discrete form of the received signal of each monitoring node in the multipath condition is: Among them, x 1 (n) represents the discrete form of the received signal of the reference node in the case of multipath, x l (n) represents the discrete form of the received signal of the non-reference node in the case of multipath, M represents the number of multipaths, and N l represents the length of the received signal, λ 1m , lm They represent x in multipath environment. 1 (n), x l (n) amplitude coefficient; τ 1m is the time delay between the mth multipath component in the reference signal and the transmitted signal; τ lm is the time delay between the mth multipath component in the lth received signal and the transmitted signal; s(n) represents the transmitted signal of the unknown radiation source, s(n-τ 1m ) represents the time delay τ 1m The emission signal of the unknown radiation source, s(n-τ lm ) represents the time delay τ lm The emission signal of the unknown radiation source; σ 1 (n), σ l (n) are the reference nodes q 1 , monitoring node q l The noise generated when receiving the signal is considered to be zero-mean additive Gaussian noise; Step 2), evenly segment the L received signals, and divide each received signal into K segments; use the received signal of the reference node as the reference signal; for each of the other L-1 received signals, perform cross-correlation between the k-th segment received data and the k-th segment received data of the reference signal, k = 1, ..., K, and obtain (L-1) × K cross-correlation functions; perform discrete Fourier transform on the obtained cross-correlation functions in turn to obtain their corresponding cross-spectra; Step 2.1), assume that there is no multipath component in the reference signal, and set x 1 (n) and x l (n) can be rewritten as Among them, s 1 (n) = λ 11 s(n-Δτ l ), μ lm =λ lm / λ 11 , Δτ lm =τ lm -τ 1 is the time delay difference between other multipath signals and the reference signal, τ 1 =||q 1 -p|| / c, c represents the propagation speed of electromagnetic waves, s 1 (n-Δτ lm ) represents the time delay Δτ lm The emission signal of the unknown radiation source; Since data segmentation of each received signal does not change the correlation between the received signals, and in order to reduce the amount of subsequent matrix operations, the length of N l The received data is evenly divided into K parts, each of which has a length of N 0 =N l / K;x 1 (n), x l The kth segment x of (n) 1k (n k ), x lk (n k ) is expressed as Among them, s 1 (n k )、σ 1 (n k )、s 1 (n k -Δτ lm )、σ l (n k ) respectively represent the k-th segment of s 1 (n), σ 1 (n), s 1 (n - Δτ lm ), σ l (n); Step 2.2), let τ be the time interval variable, then receive data x lk (n k ) and reference data x 1k (n k )'s cross-correlation function The definition of Among them, E(·) represents expectation, (·) H represents the conjugate transpose; let the signal and noise be uncorrelated, Expand, Re-write as: in, Yes 1 (n k ) is the autocorrelation function, and Also 1 (n k ) and s 1 (n k -Δτ lm )’s cross-correlation function; Step 2.3), since x 1k (n k )=s 1 (n k )+σ 1 (n k ),have in and They are 1 (n k ), x lk (n k ) and σ 1k (n k )'s power spectrum; Perform discrete Fourier transform and get x lk (n k ) and x 1k (n k ) Among them, is regarded as a noise term that follows a Gaussian distribution; Step 3), eigendecomposing the covariance matrix of the effective part of the cross spectrum to obtain its noise subspace and signal subspace, and constructing the cost function under the cross spectrum according to the orthogonality of the noise subspace and the signal subspace; Step 4), determine the positioning area and the search step, grid the positioning area, perform peak search on the cost function values ​​at all grid points, and the coordinates corresponding to the peak are the estimation results for the radiation source.

2. The direct positioning method based on cross-spectral subspace orthogonality according to claim 1, It is characterized in that The step 3) comprises the following steps: Step 3.1), select one of the cross-spectra obtained in step 2), and select the first cross-spectrum whose amplitude is greater than the preset minimum amplitude threshold TH 1 The corresponding position is recorded as the starting point, and the last one greater than TH 1 The position is recorded as the end point, and the number of points and the interval selected are determined according to the length between the start point and the end point, so that the distances between the points included in the start point and the end point are equal, and at the same time, the distances between the selected points do not exceed the preset maximum distance threshold TH 2 , thereby obtaining the effective part of the cross-spectrum uniformly selected, and applying the determined starting point, end point, number of selected points and interval under the cross-spectrum to all cross-spectra obtained in step 2) to obtain the effective part y of all cross-spectra lk (n d ),Right now Among them, y lk is the observation vector of the kth segment of the lth received data, Indicates the angular frequency corresponding to the subscript, k d is the index of the lower frequency limit, D is the interval of data selection, D>M, N d is the length of the observation vector, 1≤n d ≤N d ; y lk The vector form is in: (·) T represents transpose; Step 3.2), solve the covariance matrix of the effective part of each cross spectrum to obtain the covariance matrix R of the effective part of the cross spectrum of the kth segment of the lth received signal and the kth segment of the reference signal lk : in, In order to balance the effect of noise, Average G 11 (ω) Replacement Because A lk and a lkm Regardless of the number of segments k, they can be rewritten as Add the covariance matrices of all segments to get their average covariance matrix R l : When there is no coherent signal, P l is non-singular, R l Perform eigendecomposition and obtain N d eigenvalues ​​and corresponding eigenvectors, and obtain Among them, R l N d The eigenvalues ​​are arranged from large to small. By R l The signal subspace spanned by the eigenvectors corresponding to the first M eigenvalues ​​of By R l After N d -The noise subspace spanned by the eigenvectors corresponding to the M eigenvalues, They are R l The first M eigenvalues ​​and the last N d -A diagonal matrix of M eigenvalues; Step 3.3), in the ideal case without noise, according to the orthogonal relationship between the signal and noise subspaces, the following relationship is derived: Then there is Step 3.4), choose to estimate the radiation source position through peak search, and establish the position of any point q in the selected area. 0 The cost function in:

3. The direct positioning method based on cross-spectral subspace orthogonality according to claim 2, It is characterized in that In step 3), in order to avoid the coherence of the signal, the cross spectrum is processed by spatial smoothing to obtain a smoothed covariance matrix and then the feature decomposition is performed on it; that is, the cross spectrum in the frequency domain is spatially smoothed to obtain a smoothed covariance matrix, which is used to replace R l : definition for The d-th observation vector of is expressed as: Define Y lk is the combination of all observation vectors, and its expression is: Where D is the number of selected observation vectors; The smoothed covariance matrix is Performing feature decomposition, we get in, yes The signal subspace of yes The noise subspace of N d The eigenvalues ​​are arranged from large to small. They are The first M, the last N d -A diagonal matrix consisting of M eigenvalues; Use the obtained Replace the The cost function for the case where the signal has coherence is obtained:

Images