WO2022166477A1 - 定位方法、装置、基站、计算机设备和存储介质 - Google Patents

定位方法、装置、基站、计算机设备和存储介质 Download PDF

Info

Publication number
WO2022166477A1
WO2022166477A1 PCT/CN2021/142371 CN2021142371W WO2022166477A1 WO 2022166477 A1 WO2022166477 A1 WO 2022166477A1 CN 2021142371 W CN2021142371 W CN 2021142371W WO 2022166477 A1 WO2022166477 A1 WO 2022166477A1
Authority
WO
WIPO (PCT)
Prior art keywords
azimuth
delay
base station
frequency domain
vector
Prior art date
Application number
PCT/CN2021/142371
Other languages
English (en)
French (fr)
Inventor
潘孟冠
齐望东
刘升恒
黄永明
尤肖虎
王绍磊
徐佳
Original Assignee
网络通信与安全紫金山实验室
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Priority claimed from CN202110146480.2A external-priority patent/CN112929962B/zh
Priority claimed from CN202110144931.9A external-priority patent/CN112469119B/zh
Application filed by 网络通信与安全紫金山实验室 filed Critical 网络通信与安全紫金山实验室
Publication of WO2022166477A1 publication Critical patent/WO2022166477A1/zh

Links

Images

Classifications

    • H—ELECTRICITY
    • H04—ELECTRIC COMMUNICATION TECHNIQUE
    • H04W—WIRELESS COMMUNICATION NETWORKS
    • H04W64/00—Locating users or terminals or network equipment for network management purposes, e.g. mobility management

