A HYPER HARMONIC RESOLUTION BY USING THE DISCRETE PROLATE SPHEROIDAL WAVE FUNCTIONS

A. KHENCHAF - J. SAILLARD I RESTE

( Institut de Recherche et d'Enseignement Superieur aux techniques de l'electronique 1

'Laboratoire SZHF ( Syst&mes et Signaux Hautes Frequences 1 La Chantrerie - Cp 3003 44087 NANTES CEDEX 03

FRANCE TB1 : (33) 40.68.30.74 & (33) 40.68.30.64

ABS"

Based on harmonic analysis method, discrete spheroidal wave functions are used to achieve a very high spectral resolution.

For an N samples vector of the measured signal, the classical Fourier transform (F.T) discrimination capability is 1/N, while in the method proposed here, and for the same number of samples (NI, the discrimination capability might attain the value of (2/P) with P is the dimension of F.F.T . Additionnal comparative study of others classical windows (rectangular, Hamming and Tsengl is made.

INTRODUCTION

This paper falls in the main frame of target recognition by its backscattered complex field observed in an experimental site. Let f denote the frequency of a monochromatic plane wave illuminating the target, (f2-fl) the bandwidth to be analysed and 6f the frequency increment. For a fixed transmission-reception polarisation, discrete values of the electric field H(f), for different target bearings. are measured.

Assuming the illuminated target (fig 1) can be represented by p points scatterers, all situated rouahlv in the same olane. then the backscattered

. I

fieid -H(fl, as a function of the incident frequency, can be written as : 4'11 P j@, j ( T)6kf

H(f) = 1 l%l e e k= 1

where lakl is the amplitude of scatterer k

@k

ak scatterer k distance-projection on

is the phase of scatterer k

radar-target direction

c speed of EM wave

f frequency of the incident wave

f

direction of

o+ lakl 'i I radar 6k I i . * 1 propagation

target FIGURE 1

CH2882-919010000-0526$1 .OO 0 (1 990 IEEE)

wave

(1)

The objective is to estimate the number p and its associated parameters, the distance ak and the complex amplitude of the kth scatterer .

DISCRETE PROLATE SPHEROIDAL WAVE FUNCTIONS

The discrete prolate spheroidal wave functions are results of the solutions of the wave equation :

A++k2$=0 (2)

These functions can also be defined as the solutions of the following Fredholm integral equation :

j-)inNIT(f-f' sinn(f-f, ) )Uk(N,W;f' )df'=h (N,W)Uk(N,W;f) (3) k

(k=O, . . . , N-1) and V f E R

where U ( 0 < W < -& 1 is the local pass-band, of the frequency considered, and is of order 1, N being the number of measured samples.These functions are arranged according to their eigenvalues :

N

0 <hN-1(N,W)<AN-2(N, W)<..*<h (N,W)<hO(N,W)<l (4)

There exist N eigenvalues, and we only use the [ZNUI spheroidal functions corresponding to the first IZNWI eigenvalues which are very close to unity. In this case the U (N,W;f), and especially UO(N.W;f), contain most of the energy in the band [-W,WI. The asymptotic formula of Slepian :

1

k

l-hk(N, W) 1 E k! (8NnWIk++ e-2NnW ( 5 )

gives a fraction of the total energy of the spectral window outside the main lobe (i.e. outside [-W.Wl), and are also solutions of the differential equation :

lvalide for ksI2NWI and large N I

1 dUk(!;W;f) 1 +[ [ cos2nf - cos2nW +

ProDerties

m m The functions Uk(N,W;fl are of double orthogonality i.e they are orthogonal in [-W,Wl :

W 7 1 1 Um(N,W5f)Uk(N.W;f) df = 6mk ( 7 )

1 1 and orthonormal in [-T :

526 IEEE INTERNATIONAL RADAR CONFERENCE

discrete sequences (VF)) are also of double orthogonality: orthogonal through [ -m,+m] and orthonormal throunh [O.N-ll :

1 /2

U (N,W;f) Uk(N,W;f) df = ( 8 )

mm The Fourier transform of the discrete spheroidal wave functions are called "discrete spheroidal sequences" and they are defined by :

I,/,

for any value of n, with k =0,1,. . . ,N-l. 1 if k is even

j if k is odd and ck = {

The second orthogonality of these functions leads to a second Fourier transform defined as : V(k)(N,W) =

1 /2 Uk(N,W;f) e-J2n'n-(N-1)'21 df (10) -+dl/2

for n.k = 0,1, . . . , N-1 . The dlscrete spheroidal wave function is the

result of the inverse Fourler transform of the dlscrete spheroldal sequence :

(k= 0.1,. . . ,N-l) Wlth the same notations used above, the sequences

V(k)(N,W) satisfy the following matrix equations :

where A I s an NxN matrix defined as :

sin2n(i-J)W A I J = n(i-J)

0 5 i,J

(13)

6 N-1

where H(N,W) Is an NxN trldiagonal symetric matrix deflned as :

Li (N - i ) if J=i-1 2

H(N,WIiJ = (+ - iI2cos 2nW if j=i (15)

- 4 ( i + l)(N-1-i) if J=i+l if IJ-11 > 1

i. J=O. I,. . . ,N-1 l o

The sequences V(k)(N. W) might be calculated using the equations Uk(N.W;f).

Simllar to the spheroidal wave functions, the

i.k = 0.1 ,..., N-1 . Consequently (fig.2) show the spectrum of first fixe sequences for 14=100, W=0.04 and k=O, 1,2.3.4,5.The k discrete spheroidal sequences has k zeros through [-W,Wl and for. low window orders there is a hlgher energy concentration near the main lobe.

k=2 -60.00

-~4o.oo-l~-t-+l.ct" 0.00 0.05 0.10 0 . 1 5 0 . 9 0

-6O.00-.

-10.00-

-100.00- I

FIGURE 2 Spectral windows corresponding to the prolate spheroidal wave function data windows.

GENERALITIES

The general Cramer spectral representation for a stationary process is :

J-1/2 in which dZ(f) is a zero-mean orthogonal increment

dZ(f) is related to the spectrum S(f) by the process.

relation 2 S(f)df = E(IdZ(f)I 1. Instead of using the time domain we are going to

use the frequency domain, by first considering the

D.F.T of the samples as :

and the data are retained by the 1.F.T :

527 IEEE INTERNATIONAL RADAR CONFERENCE

1 / 2

y(f) = ( 2 0 ) sinNn(f-v) dZ(v) si en( f -u ) I,,,

(which is the basic equation for the spectrum estimation 1.

This equation is considered as Fredholm integral equation and represents projection of dZ(f) on y(f) in the frequency domain.

After using spheroidal windows and Fourier transform, we obtain the coefficients :

W Yk(f) = hk(N,W) 1 /-yU*(N,W;u)y(f+ul dv (21)

Which might be expressed as a function of dZ as

1 /z follows :

yk(f)= Uk(N,W;E) dZ(E+f) (22) i,,, As an interesting alternative, yk(fl might be

obtained using equations (18) and (21):

which is an expression corresponding to the D.F.T of x(n) weighted by the discrete spheroidal sequences (acting as a data window).

Consider a random signal x(t) composed of the sum of deterministic complex exponentials m(t) and of gaussian noise x (t) (in the simulation we consider

that xb(t)=O then :

x(t) = E(x(t)) + xb(t) = m(t) + xb(l) (24)

j2nfmt with m(t) = 1 IJ e (25)

m=O (periodic process of discrete spectrum)

from which, if M(f) is the F.T. of m(tl then it can be written as :

m(t) = E(x(t)) = eJznrt M(f) df (261

JR jZnfmt since m(t) is the sum of complex sinusoids (e then :

P

from which : &aft

dZb(fl (28)

and by extending eq.(17) E(dZ(f)) becomes : P E(dZ(f)) = 1 pm 6(f-fm) df (29)

Since eq. (22) gives the coefficient yk(f) at the m=O

frequency f at the output of the K th filter then :

1/2 P

(301

from which, and considering an isolated component

(31)

(one spectral line) of frequency f, we have :

E(Y k (f)) = p(f) Uk(N.W;f-f.) k=O, . . . , K-1 avec K=[ZNW]

then yk(f) = p(f) Uk(N.W;f-f.) + b k (f) (32)

To estimate p ( f ) least square regression method will be used.

k=O,. . . .K-1

Minimizing the expression :

we obtain

Kf Yk(f)Uk(f-f.) j ( f ) = k=O K- 1 (34)

1 ((f-f.) k=O

then for a single spectral line of frequency f.

eq. (34) becomes : Kz 1

an expression which could be written as a function of filters of the even range ( since for odd k U (N,W;O) = U ( O ) = O )

k

To have an indication about the existence probability of a spectral line of frequency f, Fisher test of 2 and 2(k-1) degreerof freedom is applied on eq. ( 3 2 ) . gq. (32) might be written in vectorial form for f=f. as

Y(f.1 = p(f.1 U(0) + b ( f , ) (36) Accordingly for any spectral line of frequency f we have a similar expression :

with Y(f) known as a vector of I R K (yk(f) k=O, . . . , K-1) observable. Y(f) is the sum of two vectors, p(f)U(O), an unknown deterministic ( = E ( Y ( f ) ) ) and random vector b(f).

Y(f) = p(f) U(0) + b(f) (37)

Then F(f1 of Fisher test can be expressed as :

IIY(f) /12/2 F ( f ) =

IIY(f)Y(f) 1I2/2(K-l)

(381

whey PU : spectral line estimated power A[PU): the difference between the actual and the

estimated power.

528

IEEE INTERNATIONAL RADAR CONFERENCE

As can been seen, probability of finding such a spectral line is directly proportional to F(f).

Application

sinusoids. defined as : Consider a signal S(n), being the sum of three

1

2 1

3 1

S(n)

(The

N=64

'k fk r k

0.1 -60 0.1016 13

60 0.1953 25

-30 0.2031 26

Table (I) shows the parameters of eq.(39) for and for the three spectrums (k = 1, 2, 3).

n'= 128

FFT HAH. FFT TSENG 1024 1024

f, 0.1016 0.1016

a 0.0979 0.0944

-56.36 -59.97

FFT SPH. V R A I E 256

0.1016 0.1016

0.1009 0 . 1

-61.80 -60

TABLE I

Samp 1 es

FIGURE 3 SIGNAL S

60'ooT (dB)

40.00

20.00

0.00 A -20.00

-40.00

-60.00 .50

normalised frequency

FIGURE 4 Signal S through rectangular wlndow

60.00

40.00

20.00

O . O O t n -20.00

-40.00

-60.00 .so

normalised frequency

FIGURE 5 Signal S through Hamming window 0.00-

-20.001

-140.00+ I , I , I , I , 0.00 0 . 1 0 0'. 20 OI.30 0'. 40 O' .sO

normalised frequency

FIGURE 6 Signal S through Tseng window

-120.00

-1 40.00

normalised frequency

FIGURE 7 Monodimensional Fisher test

529 IEEE INTERNATIONAL RADAR CONFERENCE

PARAMETERS ESTIMATION OF NEIGHBOURED SPECTRAL LINES Extension

The possibility to discrimination between two successive spectral lines. Principles

Consider the signal x generated by the equation : k xk = psk + utk+ bk

k=O,. . . , K where K is the sample size

Using the method of least-square regression , and for the estimation of parameters p and U , equation ( 4 0 ) is written as :

x = ps + ut + b (41) (41)

f x’(xk) k=O,. . . ,K-1 { s=(sk) k=O, ..., K-1 t=( ‘k) k=O,. . . ,K-I b=(bk) k=O, ..., K-1

Now we have to search for the values of p and U which minimize the difference x-ps-ut between the signal (the process) and its deterministic part with respect to s and t.

The equations of projection with respect to s and t are :

(42)

(431

(x-ps-ut, s)=(x, s ) - plls112 - u(t. s ) = 0

(x-ps-ut, t)=(x, t) - p(s. t) - ullt112= 0 which results in system that allows the following solution : - (x, s)llt112 - (x, t)(t, s )

1 1 ~ 1 1 2 llt1I2 (x, t)11s112 - (s , t)(x, s )

11s112 11t112

P = (44) - ( s * t)(t. s )

and - U = (45)

- ( s , t)(t, s )

nn Applying this on eq.(30), by considering only and f2, it becomes :

(46)

and yk(f)=pUk(f-f2) + uUk(f-fl) + bk(f) (47)

and for f=f

two spectral lines of frequency f

E(yk(f)) = pUk(f-f2) + vUk(f-fl)

k=O,. . . , K-1 2

y k 2 (f )=pUk(0) + uUk(f2-fl) + bk(f2) (48)

k=O,. . . ,K-1 Eq. (48) is equivalent to eq. (40) .

Sk = Uk(0) t = U ( f - f ) k k 2 1

with { - 1 b k = b ( f ) k 2 A

The estimators p(fl,f2) of p and u(f f ) of U are 1’ 2 given by ( 4 4 ) and (45) which are combined as :

Since the spheroidal functions form an orthogonal base, then A might be pressed as :

As in the monodimensional case we can calculate Fisher test of 2 and 2(k-2) degrees of freedom to find the joint probability of two spectral lines exlstence.

As mentioned previously, eq.(49) might take the form :

x = Y p + b (52) and Fisher test as :

II.112/2 (2-1 F(fl,f2) = (53) 11 xi 11 2/2 ( K-2 1

resulting in : F(f.,f,) =

From eq.(54) we can see that the existence probability of a couple of spectral lines (f .f is

directly related to F. This test is calculated for

Application Considering the signal S(nl defined previously,

(fig.8) shows the graphical representation of Fisher test F(f,Af) as a function of f and Af wlth f changing between L0.175 ; 0.2251 and Af between [0;0.051. In this figure we notice clearly the presence of two peaks at f = 0.1953 and Af’ = 0.0078 ( i . e f2=0.2031).

(fig. 9) shows the graphical representation of F(f f 1 Max with respect to f and Af. Searching for close spectral lines Fisher test was calculated for the normalized frequency window of LO. 175; 0.2251.

We notice that increasing the dimensions of FFT result obtained is a better separation of the peaks.

Using the rectangular Hamming and Tseng windows, only the first spectral line (of frequency 0.1016) can be disclosed while the two other lines are mixed.

1 2

f2 > fl.

1’ 2

F(f,Af)

1‘

a.175 II

d f l

FIGURE 8 13idimensionnal Fisher test

530

IEEE INTERNATIONAL RADAR CONFERENCE

normalised frequency

(FFT-128)

a

a

9,

$2

93

knpl i tuda 1 .oo

T I 1

1

1 0.9997

-60 -~~ -61.82

1.9014 0.9076

1.9014 0.9076

-56.35 -59.97

14.98 15 60 60.05

14.98 15 -30 -30.06

0.60 IIi :"t : I I 1;I j\ I

0.00 0 . 1 7 0 .18 0 . 1 9 0 . 2 0 0.21..22

normalised frequency (FFT=1024)

FIGURE 9 Max. of F as a function of f

Table 111 shows a certain error in amplitude and phase of the parameters calculated using the two windows mentioned above compared with the values calculated by the method proposed.

HAHHINC TSENC

0.1 0.1009

CONCLUSION

By classical methods, discrimination of point scatterers is rather critical. In this paper, with a harmonic-analysis based method, discrimination capability, using spheroidal windows. was improved to 2/P (P is the FFT dimension) instead of 1/N (N Is the sample size) for any dimension of F.F.T. By this method it is also possible to pinpoint the spectral lines and estimate with excellent precision all the parameters of point scatterer signal with noise. This method might also be used for high resolution spectrum analysis of output signal of any receiver when such resolution is required. A satisfying results has been obtained by off-line processing of measured data.

ACKNOWLEDGMENT The authors want to thank the C.E.L.A.R (Centre d'Electronique de 1'ARmement 35170 Bruz FRANCE) for its help during this study.

REFERENCES

111 SAILLARD J., COATANHAY J.L, GADENNE Ph. "Extraction des caracteristiques d'une cible radar ii partir d'une analyse fr.4quentielle.Extension de la m6thode"ll &me colloque GRETSI Nice, juin 1987. p. 340.

12 I THOMSON D. "Spectrum estimation and harmonic ana1ysis"Proc. IEEE, vol. 70, n" 9, 1982, p. 1055-1096.

13 J TRICOMI F. C. "Integral Equations".New-York : Wiley-Interscience, 1957.

141 SMITHIES F. "Integral Equations"New-York : Cambridge

Univ. Press, 3962.

151 TITCHMARSH E.C.

Oxford Univ. Press, 1962.

16 I BELLANGER M.

"Eigenfunction Expansions" I & 11. New-York :

"Traitement num.4rique du signal, theorie et prat ique"ed. Masson.

171 BOULEAU N. "Processus stochastique et applications" Hermann.

181 SLEPIAN D. "Prolate spheroidal wave functions, Fourier Analysis, and uncertainty. The discrete case". Bell system Technical journal, vol. 57, p. 1371-1428, 1978.

19 I FUNG I. TSENG "A novel window for harmonic analysis" IEEE transactions on acoustics,speech,and signal processing, vol. ASSP-29, No. 2, April 1981.

53 1

IEEE INTERNATIONAL RADAR CONFERENCE