Joint Doppler and DOA estimation using (Ultra-)Wideband FMCW signals

Joint Doppler and DOA estimation using (Ultra-)Wideband FMCW signals

Xu, Shengzhi; Kooij, Bert Jan; Yarovoy, Alexander

DOI

10.1016/j.sigpro.2019.107259

Publication date

2020

Document Version

Final published version

Published in

Signal Processing

Citation (APA)

Xu, S., Kooij, B. J., & Yarovoy, A. (2020). Joint Doppler and DOA estimation using (Ultra-)Wideband FMCW

signals. Signal Processing, 168, [107259]. https://doi.org/10.1016/j.sigpro.2019.107259

Important note

To cite this publication, please use the final published version (if applicable).

Please check the document version above.

Copyright

Other than for strictly personal use, it is not permitted to download, forward or distribute the text or part of it, without the consent of the author(s) and/or copyright holder(s), unless the work is under an open content license such as Creative Commons. Takedown policy

Please contact us and provide details if you believe this document breaches copyrights. We will remove access to the work immediately and investigate your claim.

This work is downloaded from Delft University of Technology.

Green Open Access added to TU Delft Institutional Repository

‘You share, we take care!’ – Taverne project

https://www.openaccess.nl/en/you-share-we-take-care

Otherwise as indicated in the copyright section: the publisher

is the copyright holder of this work and the author uses the

Dutch legislation to make this work public.

ContentslistsavailableatScienceDirect

Signal

Processing

journalhomepage:www.elsevier.com/locate/sigpro

Joint

Doppler

and

DOA

estimation

using

(Ultra-)Wideband

FMCW

signals

Shengzhi

Xu

∗

,

Bert

Jan

Kooij

,

Alexander

Yarovoy

Delft University of Technology, Delft 2628 CD, the Netherlands

a

r

t

i

c

l

e

i

n

f

o

Article history:

Received 28 January 2019 Revised 26 May 2019 Accepted 15 August 2019 Available online 4 September 2019

Keywords:

Ultra-Wideband (UWB) Direction-of-Arrival (DOA) Doppler

2 Dimensional multiple signal classification (2D MUSIC)

Lanczos algorithm Rayleigh–Ritz step

a

b

s

t

r

a

c

t

ThejointDopplerandDirection-of-Arrival(DOA)estimationofmovingtargetsusingan(Ultra-)Wideband (UWB)frequencymodulatedcontinuous-wave (FMCW)antennaarrayradarisinvestigated.Besidesthe well-knownrangemigrationproblem,anotherconcernforwidebandsignalsistheDOAestimation prob-lem.Forthefirsttime,bothproblemsareconsideredinthispapersimultaneously,wherethewideband DOAistransformedintoasecond-ordercouplingsystemsimilartotherangemigrationproblembyusing thepropertyoftheFMCWsignal.Anovelembeddedcompensationapproachtoeliminatethecoupling termscausedbyrangemigrationandwidebandDOAisproposedand2Dmultiplesignalclassification(2D MUSIC)algorithmissubsequentlyappliedwithdynamicnoisesubspacetojointestimationofDoppler and DOA.Further,to reducethe computationalload caused bymultiple eigendecompositionsof large matrices,efficientimplementationmethodsareproposedandtheirperformanceinspeed,accuracyand robustnessiscompared.Theperformanceoftheproposedmethodsisvalidatedbythenumerical simu-lationsandiscomparedwithKeystoneMUSIC.Finally,itisshownthatforasmallnumberoftargets,the Rayleigh-Ritzisthemostefficientapproachamongthem.

© 2019ElsevierB.V.Allrightsreserved.

1. Introduction

Detection and localization ofmoving targetsare importantin many fields such as automotive radar [1], ground moving target indication(GMTI)[2],underwateracousticarray[3].Themost im-portant parameters of moving targets are the range, azimuthal (and,in3Dspace,elevation)angle(orDOA)andvelocity.The tar-get range and angle together determine the location of a target. While the Doppler (along with range) velocity is determined in coherent radars by means of phase shift between chirps within the coherent processing interval (CPI), and DOA is determined from the phase shift of signals received by different antennas within antenna array. Both phase shifts can be easily measured separately using narrowband radar. Using the de-chirping tech-nique forFMCWradar [2,4],the receivedsignalsare transformed into multi-dimensionalcomplex sinusoids(whose phase depends on the fast-time - range, slow-time - Doppler velocity and ar-ray element- DOA).Then theestimationof targetsparametersis transformedintothefrequenciesestimationproblem.Byextending traditional single-frequency estimators to joint multiple frequen-cies estimators, such as matched filter, 2D-Capon [5], 2D-MUSIC

∗ _{Corresponding author.}

E-mail address: s.xu-4@tudelft.nl (S. Xu).

[6], 2D-ESPRIT [7] andsparse representation methods [8,9], joint range-Dopplerestimation algorithms have beendeveloped. These algorithmsperformwell undernarrow-band signalcondition. Tar-get movement causes, however, change in the target range dur-ingoneCPI(physically)andthecross-couplingsbetweenfast-time andslow-time (mathematically), which is called rangemigration orrangewalk in GMTI[2,10,11]. Thecross-couplingterms spread the Fourierspectrum and consequently lead to estimation errors fortheseclassicalgorithms:thelargerthesignalbandwidthorthe higherthetargetvelocity, thehighertheestimationerrorof con-ventionalmethods[10,12].

Recently,sincewidebandsignalsarewidelyusedduetothe de-mandofincreasingly higherrangeresolution,therangemigration problemhasattractedsignificantattention.Tosolvethetarget mi-gration problem,the relaxation-basedsuper-resolutionalgorithms have been proposed in [2,13,14] for multiple moving target fea-tureextractions.However,theyconsiderawide-bandapproachfor the range profile, while they assume a narrow-band approach for the steering vector. In [15], the authors present the iterative adaptivealgorithm(IAA)forjointmultipleparameters estimation, whichprovides super-resolution by iteratively calculatingthe co-variance matrix together with estimation results. In [16,17], IAA is extended to the wideband waveform case together with the rangemigrationproblem.However, IAAconsumesa hugeamount ofmemory andtime when the rawdata dimension is large and https://doi.org/10.1016/j.sigpro.2019.107259

the scanning area is divided into dense bins, which makes the algorithm impractical for real-time applications. The Keystone transformandmatchedfilterare used in[12,18] to eliminatethe rangewalkresidualandtheRadonFouriertransform(RFT)is pro-posedtoconsiderevenhigher-ordercouplingproblemsbylineor curvesearchinginthetimedomainin[10,19].Unfortunately,these approachesneedalargeamountofrawdatatodointerpolationor coherentintegration, therefore.They could not provide the same fine resolution as the super-resolution algorithms in [6,15]. Im-plementationoftheRFT alsorequiresa largeamountof comput-ingpowerforthelinesearchinginmulti-dimensionaldata.Some waveformdesignmethodsarealsoproposedtosolvetherange mi-gration problem [11], but thesealgorithms increase the system’s complexityandtheachievedresolutionisnotashighasthat ob-tainedbysuper-resolutionalgorithms.