Definitions

  • the present application relates to the field of positioning technology, and in particular, to a positioning method, apparatus, base station, computer equipment and storage medium.
  • Satellite navigation and positioning technology has the advantages of wide-area coverage and good ubiquity, but due to low signal power and weak penetration, it is mainly used for terminal positioning in outdoor open environments, and cannot provide navigation and positioning in sheltered and indoor environments. Serve.
  • the 5G cellular mobile network uses key technologies including massive MIMO, ultra-dense networking, and large-bandwidth signals. With the massive deployment of multi-antenna 5G base stations in the future, the use of 5G signals to achieve high-precision positioning will have very broad prospects for development.
  • a positioning method comprising: receiving uplink frequency domain sounding reference signals from a terminal through multiple channels of a base station; determining a channel frequency domain response vector according to the uplink frequency domain sounding reference signals of the multiple channels; The channel frequency domain response vector is converted into an overcomplete response vector representing the channel frequency domain response at a plurality of time delay-azimuth grid points in the signal range-azimuth domain of the base station; As the observation vector, the delay-azimuth two-dimensional spectral vector is determined by the signal amplitude at each of the delay-azimuth grid points to be solved, and the delay-azimuth two-dimensional spectral vector is established according to the observation vector.
  • a solution equation for solving with a dimensional spectral vector iteratively estimating the solution equation to determine a delay-azimuth spectrum formed by the signal amplitude values at each of the delay-azimuth grid points; and According to the time delay-azimuth spectrum, the position of the terminal is determined by solution.
  • a positioning apparatus includes: a sounding reference signal receiving module for receiving uplink frequency domain sounding reference signals from a terminal through multiple channels of a base station; a channel frequency domain response vector determining module for receiving according to all channels of the multiple channels.
  • the uplink frequency domain sounding reference signal is used to determine the channel frequency domain response vector;
  • the overcomplete response vector determination module is used to convert the channel frequency domain response vector into a plurality of Delay-overcomplete response vector of the frequency domain response of the channel at the azimuth grid point;
  • the solution equation establishment module is used to use the overcomplete response vector as the observation vector, and use the overcomplete response vector to be solved for each delay-
  • the signal amplitude at the azimuth grid point determines the delay-azimuth two-dimensional spectral vector, and establishes a solution equation for solving the delay-azimuth two-dimensional spectral vector according to the observation vector; delay-azimuth an angle spectrum determination module for iteratively estimating the solution equation to determine a delay-azimuth
  • a base station includes a baseband processing unit and one or more antenna processing units connected to the baseband processing unit, each of the antenna processing units includes at least one antenna array, and each of the antenna arrays includes a plurality of arrays arranged
  • the multiple array elements included in each of the antenna arrays provide corresponding multiple receiving channels; wherein, the baseband processing unit uses the multiple receiving channels provided by the antenna arrays of the respective antenna processing units The channel receives the uplink frequency domain sounding reference signal from the terminal, and executes the above positioning method to realize the positioning of the terminal.
  • a computer device includes a memory and a processor, the memory stores a computer program, and the processor implements the above positioning method when the computer program is executed.
  • Fig. 1 is the application environment diagram of the positioning method in one embodiment
  • FIG. 2 is a schematic flowchart of a positioning method in one embodiment
  • FIG. 3 is a schematic flowchart of a step of determining a channel frequency-domain response vector in one embodiment
  • Fig. 4 is a schematic flowchart of steps of converting to obtain an overcomplete response vector in one embodiment
  • FIG. 5 is a schematic flowchart of a step of determining the time delay-azimuth spectrum in one embodiment
  • FIG. 6 is a schematic flowchart of a step of determining a time delay-azimuth spectrum in another embodiment
  • FIG. 7 is a schematic flowchart of steps of determining the location of a terminal in one embodiment
  • FIG. 8 is a schematic flowchart of steps of determining the location of a terminal in one embodiment
  • 10 is a cross-sectional view of the azimuth dimension of the distance-azimuth spectrum determined by the two-dimensional space-frequency MUSIC algorithm in one embodiment
  • 11 is a cross-sectional view of the distance dimension of the distance-azimuth spectrum determined by the two-dimensional space-frequency MUSIC algorithm in one embodiment
  • 13 is a cross-sectional view of the azimuth dimension of the distance-azimuth spectrum determined by the two-dimensional space-frequency joint IAA-APES algorithm in one embodiment
  • 15 is a range-azimuth spectrum diagram determined by a two-dimensional SLIM algorithm in one embodiment
  • 16 is an azimuth estimation RMSE diagram of the two-dimensional space-frequency MUSIC algorithm and the two-dimensional SLIM algorithm in one embodiment
  • 17 is a distance estimation RMSE diagram of the two-dimensional space-frequency MUSIC algorithm and the two-dimensional SLIM algorithm in one embodiment
  • FIG. 18 is a structural block diagram of a positioning apparatus in an embodiment.
  • the positioning method provided by the present application can be applied to the application environment shown in FIG. 1 .
  • the positioning server 102 is in communication connection with one or more base stations 104 .
  • Each base station 104 may have one or more antenna arrays.
  • the base station 104 receives the uplink frequency domain Sounding Reference Signal (SRS) of the terminal 106 through multiple receiving channels of each antenna array of the base station 104, and executes this
  • SRS Sounding Reference Signal
  • the positioning methods of various embodiments are applied to locate the terminal 106 .
  • the terminal 106 may be, but is not limited to, various personal computers, notebook computers, smart phones, tablet computers, portable wearable devices, and the like.
  • the base station 104 may be a 5G base station, or may be any other type of base station suitable for positioning by the method of the present application, which is not limited in the present application.
  • the base station 104 can be a small base station with a single antenna array, or a base station with a plurality of distributed antenna arrays.
  • a group of receiving channels corresponding to the antenna array receives a group of uplink frequency-domain sounding reference signals from the terminal 106, and respectively performs the positioning method of the embodiment of the present application for each group of uplink frequency-domain sounding reference signals to determine a pair of azimuths of the direct path Then the base station 104 can send the azimuth and delay estimates corresponding to each direct path determined by the base station for each antenna array to the positioning server 102.
  • Multiple base stations 104 receive in real time the azimuth angle estimates and delay estimates corresponding to one or more direct paths to the same terminal 106, and calculate and determine based on the received azimuth estimates and time delay estimates corresponding to these direct paths. The location of the terminal 106 .
  • a positioning method is provided, and the method is applied to the base station 104 in FIG. 1 as an example for description, including the following steps S210-S260.
  • Step S210 Receive uplink frequency domain sounding reference signals from the terminal through multiple channels of the base station.
  • the base station 102 receives the uplink frequency domain SRS sent by the terminal through the multiple receiving channels of the antenna array of the base station 102. If the antenna array of the base station 102 has a total of N array elements, then the antenna array has a total of N elements.
  • the receiving channel, the number of subcarriers occupied by the SRS is M 0 , the uplink frequency domain SRS received from each receiving channel can be expressed as a vector where X m,n represents the frequency domain SRS received by the nth receiving channel and the mth subcarrier.
  • Step S220 Determine the channel frequency domain response vector according to the uplink frequency domain sounding reference signals of the multiple channels.
  • step S220 may include steps S221, S223 and S225.
  • Step S221 Determine the received signal matrix of the multiple channels according to the uplink frequency domain sounding reference signals of the multiple channels.
  • the base station 102 can express the received signal matrix of all channels as: in represents the complex space, Indicates that the received signal matrix X is an M 0 ⁇ N-dimensional complex number matrix.
  • step S221 and before step S223, the method further includes:
  • Step S222 performing antenna array correction on the received signal matrices of multiple channels
  • the channel amplitude and phase errors and the antenna amplitude and phase errors of each element of the antenna array on different subcarriers can be measured in advance, and the sum is ⁇ m,n is the total amplitude and phase error of the mth subcarrier of the nth channel, and then use these channel amplitude and phase errors and antenna amplitude and phase errors to correct the original received signal matrix, and the corrected received signal matrix X is:
  • the antenna array of the signal is an equidistant linear array (Uniform Linear Array ULA) as an example, and the element spacing of the equidistant linear array is d.
  • the delay, azimuth, and received signal amplitude of the kth path are: and in, Defined as the angle between the signal incident direction and the ULA normal direction.
  • the time delay of signal transmission may represent the distance of the signal transmission, and the time delay and the distance may be converted into each other through operations.
  • the received signal matrix X of multiple channels can be expressed as:
  • the diag( ⁇ ) operator means that each element of the vector is used as the main diagonal element to obtain the diagonal matrix.
  • is the matching vector for the signal delay domain is a matching vector function representing the delay domain of an Orthogonal Frequency Division Multiplexing (OFDM) signal, whose input is the path delay ⁇ , and T is the space composed of all possible delays ⁇ , namely in Represents the real space, and the output is an M-dimensional vector.
  • OFDM Orthogonal Frequency Division Multiplexing
  • nth element is: is a noise matrix, and the elements of the mth row and the nth column represent the noise components on the mth subcarrier and the nth receiving channel.
  • Step S223 Perform channel estimation according to the received signal matrix to obtain a channel frequency domain response matrix.
  • the channel estimation module of the base station 104 performs channel estimation according to the received signal matrix formed by the frequency domain SRS received by the multiple receiving channels, and obtains the channel frequency domain response matrix. It is assumed that the receiver uses a known SRS sequence. Perform channel estimation, then the channel frequency domain response matrix obtained by performing multi-channel frequency domain channel estimation can be expressed as:
  • the method may further include: step S224, performing dimension reduction processing on the subcarrier dimension of the channel frequency domain response matrix.
  • the base station 104 estimates the channel The sub-carrier dimension of the matrix is dimensionally reduced, that is, the sub-carrier dimension of the channel estimation matrix is extracted, and the extraction rate is v, and the channel frequency domain response matrix obtained after extraction is in, Indicates the round-down operator, then the element in the mth row of the H matrix can be expressed as:
  • H(m,:) represents all elements of the mth row of matrix H
  • H 0 (vm,:) represents all elements of the vm row of matrix H 0 .
  • Step S225 vectorize the channel frequency domain response matrix to obtain the channel frequency domain response vector.
  • the base station 104 performs vectorization on the channel frequency domain response matrix obtained in the preceding steps to obtain a channel frequency domain response vector
  • vec( ) represents the matrix vectorization operator, there are:
  • T and ⁇ are the space formed by all the delay and azimuth angles, respectively, and have:
  • Step S230 Convert the channel frequency domain response vector into an overcomplete response vector representing the channel frequency domain response at multiple delay-azimuth grid points in the signal range-azimuth domain of the base station.
  • step S230 may include steps S231-S234.
  • Step S231 obtaining the range of the working distance and the range of the receiving azimuth of the antenna array of the base station;
  • acquiring the operating distance range and the receiving azimuth angle range of the antenna array of the base station includes: determining the receiving azimuth angle range of the antenna array according to the structure information of the antenna array of the base station and the orientation information of the antenna array; The power and the sensitivity of the antenna array determine the maximum operating distance of the antenna array, and determine the operating distance range of the antenna array according to the maximum operating distance.
  • the base station 104 determines the maximum operating range of the base station antenna array according to indicators such as the transmit power of the terminal and the sensitivity of the base station receiver, thereby determining the maximum delay ⁇ max that the base station antenna array can receive the SRS, the maximum The time delay ⁇ max can represent the maximum operating distance of the antenna array.
  • the time delay range that the antenna array can receive SRS is [0, ⁇ max ]
  • the time delay range of [0, ⁇ max ] can represent the range of the antenna array.
  • Acting distance range According to the antenna array structure of the base station and the main lobe width of the array element pattern, determine the receiving azimuth range of the base station antenna array [ ⁇ min , ⁇ max ].
  • Step S232 Determine the working distance-azimuth angle domain of the antenna array based on the working distance range and the receiving azimuth angle range.
  • Step S233 using a uniform grid to divide the range-azimuth domain to determine a plurality of time-delay-azimuth grid points uniformly distributed on the range-azimuth domain.
  • the base station 104 uses a uniform grid Divide the time delay range [0, ⁇ max ] and the receiving azimuth angle range [ ⁇ min , ⁇ max ], so as to determine the time delay range [0, ⁇ max ] and the receiving azimuth angle range [ ⁇ min , ⁇ max ] are jointly determined
  • each delay-azimuth angle grid point may correspond to a delay-azimuth angle coordinate pair, for example, the delay-azimuth angle coordinate pair may be the corresponding delay-azimuth angle grid point. Coordinate pair formed by delay and azimuth.
  • Each delay-azimuth grid point can also be represented by a subscript, where p represents the grid point number in the delay domain, P represents the total number of grid points in the delay domain, and ⁇ p represents the p-th delay domain grid
  • the delay corresponding to the grid point q represents the grid point number in the azimuth domain
  • Q represents the total number of grid points in the azimuth domain
  • ⁇ q represents the azimuth angle corresponding to the qth azimuth domain grid point.
  • the delay-azimuth grid point (p,q) corresponds to the delay-azimuth coordinate pair ( ⁇ p , ⁇ q ).
  • Step S234 Convert the channel frequency domain response matrix into an overcomplete response vector of the channel frequency domain response at multiple delay-azimuth grid points.
  • the channel frequency domain response matrix can be transformed into an overcomplete response vector of the channel frequency domain response at all delay-azimuth grid points in the range-azimuth domain.
  • the number in , q k indicates that the azimuth of the kth path corresponds to the grid set in the azimuth domain
  • a ⁇ . ⁇ ( ⁇ p , ⁇ q ) represents the space-time delay domain two-dimensional matching vector when the time delay is ⁇ p and the azimuth angle is ⁇ q , in, ⁇ p is an alternative input variable, the specific form of a ⁇ ( ⁇ p ) can be similarly referred to The specific form of ; ⁇ q , As an alternative input variable, the specific form of a ⁇ ( ⁇ q ) can be similarly referred to in specific form.
  • step S240 the overcomplete response vector is used as the observation vector, and the delay-azimuth two-dimensional spectral vector is determined by the signal amplitude at each delay-azimuth grid point to be solved, and the delay-azimuth pair is established according to the observation vector.
  • the solution equation for solving the angular 2D spectral vector is used as the observation vector, and the delay-azimuth two-dimensional spectral vector is determined by the signal amplitude at each delay-azimuth grid point to be solved, and the delay-azimuth pair is established according to the observation vector.
  • the base station 104 has transformed the azimuth and distance (delay) estimation problems into The reconstruction problem of solving the delay-azimuth two-dimensional spectral vector ⁇ for the observation vector.
  • Step S250 iteratively estimate the solution equation to determine the delay-azimuth spectrum formed by the signal amplitude values at each delay-azimuth grid point.
  • step S250 A variety of different algorithms can be used to implement the above step S250.
  • two example algorithms an iterative adaptive amplitude-phase estimation algorithm and an iterative minimization sparse learning algorithm, will be provided to perform the solution equation. Iterative estimation to determine the delay-azimuth spectrum.
  • step S250 includes: using an iteratively adaptive amplitude-phase estimation algorithm to iteratively estimate the solution equation to determine a delay- Azimuth spectrum.
  • the base station 104 determines the objective function according to the solution equation, and uses the Iterative Adaptive Approach for Amplitude and Phase Estimation (IAA-APES) algorithm to This objective function is solved to obtain a solution to the above solution equation.
  • IAA-APES Iterative Adaptive Approach for Amplitude and Phase Estimation
  • 2 represents the l 2 norm of the vector.
  • R p,q represents the interference covariance matrix at the (p,q)th grid point, where the interference consists of signal components other than the current grid point (p,q), which can be expressed as:
  • APES algorithm is a kind of spectral estimation algorithm based on adaptive narrowband filter bank.
  • the adaptive filter coefficients are calculated based on the minimum variance distortion-less response (MVDR) criterion, which can guarantee On the premise that the energy of the current spectral center is not lost, the energy of other spectral positions is minimized. Therefore, on the one hand, the APES algorithm is more accurate in the estimation of the amplitude and phase of the real spectral peak position.
  • the spectral energy of the APES algorithm is relatively sharp and has a certain super-resolution ability.
  • the adaptive filter coefficients under the adaptive beamforming (Minimum variance distortionless response, MVDR) algorithm criterion need to be calculated according to the inverse of the covariance matrix.
  • the traditional APES algorithm uses multi-block beats or the moving average of molecular arrays to obtain more accurate The covariance matrix estimate of .
  • the IAA-APES algorithm is applied to the air-frequency domain joint processing of the uplink SRS, and the covariance matrix and the delay-azimuth spectrum are iteratively solved. .
  • step S250 may include steps S251 to S253.
  • Step S251 based on the overcomplete response vector, determine a matching matrix, perform two-dimensional space-frequency matched filtering, and obtain an initial estimated value of the delay-azimuth two-dimensional spectral vector;
  • the base station 104 initializes the algorithm.
  • a ⁇ and ⁇ in the solution equation determined based on the overcomplete response vector are used as matching matrices, and two-dimensional space-frequency matched filtering is performed to obtain the time delay.
  • Two-dimensional space-frequency matched filtering is performed to obtain the time delay.
  • M represents the total number of sub-carriers after decimation
  • N represents the number of elements of the antenna array
  • ( ⁇ ) H represents the operator that takes the conjugate transpose of the matrix or vector.
  • Step S252 based on the initial estimated value of the delay-azimuth two-dimensional spectral vector, iteratively update the estimated value of the power matrix, the estimated value of the covariance matrix, and the estimated value of the delay-azimuth two-dimensional spectral vector, until the current number of iterations When the change degree value between the estimated value of the delay-azimuth two-dimensional spectral vector and the estimated value of the delay-azimuth two-dimensional spectral vector of the last iteration number is less than a predetermined threshold, the iterative update is stopped;
  • the base station 104 uses the IAA-APES algorithm to iteratively update the estimation of the power matrix based on the initial estimated value ⁇ (0) of the delay-azimuth two-dimensional spectral vector ⁇ obtained by the two-dimensional space-frequency matched filtering. value, an estimate of the covariance matrix, and an estimate of the delay-azimuth two-dimensional spectral vector ⁇ .
  • i denote the ith iteration
  • P(i), R(i) and ⁇ (i) respectively denote the power matrix, covariance matrix and delay-azimuth two-dimensional spectral vector obtained from the ith iteration
  • the iterative process is as follows:
  • the base station 104 determines whether the algorithm has converged. It can be considered that when the estimated value results of the two ⁇ vectors before and after are no longer improved, that is, when the degree of change value between the estimated values of the two ⁇ vectors before and after is less than a predetermined threshold, The algorithm converges.
  • the degree of change value is a numerical value representing the degree of change between the estimated values of the ⁇ vector before and after two times, and the degree of change value can be, for example, the difference between the estimated value of the ⁇ vector of the current number and the estimated value of the ⁇ vector of the previous number of l 2
  • the value obtained by dividing the square of the norm by the square of the l2 norm estimate of the ⁇ vector of the previous number, that is, the convergence of the algorithm can be determined when the following formula is satisfied:
  • ⁇ is the set threshold value.
  • the IAA-APES algorithm has a fast convergence speed and can generally converge within 15 iterations.
  • Step S253 Obtain the signal amplitude values at each delay-azimuth grid point corresponding to the estimated value of the delay-azimuth two-dimensional spectral vector of the current iteration number when the iterative update is stopped, and based on each delay-azimuth grid The signal amplitude values at the grid points form the corresponding delay-azimuth spectrum.
  • the determined estimated value of the delay-azimuth two-dimensional spectrum contains the determined signal amplitude values at each delay-azimuth grid point, and each delay-azimuth grid point (p, q) corresponds to the delay-azimuth coordinate pair ( ⁇ p , ⁇ q ), so the base station 104 can generate the corresponding delay-azimuth based on the delay-azimuth coordinate pair ( ⁇ p , ⁇ q ) and the signal amplitude value corresponding to each delay-azimuth grid point (p, q) angle spectrum.
  • step S250 includes: using an iterative minimization sparse learning algorithm to iteratively estimate the solution equation to determine a delay-azimuth formed by the signal amplitude values at each delay-azimuth grid point angle spectrum.
  • the base station 104 determines according to the solution equation
  • the objective function of the norm-regularized least squares problem is to use the Sparse Learning via Iterative Minimization (SLIM) algorithm to iteratively solve the objective function through the loop minimization algorithm to obtain the sparseness of the above solution equation. untie.
  • SLIM Sparse Learning via Iterative Minimization
  • the objective function is expressed as:
  • the SLIM algorithm is a parameter-free sparse reconstruction algorithm. All parameters are solved in the iterative process, which avoids the influence of inaccurate parameter selection on the results, and has strong practicability.
  • step S250 may include steps S254-S257.
  • Step S254 based on the overcomplete response vector, determine a matching matrix, perform two-dimensional space-frequency matched filtering, and obtain an initial estimated value of the delay-azimuth two-dimensional spectral vector;
  • the base station 104 initializes the algorithm.
  • a ⁇ and ⁇ in the solution equation determined based on the overcomplete response vector are used as matching matrices, and two-dimensional space-frequency matched filtering is performed to obtain the time delay.
  • Two-dimensional space-frequency matched filtering is performed to obtain the time delay.
  • M represents the total number of sub-carriers after decimation
  • N represents the number of elements of the antenna array
  • ( ⁇ ) H represents the operator that takes the conjugate transpose of the matrix or vector.
  • Step S255 calculate the estimated mean square error according to the initial estimated value of the delay-azimuth two-dimensional spectral vector, as the initial estimated value of the noise power;
  • the base station 104 calculates the estimated mean square error according to the initial estimated value of the delay-azimuth two-dimensional spectral vector ⁇ , as the initial estimated value of the noise power ⁇ :
  • Step S256 based on the initial estimated value of the delay-azimuth two-dimensional spectral vector and the initial estimated value of the noise power, iteratively update the estimated value of the power matrix, the estimated value of the delay-azimuth two-dimensional spectral vector, and the estimated value of the noise power until the change degree value between the estimated value of the delay-azimuth two-dimensional spectral vector of the current iteration number and the estimated value of the delay-azimuth two-dimensional spectral vector of the previous iteration number is less than the predetermined threshold, stop the iteration renew;
  • the base station 104 is based on the initial estimated value ⁇ (0) of the delay-azimuth two-dimensional spectral vector ⁇ and the initial estimated value ⁇ (0) of the noise power ⁇ based on the two-dimensional space-frequency matched filtering
  • the SLIM algorithm iteratively updates the estimated value of the power matrix, the estimated value of the delay-azimuth two-dimensional spectral vector ⁇ , and the estimated value of the noise power ⁇ , denoting i to represent the ith iteration, correspondingly, P(i), ⁇ ( i) and ⁇ (i) represent the power matrix, delay-azimuth two-dimensional spectral vector and noise power obtained by the ith iteration, respectively.
  • the iterative processes are as follows:
  • the base station 104 determines whether the algorithm has converged. It can be considered that when the estimated value results of the two ⁇ vectors before and after are no longer improved, that is, when the degree of change value between the estimated values of the two ⁇ vectors before and after is less than a predetermined threshold, The algorithm converges.
  • the degree of change value is a numerical value representing the degree of change between the estimated values of the ⁇ vector before and after two times, and the degree of change value can be, for example, the difference between the estimated value of the ⁇ vector of the current number and the estimated value of the ⁇ vector of the previous number of l 2
  • the value obtained by dividing the square of the norm by the square of the l2 norm of the estimated value of the ⁇ vector of the previous number, that is, the convergence of the algorithm can be determined when the following formula is satisfied:
  • ⁇ is the set threshold value.
  • the SLIM algorithm has a fast convergence speed and can generally converge within 15 iterations.
  • Step S257 Obtain the signal amplitude value at each delay-azimuth grid point corresponding to the estimated value of the delay-azimuth two-dimensional spectral vector of the current iteration number when the iterative update is stopped, and based on each delay-azimuth grid The signal amplitude values at the grid points form the corresponding delay-azimuth spectrum.
  • the determined estimated value of the delay-azimuth two-dimensional spectrum contains the determined signal amplitude values at each delay-azimuth grid point, and each delay-azimuth grid point (p, q) corresponds to the delay-azimuth coordinate pair ( ⁇ p , ⁇ q ), so the base station 104 can generate the corresponding delay-azimuth based on the delay-azimuth coordinate pair ( ⁇ p , ⁇ q ) and the signal amplitude value corresponding to each delay-azimuth grid point (p, q) angle spectrum.
  • Step S260 Calculate and determine the location of the terminal according to the time delay-azimuth spectrum.
  • the position of the terminal can also be determined according to the time delay-azimuth spectrum solution through different approaches in this step S260 accordingly.
  • the determination of the position of the terminal according to the time delay-azimuth spectrum solution in this step will also be described below with two example algorithms, the iterative adaptive amplitude phase estimation algorithm and the iterative minimization sparse learning algorithm.
  • step S260 includes steps S261 to S263.
  • Step S261 performing peak detection on the delay-azimuth spectrum to obtain multiple spectrum peaks, so as to determine multiple sets of azimuth angle estimates, delay estimates and signal amplitude estimates corresponding to the multiple paths;
  • the base station 104 compares the azimuth and time delay ( ⁇ p , ⁇ q ) in the delay-azimuth coordinate pair of the delay-azimuth grid points corresponding to each spectral peak, and the spectral peak
  • the magnitude of which is determined as the estimated value of delay, estimated value of azimuth and estimated value of signal amplitude for each group of paths.
  • K max is the total number of spectral peaks
  • k' is the path number
  • the azimuth angle, estimated time delay and estimated signal amplitude of each path are the outputs of the parameter estimation/measurement module in the positioning system proposed in the present application.
  • Step S262 according to multiple sets of azimuth angle estimates, delay estimates and signal amplitude estimates corresponding to the multiple paths, determine the direct path from the multiple paths;
  • the direct path identification module of the base station 104 direct path identification module according to the input Kmax group path parameters
  • Each group includes the azimuth path delay and signal strength
  • the criterion of the strongest direct path power and the shortest delay can be used jointly at this time, or according to multiple methods.
  • the method of clustering analysis of frame data is used for identification, and finally the identification result is output, that is, the direct path is determined from multiple paths, and the estimated delay value corresponding to the direct path is determined.
  • Step S263 based on the estimated azimuth angle, the estimated time delay and the estimated value of the signal amplitude corresponding to the direct path, solve and determine the position of the terminal.
  • the base station 104 may send the azimuth estimated value, time delay estimated value and signal amplitude estimated value corresponding to the direct path determined by the base station 104 to the positioning server 102, and the positioning tracking module in the positioning server 102 will use the received direct path based on the received direct path.
  • the corresponding azimuth angle estimation value, time delay estimation value and signal amplitude estimation value are used to calculate the coordinates of the terminal 106 and its continuous positioning and tracking by using the least squares method, Kalman filter algorithm or particle filter and other algorithms.
  • the positioning server 102 may use multiple sets of azimuth estimated values, delay estimated values, and signal amplitude estimated values corresponding to multiple direct paths determined by multiple base stations 104 according to the positioning method implemented in the present application, respectively.
  • Terminal 106 performs joint positioning.
  • All U base stations that receive the uplink SRS signal of the same terminal 106 will calculate the path parameters corresponding to the direct path, including the estimated azimuth angle, the estimated time delay and the estimated value of the signal amplitude, that is, Uploaded to the positioning server 102, where u represents the u-th base station; the positioning server 102 uses the least squares method, Kalman filter or particle filter algorithm according to the path parameters of each direct path, combined with the motion model of the terminal 106, the position of the terminal 106 is continuously The sampling times t 1 , t 2 , . . . , t r are filtered, estimated and predicted, so as to realize the positioning and tracking of the terminal 106 .
  • step S260 may include steps S264-S266.
  • Step S264 extracting multiple spectral peaks in the time delay-azimuth spectrum to determine multiple pairs of azimuth angle estimates and time delay estimates corresponding to the multiple paths;
  • the base station 104 determines the azimuth and time delay ( ⁇ p , ⁇ q ) in the delay-azimuth coordinate pair of the delay-azimuth grid points corresponding to each spectral peak as each Delay estimates and azimuth estimates for paths K max is the total number of spectral peaks, k' is the path number, represents the estimated delay of the k'th path, represents the azimuth estimate of the k'th path.
  • the estimated values of the azimuth and time delay of each path are the outputs of the parameter estimation/measurement module in the positioning system proposed in the present application.
  • Step S265 according to the multiple pairs of azimuth angle estimates and time delay estimates corresponding to the multiple paths, determine the direct path from the multiple paths;
  • the direct path identification module of the base station 104 comprehensively utilizes that the direct path is shorter than the multi-path delay and the multi-frame direct path is smaller than the multi-path variance according to the azimuth angle and time delay characteristics of each path of the multi-frame. and other criteria to determine the direct path from multiple paths.
  • Step S266 based on the estimated azimuth angle and the estimated time delay corresponding to the direct path, solve and determine the position of the terminal.
  • the base station 104 may send the azimuth angle estimate value and the delay estimate value corresponding to the direct path determined by the base station 104 to the positioning server 102, and the positioning and tracking module in the positioning server 102 estimates the azimuth angle corresponding to the received direct path based on the received direct path.
  • Value and delay estimation value use least squares method, Kalman filter algorithm or particle filter and other algorithms to realize the solution of the coordinates of the terminal 106 and its continuous positioning and tracking.
  • the channel frequency domain response vector is converted into an overcomplete response vector of the channel frequency domain response at multiple delay-azimuth grid points in the signal range-azimuth domain of the base station, and the solution equation is established.
  • using an iterative adaptive amplitude-phase estimation algorithm or an iterative minimization sparse learning algorithm or other algorithms to iteratively estimate the solution equation to determine the delay formed by the signal amplitude values at each delay-azimuth grid point -Azimuth spectrum, and then solve to determine the location of the terminal.
  • the sparse distribution characteristics of the number of paths in the received signal in the signal range-azimuth angle domain of the base station can be used to simultaneously determine the estimated value of the azimuth angle and the estimated value of the delay according to the received signal, which effectively improves the positioning efficiency of the terminal. and positioning accuracy.
  • steps in the flowcharts of FIGS. 2-8 are shown in sequence according to the arrows, these steps are not necessarily executed in the sequence shown by the arrows. Unless explicitly stated herein, the execution of these steps is not strictly limited to the order, and the steps may be executed in other orders. Moreover, at least a part of the steps in FIGS. 2-8 may include multiple sub-steps or multiple stages. These sub-steps or stages are not necessarily executed and completed at the same time, but may be executed at different times. These sub-steps or stages are not necessarily completed at the same time. The order of execution of the steps is not necessarily sequential, but may be performed alternately or alternately with other steps or at least a part of sub-steps or stages of other steps.
  • MUSIC space-frequency two-dimensional Multiple Signal Classification
  • the physical layer of the multi-channel WiFi receiver is first
  • CSI channel state information
  • the azimuth-distance two-dimensional spectrum estimation is carried out by combining the joint matching vector of the air domain and the frequency domain, and then the two-dimensional spectral peak search is used to realize the simultaneous distance and azimuth angle. estimate.
  • the method for estimating azimuth, time delay and signal amplitude based on the space-frequency two-dimensional MUSIC algorithm has the following problems: (1) When there are coherent signals caused by multipath, the sampling covariance matrix in the MUSIC algorithm will appear rank deficient phenomenon , causing the signal subspace to spread to the noise subspace, and the algorithm fails; (2) In order to solve the problem of coherent signal spectrum estimation, a decoherent MUSIC algorithm based on spatial smoothing can be used.
  • the MUSIC algorithm requires the number of incident signals to be known, and this parameter is not easy to obtain in the actual environment; (4) The performance of the MUSIC algorithm depends on the accuracy of the covariance matrix estimation, so it is less fast The performance of the algorithm deteriorates significantly under the beat and low SNR; (5) The azimuth-range spectrum obtained by the space-frequency two-dimensional MUSIC algorithm is a pseudo-spectrum, that is, the amplitude of the MUSIC spectrum at each grid point has no actual physical meaning. It cannot reflect the signal energy, and in the multipath environment, the energy information of each path plays a crucial role in the identification of the direct path and even the final positioning.
  • the comparison between the above-mentioned two-dimensional space-frequency MUSIC algorithm and the iterative adaptive amplitude phase estimation (IAA-APES) algorithm proposed in this The estimation accuracy of the estimate of the signal amplitude.
  • the center frequency of the SRS signal is set to 2.565GHz
  • the subcarrier spacing is 30kHz
  • the number of subcarriers occupied by the SRS is 3264
  • the number of receiving array elements of the antenna array is 4
  • the receiving array element spacing is 5.8cm.
  • the distance-azimuth spectrum obtained by the space-frequency combined IAA-APES algorithm is as follows Figure 12.
  • Figure 9 and Figure 12 the spectral peaks are represented by x, and the real position of the target is represented by o.
  • Figures 10 and 11 are the cross-sectional views of the two-dimensional space-frequency MUSIC spectrum at the peak points in the distance and azimuth dimensions, respectively
  • Figures 13 and 14 are the cross-sectional views of the space-frequency IAA-APES spectrum at the peak points in the distance and azimuth dimensions, respectively.
  • the IAA-APES algorithm can obtain more accurate estimates of the time delay and azimuth parameters; and the spectrum obtained by the two-dimensional MUSIC algorithm is a pseudo-spectrum, which cannot reflect the real power of the signal. IAA -APES can also obtain a more accurate estimate of signal power (ie, signal amplitude).
  • the positioning method for uplink sounding reference signals based on simultaneous amplitude and phase estimation combines the broadband SRS uplink signal and the frequency domain and spatial domain information brought by the multi-channel base station, and uses the two-dimensional IAA-APES algorithm to achieve distance - Two-dimensional estimation of azimuth spectrum, accurate information of time delay, azimuth and signal power of each path can be extracted from the estimated two-dimensional spectrum.
  • using the space-frequency two-dimensional IAA-APES algorithm to estimate and determine the range-azimuth spectrum has the following advantages: (1) it can directly deal with coherent sources without smoothing operations; (2) it does not require information The number of sources is a priori; (3) In the case of a single snapshot, the accuracy of angle measurement and ranging is better than that of the two-dimensional space-frequency MUSIC algorithm; (4) Accurate estimation of the power of each path can be achieved.
  • the accuracy of the azimuth and delay estimates determined in the above-mentioned two-dimensional space-frequency MUSIC algorithm and the sparse learning iterative minimization (SLIM) algorithm proposed in the embodiments of the present application is compared through simulation experiments in a multipath environment.
  • the center frequency of the SRS signal is set to 2.565GHz
  • the subcarrier spacing is 30kHz
  • the number of subcarriers occupied by the SRS is 3264
  • the number of receiving array elements of the antenna array is 4, and the receiving array element spacing is 5.8cm.
  • RMSE ⁇ represents the RMSE of the azimuth angle ⁇
  • K is the number of paths
  • L is the number of Monte Carlo experiments
  • ⁇ l,k is the true value of the azimuth angle of the kth path in the lth experiment
  • ⁇ l,k is the true value of the azimuth angle of the kth path in the lth experiment
  • ⁇ l,k is the true value of the azimuth angle of the kth path of the lth experiment.
  • the distance RMSE is calculated as follows:
  • RMSE R represents the RMSE of the distance R
  • R l,k is the true value of the azimuth angle of the k-th path of the l-th experiment
  • the azimuth-range super-resolution estimation method based on the uplink signal proposed in the embodiment of the present application, in order to solve the problem of simultaneous super-resolution estimation of azimuth and distance, it is proposed to use the channel frequency domain response matrix signal model, the azimuth domain and the distance domain.
  • the sparsity of establish the solution equation, and use the SLIM algorithm to solve the solution equation.
  • the azimuth and distance estimation using the space-frequency two-dimensional SLIM sparse reconstruction algorithm in this application has the following advantages: (1) It can directly process coherent sources without smoothing operations; (2) There is no need for a priori on the number of sources; (3) In the case of a single snapshot, the accuracy of angle measurement and ranging is better than that of the two-dimensional space-frequency MUSIC algorithm.
  • a positioning apparatus 1800 including: a sounding reference signal receiving module 1810, a channel frequency domain response vector determination module 1820, an overcomplete response vector determination module 1830, and a solution equation establishment Module 1840, Delay-Azimuth Spectrum Determination Module 1850, and Position Determination Module 1860, wherein:
  • the sounding reference signal receiving module 1810 is configured to receive the uplink frequency domain sounding reference signal from the terminal through multiple channels of the base station;
  • a channel frequency domain response vector determination module 1820 configured to determine a channel frequency domain response vector according to uplink frequency domain sounding reference signals of multiple channels;
  • the overcomplete response vector determination module 1830 is configured to convert the channel frequency domain response vector into an overcomplete response representing the channel frequency domain response at a plurality of delay-azimuth grid points in the signal range-azimuth domain of the base station vector;
  • the solution equation establishment module 1840 is used to use the overcomplete response vector as the observation vector, determine the delay-azimuth two-dimensional spectral vector with the signal amplitude at each time delay-azimuth grid point to be solved, and establish a delay-azimuth two-dimensional spectral vector according to the observation The solution equation for the vector to solve the delay-azimuth two-dimensional spectral vector;
  • a delay-azimuth spectrum determination module 1850 configured to iteratively estimate the solution equation to determine the delay-azimuth spectrum formed by the signal amplitude values at each delay-azimuth grid point;
  • the location determination module 1860 is configured to solve and determine the location of the terminal according to the time delay-azimuth spectrum.
  • Each module in the above-mentioned positioning apparatus 1800 may be implemented in whole or in part by software, hardware and combinations thereof.
  • the above modules can be embedded in or independent of the processor in the computer device in the form of hardware, or stored in the memory in the computer device in the form of software, so that the processor can call and execute the operations corresponding to the above modules.
  • a base station in one embodiment, includes a baseband processing unit and one or more antenna processing units connected to the baseband processing unit, each antenna processing unit includes at least one antenna array, each antenna array is composed of multiple A plurality of array elements are arranged; the plurality of array element units included in each antenna array provide a plurality of corresponding receiving channels.
  • the baseband processing unit receives the uplink frequency domain sounding reference signal from the terminal through the multiple receiving channels provided by the antenna arrays of each antenna processing unit, and executes the positioning method in any of the above embodiments of the present application to locate the terminal.
  • the baseband processing unit may process the signals received from the antenna array.
  • the baseband processing unit may include an indoor baseband processing unit (Building Base Band Unit, BBU), or may include a centralized unit (Centralized Unit, CU) and a distribution unit ( Distributed Unit, DU), or any other structure capable of implementing the signal processing function required by the embodiments of the present application.
  • BBU Building Base Band Unit
  • CU Centralized Unit
  • DU Distributed Unit
  • Each antenna processing unit may include a remote radio unit (Remote Radio Unit) connected to the baseband processing unit and an antenna array connected to the remote radio unit through a feeder, or may include an active antenna unit (Active Radio Unit) connected to the baseband processing unit.
  • Antenna Unit, AAU Antenna Unit, an antenna array is integrated in the active antenna unit, or it can be any other suitable structure including an antenna array.
  • the baseband processing unit of the base station is connected in communication with the positioning server, and the baseband processing unit of the base station determines the direct path corresponding to the positioning method according to any of the above-mentioned embodiments.
  • the azimuth estimated value, delay estimated value and signal amplitude estimated value of the direct path are sent to the positioning server, and the azimuth angle of the direct path is determined by the positioning server.
  • the estimated value, the estimated time delay value and the estimated value of the signal amplitude are continuously tracked and filtered to realize the positioning of the terminal.
  • the baseband processing unit of the base station may also directly perform continuous processing on the determined azimuth angle estimation value, delay estimation value and signal amplitude estimation value of the direct path at the local end.
  • the tracking filter can realize the positioning of the terminal.
  • the baseband processing unit of the base station is communicatively connected to the positioning server, and the base determines the azimuth estimated value corresponding to the direct path and The estimated value of the time delay, and the estimated value of the azimuth angle and the estimated value of the delay corresponding to the determined direct path are sent to the positioning server. Positioning of the terminal.
  • the baseband processing unit of the base station can also directly perform continuous tracking and filtering on the determined azimuth angle estimation value and delay estimation value of the direct path at the local end, so as to realize the tracking and filtering of the terminal positioning.
  • a positioning system comprising a positioning server and one or more base stations as in any of the above embodiments, the one or more base stations being respectively connected in communication with the positioning server.
  • Each of the base stations respectively receives the uplink frequency domain sounding reference signal of the terminal, and executes the positioning methods of the foregoing embodiments of the present application to locate the terminal.
  • a computer device including a memory and a processor, a computer program is stored in the memory, and the processor implements the following steps when executing the computer program:
  • the overcomplete response vector is used as the observation vector
  • the delay-azimuth two-dimensional spectral vector is determined by the signal amplitude at each delay-azimuth grid point to be solved, and the delay-azimuth two-dimensional spectral vector is established according to the observation vector.
  • the solution equation for solving the spectral vector
  • the position of the terminal is determined by solving.
  • the processor when the processor executes the computer program, it also implements the steps of the positioning method in any of the above embodiments, and has corresponding beneficial effects.
  • a computer-readable storage medium on which a computer program is stored, and when the computer program is executed by a processor, the following steps are implemented:
  • the overcomplete response vector is used as the observation vector
  • the delay-azimuth two-dimensional spectral vector is determined by the signal amplitude at each delay-azimuth grid point to be solved, and the delay-azimuth two-dimensional spectral vector is established according to the observation vector.
  • the solution equation for solving the spectral vector
  • the position of the terminal is determined by solving.
  • the computer program when executed by the processor, it also implements the steps of the positioning method in any of the above embodiments, and has corresponding beneficial effects.
  • Nonvolatile memory may include read only memory (ROM), programmable ROM (PROM), electrically programmable ROM (EPROM), electrically erasable programmable ROM (EEPROM), or flash memory.
  • Volatile memory may include random access memory (RAM) or external cache memory.
  • RAM is available in various forms such as static RAM (SRAM), dynamic RAM (DRAM), synchronous DRAM (SDRAM), double data rate SDRAM (DDRSDRAM), enhanced SDRAM (ESDRAM), synchronous chain Road (Synchlink) DRAM (SLDRAM), memory bus (Rambus) direct RAM (RDRAM), direct memory bus dynamic RAM (DRDRAM), and memory bus dynamic RAM (RDRAM), etc.
  • SRAM static RAM
  • DRAM dynamic RAM
  • SDRAM synchronous DRAM
  • DDRSDRAM double data rate SDRAM
  • ESDRAM enhanced SDRAM
  • SLDRAM synchronous chain Road (Synchlink) DRAM
  • SLDRAM synchronous chain Road (Synchlink) DRAM
  • Rambus direct RAM
  • DRAM direct memory bus dynamic RAM
  • RDRAM memory bus dynamic RAM

