# Multi-target Position and Velocity Estimation Using OFDM Communication Signals

In this paper, we consider a passive radar system that estimates the positions and velocities of multiple moving targets by using OFDM signals transmitted by a totally un-coordinated and un-synchronizated illuminator and multiple receivers. It is assumed that data demodulation is performed separately by a communication receiver based on the direct-path signal, and the error-prone estimated data symbols are made available to the passive radar receivers, which estimate the positions and velocities of the targets in two stages. First, we formulate a problem of joint estimation of the delay-Doppler of reflectors and the demodulation errors, by exploiting two types of sparsities of the system, namely, the numbers of reflectors (i.e., targets and clutters) and demodulation errors are both small. This problem is non-convex and a conjugate gradient descent method is proposed to solve it. Then in the second stage we determine the positions and velocities of targets based on the estimated delay-Doppler in the first stage. And two methods are proposed: the first is based on numerically solving a set of nonlinear equations, while the second is based on the back propagation neural network, which is more efficient. The performance of the proposed passive OFDM radar receiver algorithm is evaluated through extensive simulations. It is seen that accurate target estimation results are obtained even with high demodulation error rate.

## Authors

• 5 publications
• 49 publications
• 3 publications
06/28/2017

### Passive Polarimetric Multistatic Radar Detection of Moving Targets

We study the exploitation of polarimetric diversity in passive multistat...
11/11/2021

### Joint Radar-Communications Processing from a Dual-Blind Deconvolution Perspective

We consider a general spectral coexistence scenario, wherein the channel...
03/20/2021

### Radar Information Theory for Joint Communication and Parameter Estimation with Passive Targets

In this paper, we derive the performance bounds for joint communication ...
03/29/2021

### MIMO-OFDM Joint Radar-Communications: Is ICI Friend or Foe?

Inter-carrier interference (ICI) poses a significant challenge for OFDM ...
08/10/2018

### Performance Tradeoff in a Unified System of Communications and Passive Radar: A Secrecy Capacity Approach

In a unified system of passive radar and communication systems of joint ...
09/19/2018

### Direct Reconstruction of Saturated Samples in Band-Limited OFDM Signals

Given a set of samples, a few of them being possibly saturated, we propo...
02/13/2020

### Cooperative Observation of Targets moving over a Planar Graph with Prediction of Positions

Consider a team with two types of agents: targets and observers. Observe...
##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## I Introduction

Passive radar systems can detect targets by utilizing readily available, non-cooperative illuminators of opportunity (IOs) [1, 2, 3, 4], and possess a number of benefits compared with active radar systems. In particular, a passive radar is smaller and less expensive because it does not need a transmitter. And many IOs (e.g., cellular base stations [5], digital television broadcasting [6], digital audio broadcasting (DAB) [7]) are available for passive sensing, as such a passive radar can operate without causing interference to existing communication systems.

A main challenge associated with passive radar is that the IOs are non-cooperative, and the transmitted signals are unknown and not under control. Hence, the conventional matched filter cannot be easily implemented. In addition, the direct-path signal is much stronger than the target reflections, making it difficult to detect and track targets. To solve those problems, a passive radar usually makes use of an additional separate channel, referred to as the reference channel, to collect the transmitted signal in order to eliminate the direct-path signal and clutters in the surveillance channels (SCs) [8, 9]. In addition, the reference signal can also be used to implement approximate matched filtering to detect targets. However, the reference signal is noisy and the target detection performance is usually significantly degraded [1, 3].

Orthogonal frequency-division multiplexing (OFDM) techniques are widely employed in many modern wireless communication systems, e.g., 4G wireless cellular [10], digital video broadcasting (DVB) [11], DAB [7] and wireless local area network (LAN) [12, 13]. For an OFDM passive radar system, demodulation can be implemented by using the reference signal [14, 15, 3]. Since demodulation provides better accuracy than directly using the reference signal, a more accurate matched filter can be implemented based on the estimated data symbols, and the performance of passive radar can be greatly improved. Moreover, the reference channel is not always necessary because data symbols can also be directly demodulated based on the received signal in SC.

Some target detection algorithms using OFDM signals have been proposed [3, 11, 15, 16, 12]. In [16], a method for detecting a moving target in the presence of multi-path reflections is proposed based on adaptive OFDM radar. In [15, 11, 12], by assuming that the demodulation is perfect, the delays and Doppler shifts of targets are estimated based on matched filtering. And in [15]

, the MUltiple SIgnal Classifier (MUSIC) and the compressed sensing (CS) techniques are employed to obtain a better target resolution