Inadditiontotherangemigration,anotherlimitationfor (Ultra-)WBsignalinthecollocatedarrayormultiple-inputand multiple-output (MIMO) application is the DOA estimation. Although the range migration problem has been studied intensively, the cur-rent algorithms jointly dealing with range migration and DOA estimationfail toprovide agood solutionto widebandDOA esti-mationbysimplifyingthesignalmodelwithnarrowband DOA as-sumption[2,14,20]. The traditional DOA estimators are based on narrowbandassumptions bythe interferometry information,such asCapon, MUSIC,etc.Toapply thetraditional narrowband super-resolutionalgorithmforwideband cases,twoof themainstreams of wideband DOA estimation are proposed, namely the incoher-entsignal subspace method (ISSM) andthe coherent signal sub-spacemethod(CSSM). ISSM solves thisproblemvia a filterbank todecompose the array output into its independent narrowband components.Thenthesubspace-basedalgorithmisappliedtoeach narrowband output, andDOA estimates can be average in some way.However,eachofnarrowbandestimatesdoesnotfullyexploit thetotalemitterpowerandsome ofthenarrowbandcomponents mayhavealowsignal-to-noiseratio(SNR), andthefinalDOA es-timatesmaybe adversely affected byfew inaccurate narrowband estimates. CSSM combines the different narrowband signal sub-spaceintoasinglesignalsubspacethat obeysthenarrowband ar-ray model.Although itis shownin [21] that the performance of CSSM is superior to ISSM, the forming of focusing matrices and universalspatialcovariancematrix(USCM)canincrease the com-putationalcomplexitysignificantly.Inaddition,theaccuracyofthe focussingmatriceshighlydependsonandissensitive tothe pre-liminaryestimateofthetrueDOAs[22].Insome other communi-cationproblems,jointtime-of-arrival(TOA)andDOAestimationin impulseradio(IR)-UWBisstudied,unfortunately,theDOAsare es-timatedbythepulsedelaywhichisdecidedbythebandwidthand notsuitableforthecompactarray.AnotherpowerfultoolforDOA estimationisthetime-frequency(TF)-MUSIC[23,24] whichisused todealwithnon-stationarysourcesanditisalsoappliedfor wide-bandDOAestimationin asimilar wayasCSSM[25].However, in FMCWradar de-chirped signals (beat frequency signals) foreach antennaelement anda singleFMCW chirp behaveas “stationary sources”,soanapplicationofTF-MUSICtothemisnothelpful.

Although both rangemigration andwidebandDOAestimation are intensively studied separately, there are few articles address bothproblemssimultaneously.Inthispaper,aMUSIC-based algo-rithmis proposed for theproblemof jointDoppler and DOA es-timationusing an (ultra)wide-bandarray-basedradar considering bothrangemigrationandwidebandDOAissues.Therange migra-tionmodelhasbeenstudiedcomprehensivelyandpresentedasthe second-ordercouplingbetweenfast-timeandslow-time.Combine thefact that the steering vector is the function ofthe frequency ofwideband DOAand the frequency is thefunction offast-time inFMCWsignal, theconventionalCSSMandISSM canbeavoided bytransformingthesteeringvector intothefunction offast-time.

Thus, the wideband DOA problemis transformed into the inter-couplingbetweenthefast-timeandtheelementindicesanalogue totherangemigrationproblem.Bythistransform,bothrange mi-grationandwidebandDOAproblempresentascouplingtermsand can be eliminated in the same way.After establishing the signal model,theclassic 2D MUSIC-basedalgorithmforjointestimation ofDoppler andDOAis presented.Unfortunately, conventional2D MUSICalgorithm cannot correctlyestimate theparameters inthe presence ofthe couplings.Toeliminate the influenceof the cou-plingtermsforaccurateparametersestimation,aphase compensa-tionmethodisproposedforbothcouplingsofrangemigrationand widebandDOA.Thecompensationmethodneeds, however, multi-plelarge-sizematrixeigendecompositionswhichare computation-ally heavy. Therefore, two efficient implementations, namely the LanczosalgorithmandRayleigh–Ritzstep,areintroduced.We com-parethetwoproposedmethodswiththeinversemethod,whichis alsoageneralMUSICacceleratingapproachpresentedin[26].The advantages of the proposed techniques are shown via numerical simulations.

Therestofthe paperisorganisedasfollows.InSection 2,we establish the signal model of multiplemoving targets for (Ultra-)WBFMCWantennaarrayradar.InSection3,theclassic2DMUSIC isappliedtojointestimationofDopplerandDOA.Then,the com-pensationalgorithmisproposed.Theefficientimplementationsare introducedandcomparedinSection4.Simulationresultsare pre-sentedinSection5.

Notationsusedinthispaperareasfollows.Scalarsaredenoted by lower-case letters, vectors andmatrices are written as lower-case and uppercase bold-face letters, respectively. (· )T_{, (· )}H _and

(_{· )}∗denotetranspose,conjugatetransposeandcomplexconjugate of a vector or matrix, respectively.  and  represent the Kro-necker and Hadamard product, respectively.



x



gives the near-estintegerlessthan orequaltox. I M denotes theM× Midentity

matrix.Someotherspecialmatricesare 1 M=[1,1,...,1]T∈RM×1,

d M=[0,1,...,M− 1]T∈RM×1.Further,the symbolsE

()

,Tr(), ()

and_O

()

representtheexpectationcomputation,thetraceofa ma-trix, the real part extraction operation and asymptotic notation, respectively.

2. UWB FMCW antenna array model 2.1. Arraysignalmodel

Inthissection,thesignalmodelusingthemono-staticantenna array withonetransmitterandLreceiversisestablished. Assume Ipoint targetswithunknown range r ₌[R1,R2,...,RI], radial

ve-locity v =[

v

1,

v

2,...,

v

I]andangle

θ

=[

θ

1,

θ

2,...,

θ

I]arelocatedin

the observed far-field, where v and

θ

are expected to be jointly estimated. The radar transmits a group of chirpsduring one CPI withthechirp duration T0 andthe pulserepetition interval(PRI)

T.Anormalizedsinglechirpsignal withthebandwidthBhasthe form s0

(

t

)

=



ej2π(f0t+0.5μt2) _t∈_[0,_T 0

)

, 0 t∈

(

−∞,0]∧t∈[T0,∞

)

, (1) where f0 denotes the startingfrequencyand

μ

=

B T0

denotes the frequencymodulation rate. Theperiodic transmittedsignal is de-composedintofast-timedomaint andchirpnumberdomainm=



_Tt



as(2),wherem=0,1,...,M− 1,andMisthetotal number ofchirpsinoneCPI.

t=t +mT t ∈[0,T0

)

. (2)

Thentheperiodictransmittedsignalhas

Consider theithscattererintheobservationdomainwiththe ve-locityviandthe rangeRi,the roundtrip time delayofthis

scat-tereris

τ

i

(

m,t

)

= 2

(

Ri+

v

i

(

t +mT

))

c =

γ

i+ 2

v

i c

(

t +mT

)

, (4) where

γ

i= 2Ri

c istheinitialroundtripdelayoftheithtarget for thefirstchirpandcisthespeedoflight.Usingthe0thelementas areference, thereceivedsignalofithscattererbythelthelement canbewrittenas

r_i(l)

(

m,t

)

=

α

iexp

(

j

ϕ

_i(l)

)

s

(

t +mT−

τ

i

(

m,t

))

=

α

iexp

(

j

ϕ

i(l)

)

s0

(

t −

τ

i

(

m,t

))

=

α

iexp

(

j

ϕ

i(l)

)

exp[j2

πφ

i

(

m,t

)

],

witht ∈[

τ

i

(

m,t0

)

,T0

)

, (5)

wheret0=

τ

i

(

m,t0

)

=

τ

i

(

m,0

)

/

(

1− 2

v

i/c

)

is the roundtrip delay

andusually t0<<T0,the superscript(l) denotesthe lthelement,

l₌0_,1_,...,L_{− 1} denotesthe indicesof theelement andLis the total numberof theelements,

α

i is theconstant complex

ampli-tude of theithscatterer, exp

(

j

ϕ

_i(l)

)