Landscapes

  • Engineering & Computer Science (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Signal Processing (AREA)
  • Mobile Radio Communication Systems (AREA)

Abstract

本申请涉及定位方法、装置、基站、计算机设备和存储介质，方法包括：通过基站的多个通道接收终端的上行频域探测参考信号；确定信道频域响应向量；将信道频域响应向量转化为代表在基站的信号作用距离-方位角域内的多个时延-方位角栅格点处的信道频域响应的超完备响应向量；建立解算方程；对解算方程进行迭代估计，以确定由各个时延-方位角栅格点处的信号幅度值形成的时延-方位角谱；以及根据时延-方位角谱，解算确定终端的位置。

Description

定位方法、装置、基站、计算机设备和存储介质
相关申请的交叉引用
本申请要求于2021年02月03日提交中国专利局、申请号为202110144931.9、名称为“定位方法、装置、计算机设备和存储介质”的中国专利申请，以及于2021年02月03日提交中国专利局、申请号为202110146480.2、名称为“定位方法、装置、计算机设备和存储介质”的中国专利申请的优先权，二者的全部内容通过引用结合在本申请中。
技术领域
本申请涉及定位技术领域，特别是涉及一种定位方法、装置、基站、计算机设备和存储介质。
背景技术
随着工业互联网、物联网和车联网的快速发展，高精度定位成为智能机器人、无人车等移动终端不可或缺的关键支撑服务。卫星导航定位技术具有广域覆盖、普适性好的优点，但是因为信号功率低、穿透力弱，主要用于室外开阔环境下的终端定位，无法在受遮蔽的环境和室内环境提供导航定位服务。
5G蜂窝移动网络运用了包括大规模MIMO、超密集组网、大带宽信号等关键技术，随着未来将大量部署的多天线5G基站，利用5G信号实现高精度定位将具有非常广阔的发展前景。
然而，目前基于5G信号进行定位的方法，在定位精度上还存在提升的空间。
发明内容
基于此，有必要提供一种能够提升定位精度的定位方法、装置、基站、计算机设备和存储介质。
一种定位方法，包括：通过基站的多个通道接收来自终端的上行频域探测参考信号；根据所述多个通道的所述上行频域探测参考信号，确定信道频域响应向量；将所述信道频域响应向量转化为代表在所述基站的信号作用距离-方位角域内的多个时延-方位角栅格点处的信道频域响应的超完备响应向量；以所述超完备响应向量作为观测向量，以待解算的各个所述时延-方位角栅格点处的信号幅度确定时延-方位角二维谱向量，建立根据所述观测向量对所述时延-方位角二维谱向量进行解算的解算方程；对所述解算方程进行迭代估计，以确定由各个所述时延-方位角栅格点处的信号幅度值形成的时延-方位角谱；以及根据所述时延-方位角谱，解算确定所述终端的位置。
一种定位装置，包括：探测参考信号接收模块，用于通过基站的多个通道接收来自终端的上行频域探测参考信号；信道频域响应向量确定模块，用于根据所述多个通道的所述上行频域探测参考信号，确定信道频域响应向量；超完备响应向量确定模块，用于将所述信道频域响应向量转化为代表在所述基站的信号作用距离-方位角域内的多个时延-方位角栅格点处的信道频域响应的超完备响应向量；解算方程建立模块，用于以所述超完备响应向量作为观测向量，以待解算的各个所述时延-方位角栅格点处的信号幅度确定时延-方位角二维谱向量，建立根据所述观测向量对所述时延-方位角二维谱向量进行解算的解算方程；时延-方位角谱确定模块，用于对所述解算方程进行迭代估计，以确定由各个所述时延-方位角栅格点处的信号幅度值形成的时延-方位角谱；以及位置确定模块，用于根据所述时延-方位角谱，解算确定所述终端的位置。
一种基站，包括基带处理单元以及与所述基带处理单元连接的一个或多个天线处理单元，每个所述天线处理单元包括至少一个天线阵，每个所述天线阵包括排列的多个阵元；每个所述天线阵包括的所述多个阵元提供对应的多个接收通道；其中，所述基带处理单元通过各个所述天线处理单元的所述天线阵提供的所述多个接收通道接收来自终端的上行频域探测参考信号，并执行上述定位方法，以实现对所述终端的定位。
一种计算机设备，包括存储器和处理器，所述存储器存储有计算机程序，所述处理器执行所述计算机程序时实现上述定位方法。
一种计算机可读存储介质，其上存储有计算机程序，所述计算机程序被处理器执行时实现上述定位方法。
本申请的一个或多个实施例的细节在下面的附图和描述中提出。本申请的其它特征、目的和优点将从说明书、附图以及权利要求书变得明显。
附图说明
为了更清楚地说明本申请实施例中的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本申请的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他实施例的附图。在附图中：
图1为一个实施例中定位方法的应用环境图；
图2为一个实施例中定位方法的流程示意图；
图3为一个实施例中确定信道频域响应向量步骤的流程示意图；
图4为一个实施例中转化得到超完备响应向量步骤的流程示意图；
图5为一个实施例中确定时延-方位角谱步骤的流程示意图；
图6为另一个实施例中确定时延-方位角谱步骤的流程示意图；
图7为一个实施例中确定终端位置步骤的流程示意图；
图8为一个实施例中确定终端位置步骤的流程示意图；
图9为一个实施例中二维空频MUSIC算法确定的距离-方位角谱图；
图10为一个实施例中二维空频MUSIC算法确定的距离-方位角谱的方位角维度截面图；
图11为一个实施例中二维空频MUSIC算法确定的距离-方位角谱的距离维度截面图；
图12为一个实施例中二维空频联合IAA-APES算法确定的距离-方位角谱图；
图13为一个实施例中二维空频联合IAA-APES算法确定的距离-方位角谱的方位角维度截面图；
图14为一个实施例中二维空频联合IAA-APES算法确定的距离-方位角谱的距离维度截面图；
图15为一个实施例中二维SLIM算法确定的距离-方位角谱图；
图16为一个实施例中二维空频MUSIC算法和二维SLIM算法的方位角估计RMSE图；
图17为一个实施例中二维空频MUSIC算法和二维SLIM算法的距离估计RMSE图；
图18为一个实施例中定位装置的结构框图。
具体实施方式
为了使本申请的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本申请进行进一步详细说明。应当理解，此处描述的具体实施例仅仅用以解释本申请，并不用于限定本申请。
本申请提供的定位方法，可以应用于如图1所示的应用环境中。其中，定位服务器102与一个或多个基站104通信连接。每个基站104可具有一个或多个天线阵，基站104通过基站104的每个天线阵具有的多个接收通道接收终端106的上行频域探测参考信号(Sounding Reference Signal，SRS)，并且执行本申请各个实施例的定位方法，以对终端106进行定位。其中，终端106可以但不限于是各种个人计算机、笔记本电脑、智能手机、平板电脑和便携式可穿戴设备等。基站104可以是5G基站，或者也可以是其他任何适用于本申请方法进行定位的类型的基站，本申请对此不作限定。
基站104可以是具有单个天线阵的小型基站，或具有分布的多个天线阵的基站，基站104的每个天线阵具有的一组阵元提供一组接收通道，基站104可以分别从其每个天线阵对应的一组接收通道接收来自终端106的一组上行频域探测参考信号，并分别针对每组上行频域探测参考信号执行本申请实施例的定位方法，以确定直达径的一对方位角估计值和时延估计值，然后基站104可以将其针对每个天线阵确定的每个直达径对应的方位角估计值和时延估计值发送至定位服务器102，定位服务器102可以从一个或多个基站104实时接收针对同一终端106的一个或多个直达径对应的方位角估计值和时延估计值，并基于接收的这些直达径对应的方位角估计值和时延估计值解算确定该终端106的位置。
在一个实施例中，如图2所示，提供了一种定位方法，以该方法应用于图1中的基站104为例进行说明，包括以下步骤S210-S260。
步骤S210，通过基站的多个通道接收来自终端的上行频域探测参考信号。
示例地，在本步骤中，基站102通过基站102的天线阵的多个接收通道接收终端发送的上行频域SRS，设基站102的该天线阵共有N个阵元，则该天线阵共有N个接收通道，SRS占用的子载波数量为M 0，则从每个接收通道接收的上行频域SRS可表示为向量
Figure PCTCN2021142371-appb-000001
其中X m,n表示第n个接收通道、第m个子载波接收的频域SRS。
步骤S220，根据多个通道的上行频域探测参考信号，确定信道频域响应向量。
在一个实施例中，如图3所示，步骤S220可以包括步骤S221、S223和S225。
步骤S221，根据多个通道的上行频域探测参考信号，确定多个通道的接收信号矩阵。
示例地，当在上一步骤S210中接收了每个接收通道的频域SRS向量
Figure PCTCN2021142371-appb-000002
后，在本步骤中，基站102可以将所有通道的接收信号矩阵表示为：
Figure PCTCN2021142371-appb-000003
其中
Figure PCTCN2021142371-appb-000004
表示复数空间，
Figure PCTCN2021142371-appb-000005
表示接收信号矩阵X为M 0×N维复数矩阵。
在一个实施例中，在步骤S221之后，步骤S223之前，方法还包括：
步骤S222，对多个通道的接收信号矩阵进行天线阵校正；
实际系统中，天线阵中各个天线阵元和各射频通道之间都存在相位不一致性。在执行步骤S222之前，可以预先测量得到天线阵的各个阵元在不同子载波上的通道幅相误差和天线幅相误差，合计为
Figure PCTCN2021142371-appb-000006
Φ m,n为第n个通道第m个子载波的总幅相误差，然后利用这些通道幅相误差和天线幅相误差对原接收信号矩阵进行校正，校正后的接收信号矩阵X为：
X＝conj(Φ)⊙X 0
其中conj(·)表示取复数共轭，⊙为Hadamard积，X 0为校正前的接收信号矩阵。通过以上校正，可以对天线阵中各个天线阵元和各射频通道之间存在的相位不一致性进行补偿，从而实现对信号的准确测量，以实现精确定位。
设SRS频域发送序列为S[m],m＝1,2,…,M 0，发送信号中心载频为f c，对应波长为λ，子载波间隔为Δf，不失一般性，以接收信号的天线阵为等距线阵(Uniform Linear Array ULA)为例，且该等距线阵的阵元间距为d。另外，假设信号传输的总路径数量为K，第k条径的时延、方位角、接收信号幅度分别为：
Figure PCTCN2021142371-appb-000007
Figure PCTCN2021142371-appb-000008
和
Figure PCTCN2021142371-appb-000009
其中，
Figure PCTCN2021142371-appb-000010
定义为信号入射方向与ULA法线方向的夹角。可以理解，在本申请的各个实施例中，信号传输的时延可代表该信号传输的距离，时延与距离之间可以通过运算而相互转化。则多个通道的接收信号矩阵X可表示为：
Figure PCTCN2021142371-appb-000011
其中，S＝diag([S[1],S[2],…,S[M 0]])，diag(·)运算符表示以向量的每个元素作为主对角线元素，获得对角矩阵。
Figure PCTCN2021142371-appb-000012
为信号时延域匹配矢量，
Figure PCTCN2021142371-appb-000013
表示正交频分复用(Orthogonal Frequency Division Multiplexing,OFDM)信号时延域的匹配向量函数，其输入为路径时延τ，T为所有可能的时延τ组成的空间，即
Figure PCTCN2021142371-appb-000014
其中
Figure PCTCN2021142371-appb-000015
表示实数空间，输出为M维向量。根据信号时延在OFDM多载波上造成的相位差异，
Figure PCTCN2021142371-appb-000016
向量的第m个元素为：
Figure PCTCN2021142371-appb-000017
Figure PCTCN2021142371-appb-000018
为空域导向矢量，
Figure PCTCN2021142371-appb-000019
表示接收阵列导向矢量函数，其中j代表虚数单位，有j 2＝-1，其输入为信号到达方位角θ，Θ为所有可能的方位角组成的空间，即
Figure PCTCN2021142371-appb-000020
输出为N维向量。
Figure PCTCN2021142371-appb-000021
的具体形式取决于阵列结构，当接收阵为ULA时，其第n个元素为：
Figure PCTCN2021142371-appb-000022
Figure PCTCN2021142371-appb-000023
为噪声矩阵，其第m行第n列的元素表示第m个子载波、第n个接收通道上的噪声分量。
步骤S223，根据接收信号矩阵，进行信道估计，得到信道频域响应矩阵。
示例地，在本步骤中，基站104的信道估计模块根据多个接收通道接收的频域SRS构成的接收信号矩阵，进行信道估计，得到信道频域响应矩阵，假设接收机利用已知的SRS序列进行信道估计，则执行多通道频域信道估计得到的信道频域响应矩阵可表示为：
Figure PCTCN2021142371-appb-000024
上式中，
Figure PCTCN2021142371-appb-000025
表示信道频域响应矩阵，其第n列为第n个接收通道估计的信道频域响应。
Figure PCTCN2021142371-appb-000026
表示信道频域响应矩阵中的噪声分量。
在一个实施例中，可选地，在步骤S223之后，在步骤S225之前，方法还可以包括：步骤S224，对信道频域响应矩阵的子载波维度进行降维处理。
当基站102为5G基站时，由于5G信道中，子载波间隔较小，子载波数量较大，而空域和频域两维联合后，信号维度较大，为降低运算量，基站104对信道估计矩阵的子载波维进行降维处理，即对信道估 计矩阵的子载波维进行抽取，记抽取率为v，抽取后得到的信道频域响应矩阵为
Figure PCTCN2021142371-appb-000027
其中，
Figure PCTCN2021142371-appb-000028
Figure PCTCN2021142371-appb-000029
表示向下取整运算符，则H矩阵第m行元素可表示为：
H(m,:)＝H 0(vm,:),m＝1,2,…,M
其中，H(m,:)表示矩阵H第m行的所有元素，H 0(vm,:)表示矩阵H 0的vm行的所有元素。
步骤S225，对信道频域响应矩阵进行向量化，得到信道频域响应向量。
示例地，在本步骤中，基站104对前述步骤得到的信道频域响应矩阵进行向量化，得到信道频域响应向量
Figure PCTCN2021142371-appb-000030
其中，vec(·)表示矩阵向量化运算符，有：
Figure PCTCN2021142371-appb-000031
上式中，
Figure PCTCN2021142371-appb-000032
表示噪声向量，
Figure PCTCN2021142371-appb-000033
表示空域和时延域两维联合匹配矢量函数，两个输入参数分别为时延τ和方位角θ，如前所述，T和Θ分别为所有时延和方位角构成的空间，且有：
Figure PCTCN2021142371-appb-000034
上式中，
Figure PCTCN2021142371-appb-000035
和
Figure PCTCN2021142371-appb-000036
分别为步骤S221中引入的信号时延域匹配矢量以及空域导向矢量。
步骤S230，将信道频域响应向量转化为代表在基站的信号作用距离-方位角域内的多个时延-方位角栅格点处的信道频域响应的超完备响应向量。
通过本步骤，可以建立对接收信号在空域-时延域的超完备表示。
在一个实施例中，如图4所示，步骤S230可以包括步骤S231-步骤S234。
步骤S231，获取基站的天线阵的作用距离范围和接收方位角范围；
在一个实施例中，获取基站的天线阵的作用距离范围和接收方位角范围包括：根据基站的天线阵的结构信息以及天线阵的朝向信息，确定天线阵的接收方位角范围；根据终端的发射功率以及天线阵的灵敏度，确定天线阵的最大作用距离，根据最大作用距离确定天线阵的作用距离范围。
示例地，在本步骤中，基站104根据终端的发射功率、基站接收机的灵敏度等指标，确定基站天线阵的最大作用范围，从而确定基站天线阵能接收SRS的最大时延τ max，该最大时延τ max即可代表天线阵的最大作用距离，如此，天线阵能接受SRS的时延范围为[0,τ max]，该时延范围为[0,τ max]即可代表天线阵的作用距离范围；根据基站的天线阵结构以及阵元方向图主瓣宽度，确定基站天线阵的接收方位角范围[θ min,θ max]。
步骤S232，基于作用距离范围和接收方位角范围，确定天线阵的作用距离-方位角域。
步骤S233，使用均匀栅格对作用距离-方位角域进行分割，以确定作用距离-方位角域上均匀分布的多个时延-方位角栅格点。
示例地，在本步骤中，基站104使用均匀栅格
Figure PCTCN2021142371-appb-000037
对时延范围[0,τ max]、接收方位角范围[θ min,θ max]进行划分，从而确定时延范围[0,τ max]和接收方位角范围[θ min,θ max]共同确定的二维作用距离-方位角域内均匀分布的多个时延-方位角栅格点。其中，每个时延-方位角栅格点可对应有一个时延-方位角坐标对，示例地，该时延-方位角坐标对可以为该对应的时延-方位角栅格点处的时延和方位角形成的坐标对。每个时延-方位角栅格点还可以用下标来表示，其中，p表示时延域栅格点编号，P表示时延域总栅格点数，τ p表示第p个时延域栅格点对应的时延，q表示方位角域栅格点编号，Q表示方位角域总栅格点数，θ q表示第q个方位角域栅格点对应的方位角。从而时延-方位角栅格点(p,q)对应有时延-方位角坐标对(τ p，θ q)。
步骤S234，将信道频域响应矩阵转化为在多个时延-方位角栅格点处的信道频域响应的超完备响应向量。
例如，可以将信道频域响应矩阵转化为在作用距离-方位角域内的所有时延-方位角栅格点处的信道频域响应的超完备响应向量。
示例地，在本步骤中，基站104假设存在p k、q k,k＝1,…,K，使得
Figure PCTCN2021142371-appb-000038
其中，k表示路径编号，K表示环境中的总路径数，包括一条直达径以及(K-1)条反射径，p k表示第k条路径时延对应在时延域栅格集合
Figure PCTCN2021142371-appb-000039
中的编号，q k表示第k条路径方位角对应在方位角域栅格集合
Figure PCTCN2021142371-appb-000040
中 的编号，且以γ p,q,p＝1,…,P,q＝1,…,Q表示第p个时延域栅格点以及第q个方位角域栅格点上的幅度，根据以上假设，有在p＝p k,q＝q k时为
Figure PCTCN2021142371-appb-000041
在其它栅格上γ p,q＝0。自此基站104得到了信道频域响应矩阵在距离-方位角二维空间上的超完备表示，
Figure PCTCN2021142371-appb-000042
可重新表示为超完备响应向量：
Figure PCTCN2021142371-appb-000043
上式中，a τ.θ(τ p,θ q)表示时延为τ p、方位角为θ q时的空域-时延域两维匹配矢量，
Figure PCTCN2021142371-appb-000044
其中，
Figure PCTCN2021142371-appb-000045
τ p作为可替代的输入变量，a τ(τ p)的具体形式可同理地参阅
Figure PCTCN2021142371-appb-000046
的具体形式；θ q、
Figure PCTCN2021142371-appb-000047
作为可替代的输入变量，a θ(θ q)的具体形式可同理地参阅
Figure PCTCN2021142371-appb-000048
的具体形式。
步骤S240，以超完备响应向量作为观测向量，以待解算的各个时延-方位角栅格点处的信号幅度确定时延-方位角二维谱向量，建立根据观测向量对时延-方位角二维谱向量进行解算的解算方程。
示例地，在本步骤中，记
Figure PCTCN2021142371-appb-000049
Figure PCTCN2021142371-appb-000050
上述A τ,θ表示空域-时延域二维匹配矩阵，γ表示幅度向量，它们的维度分别为MN×PQ以及PQ×1。因此，可将前述步骤得到的超完备响应向量
Figure PCTCN2021142371-appb-000051
转化为解算方程：
Figure PCTCN2021142371-appb-000052
至此，基站104将方位角和距离(时延)的估计问题转化为根据以超完备响应向量
Figure PCTCN2021142371-appb-000053
为观测向量对时延-方位角二维谱向量γ进行求解的重构问题。
步骤S250，对解算方程进行迭代估计，以确定由各个时延-方位角栅格点处的信号幅度值形成的时延-方位角谱。
可以采用多种不同的算法来实现上述步骤S250，在本申请的下述实施例中，将提供迭代自适应的幅度相位估计算法和迭代最小化稀疏学习算法两种示例算法来对解算方程进行迭代估计以确定时延-方位角谱。
在一个实施例中，步骤S250包括：利用迭代自适应的幅度相位估计算法，对解算方程进行迭代估计，以确定由各个时延-方位角栅格点处的信号幅度值形成的时延-方位角谱。
在本实施例中，针对步骤S240建立的解算方程，基站104根据解算方程确定目标函数，并利用迭代自适应的幅度相位估计(Iterative Adaptive Approach for Amplitude and Phase Estimation，IAA-APES)算法来对该目标函数进行求解以获得以上解算方程的解。示例地，IAA-APES算法在求解第(p,q)个栅格点的时延-方位角谱的目标函数为：
Figure PCTCN2021142371-appb-000054
其中，||·|| 2表示向量的l 2范数，如对n维向量a，有
Figure PCTCN2021142371-appb-000055
R p,q表示第(p,q)个栅格点处的干扰协方差矩阵，其中的干扰由当前栅格点(p,q)以外的信号分量构成，即可表示为：
Figure PCTCN2021142371-appb-000056
其中，
Figure PCTCN2021142371-appb-000057
为超完备响应向量
Figure PCTCN2021142371-appb-000058
的协方差矩阵，E[·]表示取期望。IAA-APES的目标函数中，R为未知量，因此IAA-APES采用迭代的方式，对R和谱值γ进行交替求解。
APES算法是一类基于自适应窄带滤波器组的谱估计算法，为提高分辨能力，自适应滤波器系数以最小方差无失真响应(Minimum Variance Distortion-less Response,MVDR)为准则计算得到，能够保证当前谱中心能量无损失的前提下，其它谱位置的能量最小化，因此，一方面，APES算法在真实的谱峰位置，幅度和相位的估计较为准确，另一方面，由于抑制了其它位置的谱能量，APES算法的谱峰较为尖锐，具有一定的超分辨能力。
自适应波束形成(Minimum variance distortionless response，MVDR)算法准则下的自适应滤波器系数需要根据协方差矩阵的逆来计算，传统APES算法利用多块拍，或分子阵滑动平均的方式来获得较为准确的协方差矩阵估值。为提高定位系统的实时性、降低存储量，本申请实施例中，将IAA-APES算法应用于上行SRS的空域-频域联合处理中，迭代地对协方差矩阵和时延-方位角谱求解。
在一个实施例中，如图5所示，对于迭代自适应的幅度相位估计算法，步骤S250可以包括步骤S251-步骤S253。
步骤S251，基于超完备响应向量，确定匹配矩阵，进行二维空频匹配滤波，得到时延-方位角二维谱向量的初始估计值；
示例地，在本步骤中，基站104对算法进行初始化，初始化时，以基于超完备响应向量确定的解算方程中的A τ,θ为匹配矩阵，进行二维空频匹配滤波，得到时延-方位角二维谱向量γ的初始估计值：
Figure PCTCN2021142371-appb-000059
上式中，M表示抽取后子载波总数，N表示天线阵的阵元数，(·) H表示对矩阵或向量取共轭转置的运算符。
步骤S252，基于时延-方位角二维谱向量的初始估计值，迭代更新功率矩阵的估计值、协方差矩阵的估计值以及时延-方位角二维谱向量的估计值，直至当前迭代次数的时延-方位角二维谱向量的估计值与上一迭代次数的时延-方位角二维谱向量的估计值之间的变化度值小于预定阈值时，停止迭代更新；
示例地，在本步骤中，基站104基于二维空频匹配滤波得到的时延-方位角二维谱向量γ的初始估计值γ(0)，利用IAA-APES算法迭代地更新功率矩阵的估计值、协方差矩阵的估计值以及时延-方位角二维谱向量γ的估计值。记i表示第i次迭代，相应地，P(i)、R(i)以及γ(i)分别表示第i次迭代得到的功率矩阵、协方差矩阵以及时延-方位角二维谱向量，迭代过程分别如下：
更新功率矩阵：
P(i)＝diag(|γ 1,1| 2,…,|γ P,Q| 2)
更新协方差矩阵：
Figure PCTCN2021142371-appb-000060
更新γ：
Figure PCTCN2021142371-appb-000061
在每执行一次更新后，基站104判断算法是否收敛，可以认为前后两次γ向量的估计值结果不再改善时，即前后两次γ向量的估计值之间的变化度值小于预定阈值时，算法收敛。变化度值是表征前后两次γ向量的估计值之间的变化程度的数值，变化度值例如可以是当前次数的γ向量的估计值与前一次数的γ向量的估计值之差的l 2范数平方除以前一次数的γ向量的l 2范数估计值平方得到的值，即满足下式时可判定算法收敛：
Figure PCTCN2021142371-appb-000062
ε为设定的门限值。当算法收敛时，得到时延-方位角二维谱的估计值：
Figure PCTCN2021142371-appb-000063
IAA-APES算法收敛速度较快，一般能够在15次迭代以内收敛。
步骤S253，获取停止迭代更新时的当前迭代次数的时延-方位角二维谱向量的估计值对应的各个时延-方位角栅格点处的信号幅度值，基于各个时延-方位角栅格点处的信号幅度值，形成对应的时延-方位角谱。
示例地，在本步骤中，在算法收敛时，即停止迭代更新时，确定的时延-方位角二维谱的估计值
Figure PCTCN2021142371-appb-000064
中，包含有确定的各个时延-方位角栅格点处的信号幅度值，而每个时延-方位角栅格点(p,q)对应有时延-方位角坐标对(τ p，θ q)，从而基站104可以基于各个时延-方位角栅格点(p,q)对应的时延-方位角坐标对(τ p，θ q)和信号幅度值，生成对应的时延-方位角谱。
在另一个实施例中，步骤S250包括：利用迭代最小化稀疏学习算法，对解算方程进行迭代估计，以确定由各个时延-方位角栅格点处的信号幅度值形成的时延-方位角谱。
在本实施例中，针对步骤S240建立的解算方程，基站104根据解算方程确定
Figure PCTCN2021142371-appb-000065
范数正则化最小二乘问题的目标函数，利用迭代最小化稀疏学习(Sparse Learning via Iterative Minimization，SLIM)算法，通过循环最小化算法迭代地对该目标函数进行求解来获得上述解算方程的稀疏解。示例地，该目标函数表示为：
Figure PCTCN2021142371-appb-000066
其中，η为噪声功率，
Figure PCTCN2021142371-appb-000067
为向量γ的l 1范数。
Figure PCTCN2021142371-appb-000068
表示向量
Figure PCTCN2021142371-appb-000069
的l 2范数。
SLIM算法是一种无参数的稀疏重构类算法，所有参数都在迭代过程中求解，避免了参数选取不准确对结果的影响，实用性较强。
在一个实施例中，如图6所示，对于迭代最小化稀疏学习算法，步骤S250可以包括步骤S254-步骤S257。
步骤S254，基于超完备响应向量，确定匹配矩阵，进行二维空频匹配滤波，得到时延-方位角二维谱向量的初始估计值；
示例地，在本步骤中，基站104对算法进行初始化，初始化时，以基于超完备响应向量确定的解算方程中的A τ,θ为匹配矩阵，进行二维空频匹配滤波，得到时延-方位角二维谱向量γ的初始估计值：
Figure PCTCN2021142371-appb-000070
上式中，M表示抽取后子载波总数，N表示天线阵的阵元数，(·) H表示对矩阵或向量取共轭转置的运算符。
步骤S255，根据时延-方位角二维谱向量的初始估计值，计算估计均方误差，作为噪声功率的初始估计值；
示例地，在本步骤中，基站104根据时延-方位角二维谱向量γ的初始估计值计算估计均方误差，作为噪声功率η的初始估计值：
Figure PCTCN2021142371-appb-000071
步骤S256，基于时延-方位角二维谱向量的初始估计值和噪声功率的初始估计值，迭代更新功率矩阵的估计值、时延-方位角二维谱向量的估计值以及噪声功率的估计值，直至当前迭代次数的时延-方位角二维谱向量的估计值与上一迭代次数的时延-方位角二维谱向量的估计值之间的变化度值小于预定阈值时，停止迭代更新；
示例地，在本步骤中，基站104基于二维空频匹配滤波得到的时延-方位角二维谱向量γ的初始估计值γ(0)以及噪声功率η的初始估计值η(0)，SLIM算法迭代地更新功率矩阵的估计值、时延-方位角二维谱向量γ的估计值和噪声功率η的估计值，记i表示第i次迭代，相应地，P(i)、γ(i)和η(i)分别表示第i次迭代得到的功率矩阵、时延-方位角二维谱向量以及噪声功率，迭代过程分别如下：
更新功率矩阵：
P(i)＝diag(|γ 1,1| 2,…,|γ P,Q| 2)
更新γ：
Figure PCTCN2021142371-appb-000072
上式中，
Figure PCTCN2021142371-appb-000073
表示MN维单位矩阵。
更新η：
Figure PCTCN2021142371-appb-000074
在每执行一次更新后，基站104判断算法是否收敛，可以认为前后两次γ向量的估计值结果不再改善时，即前后两次γ向量的估计值之间的变化度值小于预定阈值时，算法收敛。变化度值是表征前后两次γ向量的估计值之间的变化程度的数值，变化度值例如可以是当前次数的γ向量的估计值与前一次数的γ向量的估计值之差的l 2范数平方除以前一次数的γ向量的估计值的l 2范数的平方得到的值，即满足下式时可判定算法收敛：
Figure PCTCN2021142371-appb-000075
ε为设定的门限值。当算法收敛时，得到时延-方位角二维谱的估计值：
Figure PCTCN2021142371-appb-000076
SLIM算法收敛速度较快，一般能够在15次迭代以内收敛。
步骤S257，获取停止迭代更新时的当前迭代次数的时延-方位角二维谱向量的估计值对应的各个时延-方位角栅格点处的信号幅度值，基于各个时延-方位角栅格点处的信号幅度值，形成对应的时延-方位角谱。
示例地，在本步骤中，在算法收敛时，即停止迭代更新时，确定的时延-方位角二维谱的估计值
Figure PCTCN2021142371-appb-000077
中，包含有确定的各个时延-方位角栅格点处的信号幅度值，而每个时延-方位角栅格点(p,q)对应有时延-方位角坐标对(τ p，θ q)，从而基站104可以基于各个时延-方位角栅格点(p,q)对应的时延- 方位角坐标对(τ p，θ q)和信号幅度值，生成对应的时延-方位角谱。
步骤S260，根据时延-方位角谱，解算确定终端的位置。
根据上述步骤S250中所采用的算法的不同，本步骤S260也可以相应地通过不同的途径来根据时延-方位角谱解算确定终端的位置。以下也将以迭代自适应的幅度相位估计算法和迭代最小化稀疏学习算法两种示例算法来对本步骤中根据时延-方位角谱解算确定终端的位置进行说明。
在一个实施例中，如图7所示，对于迭代自适应的幅度相位估计算法，步骤S260包括步骤S261-步骤S263。
步骤S261，对时延-方位角谱进行峰值检测，得到多个谱峰，以确定多个路径对应的多组方位角估计值、时延估计值和信号幅度估计值；
示例地，在本步骤中，基站104将各个谱峰对应的时延-方位角栅格点的时延-方位角坐标对中的方位角和时延(τ p，θ q)，以及谱峰的幅度，确定为各路径的各组时延估计值、方位角估计值和信号幅度估计值
Figure PCTCN2021142371-appb-000078
K max为谱峰总数，k’表示路径编号，
Figure PCTCN2021142371-appb-000079
表示第k’条路径的时延估计值，
Figure PCTCN2021142371-appb-000080
表示第k’条路径的方位角估计值，
Figure PCTCN2021142371-appb-000081
表示第k’条路径的信号幅度估计值。各路径的方位角、时延估计值和信号幅度估计值为本申请所提出的定位系统中的参数估计/测量模块的输出。
步骤S262，根据多个路径对应的多组方位角估计值、时延估计值和信号幅度估计值，从多个路径中确定直达径；
示例地，在本步骤中，基站104的直达径识别模块直达径识别模块根据输入的K max组路径参数
Figure PCTCN2021142371-appb-000082
每组包括方位角
Figure PCTCN2021142371-appb-000083
路径时延
Figure PCTCN2021142371-appb-000084
以及信号强度
Figure PCTCN2021142371-appb-000085
对直达径(Line-of-Sight，LOS)进行识别，由于APES类算法能够较准确地恢复各路径的幅度，因此此时可联合利用直达径功率最强和时延最短的准则，或根据多帧数据进行聚类分析的方法进行识别，最后输出识别结果，即从多个路径中确定直达径，并确定直达径对应的时延估计值
Figure PCTCN2021142371-appb-000086
方位角估计值
Figure PCTCN2021142371-appb-000087
以及信号幅度估计值
Figure PCTCN2021142371-appb-000088
步骤S263，基于直达径对应的方位角估计值、时延估计值和信号幅度估计值，解算确定终端的位置。
在本步骤中，基站104可以将其确定的直达径对应的方位角估计值、时延估计值和信号幅度估计值发送至定位服务器102，由定位服务器102中的定位跟踪模块基于接收的直达径对应的方位角估计值、时延估计值和信号幅度估计值，使用最小二乘方法、卡尔曼滤波算法或粒子滤波等算法，实现对终端106坐标的解算以及其连续定位和跟踪。
在本申请实施例中，定位服务器102可以利用多个基站104分别根据本申请实施的定位方法确定的多个直达径对应的多组方位角估计值、时延估计值和信号幅度估计值来对终端106进行联合定位。所有接收到同一终端106上行SRS信号的U个基站将各自解算的直达径对应的路径参数，包括方位角估计值、时延估计值和信号幅度估计值，即
Figure PCTCN2021142371-appb-000089
上传至定位服务器102，其中u代表第u个基站；定位服务器102使用最小二乘法、卡尔曼滤波或粒子滤波算法根据各个直达径的路径参数，结合终端106运动模型，对终端106的位置在连续的采样时刻t 1,t 2,…,t r进行滤波、估计和预测，从而实现对终端106的定位与跟踪。
在另一个实施例中，如图8所示，对于迭代最小化稀疏学习算法，步骤S260可以包括步骤S264-步骤S266。
步骤S264，提取时延-方位角谱中的多个谱峰，以确定多个路径对应的多对方位角估计值和时延估计值；
示例地，在本步骤中，基站104将各个谱峰对应的时延-方位角栅格点的时延-方位角坐标对中的方位角和时延(τ p，θ q)，确定为各路径的时延估计值和方位角估计值
Figure PCTCN2021142371-appb-000090
K max为谱峰总数，k’表示路径编号，
Figure PCTCN2021142371-appb-000091
表示第k’条路径的时延估计值，
Figure PCTCN2021142371-appb-000092
表示第k’条路径的方位角估计值。各路径的方位角和时延估计值为本申请所提出的定位系统中的参数估计/测量模块的输出。
步骤S265，根据多个路径对应的多对方位角估计值和时延估计值，从多个路径中确定直达径；
示例地，在本步骤中，基站104的直达径识别模块根据多帧的各路径的方位角和时延特性，综合利用直达径比多径时延短、多帧间直达径比多径方差小等准则，从多个路径中确定直达径。
步骤S266，基于直达径对应的方位角估计值和时延估计值，解算确定终端的位置。
在本步骤中，基站104可以将其确定的直达径对应的方位角估计值和时延估计值发送至定位服务器102， 由定位服务器102中的定位跟踪模块基于接收的直达径对应的方位角估计值和时延估计值，使用最小二乘方法、卡尔曼滤波算法或粒子滤波等算法，实现终端106坐标的解算以及其连续定位和跟踪。
上述定位方法中，将信道频域响应向量转化为在基站的信号作用距离-方位角域内的多个时延-方位角栅格点处的信道频域响应的超完备响应向量，建立解算方程，利用迭代自适应的幅度相位估计算法或迭代最小化稀疏学习算法或其他算法，对解算方程进行迭代估计，以确定由各个时延-方位角栅格点处的信号幅度值形成的时延-方位角谱，进而解算确定终端的位置。从而能够利用在基站的信号作用距离-方位角域内，接收信号中路径数的稀疏分布特性，根据接收的信号同时确定方位角的估计值和时延的估计值，有效提升了对终端的定位效率和定位精度。
应该理解的是，虽然图2-8的流程图中的各个步骤按照箭头的指示依次显示，但是这些步骤并不是必然按照箭头指示的顺序依次执行。除非本文中有明确的说明，这些步骤的执行并没有严格的顺序限制，这些步骤可以以其它的顺序执行。而且，图2-8中的至少一部分步骤可以包括多个子步骤或者多个阶段，这些子步骤或者阶段并不必然是在同一时刻执行完成，而是可以在不同的时刻执行，这些子步骤或者阶段的执行顺序也不必然是依次进行，而是可以与其它步骤或者其它步骤的子步骤或者阶段的至少一部分轮流或者交替地执行。
仿真实验
存在一种基于WiFi信号的空频二维多重信号分类(Multiple Signal Classification，MUSIC)算法对方位角、时延和信号幅度进行估计的方法，在该方法中，首先对多通道WiFi接收机物理层反馈的信道状态信息(Channel State Information,CSI)进行重排，并结合空域和频域的联合匹配矢量进行方位角-距离二维谱估计，随后通过二维谱峰搜索实现距离和方位角的同时估计。
基于空频二维MUSIC算法对方位角、时延和信号幅度进行估计的方法存在如下问题：(1)在存在多径引起的相干信号时，MUSIC算法中的采样协方差阵会出现秩亏现象，导致信号子空间扩散至噪声子空间，算法失效；(2)为解决相干信号谱估计问题，可采用基于空间平滑的解相干MUSIC算法，然而空间平滑一方面会导致孔径损失，另一方面也增大了运算量；(3)MUSIC算法要求已知入射信号个数，在实际环境中，这个参数不容易获得；(4)MUSIC算法性能依赖于协方差矩阵估计的准确度，因此在少快拍和低SNR下，算法性能显著恶化；(5)空频二维MUSIC算法获得的方位角-距离谱是一种伪谱，即MUSIC谱在各栅格点上的幅值没有实际物理含义，无法反映信号能量，而在多径环境中，各路径的能量信息对于直达径的识别乃至最终的定位都有着至关重要的作用。
通过仿真实验对比多径环境下，上述二维空频MUSIC算法和本申请实施例所提出的迭代自适应的幅度相位估计(IAA-APES)算法对方位角、时延(可换算成距离)和信号幅度的估计值的估计精度。仿真实验中，设置SRS信号中心频率2.565GHz，子载波间隔为30kHz，SRS占据子载波个数为3264，天线阵的接收阵元数为4，接收阵元间距为5.8cm。设置6条相干路径，6条路径的真实方位角和距离分别为：[34.6196°,14.5595°,-24.3652°,50.3258°,19.6609°,26.7446°]和[27.6448m,7.4385m,18.6720m,29.1222m,23.9435m,14.8057m]。
当信噪比(Signal to Noise Ratio，SNR)＝0dB时，平滑后二维空频MUSIC的距离-方位角谱如图9所示，空频联合IAA-APES算法得到的距离-方位角谱如图12所示。图9与图12中，谱峰用x号表示，目标真实位置用o表示。图10和图11分别为二维空频MUSIC谱在距离维和方位维峰值点处的截面图，图13和图14则分别为空频IAA-APES谱在距离维和方位维峰值点的截面图。从这些图中可以看出，整体上，IAA-APES算法能够获得更准确的时延、方位角参数的估计值；且二维MUSIC算法得到的谱是伪谱，无法反映信号的真实功率，IAA-APES则还能够得到信号功率(即信号幅度)较为准确的估计值。
本申请实施例所提出的基于幅度相位同时估计的上行探测参考信号的定位方法，联合了宽带SRS上行信号以及多通道基站带来的频域和空域信息，并使用二维IAA-APES算法实现距离-方位角谱的二维估计，从估计确定的二维谱中可以提取准确的各路径的时延、方位角以及信号功率的信息。
相比于二维空频MUSIC算法，使用空频二维IAA-APES算法估计确定距离-方位角谱具有以下优点：(1)能够直接处理相干信源，无需进行平滑操作；(2)无需信源个数先验；(3)单快拍情况下，测角和测距精度均优于二维空频MUSIC算法；(4)能够实现各路径功率的准确估计。
通过仿真实验对比多径环境下，上述二维空频MUSIC算法和本申请实施例所提出的迭代最小化稀疏学习(SLIM)算法中确定的方位角和时延估计值的精度。仿真实验中，设置SRS信号中心频率2.565GHz，子载波间隔为30kHz，SRS占据子载波个数为3264，天线阵的接收阵元数为4，接收阵元间距为5.8cm。
设置6条相干路径，在每次蒙特卡洛实验中，路径的方位角和传播距离分别在[-60°,60°]和[5m,30m]的区间内随机选取，每次蒙特卡洛实验中考虑5组信噪比(SIGNAL-to-NOISE RATIO，SNR)，分别为-10dB、-5dB、0dB、5dB、10dB。
在对实验数据的处理过程中，在二维MUSIC处理前增加了空频平滑处理，以使其能够处理相干信号。以某次蒙特卡洛实验为例，6条路径的真实方位角和距离分别为：[34.6196°,14.5595°,-24.3652°,50.3258°,19.6609°,26.7446°]和[27.6448m,7.4385m,18.6720m,29.1222m,23.9435m,14.8057m]，当SNR＝0dB时，平滑后二维空频MUSIC的距离-方位角谱如图9所示，二维SLIM算法得到的相应结果如图15所示。图9与图15中，谱峰用x号表示，目标真实位置用o表示。从图中能够看到，整体上，使用二维SLIM稀疏重构算法能够获得更准确的路径参数估计。
最后，统计了400次蒙特卡洛实验中，空频二维MUSIC算法和空频二维SLIM稀疏重构算法对所有路径的方位角和距离的估计精度，并分别以均方根误差(Root Mean Square Error,RMSE)来表达，分别如图16和图17所示。其中，方位角RMSE的计算方法如下：
Figure PCTCN2021142371-appb-000093
其中，RMSE θ表示方位角θ的RMSE，K为路径个数，L为蒙特卡洛实验次数，θ l,k为第l次实验第k个路径的方位角真实值，
Figure PCTCN2021142371-appb-000094
为第l次实验第k个路径的方位角估计值。距离RMSE的计算方法如下：
Figure PCTCN2021142371-appb-000095
其中，RMSE R表示距离R的RMSE，R l,k为第l次实验第k个路径的方位角真实值，
Figure PCTCN2021142371-appb-000096
为第l次实验第k个路径的方位角估计值。
从统计结果能够看到，空频二维SLIM稀疏重构算法的方位角和距离估计精度均优于空频二维MUSIC算法。
本申请实施例提出的基于上行信号的方位角-距离超分辨估计方法中，为解决其中的方位角、距离同时超分辨估计问题，提出利用信道频域响应矩阵信号模型中，方位域和距离域的稀疏性，建立解算方程，并使用SLIM算法对解算方程进行求解。相比于二维空频MUSIC算法，本申请使用空频二维SLIM稀疏重构算法的方位角和距离估计具有以下优点：(1)能够直接处理相干信源，无需进行平滑操作；(2)无需信源个数先验；(3)单快拍情况下，测角和测距精度均优于二维空频MUSIC算法。
在一个实施例中，如图18所示，提供了一种定位装置1800，包括：探测参考信号接收模块1810、信道频域响应向量确定模块1820、超完备响应向量确定模块1830、解算方程建立模块1840、时延-方位角谱确定模块1850和位置确定模块1860，其中：
探测参考信号接收模块1810，用于通过基站的多个通道接收来自终端的上行频域探测参考信号；
信道频域响应向量确定模块1820，用于根据多个通道的上行频域探测参考信号，确定信道频域响应向量；
超完备响应向量确定模块1830，用于将信道频域响应向量转化为代表在基站的信号作用距离-方位角域内的多个时延-方位角栅格点处的信道频域响应的超完备响应向量；
解算方程建立模块1840，用于以超完备响应向量作为观测向量，以待解算的各个时延-方位角栅格点处的信号幅度确定时延-方位角二维谱向量，建立根据观测向量对时延-方位角二维谱向量进行解算的解算方程；
时延-方位角谱确定模块1850，用于对解算方程进行迭代估计，以确定由各个时延-方位角栅格点处的信号幅度值形成的时延-方位角谱；
位置确定模块1860，用于根据时延-方位角谱，解算确定终端的位置。
关于定位装置1800的具体限定可以参见上文中对于定位方法的限定，在此不再赘述。上述定位装置1800中的各个模块可全部或部分通过软件、硬件及其组合来实现。上述各模块可以硬件形式内嵌于或独立于计算机设备中的处理器中，也可以以软件形式存储于计算机设备中的存储器中，以便于处理器调用执行以上各个模块对应的操作。
在一个实施例中，提供了一种基站，该基站包括基带处理单元以及一个或多个与基带处理单元连接的天线处理单元，每个天线处理单元包括至少一个天线阵，每个天线阵由多个阵元排列而成；每个天线阵包 括的多个阵元单元提供对应的多个接收通道。其中，基带处理单元通过各个天线处理单元的天线阵提供的多个接收通道接收来自终端的上行频域探测参考信号，并执行本申请上述任意实施例的定位方法，以对终端进行定位。
基带处理单元可以对从天线阵接收的信号进行处理，示例地，基带处理单元可以包括室内基带处理单元(Building Base Band Unit，BBU)，或者可以包括集中单元(Centralized Unit，CU)和分布单元(Distributed Unit，DU)，或者其他能够实现本申请实施例所需的信号处理功能的任意结构。
每个天线处理单元可以包括与基带处理单元连接的射频拉远单元(Remote Radio Unit)以及通过馈线与射频拉远单元连接的天线阵，或者可以包括与基带处理单元连接的有源天线单元(Active Antenna Unit，AAU)，该有源天线单元内集成有天线阵，或者也可以是包括天线阵的其他任意合适的结构。
进一步地，对于迭代自适应的幅度相位估计算法，在一个实施例中，该基站的基带处理单元与定位服务器通信连接，该基站的基带处理单元在通过上述任意实施例的定位方法确定直达径对应的方位角估计值、时延估计值和信号幅度估计值，并将直达径对应的方位角估计值、时延估计值和信号幅度估计值发送至定位服务器，由定位服务器对直达径的方位角估计值、时延估计值和信号幅度估计值进行连续的跟踪滤波，实现对终端的定位。
对于迭代自适应的幅度相位估计算法，在另一个实施例中，基站的基带处理单元也可以直接在本端对确定的直达径的方位角估计值、时延估计值和信号幅度估计值进行连续的跟踪滤波，实现对终端的定位。
进一步地，对于迭代最小化稀疏学习算法，在一个实施例中，该基站的基带处理单元与定位服务器通信连接，该基在通过上述任意实施例的定位方法确定直达径对应的方位角估计值和时延估计值，并将确定的直达径对应的方位角估计值和时延估计值发送至定位服务器，由定位服务器对直达径的方位角估计值和时延估计值进行连续的跟踪滤波，实现对终端的定位。
对于迭代最小化稀疏学习算法，在另一个实施例中，基站的基带处理单元也可以直接在本端对确定的直达径的方位角估计值和时延估计值进行连续的跟踪滤波，实现对终端的定位。
在一个实施例中，提供了一种定位系统，该定位系统包括定位服务器以及一个或多个如上任意实施例的基站，该一个或多个基站分别与定位服务器通信连接。其中每个基站分别接收终端的上行频域探测参考信号，并且执行本申请上述各个实施例的定位方法，以对终端进行定位。
在一个实施例中，提供了一种计算机设备，包括存储器和处理器，存储器中存储有计算机程序，该处理器执行计算机程序时实现以下步骤：
通过基站的多个通道接收来自终端的上行频域探测参考信号；
根据多个通道的上行频域探测参考信号，确定信道频域响应向量；
将信道频域响应向量转化为代表在基站的信号作用距离-方位角域内的多个时延-方位角栅格点处的信道频域响应的超完备响应向量；
以超完备响应向量作为观测向量，以待解算的各个时延-方位角栅格点处的信号幅度确定时延-方位角二维谱向量，建立根据观测向量对时延-方位角二维谱向量进行解算的解算方程；
对解算方程进行迭代估计，以确定由各个时延-方位角栅格点处的信号幅度值形成的时延-方位角谱；
根据时延-方位角谱，解算确定终端的位置。
在其他实施例中，处理器执行计算机程序时还实现如上任意实施例的定位方法的步骤，并具有相应的有益效果。
在一个实施例中，提供了一种计算机可读存储介质，其上存储有计算机程序，计算机程序被处理器执行时实现以下步骤：
通过基站的多个通道接收来自终端的上行频域探测参考信号；
根据多个通道的上行频域探测参考信号，确定信道频域响应向量；
将信道频域响应向量转化为代表在基站的信号作用距离-方位角域内的多个时延-方位角栅格点处的信道频域响应的超完备响应向量；
以超完备响应向量作为观测向量，以待解算的各个时延-方位角栅格点处的信号幅度确定时延-方位角二维谱向量，建立根据观测向量对时延-方位角二维谱向量进行解算的解算方程；
对解算方程进行迭代估计，以确定由各个时延-方位角栅格点处的信号幅度值形成的时延-方位角谱；
根据时延-方位角谱，解算确定终端的位置。
在其他实施例中，计算机程序被处理器执行时还实现如上任意实施例的定位方法的步骤，并具有相应的有益效果。
本领域普通技术人员可以理解实现上述实施例方法中的全部或部分流程，是可以通过计算机程序来指令相关的硬件来完成，所述的计算机程序可存储于一非易失性计算机可读取存储介质中，该计算机程序在执行时，可包括如上述各方法的实施例的流程。其中，本申请所提供的各实施例中所使用的对存储器、存储、数据库或其它介质的任何引用，均可包括非易失性和/或易失性存储器。非易失性存储器可包括只读存储器(ROM)、可编程ROM(PROM)、电可编程ROM(EPROM)、电可擦除可编程ROM(EEPROM)或闪存。易失性存储器可包括随机存取存储器(RAM)或者外部高速缓冲存储器。作为说明而非局限，RAM以多种形式可得，诸如静态RAM(SRAM)、动态RAM(DRAM)、同步DRAM(SDRAM)、双数据率SDRAM(DDRSDRAM)、增强型SDRAM(ESDRAM)、同步链路(Synchlink)DRAM(SLDRAM)、存储器总线(Rambus)直接RAM(RDRAM)、直接存储器总线动态RAM(DRDRAM)、以及存储器总线动态RAM(RDRAM)等。
以上实施例的各技术特征可以进行任意的组合，为使描述简洁，未对上述实施例中的各个技术特征所有可能的组合都进行描述，然而，只要这些技术特征的组合不存在矛盾，都应当认为是本说明书记载的范围。
以上所述实施例仅表达了本申请的几种实施方式，其描述较为具体和详细，但并不能因此而理解为对发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本申请构思的前提下，还可以做出若干变形和改进，这些都属于本申请的保护范围。因此，本申请专利的保护范围应以所附权利要求为准。

Claims (20)

  1. 一种定位方法，所述方法包括：
    通过基站的多个通道接收来自终端的上行频域探测参考信号；
    根据所述多个通道的所述上行频域探测参考信号，确定信道频域响应向量；
    将所述信道频域响应向量转化为代表在所述基站的信号作用距离-方位角域内的多个时延-方位角栅格点处的信道频域响应的超完备响应向量；
    以所述超完备响应向量作为观测向量，以待解算的各个所述时延-方位角栅格点处的信号幅度确定时延-方位角二维谱向量，建立根据所述观测向量对所述时延-方位角二维谱向量进行解算的解算方程；
    对所述解算方程进行迭代估计，以确定由各个所述时延-方位角栅格点处的信号幅度值形成的时延-方位角谱；以及
    根据所述时延-方位角谱，解算确定所述终端的位置。
  2. 根据权利要求1所述的方法，其中，所述根据所述多个通道的所述上行频域探测参考信号，确定信道频域响应向量包括：
    根据所述多个通道的上行频域探测参考信号，确定所述多个通道的接收信号矩阵；
    根据所述接收信号矩阵，进行信道估计，得到信道频域响应矩阵；以及
    对所述信道频域响应矩阵进行向量化，得到信道频域响应向量。
  3. 根据权利要求2所述的方法，其中，在所述确定所述多个通道的接收信号矩阵之后，所述根据所述接收信号矩阵，进行信道估计，得到信道频域响应矩阵之前，还包括：
    对所述多个通道的接收信号矩阵进行天线阵校正。
  4. 根据权利要求2所述的方法，其中，在所述得到信道频域响应矩阵之后，在所述对所述信道频域响应矩阵进行向量化之前，还包括：
    对所述信道频域响应矩阵的子载波维度进行降维处理。
  5. 根据权利要求1所述的方法，其中，所述将所述信道频域响应向量转化为代表在所述基站的信号作用距离-方位角域内的多个时延-方位角栅格点处的信道频域响应的超完备响应向量，包括：
    获取所述基站的天线阵的作用距离范围和接收方位角范围；
    基于所述作用距离范围和所述接收方位角范围，确定所述天线阵的作用距离-方位角域；
    使用均匀栅格对所述作用距离-方位角域进行分割，以确定所述作用距离-方位角域上均匀分布的多个时延-方位角栅格点；以及
    将所述信道频域响应矩阵转化为在所述多个时延-方位角栅格点处的信道频域响应的超完备响应向量。
  6. 根据权利要求5所述的方法，其中，所述获取所述基站的天线阵的作用距离范围和接收方位角范围包括：
    根据所述基站的天线阵的结构信息以及所述天线阵的朝向信息，确定所述天线阵的接收方位角范围；以及
    根据所述终端的发射功率以及所述天线阵的灵敏度，确定所述天线阵的最大作用距离，根据所述最大作用距离确定所述天线阵的作用距离范围。
  7. 根据权利要求1至6中任一项所述的方法，其中，所述对所述解算方程进行迭代估计，以确定由各个所述时延-方位角栅格点处的信号幅度值形成的时延-方位角谱，包括：
    利用迭代自适应的幅度相位估计算法，对所述解算方程进行迭代估计，以确定由各个所述时延-方位角栅格点处的信号幅度值形成的时延-方位角谱。
  8. 根据权利要求7所述的方法，其中，所述利用迭代自适应的幅度相位估计算法，对所述解算方程进行迭代估计，以确定由各个所述时延-方位角栅格点处的信号幅度值形成的时延-方位角谱，包括：
    基于所述超完备响应向量，确定匹配矩阵，进行二维空频匹配滤波，得到所述时延-方位角二维谱向量的初始估计值；
    基于所述时延-方位角二维谱向量的所述初始估计值，迭代更新功率矩阵的估计值、协方差矩阵的估计值以及所述时延-方位角二维谱向量的估计值，直至当前迭代次数的时延-方位角二维谱向量的估计值与上一迭代次数的时延-方位角二维谱向量的估计值之间的变化度值小于预定阈值时，停止所述迭代更新；以及
    获取停止所述迭代更新时的所述当前迭代次数的时延-方位角二维谱向量的估计值对应的各个所述时延-方位角栅格点处的信号幅度值，基于各个所述时延-方位角栅格点处的所述信号幅度值，形成对应的时延-方位角谱。
  9. 根据权利要求8所述的方法，其中，所述根据所述时延-方位角谱，解算确定所述终端的位置，包括：
    对所述时延-方位角谱进行峰值检测，得到多个谱峰，以确定多个路径对应的多组方位角估计值、时延 估计值和信号幅度估计值；
    根据所述多个路径对应的所述多组方位角估计值、时延估计值和信号幅度估计值，从所述多个路径中确定直达径；以及
    基于所述直达径对应的所述方位角估计值、所述时延估计值和所述信号幅度估计值，解算确定所述终端的位置。
  10. 根据权利要求1至6中任一项所述的方法，其中，所述对所述解算方程进行迭代估计，以确定由各个所述时延-方位角栅格点处的信号幅度值形成的时延-方位角谱，包括：
    利用迭代最小化稀疏学习算法，对所述解算方程进行迭代估计，以确定由各个所述时延-方位角栅格点处的信号幅度值形成的时延-方位角谱。
  11. 根据权利要求10所述的方法，其中，所述利用迭代最小化稀疏学习算法，对所述解算方程进行迭代估计，以确定由各个所述时延-方位角栅格点处的信号幅度值形成的时延-方位角谱，包括：
    基于所述超完备响应向量，确定匹配矩阵，进行二维空频匹配滤波，得到所述时延-方位角二维谱向量的初始估计值；
    根据所述时延-方位角二维谱向量的所述初始估计值，计算估计均方误差，作为噪声功率的初始估计值；
    基于所述时延-方位角二维谱向量的所述初始估计值和所述噪声功率的所述初始估计值，迭代更新功率矩阵的估计值、所述时延-方位角二维谱向量的估计值以及所述噪声功率的估计值，直至当前迭代次数的时延-方位角二维谱向量的估计值与上一迭代次数的时延-方位角二维谱向量的估计值之间的变化度值小于预定阈值时，停止所述迭代更新；以及
    获取停止所述迭代更新时的所述当前迭代次数的时延-方位角二维谱向量的估计值对应的各个所述时延-方位角栅格点处的信号幅度值，基于各个所述时延-方位角栅格点处的所述信号幅度值，形成对应的时延-方位角谱。
  12. 根据权利要求11所述的方法，其中，所述根据所述时延-方位角谱，解算确定所述终端的位置，包括：
    提取所述时延-方位角谱中的多个谱峰，以确定多个路径对应的多对方位角估计值和时延估计值；
    根据所述多个路径对应的所述多对方位角估计值和时延估计值，从所述多个路径中确定直达径；以及
    基于所述直达径对应的所述方位角估计值和所述时延估计值，解算确定所述终端的位置。
  13. 一种定位装置，其中，所述装置包括：
    探测参考信号接收模块，用于通过基站的多个通道接收来自终端的上行频域探测参考信号；
    信道频域响应向量确定模块，用于根据所述多个通道的所述上行频域探测参考信号，确定信道频域响应向量；
    超完备响应向量确定模块，用于将所述信道频域响应向量转化为代表在所述基站的信号作用距离-方位角域内的多个时延-方位角栅格点处的信道频域响应的超完备响应向量；
    解算方程建立模块，用于以所述超完备响应向量作为观测向量，以待解算的各个所述时延-方位角栅格点处的信号幅度确定时延-方位角二维谱向量，建立根据所述观测向量对所述时延-方位角二维谱向量进行解算的解算方程；
    时延-方位角谱确定模块，用于对所述解算方程进行迭代估计，以确定由各个所述时延-方位角栅格点处的信号幅度值形成的时延-方位角谱；以及
    位置确定模块，用于根据所述时延-方位角谱，解算确定所述终端的位置。
  14. 一种基站，包括基带处理单元以及与所述基带处理单元连接的一个或多个天线处理单元，每个所述天线处理单元包括至少一个天线阵，每个所述天线阵包括排列的多个阵元；每个所述天线阵包括的所述多个阵元提供对应的多个接收通道；其中，所述基带处理单元通过各个所述天线处理单元的所述天线阵提供的所述多个接收通道接收来自终端的上行频域探测参考信号，并执行根据权利要求1至12中任一项所述的定位方法，以实现对所述终端的定位。
  15. 根据权利要求14所述的基站，其中，所述基站的所述基带处理单元与定位服务器通信连接，所述基站的所述基带处理单元在执行根据权利要求9所述的定位方法时，确定所述直达径对应的所述方位角估计值、所述时延估计值和所述信号幅度估计值，并将所述直达径对应的所述方位角估计值、所述时延估计值和所述信号幅度估计值发送至所述定位服务器，由所述定位服务器对所述直达径的所述方位角估计值、所述时延估计值和所述信号幅度估计值进行连续的跟踪滤波，以实现对所述终端的定位。
  16. 根据权利要求14所述的基站，其中，所述基站在执行根据权利要求9所述的定位方法时，确定所述直达径对应的所述方位角估计值、所述时延估计值和所述信号幅度估计值，并由所述基站的所述基带处理单元对所述直达径的所述方位角估计值、所述时延估计值和所述信号幅度估计值进行连续的跟踪滤波，以实现对所述终端的定位。
  17. 根据权利要求14所述的基站，其中，所述基站的所述基带处理单元与定位服务器通信连接，所述基站的所述基带处理单元在执行根据权利要求12所述的定位方法时，确定所述直达径对应的所述方位角估计值和所述时延估计值，并将所述直达径对应的所述方位角估计值和所述时延估计值发送至所述定位服务器，由所述定位服务器对所述直达径的所述方位角估计值和所述时延估计值进行连续的跟踪滤波，以实现对所述终端的定位。
  18. 根据权利要求14所述的基站，其中，所述基站在执行根据权利要求12所述的定位方法时，确定所述直达径对应的所述方位角估计值和所述时延估计值，并由所述基站的所述基带处理单元对所述直达径的方位角估计值和所述时延估计值进行连续的跟踪滤波，以实现对所述终端的定位。
  19. 一种计算机设备，包括存储器和处理器，所述存储器存储有计算机程序，其中，所述处理器执行所述计算机程序时实现权利要求1至12中任一项所述方法的步骤。
  20. 一种计算机可读存储介质，其上存储有计算机程序，其中，所述计算机程序被处理器执行时实现权利要求1至12中任一项所述的方法的步骤。
PCT/CN2021/142371 2021-02-03 2021-12-29 定位方法、装置、基站、计算机设备和存储介质 WO2022166477A1 (zh)

Applications Claiming Priority (4)

Application Number Priority Date Filing Date Title
CN202110146480.2A CN112929962B (zh) 2021-02-03 2021-02-03 定位方法、装置、计算机设备和存储介质
CN202110144931.9A CN112469119B (zh) 2021-02-03 2021-02-03 定位方法、装置、计算机设备和存储介质
CN202110144931.9 2021-02-03
CN202110146480.2 2021-02-03

Publications (1)

Publication Number Publication Date
WO2022166477A1 true WO2022166477A1 (zh) 2022-08-11

Family

ID=82740840

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2021/142371 WO2022166477A1 (zh) 2021-02-03 2021-12-29 定位方法、装置、基站、计算机设备和存储介质

Country Status (1)

Country Link
WO (1) WO2022166477A1 (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2024067730A1 (zh) * 2022-09-29 2024-04-04 中兴通讯股份有限公司 多目标位置分类无线感知方法、设备、介质

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104181499A (zh) * 2014-08-23 2014-12-03 电子科技大学 方位角先验条件下基于稀布线阵的测距无源定位方法
CN106559118A (zh) * 2015-09-24 2017-04-05 中国电信股份有限公司 用于大规模天线下用户终端方位角估计的方法和装置
CN107466043A (zh) * 2016-06-03 2017-12-12 中国移动通信集团河北有限公司 一种确定基站天线的方位角的方法和设备
CN108293172A (zh) * 2015-10-08 2018-07-17 波尔特公司 用于追踪对象的到达角度定位系统
US20180358688A1 (en) * 2016-12-16 2018-12-13 Futurewei Technologies, Inc. High dimensional (hidi) radio environment characterization and representation
CN112469119A (zh) * 2021-02-03 2021-03-09 网络通信与安全紫金山实验室 定位方法、装置、计算机设备和存储介质
CN112929962A (zh) * 2021-02-03 2021-06-08 网络通信与安全紫金山实验室 定位方法、装置、计算机设备和存储介质

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104181499A (zh) * 2014-08-23 2014-12-03 电子科技大学 方位角先验条件下基于稀布线阵的测距无源定位方法
CN106559118A (zh) * 2015-09-24 2017-04-05 中国电信股份有限公司 用于大规模天线下用户终端方位角估计的方法和装置
CN108293172A (zh) * 2015-10-08 2018-07-17 波尔特公司 用于追踪对象的到达角度定位系统
CN107466043A (zh) * 2016-06-03 2017-12-12 中国移动通信集团河北有限公司 一种确定基站天线的方位角的方法和设备
US20180358688A1 (en) * 2016-12-16 2018-12-13 Futurewei Technologies, Inc. High dimensional (hidi) radio environment characterization and representation
CN112469119A (zh) * 2021-02-03 2021-03-09 网络通信与安全紫金山实验室 定位方法、装置、计算机设备和存储介质
CN112929962A (zh) * 2021-02-03 2021-06-08 网络通信与安全紫金山实验室 定位方法、装置、计算机设备和存储介质

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
HUAWEI, HISILICON: "Discussion on CSI enhancement", 3GPP DRAFT; R1-1812242, 3RD GENERATION PARTNERSHIP PROJECT (3GPP), MOBILE COMPETENCE CENTRE ; 650, ROUTE DES LUCIOLES ; F-06921 SOPHIA-ANTIPOLIS CEDEX ; FRANCE, vol. RAN WG1, no. Spokane, USA; 20181112 - 20181116, 11 November 2018 (2018-11-11), Mobile Competence Centre ; 650, route des Lucioles ; F-06921 Sophia-Antipolis Cedex ; France , XP051554124 *

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2024067730A1 (zh) * 2022-09-29 2024-04-04 中兴通讯股份有限公司 多目标位置分类无线感知方法、设备、介质

Similar Documents

Publication Publication Date Title
CN112929962B (zh) 定位方法、装置、计算机设备和存储介质
Lee et al. Low-complexity estimation of 2D DOA for coherently distributed sources
CN112469119B (zh) 定位方法、装置、计算机设备和存储介质
JP4339801B2 (ja) 固有値分解を利用しない信号到来方向推定手法および受信ビーム形成装置
Liu et al. Joint estimation of DOA and TDOA of multiple reflections in mobile communications
US20060007043A1 (en) Method and device for tracking the directions-of-arrival of radio waves
CN104977558A (zh) 一种基于贝叶斯压缩感知的分布源中心波达方向估计方法
Zhang et al. High resolution 3-D angle of arrival determination for indoor UWB multipath propagation
Abdalla et al. Performance evaluation of direction of arrival estimation using MUSIC and ESPRIT algorithms for mobile communication systems
Akbari et al. MUSIC and MVDR DOA estimation algorithms with higher resolution and accuracy
WO2023000614A1 (zh) 无线定位参数估计方法、装置、系统、计算机设备及存储介质
Chen et al. Millidegree-level direction-of-arrival (DoA) estimation and tracking for terahertz wireless communications
CN112255629B (zh) 基于联合uca的序贯esprit二维不相干分布源参数估计方法
Qi et al. Time-frequency DOA estimation of chirp signals based on multi-subarray
Feng et al. Wireless channel parameter estimation algorithms: Recent advances and future challenges
Lin et al. 3D wideband mmWave localization for 5G massive MIMO systems
WO2022166477A1 (zh) 定位方法、装置、基站、计算机设备和存储介质
Molaei et al. Efficient clustering of non-coherent and coherent components regardless of sources’ powers for 2D DOA estimation
CN114884841A (zh) 基于高阶统计和非均匀阵列的欠定参数联合估计方法
Wen et al. HOSVD based multidimensional parameter estimation for massive MIMO system from incomplete channel measurements
Li et al. A novel method for DOA and time delay joint estimation in multipath OFDM environment
Wang et al. DOA estimation of smart antenna signal based on MUSIC algorithm
Pan et al. A low-complexity joint AoA and ToA estimation method for positioning with 5G signals
Qianlin et al. Fast parallel factor decomposition technique for coherently distributed source localization
CN113740797B (zh) 一种透镜阵列下的高精度单快拍目标到达角度估计方法

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 21924470

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 21924470

Country of ref document: EP

Kind code of ref document: A1