[17] and clutter removal performance. In [3], by using the received OFDM signal from an un-coordinated but synchronizated illuminator, a delay and Doppler shift estimation algorithm is proposed taking into account the demodulation error. The atomic norm (AN) is used to enforce the signal sparsity in the delay-Doppler plane and the -norm is used to enforce the sparsity of the demodulation error signal. Then, a convex semidefinite program (SDP) is solved to obtain the estimate of the target delays and Doppler shifts.

The present contribution is aimed at extending the results of [3] by using multiple receivers to achieve the target position and velocity estimation based on the OFDM signal emitted by a totally un-coordinated and un-synchronizated illuminator. Assuming that data demodulation is performed separately by a communication receiver based on the direct-path signal, and the error-prone estimated data symbols are made available to the passive radar receivers. Then, a two-stage procedure is proposed to estimate the positions and velocities of targets.

The first stage is aimed at estimating the delay-Doppler shift and demodulation error by exploiting two types of sparsity: on one hand, as targets and clutters are sparsely distributed in space, the reflected signals hitting the radar receivers are sparse; on the other, the demodulation error rate of a communication system is typically low under normal operating conditions and hence the demodulation error signal is also sparse. Since the delays and Doppler shifts of the targets are continuous parameters, conventional CS tools [18] may lead to unsatisfactory performance [19] when the signals cannot be sparsely represented by a finite discrete dictionary [20, 21, 22]. We make use of the recently developed mathematical theory of continuous sparse recovery for super-resolution [23, 24, 25], and especially the AN minimization techniques which have been successfully applied for continuous frequency recovery, line spectral estimation and direction-of-arrival estimation [25, 26, 27]. Note that, unlike the convex problem of the delay-Doppler estimation in [3] for one receiver, in our model, different receivers share the same estimated data symbols and impose the same constraint, which yields a non-convex problem due to existence of the produce term of decision variables. Hence, we use non-convex factorization (NF) to transform the problem to a smooth unconstrained optimization problem, which is then solved by a conjugate gradient descent (CGD) algorithm.

The second stage is aimed at determining the target positions and velocities based on the estimates in the first stage. Since the illuminator and receivers are un-synchronizated, we utilize the delay differences between different receivers to calculate each target position. The first method numerically solves a set of nonlinear equations, and the second method utilizes the back propagation (BP) neural network [28, 29, 30] to estimate the target position, which is more computationally efficient. The corresponding target velocity can then be determined based on the estimated target position and Doppler shift. Extensive simulation results are provided to illustrate that the proposed methodology can estimate the target positions and velocities accurately.

The remainder of the paper is organized as follows. In Section II, we present the signal model of the OFDM passive radar and set up the problem. In Section III, we develop a delay-Doppler estimator based on conjugate gradient descent. In Section IV, we discuss methods for estimating the locations and velocities. Simulation results are presented in Section V. Section VI concludes the paper.

## Ii System Descriptions & Problem Formulation

### Ii-a System Descriptions

As shown in Fig. 1, we consider a passive radar system consisting of () receivers and one non-cooperative illuminator, that aims to estimate the locations and velocities of multiple targets in a three-dimensional cartesian coordinate system. Suppose that the coordinates of the illuminator and receiver are and , respectively. Assume that there are reflectors in the surveillance area, include targets and clutters. Note that we consider clutters as zero-velocity targets. Let and be the location and velocity of the -th reflector, respectively. Then, the traveling time from the illuminator to the -th receiving antenna due to the -th reflector is

 ¯τℓ,m=1c(∥p0−xℓ∥2+∥pm−xℓ∥2), (1)

where is the speed of light in free-space; denotes the -norm. And the corresponding Doppler shift is given by [31]

 ¯fℓ,m= vTℓ(p0−xℓ)λ∥p0−xℓ∥2+vTℓ(pm−xℓ)λ∥pm−xℓ∥2, (2)

where denotes the wavelength of the carrier.

Due to the reflections of targets and clutters, the received signal at the -th receiver is given by

 ym(t)= ydm(t)+s(t)∗L∑ℓ=1cℓ,mei2π¯fℓ,mtδ(t−¯τℓ,m+Δτ)+wm(t), (3)

where denotes the convolution operator; is the unit impulse function; is the -th path’s complex gain at the -th receiving antenna; is the unknown communication signal; is the direct-path (illuminator-to-receiver) signal111Note that the direct-path signal can be suppressed by the spatial filtering method in [8] or by using a reference channel to collect the direct-path signal [1, 2, 32], i.e., using a narrow beam antenna towards the transmitter to receive the direct-path signal [32], or using the side-lobe of the receiving antenna to receive the direct-path signal [2].; is the synchronization error between the transmitter and receivers; and is a white, complex circularly symmetric Gaussian process.