denotes thephase delay rel-ativeto the0thelement, and2

πφ

_i(m,t ) is thephase ofthe re-ceived signal of the0th element, whichaccording to(1) hasthe form

φ

i

(

m,t

)

= fc

(

t −

τ

i

(

m,t

))

+0.5

μ

(

t −

τ

i

(

m,t

))

2, (6)

witht ∈[

τ

i

(

m,t0

)

,T0

)

.

From the phase of the received signal, the instantaneous fre-quencyofthereceivedsignalisextractedas

fi

(

m,t

)

= d

φ

i

(

m,t

)

dt = fc



1−d

τ

i

(

m,t

)

dt



+

μ(

t −

τ

i

(

m,t

))



1−d

τ

i

(

m,t

)

dt



≈ fc+

μ

t . (7)

Now the terms of the time delay are neglected since

τ

i(m,

t )<<T0 andvi<<c.Thenthephase delayofthelthelement is

givenby

ϕ

(l) i =2

π

fi

(

m,t

)

dl c sin

θ

i=2

π

(

fc+

μ

t

₎

dl c sin

θ

i, (8)

where

θ

i denotes the angle of the ith scatterer, fi(m, t ) denotes

the instantaneous frequencyandisobtainedfrom(7), andd_l de-notesthedistancebetweenthelthelementandthereference ele-ment,respectively.Itisseenin(8) thatthephasedelayisnotonly related to thelthelement, but alsothe fast-time t . In fact, it is verystraightforwardbecausethesteeringvectoristhefunctionof frequencyforwidebandDOAandthefrequencyisthefunctionof timeforFMCWsignal,therefore,thesteeringvectorisalsoa func-tionoftime.Then,thewidebandDOAcanbetransformedintoan additionalsecond-ordercouplingbetweentheindicesofelements andtheindicesoffast-time.Mostapproachesforjointestimation DOAandDoppler havefailedtoprovidea solutionforthis wide-bandDOAestimation[14,26,27].Besides,theconventionalISSMor CSSMforwidebandDOAcanbeavoidedbysolvingtheproblemof couplingterms.

In this paper, the targets are located in the far-field and the observation time in one CPI is very short, thus, the angles are assumed not changing in one CPI. Without losing generality, the uniformlydistributedlineararray(ULA)withomnidirectional ele-mentsisusedtoestablishthesignalmodelinthefollowing,where dl=ld anddistheinterspacebetweentheneighbouringelements.

The receivedsignal isthencross-correlatedwiththe transmit-ted signal andthe de-chirped signal ofthe ithscatterer received

bylthelementcanbewrittenas

r(_il)

(

m,t

)

s∗

(

m,t

)

=

α

iexp

(

j

ϕ

i(l)

)

exp



j2

πφ

i

(

m,t

)



s∗

(

m,t

)

=

α

iexp

(

j

ϕ

i(l)

)

× exp



− j2

π

(

fc

τ

i

(

m,t

)

− 0.5

μτ

i2

(

m,t

)

+

μ

t

τ

i

(

m,t

))



≈

α

iexp

(

j

ϕ

i(l)

)

exp



− j2

π



fc

τ

i

(

m,t

)

+

μ

t

τ

i

(

m,t

)



. (9)

SubstitutionofEq.(4) intheresultofEq.(9) yieldsthede-chirped dataz(_il)

(

m,t

)

givenby,

z_i(l)

(

m,t

)

=

α

iexp

(

j

ϕ

i(l)

)

× exp

− j2

π

fc

(

γ

i+ 2

v

i c

(

t +mT

))

+

μ

t

(

γ

i+ 2

v

i c

(

t +mT

))



≈

α

iexp

(

j

ϕ

( l) i

)

exp



− j2

π



μγ

it +fd,iTm +

μ

2

v

i c Tmt



, (10)

where fd,i=2

v

_cifc represents the Doppler frequency of the ith

scatterer andexp(j2

π

fc

γ

i) isabsorbed into

α

i (for simplicity,

α

i

isstillusedinthefollowing),andthetermsfd,it and

μ

2

v

i

c t

2 _are

omittedastheyareverysmall.Aswecanseefrom(10),thesignal modelisa2Dcomplexsinusoidwithacouplingtermbetweenthe fast-timet andtheslow-timenumberm.Ifthescatterersmigrate oneorseveralrangebinsinoneCPIas

v

i× MT

c

2B, (11)

where MT is one CPI formed by M chirps and vi is the velocity

of the ith targets.It is called range migration in GMTI and this couplingtermwilldecreasetheperformance ofestimation apply-ingclassicFouriertransform.Itisworthnotingthat,thedefinition ofconventional rangemigration,whichshownasacouplingterm inthe signal model,is basedon the Rayleigh criterion.However, theresolutionofsubspace-basedmethodshasbrokentheRayleigh criterionandthe couplingtermswill always decreasethe perfor-manceofthesemethodsmoreor lesseven ifthetargetsmigrate lessthanonerangeresolutioncellinoneCPI.

Afterthe de-chirping,thedatain Eq.(10) issampledwith re-specttofast-timewithfrequencyfsandthediscretizeddatazˆ(_il)in

thetimedomainisobtainedas ˆ z_i(l)

(

m,k

)

=

α

ia(l)

(

θ

i

)

× exp

− j2

π



−

μ

kld c fssin

θ

i+

μγ

i k fs (12) + f_d,iTm+

μ

2

v

i c fs Tmk



, where a(l)

(

θ

i

)

=exp



j2

π

ld

λ

sin

θ

i



, (13)

and

λ

=c/ fc is the wavelength of the starting frequency, k

de-notesthesamplingindex(k₌0_,1_...,K_{− 1,}inwhichKisthetotal numberofsnapshots inone chirp).Bystacking all the zˆ(l)

₍

_m,k

₎

_,

the raw data ofthe lthelement is written in the matrixformat Z (_il)_∈CM×K_conform

Z_i(l)=

α

ia(l)

(

θ

i

)(

1M[h_i(l)]T

)



(

fd,ifTr,i

)



i

, (14)

where h (_il)∈CK×1_{, f r}

h(_il)=

1,ej2πμc fldssinθi,...,ej2πμ (K−1)ld c fs sinθi



T , f_r,i=

1,e− j2πμγfsi,...,e− j2πμ γi fs(K−1)



T , fd,i=



1,e− j2πfd,iT,...,_e− j2πfd,iT(M−1)



T, ₍₁₅₎ and



i∈CM×Kisdenotedby

i

=

⎛

⎜

⎜

⎝

gT i

(

0

)

gT i

(

1

)

. . . gT i

(

M− 1

)

⎞

⎟

⎟

⎠

, (16) inwhich g i

(

m

)

∈CK×1isgivenby gi

(

m

)

=[1,e− j2πμ 2vi c fsTm,...,e− j2πμ 2vi c fsTm(K−1)]T. (17)

The coupling term 1 M[h i(l)]T in Eq. (14) betweenelement

in-dices and fast-time sampling indices is introduced by wideband modulatedfrequencies. Ifwe consider a narrowband signal, then thesignalmodelreducesto

Z(_il)=

α

ia(l)

(

θ

i

)

fd,ifTr,i. (18)

However, the coupling terms are in general too large to be ne-glectedfor(Ultra-)WBsignals,whichwe willconsider inthe fol-lowing sections.The received signal model X (l)_∈CM×K _for

mul-tiplescatterersin the presenceof whiteGaussian noise is repre-sentedby X(l)₌ I  i=1 Z(_il)+N(l)_{, (19)}

where N (l)_∈CM×K _{denotes the Gaussian noise with distribution}

CN

(

0 ,

σ

2_I

₎

_.

2.2.Unambiguousangleandvelocity

According totheNyquistsamplingcriteria,inthesampled de-chirpeddatainEq.(13) thefollowingparametersareboundedby

2

π

Tfd,i∈

(

−

π

,

π

)

,

2

π

fc

d

csin

θ

i∈

(

−

π

,

π

)

, (20)

sotheunambiguousvelocityandangleareobtainedas

v

m= c 4Tfc, sin

(

θ

m

)

=min



_c 2fcd,1



. (21)

Ifthevelocitiesoranglesaregreaterthantheunambiguousones, theywillbefoldedintotheunambiguousdomain.

3. 2D MUSIC algorithm and compensation method

Inthissection,theclassic2DMUSICalgorithmforjoint estima-tionofDOA andDoppler ispresented atfirst.However, the cou-plingterms decline the performance of theclassic 2D MUSIC, in ordertocircumventthisanovelcompensationmethodisproposed intheMUSICalgorithm toremove suchinterference.The estima-tionofthemodelorderisdiscussedatlastinthissection. 3.1.2DMUSICalgorithm

With thethree-dimensional signal model,it ispossible to ap-plythe MUSICalgorithm forjointparameterestimation ifwe ig-norethecouplingterms.Usingonedimensionofsinusoidaldataas

reference, the2D MUSICalgorithm canbe implemented forjoint two-dimensional parameter estimation. The 3D MUSIC algorithm canbefurtherusedforjointthree-dimensionalparameter estima-tion forDOA, rangeand Doppler.The noise subspace can be ex-tractedbyapplyingaspatialsmoothingtechniquetoeliminate co-herencebetweenthesources[28] orbyapplyingahighorder sin-gularvaluedecomposition(HOSVD)[29].However,itisboth time-andmemory-consumingtodirectlyapplythe3DMUSICalgorithm. Thus,the2DMUSICalgorithmisappliedhere,forinstance,to es-timateDoppler andDOA jointly.It isworth notingthat, the pro-posedmethodscanalsobeappliedforjointestimationofDoppler andrangeorDOAandrange.Toapplythe2DMUSICalgorithm,the rawdatahavetobereshapedfromthe3-dimensionaltensorform to the2-dimensionalmatrix form Y ∈CLM×K _{by stacking element}

andslow-timedimensionstogetheras

Y=

⎛

⎜

⎜

⎜

⎝

X(0) X(1) . . . X(L−1)

⎞

⎟

⎟

⎟

⎠

. (22)

Forsimplicity, Y isrewritteninmatrixnotationas:

Y= I  i

α

i

(

a_θifd,i

)(

fr,i

)

T_

_dr,i

_

_

θr,i+N, (23) where a _θ i∈C L×1_,

_

dr,i∈CLM×Kand



θr,i∈CLM×K aregivenby,

a_θi=[1,a(1)

₍

_θ

i

)

,...,a(L−1)

(

θ

i

)

]T,

dr

,i=1L

i

,



θr,i=

⎛

⎜

⎜

⎜

⎜

⎝

(

h_i(0)

)

T

(

h_i(1)

)

T . . .

(

h(_iL−1)

)

T

⎞

⎟

⎟

⎟

⎟

⎠

1M. (24)

Now,theclassic2D MUSICalgorithmisapplieddirectlyby ignor-ingthe couplingterms.First, thecovariancematrix R _∈CLM×LM _is

computedaccordingto

R=E

(

YYH

₎

_{. (25)}

Theeigendecompositionisappliedtosplitthedataspaceintothe noisesubspaceassociatedtothenoiseeigenvectors U nandthe

sig-nalsubspaceassociatedtothesignaleigenvectors U s.

R=U



UH,

U=[UsUn]. (26)

Toextractthenoisesubspaceweassumethatthenumberof scat-terersisknown.Theestimationofthenumberofscattererswillbe discussedlater. Thematchedsteering vector

α

(

v

p,

θ

q

)

∈CLM×1 for

thevelocityvpandtheangle

θ

qisformulatedas:

α

(

v

p,

θ

q

)

=a_θqfd,p. (27)

Afterthat,theMUSICspectrumatthepoint(vp,

θ

q)canbe

calcu-latedby P

(

v

p,

θ

q

)

= 1

α

H

(

_v

p,

θ

q

)

UnUnH

α

(

v

p,

θ

q

)

. (28)

3.2. Compensationforcouplingterms

Directly applying the classic MUSIC algorithm without any phasecompensationyieldsan estimationperformancethat is sig-nificantlylessaccurateduetotheinfluenceofthecouplingterms. Hence,phaseadjustmentisneededbeforetheMUSICalgorithmis

applied.AlthoughtheKeystonetransformisthemostcommon ap-proachforthecouplingtermadjustment,theinterpolationofthe Keystonetransformleadstosignificantphaseerrorswhenthedata size is small [18,30]. Despite this drawback, the performance of Keystonetransformwillbediscussedandcomparedwiththe pro-posedalgorithminSection5.Fortunately,sincethecouplingterms are functionsofvfor



_dr,iand

θ

for



_θr,i,weare abletoremove thecouplingtermsineachscanninggrid.Thecompensationterm forthegridintermsof(vp,

θ

q)isformulatedas

C=

(

dr,p



θr,q

)

∗. (29)

ThenbyHadamardproductwiththerawdatamatrix Y yields

ˆ

Y=YC. (30)

Since the compensationtermis justaphase shift,it will not in-crease thenoise power. The couplingterms ofthe new obtained dataareremovedforthegridintermsof(vp,

θ

q).Withthis

com-pensation, not only the phase is adjusted to improve the accu-racy, butalso the orthogonality betweenthe steering vector and thenoise subspaceisenhancedwhichhelpstoimprovethe reso-lution.Thecovariancematrixiscalculatedusingtheimproveddata ˆ

Y according to

ˆ

R=E

(

YˆYˆH

₎

_{. (31)}

Finally, the 2D MUSIC algorithm (28) can be applied to the improved covariance matrix. The algorithm is concluded in

Algorithm1.

Algorithm 1 2DMUSICwithCompensation. 1: Reshapetherawdataas(22)

2: for

v

pin[−

v

m,

v

m] do 3: for

θ

qin[−

θ

m,

θ

m] do 4: C :=

(

dr,p



θr,q

)

∗ 5: Y ˆ := Y  C 6: R ˆ :=E

(

_Y ˆ _Y ˆH

₎

7: R ˆ _=:U



U −1 #Eigendecomposition 8: U n:=U [:,I:end] 9:

α

(v

p,

θ

q

)

:=a _θq fd,p 10: P

(

v

p,

θ

q

)

:= 1

α

H

₍

_v

p,

θ

q

)

U nU nH

α

(

v

p,

θ

q

)

11: endfor 12: endfor

3.3.Estimationofthetargetnumber

Before implementingthe 2D MUSIC algorithm, thenumber of targetshas tobe estimatedto correctlysplit the noise subspace. Apparently, the coupling terms in the UWB signal model bring difficultiestoestimatethenumberofthetargets,sincethe eigen-values decrease more smoothly than that of narrowband data. Therefore,amethodisproposedfortheestimationofthenumber oftargetsinvolved.Usually,thenumberoftargetsismuchsmaller

Fig. 1. Time consumption comparison for (a) calculating the different number of eigenvectors with the same dimension of the Hermitian matrix of 768 × 768 and (b) calculating 10 eigenvectors with the different dimension of the Hermitian matrix.

Fig. 3. Angle-Velocity maps of B = 1 GHz and SNR = 3 dB for (a) 2D MUSIC without phase compensation, (b) 2D MUSIC after Keystone transform, (c) 2D MUSIC with phase compensation and accelerated by Lanczos algorithm, (d) 2D MUSIC with phase compensation and accelerated by Rayleigh–Ritz step and (e) 2D MUSIC with phase compen- sation and accelerated by inverse algorithm.

thanthedimensionofthecovariancematrix,sothedimensionof thesignalsubspaceisallowedtobeslightlyoverestimated. Accord-ingtothisproperty,alargermodelorderthanthetrueonecanbe selectedatfirsttoimagetheMUSICpseudo-spectrum.Thenby us-ingpeaksdetectionmethods,weestimate thenumberofthe tar-getsfrom theMUSIC pseudo-spectrum.Although it isallowed to assumethelargerdimensionofthesignalsubspace thanthetrue one, it provides better imagingresults by using the whole noise subspace.After weobtainthenumberoftargets,theproposed al-gorithmscanbe appliedtoobtainbetterestimations.The simula-tionsofsuchamethodwillbeshowninSection5.3.

4. Efficient implementation

In this section, the efficient implementations of the proposed methodare proposed by accelerating the extractionof the noise subspaceandparallelimplementation.

4.1. Efficientimplementationforthenoisesubspaceextraction As for 2D MUSIC, the dimension of the covariance matrix LM_{× LM} is usually very large. Thus, it is a heavy computational burdentoperform allthe eigendecompositionsforeach scanning gridandlimitstheproposedmethods forrealapplications. Fortu-nately, as the covariance is a Hermitian matrix, some properties of the algorithm allow opportunities to accelerate the algorithm. Thefirstoneisthatthenumberoftargetsisusuallymuchsmaller than thedimension ofthe covariancematrix.Instead of calculat-ing all eigenvectors, onecan only calculatethe needed eigenvec-torsinthesignalsubspace,whilethenoisesubspacecanbeeasily obtainedfromtheorthogonal complementsubspaceofthesignal subspaceaccordingto,

Fig. 4. Comparison of RMSEs with CRLB as a function of SNR at B = 4 GHz for one target at angle 40 ◦_{and velocity 8 m/s. (a) angle estimation and (b) velocity estimation.}

Another interesting property is that the compensation term is just a minorphase shift. Therefore, thesignal subspaces in each adjacentgridareclosetoeachother inthenorm.Based onthese properties, two acceleration methods are introduced, namely the LanczosalgorithmandtheRayleigh-Ritzstep.Theinverseapproach in[26] isalsomentionedforcomparison.

4.1.1. Lanczosalgorithm

TheLanczosalgorithmisaniterativemethodforcalculatingthe eigendecomposition of large Hermitian/symmetric matrix [31]. It saves a lot of computationby only computing the largest eigen-valuesandtheir corresponding eigenvectors.Thus, itcanbe used inscenarioswhereonlythesignalsubspaceisneededandthe di-mensionof thesignal subspace ismuch smallerthan that ofthe covariancematrix.

4.1.2. Rayleigh–Ritzstep

Lanczos ismuchfastertoextract thesignal subspacethanthe default eig function if the dimension of the signal subspace is small.However, it isnot fastenough andwe do nottake advan-tageofthefactthattheadjacentsignalsubspacesareclosetoeach other.Thesignalsubspacesforneighbouringgridsareclosetoeach otherinthenormsotheprevioussignalsubspaceprovidesagood initialguesstocalculatethenextone.Thus,theRayleigh–Ritzstep method[32] isadoptedtousetheprevioussignalsubspaceasthe initialguesstoapproachthecurrentsignalsubspace.Accordingto thesimulation,justonestepisneededtoobtainasufficientlygood eigenvectorapproximation.

4.1.3. InversealgorithmwithoutEVD

In [26], the authors propose to usethe inverse ofthe covari-ance matrix to replace the noise subspace. Certainly, calculating the inverse or pseudo-inverse of a huge matrix will save a lot oftimecomparedtocalculatingtheeigendecomposition.It shows comparableresults withMUSICbutis fasterthan eigendecompo-sition MUSIC. However, this method can only work in high SNR condition. If the SNR is low, the approximation of this method no longer valid. Besides, the convergence performance of matrix inverse/pseudo-inverse is not monotonically increasing with the snapshot/SNR [33]. Thus this algorithm is not stable and robust. Despiteitsdisadvantages,fromthecomputationalburdenand es-timationperformanceperspectives,weusethisalgorithmasa ref-erencetocompareitwithourproposedalgorithms.

4.1.4. Comparison

To compare the priority of selected methods on time con-sumption,we firstfixthedimensionofthecovariancematrixand

simulate the time consumption with the different number of eigenvectors associated withthe largest eigenvalues. The simula-tionresultswithameantimeof100repeatsusingPython3.5with SciPy0.19 underInter(R) Core i5-6500 @ 3.20GHz are shown in

Fig.1(a).(Itisworth notingthat theresultsusing MATLABcould bedifferent.)

Thetimeconsumptionisnotonlyinfluencedbythenumberof eigenvectors associated with the largest eigenvalues, but also by thedimension of thecovariance matrix.Thus, the comparisonof timeconsumptionwiththedifferentdimensionsofthecovariance matrixisshowninFig.1(b).Aswecansee,thetimeconsumption ofthedefaulteig functionincreasessignificantly withthe dimen-sionofthecovariancematrixasitneedsO

(

n3

₎

_{flops(floatnumber}

operations).TheRayleigh-Ritzstepshowstobethemostefficient methodamongallofthem.

The computational complexity is shown in Table 1, where n represents the dimension of the Hermitian matrix and

β

repre-sents the number of eigenvectors associated with the

β

largest eigenvalues.

Accordingtotheaboveanalysis,usingtheRayleigh-Ritzstepas anexample,thealgorithmcanbeillustratedasAlgorithm2.

Algorithm 2 CompensationalgorithmwithRayleigh–Ritzstep. 1: Reshapetherawdata Y as(22)

2: R :=E

(

YY H

₎

3: R =:U



U H #Eigendecomposition 4: U s:= U [:,0:I− 1] 5: for

v

pin[−

v

m,

v

m] do 6: Y ˆ :₌ Y _

(

_dr,p

)

∗ 7: for

θ

qin[−

θ

m,

θ

m] do 8: Y ˆ :₌Y ˆ _

(

_θr,q

)

∗ 9: R ˆ :=E

(

Y ˆ Y ˆ H

₎

10: Z :₌R ˆ U s 11: Z =:QP #QRdecomposition 12: H :=Q H_{R ˆ Q} 13: H =:F



F H _{#Eigendecomposition} 14: U s:= QF 15:

α

(v

p,

θ

q

)

:=a _θq fd,p 16: P

(v

p,

θ

q

)

:= 1

α

H

_(v

p,

θ

q

)(

I − UsU sH

)

α

(v

p,

θ

q

)

17: endfor 18: endfor

Fig. 5. Angle-Velocity maps of B = 4 GHz and SNR = 30 dB for (a) 2D MUSIC without phase compensation, (b) 2D MUSIC after Keystone transform, (c) 2D MUSIC with phase compensation and accelerated by Lanczos algorithm, (d) 2D MUSIC with phase compensation and accelerated by Rayleigh–Ritz step and (e) 2D MUSIC with phase compensation and accelerated by inverse algorithm.

4.2.Parallelprocessing

As the 2D MUSIC algorithm is a scanning process, it is pos-sible to divide the scanning domain into several parts related to CPU cores for parallel processing. We can process each part parallelly to fully utilize the hardware. Here, by using a thread-ing package in python3.5, we divide the scanning domain into 4 parts with identical size and by using covariance matrix size of 256× 256, the computational time with and without paral-lel computing are 25.8 s and 41.9 s, respectively. In the simu-lation,62% of the computational time is saved by using parallel processing.

5. Simulations

Inthissection,theperformanceoftheproposedmethodsis dis-cussed. As the coupling terms are related to the bandwidth, the performance oftheproposedmethods withdifferentbandwidths, i.e. 1GHz and 4GHz, will be simulated. The system parameters usedforsimulationareshowninTable2.

5.1. Bandwidth1GHz

We start with considering a case with a bandwidth of B₌ 1GHz, where the relative bandwidth is 1.3%. Fig. 2 shows the

Fig. 6. Angle-Velocity maps of B = 4 GHz and SNR = 20 dB for (a) 2D MUSIC without phase compensation, (b) 2D MUSIC after Keystone transform, (c) 2D MUSIC with phase compensation and accelerated by Lanczos algorithm, (d) 2D MUSIC with phase compensation and accelerated by Rayleigh–Ritz step and (e) 2D MUSIC with phase compensation and accelerated by inverse algorithm.

Table 1

Computational complexity.

Algorithm Computational complexity Default eig Func O(n3₎

Inverse O(1 3 n 3)

Lanczos O(βn2₎

Rayleigh-Ritz O(β2 n + _β3₎

root-mean-squareerrors(RMSEs)ofestimatesofDOAandDoppler ofasinglepointscattererwiththeradialvelocity8m/s,angle40◦ andrange80masafunctionoftheSNR.Theyarecomparedwith thecorrespondingCramér-Raolowerbounds(CRLB)(seeAppendix forCRLBderivation).TheRMSEsareobtainedfrom40MonteCarlo

Table 2

Parameters of system.

Parameters Value

Number of Chirps in one CPI 16 Number of Snapshots in one Chirp 32 Number of Antenna Elements 8 Starting Frequency 77 GHz Inter-element Distance 1.899 mm Chirp Repetition Interval 0.1 ms

Chirp Duration 0.09 ms

trials.As classic MUSICis abiasedestimator, theRMSEs willnot decreasewith an increase ofSNR. One can alsoobserve that the RMSEsoftheKeystoneMUSICdonotalways decreasealongwith

Fig. 7. (a) 2D MUSIC without phase compensation and (b) 2D MUSIC after Keystone transform (c) 2D MUSIC with phase compensation and accelerated by Rayleigh–Ritz step.

Table 3

Comparison of computational time. Algorithm Computational Time Default eig Func 152.68 s

Lanczos 18.75 s Rayleigh-Ritz 5.48 s Inverse 32.66 s Keystone MUSIC 1.01 s

theincreaseoftheSNR sincetheerrorintroducedbythe interpo-lationwilldominantlydecreasethe accuracyofthe estimationat highSNRcondition.

We further focus on the performance of the proposed meth-odsina challengingscenario.Tothisend, fourcloselypositioned pointscattererswiththesameangles

θ

=45◦andthesame ampli-tudes0dB,butcloseradialvelocities

v

=[4.6,5.68,6.86,7.91]m/s andrandomrangefrom100mto200mareset.TheSNRissetto 3dB.ThenormalizedresultsareshowninFig.3.Whiletheinverse methodandtheclassic MUSIC algorithmwithout phase compen-sationare not abletoseparate scatterersfromeach other,MUSIC withphase compensation(in both Lanczosand Rayleigh–Ritz ac-celerationimplementations) achieves at least-3 dB isolation be-tweenthescatterers. TheKeystone-MUSIC(Fig.3(b))showssome separationofthescattererswithlowspeed,butfast-moving scat-terersare notseparated.Table3 showsthetime consumptionsof thefourmethodsandthecompensatedMUSICalgorithmwith de-fault eig function, where the observation domain is divided into a100× 100grid. It isworth notingthat the superioritywith re-specttothe computationaltime ofthe proposed algorithmcould be more significant if the dimension of the covariance matrix is larger.

5.2. Bandwidth4GHz

Next,we increase the bandwidth from1GHz to4GHz, where therelative bandwidthis5.1%.Accordingto themodelfrom(23), the estimation accuracy will deteriorate as the bandwidth in-creases. The same point scatterer with radial velocity 8m/s, an-gle 40◦ and range80mis setforRMSEsimulation. Theobtained RMSEs results ofthe proposed methods and the competitors are compared with CRLB in Fig. 4. The error of the no compensa-tionmethodismuchlargerthanwhat weobtainedinthecaseof B=1GHz.

Thenextsimulationistotesttheabilityoftheproposed algo-rithm to detect a relatively weak target. According to the analy-sis, ifthereis astrong migratedtarget presentinan observation domain, the energy of this target will dominantly spread into several eigenvectors. Thus, the subspace corresponding to the relativelyweaktargetwillbeallocatedtonoisesubspace.Two tar-gets, one withrange 100m, angle20◦, velocity 8m/sand ampli-tude

α

=−10dB andanotherone withrange80m,angle40◦, ra-dialvelocity5m/sandamplitude

α

₌0dB_,areset.TheSNR isset to30dB.TheresultsareshowninFig.5,whereforimproved visi-bilitythe results are normalized.From Fig. 5(a),we can see that the weak target is missing in the classic MUSIC result without phase compensation.Bothtargetsare seenin theMUSIC pseudo-spectrum obtained by Keystone MUSIC algorithms, however, the peakscorrespondingto thetargetsarewiderthan inthe conven-tional MUSIC or proposed algorithms and their relative contrast withthebackgroundismuchsmallerinmagnitude.

Then 11 point scatterers with random angles from 0◦ to 50◦, random radial velocitiesfrom 0m/s to 9m/sand random ranges from 100m to 200m random

α

from−3dB to 0dB are set and

Fig. 8. Angle-Velocity maps of B = 4 GHz , SNR = 20 dB and 128 snapshots for (a) 2D MUSIC without phase compensation, (b) 2D MUSIC after Keystone transform, (c) 2D MUSIC with phase compensation and accelerated by Lanczos algorithm, (d) 2D MUSIC with phase compensation and accelerated by Rayleigh-Ritz step and (e) 2D MUSIC with phase compensation and accelerated by inverse algorithm.

the SNR is set to 20dB. The normalized results are presentedin

Fig. 6, where the dynamic range is limited to 20dB. From the angle-velocity map, one can conclude that the peaks of estima-tionwithoutphasecompensationarebiasedtowardshigher veloc-ities and widened (especially in the azimuthal domain) in com-parison withthat oftheproposed methods.Althoughthe accura-cies ofestimation are slightlybetter, the KeystoneMUSIC suffers fromapoorresolutionofcloselyspacedtargets(especiallyin az-imuthal domain). Atthe sametime, all threecompensation algo-rithmsdemonstratetheclearseparationofalltargetsandaccurate estimation oftheir parameters. Toshow the improvement ofthe resolution oftheproposed compensationmethod, an extra simu-lation usingthe same systemparameters of 6 random point tar-gets withlarge angles and velocitiesis implemented.The results without compensation, phase adjustment by Keystone transform and phase compensation by the proposed method are presented

inFig. 7.The threecloselypositioned targets are hardlyresolved fromFig.7(a)(b),whiletheyareclearlyresolvedinFig.7(c).

In the next simulation, we keep the same parameters as be-forewhile increase the snapshots from32 to 128.The results of the11pointtargetsfromtheprevioussimulationarepresentedin

Fig.8. The imagingperformance of the inversealgorithm signifi-cantlydegradeswhenwe increasesnapshotsto 128anda strong ghost target appears at the position

(v

=0,

θ

=0

)

, while the Rayleigh-RitzandLanczos algorithmsrevealsharperpeaksrelated to targets.It is notedthat the results agree with the simulation in[33].

5.3.Numberofthetargetsestimation

The11point scatterersfromthe previoussimulation are used againinthissimulationandtheSNR issetto20dB.Fig.9 shows

Fig. 9. The number of targets estimation from MUSIC pseudo-spectrum using dif- ferent dimensions of the signal subspace assumptions.

thenumberoftarget estimationsusingthedifferentdimensionof thesignal subspaceassumption.The connectedregion label algo-rithmisusedheretoobtainthenumberoftargetsfromthe bina-rizedMUSICspectrum.Herethethresholdforbinarizationissetto −7dB,whichisathirdofthemeanvalueofthenormalizedMUSIC spectrumofalargenumberoftargetassumptions.

Fig.10 showstheMUSICspectrumofRayleigh-Ritzmethodand classic MUSIC without compensationfor the incorrect dimension ofsignalsubspaceassumptions.Fromtheresults,bothclassic

MU-SIC and compensation MUSIC will miss targets ifthe dimension ofthesignalsubspaceisunderestimated.However,comparedwith

Fig. 6,our proposed compensationMUSIC hasa higher tolerance foroverestimatingthe dimension ofsignal subspacethan that of classicMUSIC.

6. Conclusions

In thispaper, a noveljoint Doppler-DOA estimationusing the UWB FMCW array-based radar for moving targets is proposed. Wepresentanaccuratesignal modelforradarreturns from mov-ing targets which at first takes both range migration and wide-band DOA into account. Using this model we introduce a mod-ified MUSIC algorithm to eliminate the influence of the inter-couplingterms.Themethod adjuststhephase oftherawdatain each scanninggridbefore eigendecompositionto improvethe ac-curacyoftheDopplerandDOAestimation.Moreover,we propose two efficientimplementations,namelya Lanczosalgorithm anda Rayleigh-Ritzstep,toreducethecomputationalburdenspecifically fortheproposedmethod.

By comparing RMSEs and CRLB of conventional MUSIC, Key-stoneMUSICandproposed algorithm(forthebandwidthof1GHz and 4GHz) via numerical simulations, we demonstrate that the phase compensation algorithm improves the accuracies of both Doppler andDOAestimation over theconventional andKeystone MUSICandtheaccuraciesoftheproposedalgorithmimprovewith SNR. For example,the accuracies of both Doppler and DOA esti-mationsareimprovedmorethan20dB forSNR ₌20dB inFig.4. AlthoughforSNRbelow-10dBKeystoneMUSIChasaccuracy sim-ilar to the proposed method, the resolution and overall contrast

Fig. 10. Angle-Velocity maps of B = 4 GHz , SNR = 20 dB and 32 snapshots for (a) 2D MUSIC without compensation assume 10 targets, (b) Rayleigh–Ritz step assume 10 targets, (c) 2D MUSIC without compensation assume 16 targets and (d) Rayleigh-Ritz step assume 16 targets.

of the MUSIC pseudo-spectrumare worse than by the algorithm proposed.Duetothephasecompensation,thealgorithmproposed alsoresolvestargetscloselyspacedinthevelocity-angulardomain, whicharenotresolvablebothwithconventionalandKeystone MU-SIC algorithms.Further,weshow thatthe proposedLanczos algo-rithm and Rayleigh–Ritz are more robust than the inverse algo-rithminoursimulations.Inaddition,theRayleigh–Ritzstepshows superioritywithrespecttocomputational timewhenthenumber of targets is much smaller than the dimension of the signal co-variance matrix and has a high tolerance foroverestimating the dimensionofthesignalsubspace.

Declaration of Competing Interest

Theauthorsdeclarethattheyhavenoknowncompeting finan-cialinterestsorpersonalrelationshipsthatcouldhaveappearedto influencetheworkreportedinthispaper.

Acknowledgment

TheauthorswouldliketothanktheChinaScholarshipCouncil (CSC) forthe funding support. The authors thank the editor and the anonymous referees fortheir valuable suggestions and com-ments.

Appendix: CRLB derivation

ToformulatetheCRBmatrix,wefirstreshapetherawdatainto thevectorform y ∈CLMK×1_as

y= I  i=1

α

i

(

aθifd,ifr,i

)



ω

dr,i

ω

θr,i+n, (33) where

ω

dr,i∈CLMK×1,

ω

θr,i∈CLMK×1 are

ω

dr,i=1L

⎛

⎜

⎜

⎝

gi

(

0

)

gi

(

1

)

. . . gi

(

M− 1

)

⎞

⎟

⎟

⎠

,

ω

θr,i=

⎛

⎜

⎜

⎜

⎝

1Mh(i0) 1Mh(i1) . . . 1Mh(iL−1)

⎞

⎟

⎟

⎟

⎠

. (34)

Let Q ∈CLMK×1_{bethenoisecovariancematrix,whichis}

Q=E

(

nnH

₎

₌

_σ

2_I

MLK, (35)

with

σ

2 _{being the variances of the noise. According to the}

ex-tended Slepian-Bangs’formula[34],theijthelement ofthefisher informationmatrix(FIM)hastheform:

{

FIM

}

i j=Tr



Q−1

∂

Q

∂

η

i Q−1

∂

Q

∂

η

j



+2



_

∂

y

∂

η

i



H Q−1



∂

y

∂

η

j



, (36) where

η

=[

θ

,v]T, (37) with

θ

and v beingthevectorsconsistingoftheDOAsandDoppler frequencies, respectively.

∂

y /

∂

ηi

denotes the derivative of y with respecttotheithparameterof

η

.NotethattheFIMisblock diag-onalsincetheparametersin Q areindependentofthosein

μ

and viceversa.Thus,theCRBmatrixforthemotionparameterscanbe calculatedfromthesecondtermontherightsideof(36).Let

χ

θ i = j

ξ

1,i

(

dLa_θi

)

fd,ifr,i

ω

dr,i

ω

θr,i + j

ξ

2,iaθifd,ifr,i

ω

dr,i[

ω

θr,i

(

dL1MdK

)

], (38)

χ

v i = j

ζ

1,iaθi

(

dMfd,i

)

fr,i

ω

dr,i

ω

θr,i + j

ζ

2,iaθifd,ifr,i[

(

1LdMdK

)



ω

dr,i]

ω

θr,i, (39) where

ξ

1,i=2

πα

ifcd ccos

θ

i,

ξ

2,i=2

πα

i

μ

d c fs cos

θ

i,

ζ

1= −4

πα

iT fc c,

ζ

2=−4

πα

iT

μ

c fs.Let G=[

χ

θ₁ ...

χ

θ I,

χ

v1 ...

χ

vI], (40)

thentheCRBmatrixfortheparametervector

η

isgivenby

CRB

(

η

)

=[2

(

GH_Q−1_G

₎

_]−1_{. (41)}

References

[1] V. Winkler, Range doppler detection for automotive FMCW radars, in: Proceed- ings of the European Microwave Conference, 2007, pp. 1445–1448, doi: 10.1109/ EUMC.2007.4405477 .

[2] N. Jiang, R. Wu, J. Li, Super resolution feature extraction of moving tar- gets, IEEE Trans. Aerosp. Electron. Syst. 37 (3) (2001) 781–793, doi: 10.1109/ 7.953236 .

[3] B. Xerri , J.-F. Cavassilas , B. Borloz , Passive tracking in underwater acoustic, Sig- nal Process. 82 (8) (2002) 1067–1085 .

[4] S. Liu , Z. Zeng , Y.D. Zhang , T. Fan , T. Shan , R. Tao , Automatic human fall detec- tion in fractional fourier domain for assisted living, in: Proceedings of the IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), IEEE, 2016, pp. 799–803 .

[5] K. Yamaguchi, M. Saito, K. Miyasaka, H. Matsue, Design and performance of a 24 GHz band FM-CW radar system and its application, in: Proceedings of the IEEE Asia Pacific Conference on Wireless and Mobile, 2014, pp. 226–231, doi: 10.1109/APWiMob.2014.6920270 .

[6] F. Belfiori , W. van Rossum , P. Hoogeboom , Application of 2D MUSIC algorithm to range-azimuth FMCW radar data, in: Proceedings of the Ninth European Radar Conference (EuRAD), 2012, pp. 242–245 .

[7] M.D. Zoltowski , K.T. Wong , Esprit-based 2-D direction finding with a sparse uniform array of electromagnetic vector sensors, IEEE Trans. Signal Process. 48 (8) (20 0 0) 2195–2204 .

[8] H. Al-Tous, I. Barhumi, N. Al-Dhahir, Atomic-norm for joint data recovery and narrow-band interference mitigation in OFDM systems, in: Proceedings of the IEEE Twenty-Seventh Annual International Symposium on Personal, Indoor, and Mobile Radio Communications (PIMRC), 2016, pp. 1–5, doi: 10.1109/PIMRC. 2016.7794716 .

[9] L. Zheng, X. Wang, Super-resolution delay-doppler estimation for OFDM pas- sive radar, IEEE Trans. Signal Process. 65 (9) (2017) 2197–2210, doi: 10.1109/ TSP.2017.2659650 .

[10] J. Yu, J. Xu, Y.N. Peng, X.G. Xia, Radon-fourier transform for radar target de- tection (III): optimality and fast implementations, IEEE Trans. Aerosp. Electron. Syst. 48 (2) (2012) 991–1004, doi: 10.1109/TAES.2012.6178044 .

[11] F. Roos, D. Ellenrieder, N. Appenrodt, J. Dickmann, C. Waldschmidt, Range migration compensation for chirp-sequence based radar, in: Proceedings of the German Microwave Conference (GeMiC), 2016, pp. 317–320, doi: 10.1109/ GEMIC.2016.7461620 .

[12] Z. Sun, X. Li, W. Yi, G. Cui, L. Kong, Detection of weak maneuvering target based on keystone transform and matched filtering process, Signal Process. 140 (2017) 127–138, doi: 10.1016/j.sigpro.2017.05.013 . http://www.sciencedirect. com/science/article/pii/S0165168417301846

[13] N. Jiang, J. Li, Multiple moving target feature extraction for airborne HRR radar, IEEE Trans. Aerosp. Electron. Syst. 37 (4) (2001) 1254–1266, doi: 10.1109/ 7.976963 .

[14] J. Li, G. Liu, N. Jiang, P. Stoica, Moving target feature extraction for airborne high-range resolution phased-array radar, IEEE Trans. Signal Process. 49 (2) (2001) 277–289, doi: 10.1109/78.902110 .

[15] W. Roberts, P. Stoica, J. Li, T. Yardibi, F.A. Sadjadi, Iterative adaptive approaches to MIMO radar imaging, IEEE J. Sel. Top. Signal Process. 4 (1) (2010) 5–20, doi: 10.1109/JSTSP.2009.2038964 .

[16] N. Petrov, F.L. Chevalier, Iterative adaptive approach for unambiguous wide- band radar target detection, in: Proceedings of the European Radar Conference (EuRAD), 2015, pp. 45–48, doi: 10.1109/EuRAD.2015.7346233 .

[17] N. Petrov, F.L. Chevalier, Fast implementation of iterative adaptive approach for wideband unambiguous radar detection, in: Proceedings of the Twenty-Third European Signal Processing Conference (EUSIPCO), 2015, pp. 1207–1211, doi: 10. 1109/EUSIPCO.2015.7362575 .

[18] T. Shan , S. Liu , Y.D. Zhang , M.G. Amin , R. Tao , Y. Feng , Efficient architecture and hardware implementation of coherent integration processor for digital video broadcast-based passive bistatic radar, IET Radar Sonar Navig. 10 (1) (2016) 97–106 .

[19] J. Xu, J. Yu, Y.N. Peng, X.G. Xia, Radon-fourier transform for radar target detec- tion (II): blind speed sidelobe suppression, IEEE Trans. Aerosp. Electron. Syst. 47 (4) (2011) 2473–2489, doi: 10.1109/TAES.2011.6034645 .

[20] D. Zoeke, A. Ziroff, Phase migration effects in moving target localization us- ing switched MIMO arrays, in: Proceedings of the European Radar Conference (EuRAD), 2015, pp. 85–88, doi: 10.1109/EuRAD.2015.7346243 .

[21] H. Wang, M. Kaveh, Coherent signal-subspace processing for the detection and estimation of angles of arrival of multiple wide-band sources, IEEE Trans. Acoust. Speech Signal Process. 33 (4) (1985) 823–831, doi: 10.1109/TASSP.1985. 1164667 .

[22] H. Hung, M. Kaveh, Focussing matrices for coherent signal-subspace pro- cessing, IEEE Trans. Acoust. Speech Signal Process. 36 (8) (1988) 1272–1281, doi: 10.1109/29.1655 .

[23] A. Belouchrani , M.G. Amin , Time-frequency MUSIC, IEEE Signal Process. Lett. 6 (5) (1999) 109–110 .

[24] Y. Zhang , M.G. Amin , Blind separation of nonstationary sources based on spa- tial time-frequency distributions, EURASIP J. Appl. Signal Process. 2006 (2006) 198–198 .

[25] A.B. Gershman , M.G. Amin , Wideband direction-of-arrival estimation of mul- tiple chirp signals using spatial time-frequency distributions, IEEE Signal Pro- cess. Lett. 7 (6) (20 0 0) 152–155 .

[26] D. Oh, J.H. Lee, Low-complexity range-azimuth FMCW radar sensor using joint angle and delay estimation without SVD and EVD, IEEE Sens. J. 15 (9) (2015) 4799–4811, doi: 10.1109/JSEN.2015.2428814 .

[27] S. Kim, D. Oh, J. Lee, Joint DFT-ESPRIT estimation for TOA and DOA in vehicle FMCW radars, IEEE Antennas Wirel. Propag. Lett. 14 (2015) 1710–1713, doi: 10. 1109/LAWP.2015.2420579 .

[28] S. Hong, X. Wan, H. Ke, Spatial difference smoothing for coherent sources location in MIMO radar, Signal Process. 109 (2015) 69–83, doi: 10.1016/j.sigpro.2014.09.037 . http://www.sciencedirect.com/science/article/ pii/S0165168414004642 .

[29] M. Haardt , F. Roemer , G. Del Galdo , Higher-order SVD-based subspace estima- tion to improve the parameter estimation accuracy in multidimensional har- monic retrieval problems, IEEE Trans. Signal Process. 56 (7) (2008) 3198–3213 . [30] M.A. Richards , The keystone transformation for correcting range migration in

range-doppler processing, Pulse 10 0 0 (2014) 1 .

[31] Y. Saad , Iterative Methods for Sparse Linear Systems, volume 82, SIAM, 2003 . [32] B.N. Parlett , The Symmetric Eigenvalue Problem, SIAM, 1998 .

[33] M. Steiner, K. Gerlach, Fast converging adaptive processor or a structured co- variance matrix, IEEE Trans. Aerosp. Electron. Syst. 36 (4) (20 0 0) 1115–1126, doi: 10.1109/7.892662 .

[34] P. Stoica , R.L. Moses , Spectral Analysis of Signals, Pearson/Prentice Hall Upper Saddle River, NJ, 2005 .