### Ii-B OFDM-based Passive Radar Signal Model

In this paper, we assume that the signal is the OFDM signal that is widely adopted in contemporary wireless communication systems. The OFDM system consists of sub-carriers, with data sub-carriers and cyclic prefix (CP) sub-carriers. The duration of an OFDM block is with being the “sub-pulse duration.” Then, the transmitted baseband OFDM signal over blocks is given by

 s(t)=Nb−1∑n=0Nd−1∑k=0bn(k)ei2πktNdTu(t−nNT), (4)

where is the -th normalized data symbol block, such that with denoting the complex conjugate operator; and

 u(t)={1,  t∈[−NpT,NdT],0,  otherwise. (5)

At each radar receiver , suppose that the direct-path signal is first removed, and we only refer to the baseband signals by assuming that down-conversion has been performed. The CP is removed assuming that its length is no less than the maximum path delay, i.e., . Note that the data symbols are unknown, but can be estimated by demodulation using the direct-path signals  [1]. However, demodulation may be error-prone. Hence in the following we assume that an estimate of the data symbols, , is available such that

 bn(k)=^bn(k)+en(k),k=0,...,Nd−1, (6)

where denotes the corresponding demodulation error. Furthermore, we assume that the velocity of the target is low, such that . Hence the phase rotation due to the Doppler shift can be approximated as constant over an OFDM block duration , i.e., [15, 3]

 ei2π¯fℓ,mt≈ei2π¯fℓ,mnNT,t∈[nNT,(n+1)NT]. (7)

At each receiver , in the -th OFDM block matched filtering is performed to obtain, for ,

 ¯yn,m(k)= 1NdT∫nNT+NdTnNT(ym(t)−ydm(t))e−i2πktNdTdt+¯wn,m(k) = L∑ℓ=1cℓ,mNd−1∑q=0bn(q)1NdT∫nNT+NdTnNTei2π¯fℓ,mt≈ei2π¯fℓ,mnNTei2πqt−¯τℓ,m+ΔτNdTe−i2πktNdTdt+¯wn,m(k) ≈ L∑ℓ=1cℓ,mei2π¯fℓ,mnNTNd−1∑q=0bn(q)e−i2πq¯τℓ,m−ΔτNdT1NdT∫nNT+NdTnNTei2π(q−k)tNdTdtNdT⋅δ(q−k)+¯wn,m(k) (8) = (^bn(k)+en(k))L∑ℓ=1cℓ,mei2πnfℓ,me−i2πkτℓ,m+¯wn,m(k), (9)

where

 τℓ,m= ¯τℓ,m−ΔτNdT∈[0,1), (10) fℓ,m= ¯fℓ,mNT∈[0,1), (11) ¯wn,m(k)= 1NdT∫nNT+NdTnNTwm(t)e−i2πktNdTdt. (12)

### Ii-C Problem Formulation

Let us now define

 cm= [c1,m,c2,m,...,cL,m]T∈CL×1, (13) fm= [f1,m,f2,m,...,fL,m]T∈CL×1, (14) τm= [τ1,m,τ2,m,...,τL,m]T∈CL×1, (15)

and the steering vectors

 s(f)= [1,ei2πf,...,ei2π(Nb−1)f]T∈CNb×1, (16) d(τ)= [1,ei2πτ,...,ei2π(Nd−1)τ]T∈CNd×1. (17)

Correspondingly, the response matrices are defined as

 S(fm)= [s(f1,m),s(f2,m),...,s(fL,m)]∈CNb×L, (18) D(τm)= [d(τ1,m),d(τ2,m),...,d(τL,m)]∈CNd×L. (19)

Then (9) can be written as the following matrix form

 ¯Ym=(^B+E)⊙(S(fm)diag(cm)D(τm)H)+¯Wm, (20)

where denotes the Hadamard product; denotes the diagonal matrix whose diagonal entries are ; , , and are matrices whose -th element are , , and , respectively.

Further denote , , and

 ϕm=L∑ℓ=1cℓ,ma(τℓ,m,fℓ,m)∈CNbNd×1, (21)

with

 a(τ,f)=d(τ)∗⊗s(f)∈CNbNd×1, (22)

and being the Kronecker product. Then we vectorize in (20) to obtain

 ¯ym= vec(¯Ym), = diag(^b+e)(D(τm)∗∘S(fm))cm+¯wm, = diag(^b+e)ϕm+¯wm, (23)

where is the Khatri-Rao product; is a matrix whose -th column has the form of . Finally, (23) can be rewritten in the following matrix form

 Y=diag(^b+e)Φ+W, (24)

where the -th columns of , and are , and , respectively.

In this paper, we first estimate the delays and Doppler shifts contained in from the received signals . Then based on these estimates, we further estimate the locations and velocities of the reflectors , and those with are considered clutters.

## Iii Stage 1: Delay-Doppler Estimation

In this section, we propose a CGD method to estimate the delays and Doppler shifts in (24). We first formulate a non-convex optimization problem by exploiting two types of sparsity. Then we relax the non-convex optimization problem to a smooth unconstrained form. The smoothed problem can then be solved via CGD.

### Iii-a Non-convex Optimization Problem Setup

We will exploit the following two types of sparsity: firstly, the number of reflectors in (21); secondly, assuming that the demodulation error rate is low, then has a small number of non-zero entries, i.e., with being the -norm. Since the delays and the Doppler shifts take continuous values, the atomic norm [25, 33] is used to exploit the first type of sparsity. Let be the set of atoms, where is defined in (22). Then the 2D atomic norm [25] associated to for is defined as

 ∥ϕm∥A = inf{χ>0:ϕm∈χconv(A)} (25) = infcℓ,m∈C,fℓ,m∈[0,1),τℓ,m∈[0,1){∑ℓ|cℓ,m|:ϕm=∑ℓcℓ,ma(τℓ,m,fℓ,m)}.

The atomic norm can enforce sparsity in the atom set . Note that columns in are independent with their own sparsities. On this basis, our delay-Doppler estimation problem can be formulated according to (24) as:

 (^Φ,^e)=argminΦ∈CNbNd×M,e∈CNbNd×112∥Y−diag(^b+e)Φ∥2F+γM∑m=1∥ϕm∥A+η∥e∥1, (26)

where denotes the -norm, and are the weight factors.

However, finding the harmonic components via atomic norm is an infinite programming problem over all feasible and . For the convenience of calculation, we use the following equivalent form of (25) for  [25, 34]

 ∥ϕm∥A=infQm∈C(2Nb−1)×(2Nd−1),νm∈R⎧⎪ ⎪ ⎪⎨⎪ ⎪ ⎪⎩12NbNdTr(T(Qm))+νm2,s.t.[T(Qm)ϕmϕHmνm]⪰0⎫⎪ ⎪ ⎪⎬⎪ ⎪ ⎪⎭, (27)

where denotes the trace operator, stands for a positive semidefinite matrix, and takes as input a matrix

 Qm=[qm,−Nd+1,qm,−Nd+2,...,qm,Nd−1]∈C(2Nb−1)×(2Nd−1), (28)

with

 qm,n=[qm,n(−Nb+1),qm,n(−Nb+2),...,qm,n(Nb−1)]T∈C(2Nb−1)×1, (29) n=−Nd+1,−Nd+2,...,Nd−1,

and outputs an block Toeplitz matrix

 T(Qm)=⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣Toep(qm,0)Toep(qm,−1)⋯Toep(qm,−Nd+1)Toep(qm,1)Toep(qm,0)⋯Toep(qm,−Nd+2)⋮⋮⋱⋮Toep(qm,Nd−1)Toep(qm,Nd−2)⋯Toep(qm,0)⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦∈CNbNd×NbNd, (30)

where denotes the Toeplitz matrix whose first column is the last elements of the input vector. More specifically, we have

 Toep(qm,n)=⎡⎢ ⎢ ⎢ ⎢ ⎢⎣qm,n(0)qm,n(−1)⋯qm,n(−Nb+1)qm,n(1)qm,n(0)⋯qm,n(−Nb+2)⋮⋮⋱⋮qm,n(Nb−1)qm,n(Nb−2)⋯qm,n(0)⎤⎥ ⎥ ⎥ ⎥ ⎥⎦∈CNb×Nb, (31) n=−Nd+1,−Nd+2,...,Nd−1.

Equations (25) and (27) are related when achieving the optimum through the relationship

 T(Qm)= ∑ℓ,m|cℓ,m|a(τℓ,m,fℓ,m)a(τℓ,m,fℓ,m)H, (32) νm= ∑ℓ,m|cℓ,m|. (33)

By using (27), (26) can be transformed to the following optimization problem:

 (^Φ,^e)= argminΦ∈CNbNd×M,Qm∈C(2Nb−1)×(2Nd−1),e∈CNbNd×112∥Y−diag(^b+e)Φ∥2F+γ2NbNdM∑m=1Tr(T(Qm))+γ2M∑m=1νm+η∥e∥1, (34) s.t. [T(Qm)ϕmϕHmνm]⪰0, m=1,...,M.

Note that the above problem is non-convex, since it involves the product term of and . In the following subsection, we will introduce a CGD method to solve the non-convex optimization problem (34).

### Iii-B Conjugate Gradient Descent Algorithm

Define and

 Θm=[UmϕmϕHmνm], m=1,...,M. (35)

Then problem (34) is rewritten as

 (^Φ,^e)= argminΘm∈C(NbNd+1)×(NbNd+1)e∈CNbNd×112∥Y−diag(^b+e)Φ∥2F+γ2NbNdM∑m=1Tr(Um)+γ2M∑m=1νm+η∥e∥1, (36) s.t. T(P(Um))=Um, Θm⪰0, m=1,...,M,

where denotes an inverse operation on the input block Toeplitz matrix, and outputs a matrix. In particular, if we partition the block Toeplitz matrix into blocks, i.e.,

 Um=⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣¯Um,1,1¯Um,1,2⋯¯Um,1,Nd¯Um,2,1¯Um,2,2⋯¯Um,2,Nd⋮⋮⋱⋮¯Um,N,1¯Um,N,2⋯¯Um,Nd,Nd⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦∈CNbNd×NbNd, (37)

then the -th element of is given by

 P(Um)(i,j)=1βi,jb1−b2=j∑d1−d2=i¯Um,d1,d2(b1,b2), b1,b2=1,2,...,Nb; d1,d2=1,2,...,Nd, (38)

where

 βi,j= (Nb−|j|)(Nd−|i|), (39) j=−Nb+1,−Nb+2,...,Nb−1, i=−Nd+1,−Nd+2,...,Nd−1.

To solve (36) via the CGD algorithm, we need relax it to a smooth unconstrained form. Hence, we first replace the constraint by the penalty term , and approximate the -norm by a twice continuously differentiable function [35]

 ∥e∥1≈ ϕϖ(e) = ϖNdNb∑n=1log(exp(|en|/ϖ)+exp(−|en|/ϖ)2) = ϖNdNb∑n=1logcosh(|en|ϖ), (40)

where denotes the -th element in and is a weight parameter, which controls the smoothing level.

Furthermore, we remove the constraints by setting with , such that is chosen minimally according to the ranks of . In particular, we have the following lemma. The proof is given in Appendix A.

###### Lemma 1.

Suppose is the solution to (26), where

 ^ϕm=L∑ℓ=1^cℓ,ma(^τℓ,m,^fℓ,m), m=1,...,M. (41)

Then each in the solution to (36) is rank- if 222Actually, the condition is a technical requirement that used in Theorem 1 in [34]. Through simulations we found that the result still holds without such condition. and hold.

The above Lemma 1 shows that each in the solution to (36) should be rank-. It is mentioned in [36] that an algorithm can be accelerated through non-convex factorization if the solution is low-rank, since the size of the optimization variables is significantly reduced (see Fig. 2). Hence, let be the upper bound on the number of reflectors. We can then let such that the constraints and are both satisfied. Then, (36) can be rewritten as the following smooth unconstrained optimization problem

 minZm∈C(NbNd+1)×¯Le∈CNbNd×1ζ({ZmZHm}Mm=1,e), (42)

where

 ζ({ZmZHm}Mm=1,e)= 12∥Y−diag(^b+e)Φ∥2F+γ2NbNdM∑m=1Tr(Um) +γ2M∑m=1νm+ηϕϖ(e)+M∑m=1ρ2∥T(P(Um))−Um∥2F. (43)

The CGD algorithm for solving (42) performs the following iterations

 Zim= Zi−1m+μiGim, m=1,...,M, (44) ei= ei−1+μigi, (45)

where is the step size, which is chosen according to the backtracking line search [37], given in Appendix B, to guarantee that the objective function does not increase with ; and are the search directions

 Gim= −∇iZmζ+¯μiGi−1m, (46) gi= −∇ieζ+¯μigi−1, (47)

with and , where and are the gradients, which are derived in Appendix C; and

 ¯μi=M∑m=1⟨∇iZmζ,∇iZmζ−∇i−1Zmζ⟩+⟨∇ieζ,∇ieζ−∇i−1eζ⟩M∑m=1⟨Gi−1m,∇iZmζ−∇i−1Zmζ⟩+⟨gi−1,∇ieζ−∇i−1eζ⟩, (48)

where