Method and device to determine the frequency spectrum of a signal

The present invention pertains to the analysis of the frequency spectrum of a signal by a high-resolution method that is more precise than the Fourier transform when all that is available is a small number of samples of the signal. It pertains more particularly to a method for the autoregressive modelling frequency spectral analysis of a signal, known as Burg's maximum entropy method, and consists of a regularization of this method by means of a criterion similar to that already known to have been used in order to regularize the method of autoregressive spectral analysis known as the least error squares method, while at the same time keeping the advantages of Burg's maximum entropy method which are that it can be used in real time with a lattice type computation structure, that it can be extended to multisegment configurations and that it is robust with respect to computation noise, errors of quantization and rounding-off operations.

BACKGROUND OF THE INVENTION 
1. Field of the Invention 
The present invention relates to the analysis of the frequency spectrum of 
a signal by a high resolution method that is more precise than the Fourier 
transform when all that is available is a small number of samples of the 
signal. It relates more particularly to a method for the autoregressive 
modelling analysis of the frequency spectrum of a signal by the method 
known as the Burg's maximum entropy method which has been the subject of 
articles, among them: 1! Burg J. P., "Maximum entropy spectral analysis", 
in Proc. 37th Annual Inter. Meeting Society of Exploration Geophysicists 
(Oklahoma City, Okla.), Oct. 31, 1967. This method can be used in a wide 
variety of fields, especially in radar for Doppler filtering and the 
detection of a rupture, namely an event that does not follow the 
modelling. 
2. Description of the Prior Art 
In the method of autoregressive modelling analysis a complex signal x(t) 
represented by a sequence of complex samples {x.sub.n } is modelled by 
means of a law of prediction that can be used to deduce the value of a 
sample x.sub.n from a linear combination of the values of the n-p 
preceding samples of the sequence. This entails searching for a 
relationship of forward prediction having the form: 
##EQU1## 
where p is the order of the model, namely the number of previous samples 
of the sequence taken into account in the prediction, {a.sub.p,k } is a 
sequence of p complex coefficients defining the model and 
.epsilon..sub.f,n,p is a forward prediction error, by the p order model, 
of the sample x.sub.n, which it should be possible to liken to a white 
noise if the model is faithful. 
The modelling is said to be autoregressive because, when there is no noise 
and there are no errors of measurements on the samples, the higher the 
order p of the model, the more faithful is the modelling. In practice, the 
noises and measurement errors cause the modelling to diverge so that it is 
sought to limit the p order modelling to the lowest value for which the 
modelling error can be likened to a white noise. 
The modelling can be used to obtain knowledge of the frequency spectrum of 
a regularly sampled signal x(t). Indeed, it is possible, from the 
relationship (1), to deduce the relationship: 
##EQU2## 
with a.sub.p,0 =1 which by Fourier transform of its two members, becomes: 
##EQU3## 
with: 
##EQU4## 
We deduce therefrom: 
##EQU5## 
Assuming that the model is reliable, namely that .epsilon..sub.f,n,p is a 
white noise, it can be assumed that: 
EQU .vertline.E.sub.f,n,p (f).vertline..sup.2 =P.sub.f,n,p 
where P.sub.f,n,p is the energy of the forward prediction error of the p 
order model for the sample x.sub.n. The following is then got: 
##EQU6## 
which is the expression of the energy of the frequency spectrum of the 
sampled signal x(t) and which shows that this energy is deduced from the 
values of the sequence of the coefficients {a.sub.p,k } of the model. 
The autoregressive model defined by the sequence of the coefficients 
{a.sub.p,k } of the law of forward prediction of the sample x.sub.n as a 
function of the n-p previous samples also has the property, which shall be 
used hereinafter, of corresponding to a law of backward prediction of the 
sample x.sub.n-p as a function of the following n-p. samples. For, we also 
have: 
##EQU7## 
where * designates the conjugate operator and .epsilon..sub.b,n,p 
designates a error of backward prediction of the p order model on the 
value of the sample x.sub.n-p. 
The methods of autoregressive modelling spectral analysis raise the problem 
of determining the sequence of the complex coefficients {a.sub.p,k } 
defining the model adopted. 
A known approach is that of the least error squares method which is based 
on a minimizing of the sum of the squares of the forward prediction errors 
.epsilon..sub.f,n,p for n varying from 0 to N-p+1, the samples x.sub.n 
being assumed to be zero beyond the interval 1, N!. According to the 
relationship (2), we have: 
##EQU8## 
To minimize this sum of error squares, it is necessary to cancel its 
derivatives with respect to the sequence {a.sub.p,k } of the coefficients 
of the model. This is expressed by the system of linear equations: 
##EQU9## 
where again since a.sub.p,0 is equal to one by definition: 
##EQU10## 
Assuming: 
##EQU11## 
T being the transposition operator, and: 
##EQU12## 
it is possible to rewrite the system of equations in matrix form and 
arrive at the Yule-Walker equation: 
EQU .OMEGA..sub.p A.sub.p =-C.sub.p ( 7) 
whence we deduce the value of the vector A.sub.p of the coefficients of the 
p order model: 
EQU A.sub.p =-.OMEGA..sub.p.sup.-1 C.sub.p 
-1 being the inversion operator. 
It can be seen that the least error squares method necessitates an 
estimation of a p order autocorrelation matrix .OMEGA..sub.p and its 
inversion. This amounts to a large number of computations. 
To reduce the number of computations, Levinson has sought a recursive 
formula on the order of the model. To do so, he has made use of the 
following operator: 
EQU V.sup.(-) =J.V* 
where J is the antidiagonal unit matrix. With this operator, the p order 
autocorrelation matrix .OMEGA..sub.p can be rewritten as a function of the 
p-1 order autocorrelation matrix .OMEGA..sub.p-1 as follows: 
##EQU13## 
where + is the transconjugate operator. Assigning the notation 
A.sub.p.sup.m to the m first coefficients of the p order model, it is 
possible to write: 
EQU A.sub.p =A.sub.p.sup.p-1,a.sub.p,p ! 
thus making it possible to put the Yule-Walker equation in the form: 
##EQU14## 
From this, a first matrix equation is derived: 
EQU .OMEGA..sub.p-1.A.sub.p.sup.p-1 +a.sub.p,p.C.sub.p-1.sup.(-) =-C.sub.p-1 
This first matrix equation is transformed into: 
EQU A.sub.p.sup.p-1 +a.sub.p,p..OMEGA..sub.p-1.sup.-1.C.sub.p-1.sup.(-) 
=-.OMEGA..sub.p-1.sup.-1.C.sub.p-1 
Observing, according to the relationship (7) and the definition of the 
operator (-) that: 
EQU A.sub.p-1 =-.OMEGA..sub.p-1.sup.-1.C.sub.p-1 
and that: 
EQU A.sub.p-1.sup.(-) =-.OMEGA..sub.p-1.sup.-1.C.sub.p-1.sup.(-) 
we get: 
EQU A.sub.p.sup.p-1 =A.sub.p-1 +a.sub.p,p.A.sub.p-1.sup.(-) 
Whence a recursive relationship on the coefficients of the models as a 
function of their order: 
EQU a.sub.p,k =a.sub.p-1,k +a.sub.p,p.a*.sub.p-1,p-k 
which makes it possible to determine the coefficient of a p order model 
from the coefficients of the immediately lower p-1 order model, provided 
that the last coefficient a.sub.p,p has been determined beforehand and, 
therefore, provided that a structure of recursive computation is used on 
the order for the computations of all the coefficients of a model except 
for the last one. Unfortunately, the determining of the last coefficient 
a.sub.p,p is a difficult task for it makes use of the remaining equations 
of the system of Yule-Walker equations. 
Burg's method of autoregressive modelling spectral analysis takes up 
Levinson's approach for determining the p-1 first coefficients of a p 
order prediction model. But it replaces the determining of the last 
coefficient a.sub.p,p by an estimation under a least error squares 
constraint using the sum of the energies of the forward and backward 
linear prediction errors. These energy values are given by: 
##EQU15## 
which, taking account of Levinson's recursive relationship: 
EQU a.sub.p,k =a.sub.p-1,k +.mu..sub.p.a*.sub.p-1,p-k 
where .mu..sub.p is an estimation of the coefficient a.sub.p,p called a 
reflection coefficient, leads to the recursive relationships: 
##EQU16## 
The mean of the sum of the error energy values is equal to: 
##EQU17## 
The coefficient .mu..sub.p+1 which minimizes this mean of the sum of 
energies is such that: 
##EQU18## 
which, by substituting the expressions of the forward and backward 
prediction errors at the order immediately below p in the formula of 
P.sub.p+1, leads to the relationship: 
##EQU19## 
To compute the forward and backward prediction error terms in the numerator 
and denominator of the previous relationship (9) for defining a reflection 
coefficient, the relationships (8) are applied between the prediction 
errors bringing into play a recurrence on the order of the model starting 
with the sample of the signal x.sub.n. This is done by means of a computer 
structure in the form of a lattice filter. In order to prevent the errors 
from diverging as the order position of the model increases, this lattice 
filter structure requires that the inverse filter should be stable. This 
means that the modulus of the reflection coefficient should be smaller 
than one. 
The choice of a lattice structure of computer rather than a transversal 
structure has the threefold advantage of simplified hardware 
implementation, greater robustness with regard to truncation noises and 
greater robustness with regard to measurement noises. 
A major application, in radar, of autoregressive model spectral analysis by 
Burg's maximum entropy method is in the discrimination of the spectral 
signal of clutter from the spectral signature of targets. This application 
requires the ability to compute the mean spatially on several range gates 
and several azimuths for the clutter is differentiated from the targets by 
a wide spatio-Doppler correlation. A segment of radar data elements is 
defined as a set of coherent recurrences in Doppler form coming from one 
and the same resolution cell of the radar. A spatial correlation of the 
data elements implies the use of data segments coming from resolution 
cells that are adjacent in distance and azimuth. The difficulty in the 
combination of these data segments lies in the temporal discontinuity 
between the different segments. Because of this discontinuity, samples 
cannot be simply recombined to give a single string of data elements taken 
as an input signal for Burg's algorithm. To resolve this problem, Haykin 
has proposed a multisegment version of Burg's algorithm in an article: 2! 
Haykin S., Currie B. W., Kesler S. B., "Maximum entropy spectral analysis 
of radar clutter", Proceedings of the IEEE, Vol. 70, No. 9, pp. 953-962, 
September 1982. This consists in modifying Burg's algorithm to define as 
many lattice filters as there are segments and, at each step, to define a 
common reflection coefficient computed to minimize the sum of the errors 
of each filter. Then, the procedure of estimation of the reflection 
coefficient is performed under a least error squares constraint on the 
mean of the energy values of the prediction errors on all the segments. 
The reflection coefficient is then reinjected into the different lattice 
structure filters. Let it be assumed that there are 1 disjointed segments 
of samples available for analysis: 
##EQU20## 
1 lattice filters are defined, having forward and backward linear 
prediction errors: 
##EQU21## 
with n=p+2, p+3, . . . , N.sub.i and i=1, 2, . . . , I. The 
spatial-arithmetic mean of the sum of the energies of error on the 1 
segments gives: 
##EQU22## 
By substituting the expressions of the forward and backward prediction 
errors at the immediately lower order p in the above formula and by 
writing that the sought value of the reflection coefficient .mu..sub.I,p+1 
cancels the expression: 
##EQU23## 
we get: 
##EQU24## 
In an autoregressive model spectral analysis using Burg's maximum entropy 
method, as in the least error squares method, the choice of the order of 
prediction model is a crucial problem that has direct consequences for the 
resolution of the spectral analysis. A low order will deliver a smooth 
frequency spectrum erasing existing peaks while an excessively high order 
will induce spurious peaks in the frequency spectrum. Many statistical 
criteria have been proposed to estimate the optimum order at which it is 
preferable to stop: 
criteria based on the energy of the prediction error (Akaike's, Hannan's 
and Rissanen's criteria), 
criteria based on the whiteness of the residue (statistical coat-stand 
test, measurement of the flattening of the spectrum), 
criteria based on an optimization of a compromise between a bias term and a 
variance term. 
These statistical criteria which are sometimes costly in terms of 
computations are limited to the case where the horizon of analysis is 
large. Indeed, they presuppose a normal distribution of the errors which, 
when it is not verified, especially in the case of a short horizon of 
analysis, makes the estimation underlying the order of the autoregressive 
model inconsistent. 
In certain applications, such as those related to radar, the uncertainty as 
regards the determining of the order of the autoregressive model may be 
great owing to a horizon of observation that is extremely short. 
To resolve this problem of estimation of the order of the prediction model 
of a method of autoregressive spectral analysis, Kitagawa and Gersch, in 
their articles, 
3! Kitagawa G., Gersch W., "A smoothness prior time-varying AR coefficient 
modelling of nonstationary covariance time series", IEEE Trans. AC-30, No. 
1, pp. 48-56, January 1985, 
4! Kitagawa G., Gersch W., "A smoothness prior long AR model method for 
spectral estimation", IEEE Trans. AC-30, No. 1, pp. 57-65, January 1985, 
have proposed an approach within the framework of the least error squares 
method known as the regularized least error squares method. This approach 
consists in adding a regularizing function to the criteria of the least 
error squares, tending to dictate a smooth variation on the energy 
spectrum deduced from the coefficients of the prediction model in order to 
enable the maximum order of the model to be reached without allowing 
spurious peaks to develop in the energy spectrum. Instead of minimizing 
the simple criterion of the least error squares J(A.sub.p) whose 
expression resumed from the relationship (3) is: 
##EQU25## 
a criterion is minimized: 
EQU I(A.sub.p)=J(A.sub.p)+R(A.sub.p) 
where R(A.sub.p) is a regularizing function deduced from the transfer 
function H.sub.p (f) of the filter corresponding to the law of prediction 
of a p order model: 
##EQU26## 
by means of the following relationship: 
##EQU27## 
Taking only the first two orders, we get: 
##EQU28## 
or again: 
##EQU29## 
The least error squares method applied to this criterion I(A.sub.p) gives 
the coefficients a.sub.p,k of the p order model as the solution of the 
system of p linear equations: 
##EQU30## 
which can be put in the form: 
##EQU31## 
Assigning the notation D to the diagonal matrix defined by: 
##EQU32## 
and resuming the vector and matrix notations of the relationships (4, 5 
and 6), it is possible to rewrite the system of equations giving the 
coefficients of the p order prediction model in the matrix form related to 
the Yule-Walker equation: 
EQU (.OMEGA..sub.p +D.sup.T.D)A.sub.p =-C.sub.p 
Whence: 
EQU A.sub.p =-(.OMEGA..sub.p +D.sup.T.D).sup.-1.C.sub.p 
It can be seen that the regularization is equivalent to adding 
regularization terms to the diagonal terms of the measured p order 
autocorrelation matrix .OMEGA..sub.p, before inverting it. This addition 
makes it possible to give this positive autocorrelation matrix a 
semi-defined nature so as to make it preserve a physical character that it 
may have lost owing to errors of measurement. 
To optimize the coefficients .lambda..sub.0 and .lambda..sub.1 in terms of 
maximum likelihood, it can be seen that the criteria of the regularized 
least error squares leads to the maximizing of the term: 
##EQU33## 
with an a posteriori law defined by the conditional distribution of the 
data elements: 
##EQU34## 
and an a priori law: 
##EQU35## 
The relationships being Gaussian, the likelihood of these parameters 
.lambda..sub.0 and .lambda..sub.1 is also Gaussian and is obtained by 
integration: 
EQU V(P.sub.p,.lambda..sub.0,.lambda..sub.1 
/x)=.function.f(x/A.sub.p,P.sub.p).f(A.sub.p / 
.lambda..sub.0,.lambda..sub.1,P.sub.p).dA.sub.p 
and the maximizing of LogV(P.sub.p,.lambda..sub.0,.lambda..sub.1 /x)! 
gives the optimum values of the coefficients .lambda..sub.0 and 
.lambda..sub.1. 
These parameters .lambda..sub.0 and .lambda..sub.1 are chosen as a function 
of the type of application for a given measurement noise but are 
insensitive to errors of a factor of 10. For example, for an application 
involving the determination of a Doppler spectrum in radar technology, the 
parameter .lambda..sub.0 may be taken to be a first order parameter and 
the parameter .lambda..sub.1 may be taken to be a 0.001 order parameter. 
The regularized least error squares method developed by Kitagawa and Gersch 
enables the choice of a high prediction model order with the constraint of 
spectral smoothness. However, it is not appropriate when the number of 
signal samples available is small. In this case, Burg's maximum entropy 
method is preferable for it is simpler to implement (with no matrix 
inversion but solely linear computations) and more robust in the face of 
computation noises and errors of quantization and rounding-off operations. 
To regularize Burg's maximum entropy method, Silverstein and Pimbley have 
proposed a minimum free energy approach in the article: 
5! Silverstein S. D, Pimbley J. M. "The minimum free energy regularization 
connection: Linear-MFE", Marple Press, 23 ACSSC, pp. 365-370, 1989. 
It has been seen earlier that the coefficients a.sub.p,k of a p order 
autoregressive prediction model and the coefficients a.sub.p-1,k of the 
immediately lower p-1 order autoregressive prediction model are related by 
Levinson's condition: 
##EQU36## 
with: 
##EQU37## 
For the transfer functions H.sub.p (p) and H.sub.p-1 (f) of the predictive 
filters associated with the models: 
##EQU38## 
this leads to a recursive relationship having the form: 
EQU H.sub.p (f)=H.sub.p-1 (f)+.mu..sub.p.e.sup.-j2 .pi.pf.H*.sub.p-1 (f) 
The idea is to modify the criterion minimized by Burg which is the mean of 
the forward and backward prediction error energies: 
##EQU39## 
by a corrective term L.sub.p : 
##EQU40## 
so that there is a new criterion to be minimized having the form: 
EQU E.sub.p =U.sub.p -.alpha.L.sub.p 
called free energy by reference to a physical system for which U.sub.p 
would represent the energy, L.sub.p would represent entropy and .alpha. 
would represent the effective temperature. The choice of the reflection 
coefficient .mu..sub.p minimizing this new criterion is such that we have: 
##EQU41## 
which leads to the resolution of a third degree equation with real 
coefficients in .gamma..sub.p having the form: 
EQU (1-.gamma..sub.p.sup.2).(.gamma..sub.p.G.sub.p +.vertline.D.sub.p 
.vertline.)=-2.alpha..gamma..sub.p ( 10) 
with 
##EQU42## 
It can be seen that for .alpha. equals zero, Burg's classic solution is 
found again. Indeed, the polynomial then has three roots: +1, -1 and 
-.vertline.D.sub.p .vertline./G.sub.p so that there are three possible 
values for the reflection coefficient .mu..sub.p 
##EQU43## 
Since it is desired that the modulus of .mu..sub.p should be strictly 
smaller than 1, only the last value is appropriate: 
##EQU44## 
Should .alpha..noteq.0, it can be shown that there is only one approach 
giving a reflection coefficient whose modulus is strictly smaller than 1. 
The regularizing of Burg's algorithm by the minimum free energy approach 
developed by Silverstein and Pimbley therefore consists in: computing the 
terms D.sub.p and G.sub.p on the basis of their definitions of the 
relationships (11), setting up the third degree equation in .gamma..sub.m 
of the relationship (10), computing the three roots of this equation, 
selecting the root .gamma..sub.m0 whose modulus is strictly smaller than 
one, computing the reflection coefficient .mu..sub.p on the basis of the 
relationship: 
##EQU45## 
and computing the coefficients a.sub.p,k of the model by the usual 
recursive relationship: 
EQU a.sub.p,k =a.sub.p-1,k +.mu..sub.p.a*.sub.p-1,p-k 
As compared with the non-regularized Burg's method, it requires an 
additional calculation consisting of a search for the roots of the third 
degree equation in .gamma..sub.m. 
The present invention is aimed at obtaining another method of regularizing 
Burg's method which is less restrictive from the computation point of 
view. For it has been seen that, had there not been this need for 
regularization, Burg's method as compared with the least errors squares 
method has a threefold advantage in that it can be implanted in real time 
at low computation cost, with a lattice structure, that it can be extended 
to the multisegment case introduced by Haykin for frequency filtering or 
for the temporal smoothing of spectral analysis and that it is robust with 
respect to computation noise, quantization errors and rounding-off 
operations. 
SUMMARY OF THE INVENTION 
An object of the present invention is a first method for determining the 
frequency spectrum of a sampled signal {x.sub.n } consisting in assessing 
the coefficients {a.sub.p,k } of a p order prediction model that 
correspond to a law of prediction of the nth sample x.sub.n from the n-p 
previous samples {x.sub.n-k }, with k varying from 1 to p, having the 
form: 
##EQU46## 
where .epsilon..sub.f,n,p is a forward prediction error, by the p order 
model, of the sample x.sub.n, said coefficients {a.sub.p,k } giving an 
estimation of the energy of the frequency spectrum of the signal sampled 
in the form: 
##EQU47## 
where a.sub.p,0 is equal to 1 and where P.sub.p is the energy of the 
forward prediction error of the p order model for the sample x.sub.n when 
the model is reliable, namely when said error is a white noise, on the 
basis of a recursive law on the p order of the model having the form: 
EQU a.sub.p,k =a.sub.p-1,k +.mu..sub.p.a*.sub.p-1,p-k 
* designating the conjugate operator, with: 
EQU a.sub.p,0 =1 
EQU a.sub.p-1,0 -1 
EQU k.epsilon.1, . . . p-1! 
.mu..sub.p being a reflection coefficient determined from the terms: 
##EQU48## 
where .epsilon..sub.f,n,p-1 designates the forward prediction error by the 
p-1 order prediction model on the sample x.sub.n and 
.epsilon..sub.b,n-1,p-1 designates the backward prediction error by the 
p-1 order model on the sample x.sub.n-p : 
##EQU49## 
This method is noteworthy in that the reflection coefficient .mu..sub.p is 
determined by a relationship having the form: 
##EQU50## 
where .lambda..sub.0 and .lambda..sub.1 are positive real coefficients. 
An object of the invention is also a second method for determining the 
frequency spectrum of the mean of a family of 1 sampled signals {x.sub.n 
}.sub.i consisting in assessing the coefficients {a.sub.p,k } of a p order 
prediction model that correspond to a law of prediction of an nth sample 
x.sub.i,n on the basis of the n-p previous samples {x.sub.i,n-k }, k 
varying from 1 to p, having the form: 
##EQU51## 
where .epsilon..sub.f,i,n,p is a forward prediction error, by the p order 
model, of the sample x.sub.i,n of the signal i, said coefficients 
{a.sub.p,k } giving an estimation of the mean of the energies of the 
frequency spectra of the sampled signals in the form: 
##EQU52## 
where a.sub.p,0 is equal to 1 and where P.sub.I,p is the mean of the 
energies of the forward prediction errors of the p order model for the nth 
samples x.sub.i,n of the 1 signals when the model is reliable, namely when 
said errors are white noises, on the basis of a recursive law on the p 
order of the model having the form: 
EQU a.sub.p,k =a.sub.p-1,k +.mu..sub.I,p.a*.sub.p-1,p-k 
* designating the conjugate operator: 
with: 
EQU a.sub.p,0 =1 
EQU a.sub.p-1,0 =1 
EQU k.epsilon.1, . . . p-1! 
.mu..sub.p being a reflection coefficient determined on the basis of the 
terms: 
##EQU53## 
where .epsilon..sub.f,n,p-1 designates the forward prediction error by the 
p-1 order prediction model on the nth sample of the i.sup.th signal and 
.epsilon..sub.b,n-1,p-1 designates the backward prediction error, by the 
p-1 order model, on the n-p.sup.th sample of the i.sup.th signal: 
##EQU54## 
This method is noteworthy in that the reflection coefficient .mu..sub.p is 
determined by a relationship having the form: 
##EQU55## 
where .lambda..sub.0 and .lambda..sub.1 are positive real coefficients. 
An object of the invention is also a device for implementing the first 
method referred to here above.

As recalled here above, the spectral analysis of a sampled signal {x.sub.n 
} sampled by means of an autoregressive model using Burg's maximum entropy 
method consists in determining the coefficient {a.sub.p,1, . . . a.sub.p,p 
} of a p order prediction model that corresponds to a law of prediction of 
the nth sample x.sub.n on the basis of the n-p previous samples {x.sub.n-k 
}, with k varying from 1 to p, having the form: 
##EQU56## 
where .epsilon..sub.f,n,p is a positively oriented prediction error by the 
p order model, of the sample x.sub.n on the basis of a recursive law on 
the p order of the model having the form: 
EQU a.sub.p,k =a.sub.p-1,k +.mu..sub.p.a*.sub.p-1,p-k (12) 
designating the conjugate operator with: 
##EQU57## 
.mu..sub.p being a reflection coefficient determined on the basis of the 
relationship: 
##EQU58## 
with: 
##EQU59## 
where .epsilon..sub.f,n,p -1 designates the forward prediction error, by 
the p-1 order prediction model, on the sample x.sub.n and 
.epsilon..sub.b,n-1,p-1 the backward prediction error, by the p-1 order 
model, on the sample x.sub.n-p : 
##EQU60## 
This method of determining the reflection coefficient .mu..sub.p is 
equivalent to seeking a maximum entropy, namely the minimizing of the 
term: 
##EQU61## 
The drawback of this method is that its reliability is highly sensitive to 
the choice of the order p of the prediction model: an excessively low 
order gives an estimation of the frequency spectrum that is far too rough 
and does not reproduce all the existing peaks while an excessively high 
order induces spurious peaks in this estimation of the frequency spectrum. 
To make it possible to obtain a maximum order without allowing the 
development of spurious peaks in the estimation of the spectrum, it is 
proposed to add, to the criterion to be minimized, the regularization 
functions used by Kitagawa and Gersch in the framework of the spectral 
analysis by the least errors square method. 
The criterion to be minimized becomes: 
##EQU62## 
wherein: 
##EQU63## 
where H.sub.p (f) is the transfer function of the predictive filter 
associated with the p order predictive model 
##EQU64## 
Hereinafter, the description shall be limited to a first order regularizing 
function namely, the criterion to be minimized that will be adopted will 
be the term: 
EQU E.sub.p =U.sub.p +.lambda..sub.0 R.sub.p,0 +.lambda..sub.1 R.sub.p,1 
The minimizing consists in cancelling the gradient of this term with 
reference to .mu..sub.p : 
##EQU65## 
which is also written as: 
EQU .gradient..sub..mu..sbsb.p E.sub.p =.gradient..sub..mu..sbsb.p U.sub.p 
+.lambda..sub.0 .gradient..sub..mu..sbsb.p R.sub.p,0 +.lambda..sub.1 
.gradient..sub..mu..sbsb.p R.sub.p,1 =0 (18) 
To explain this term, first of all the term .gradient..sub..mu.p U.sub.p is 
computed. According to the relationship (12) we have: 
##EQU66## 
Taking account of the recursive relationship on the forward and backward 
prediction errors resulting from the recursive relationship (12) on the 
coefficients of the prediction model: 
##EQU67## 
and using the definition of the following complex derivation: 
##EQU68## 
it is shown that: 
##EQU69## 
We then obtain: 
EQU .gradient..sub..mu..sbsb.up P.sub.p =.mu..sub.p.G.sub.p +D*.sub.p 
with: 
##EQU70## 
We then go to the computation of .gradient..sub..mu..sbsb.p R.sub.p,0. By 
assumption, there is: 
##EQU71## 
so that: 
##EQU72## 
Owing to the recursive relationship (12) on the coefficients of the 
prediction model and their relationship of definition (16), there is also 
a recursive relationship between the transfer functions of the 
corresponding predictive filters as a function of the order. This 
relationship is as follows: 
EQU H.sub.p (f)=H.sub.p-1 (f)+.mu..sub.p.e.sup.-j2 .pi.fp.H*.sub.p-1 (f)(20) 
Taking account of the fact that: 
EQU .vertline.H.sub.p (f).vertline..sup.2 =H.sub.p (f).H*.sub.p (f) 
it is possible to write: 
EQU .vertline.H.sub.p (f).vertline..sup.2 =H.sub.p-1 (f)+.mu..sub.p.e.sup.-j2 
.pi.fp.H*.sub.p-1 (f)!.H*.sub.p-1 (f)+.mu.*.sub.p.e.sup.j2 
.pi.fp.H.sub.p-1 (f)! 
which can be put in the form: 
EQU .vertline.H.sub.p (f).vertline..sup.2 =(1+.vertline..mu..sub.p 
.vertline..sup.2)..vertline.H.sub.p-1 (f).vertline..sup.2 
+(.mu..sub.p.Y*.sub.p +.mu.*.sub.p.Y.sub.p) 
with: 
EQU Y.sub.p =e.sup.j2 .pi.fp.H.sub.p-1.sup.2 (f) 
From this it is deduced that: 
EQU .gradient..sub..mu.p .vertline.H.sub.p (f).vertline..sup.2 
=2..mu..sub.p..vertline.H.sub.p-1 (f).vertline..sup.+2 2.e.sup.j2 
.pi.fp.H.sub.p-1.sup.2 (f) 
Observing that: 
EQU .gradient..sub..mu.d.sbsb.p H*.sub.p (f)=2.e.sup.j2 .pi.fp.H.sub.p-1 (f) 
a first relationship is obtained: 
EQU .gradient..sub..mu..sbsb.p .vertline.H.sub.p (f).vertline..sup.2 
=2..mu..sub.p..vertline.H.sub.p-1 (f).vertline..sup.2 
+.gradient..sub..mu..sbsb.p H*.sub.p (f).H.sub.p-1(f) (21) 
This relationship can be used to make a trace back to the expression of the 
gradient of R.sub.p,0. Indeed, by definition owing to the relationship 
(16), we have: 
##EQU73## 
Taking account of the relationship (19) for the definition of the gradient 
R.sub.p,0, the following is obtained by integration of the relationship 
(21): 
##EQU74## 
Now: 
##EQU75## 
so that, from the above equation, the following relationship is deduced: 
##EQU76## 
The term R.sub.p-1,0 may be developed as follows: 
##EQU77## 
again: 
##EQU78## 
whence: 
##EQU79## 
Using the latter expression in the relationship (22) we obtain: 
##EQU80## 
To continue the development of the relationship (18), the following term 
still has to be computed: 
##EQU81## 
Owing to the recursive relationship (20) on the transfer function of the 
predictive filter associated with the model: 
EQU H.sub.p (f)=H.sub.p-1 (f)+.mu..sub.p.e.sup.-j2 .pi.fp.H*.sub.p-1 (f) 
which can be noted as: 
EQU H.sub.p (f)=H.sub.p-1 (f)+.mu..sub.p.T.sub.p (f) 
with: 
EQU T.sub.p (f)=e.sup.-j2 .pi.fp.H*.sub.p-1 (f) (26) 
it can be written that: 
##EQU82## 
is also equal to: 
##EQU83## 
This is also written as: 
##EQU84## 
in assuming: 
##EQU85## 
Passing to the gradient, we obtain: 
##EQU86## 
To assess M.sub.p (f), it is observed that, according to the relationships 
of definition (27) and (26), we have: 
##EQU87## 
To assesss: 
##EQU88## 
it is observed that: 
##EQU89## 
for: 
##EQU90## 
so that it is possible to write: 
##EQU91## 
with: 
##EQU92## 
To compute the expression (25), it is necessary to integrate each of the 
terms of the relationship (28). For the first term, there is obtained: 
##EQU93## 
which is equal to: 
##EQU94## 
The integration of the other term of the relationship (28) gives: 
##EQU95## 
which is also written as: 
##EQU96## 
The synthetic expression of the gradient deduced therefrom: 
##EQU97## 
The terms R.sub.p-1,0 and R.sub.p-1,1 remain to be expressed. It has been 
seen, with the relationship (23) that: 
##EQU98## 
For the term R.sub.p-1,1 the operation is started again from its 
definition: 
##EQU99## 
By explaining H.sub.p-1 (f), the following is got: 
##EQU100## 
or again: 
##EQU101## 
The gradient of R.sub.p,1 then takes the form: 
##EQU102## 
Whence the final expression of the energy gradient to be minimized as a 
function of the coefficients of the p-1 order prediction model and of the 
p order reflection coefficient: 
##EQU103## 
Returning to the relationship (17) which expresses the fact that a 
reflection coefficient .mu..sub.p is sought such that: 
EQU .gradient..sub..mu..sbsb.p E.sub.p =0 
a new expression of the reflection coefficient of the regularized Burg's 
spectral analysis is deduced therefrom: 
##EQU104## 
which can be written in abridged form: 
##EQU105## 
with the definitions of the previous relationships (14): 
##EQU106## 
It is noted that the regularization consists in modifying the definition in 
the form of a fraction of the reflection coefficient: 
##EQU107## 
of Burg's spectral analysis method by adding the following expression to 
the numerator: 
##EQU108## 
and the following expression to the denominator: 
##EQU109## 
Since the multisegment version of Burg's method amounts only to replacing 
the terms D.sub.p and G.sub.p of the fractional expression of the 
reflection coefficient by the terms D.sub.I,p and G.sub.I,p : 
##EQU110## 
with: 
##EQU111## 
it is deduced therefrom that the same regularization applied to the 
multisegment version of Burg's method leads to the adoption of a 
reflection coefficient having the form: 
##EQU112## 
The multisegment version of Burg's algorithm in the case of filtering can 
be made more robust with respect to the possible presence of one or more 
targets in the estimation window. indeed, the estimation of the mean 
autoregressive filter given by Burg's multisegment algorithm is 
"desensitized" since a target of high intensity of several lower 
intensities occupies a relative non-negligible number of cells of the 
estimation window. Since the mean filter is mismatched with the clutter 
environment, the detection of the target at output of the filtering is 
deteriorated or even absent. 
To overcome this drawback, it is possible to develop a technique of 
statistical weighting of data elements coming from the cells as a function 
of a probability of a posteriori membership of these cells in the majority 
population, the mean autoregressive model of which is determined by the 
multisegment Burg algorithm. This technique leads to the following 
modification of the expression of the regularized multisegment reflection 
coefficient: 
##EQU113## 
with 
##EQU114## 
In the above expression, the weighting K(X.sub.i,p-1.sup.2) is expressed as 
follows: 
##EQU115## 
s.sub.j being equal to 1 if the cell i of the estimation window contains 
only clutter and 0 if the cell i of the estimation window contains a 
target. 
Furthermore, whatever the value of i, we have: 
EQU P(s.sub.i =1)=p.sub.1 et P(s.sub.i =0)=p.sub.0 
with 
EQU (p.sub.0 +p.sub.1)=1 
p.sub.1 being the a priori probability that a cell of the estimation window 
contains only clutter and p.sub.0 being the a priori probability that a 
cell of the estimation window contains a target. 
The variable X.sub.i,p-1.sup.2 depends on the variables of the p-1 order 
algorithm, either rough data elements if p=1 or forward and backward 
prediction errors. This variable follows a relationship of X.sup.2. 
The expression of the regularized statistically weighted multisegment 
reflection coefficient (30) is obtained by minimizing the following 
regularized energy: 
##EQU116## 
where: 
##EQU117## 
with 
##EQU118## 
and where, as seen here above: 
##EQU119## 
the minimizing consisting in cancelling the gradient of this regularized 
error energy term with reference to .mu..sub.I,p. 
##EQU120## 
The expression of U.sub.p.sup.stat represents the statistical expectation 
of the total energy (and no longer mean energy) of the forward and 
backward prediction errors of the majority population of the estimation 
window, conditional to the knowledge of the p-1 order variable 
X.sub.p-1.sup.2 =(X.sub.1,p-1.sup.2, . . . ,X.sub.I,p-1): 
EQU U.sub.p.sup.stat =E(P.sub.p /X.sub.p-1.sup.2) 
This expression is developed by applying the formula of the total 
probabilities: 
##EQU121## 
It is observed that s.sub.i,k corresponds to the cell k of the estimation 
window for the event C.sub.i and is shown that it is possible to reduce 
the operation to the following expression: 
##EQU122## 
Baye's relationship gives: 
EQU P(s.sub.j =1/X.sub.j,p-1.sup.2)=K(X.sub.j,p-1.sup.2) 
The formula of the new multisegment reflection coefficient given by the 
relationship (30) may be explained by developing the following 
considerations on the random variables conditioned by events: 
Let us assume a partition of an event C.sub.i of a probabilized space E and 
a real random variable X. The formula of the total probabilities applied 
to the event {X&lt;x} gives: 
##EQU123## 
The last-named relationship can be used to compute the mean estimation 
error energy and deduce therefrom an expression of the multisegment 
reflection coefficient that possesses a statistical weighting of the data 
elements coming from each cell as a function of the probability that this 
cell belongs to the clutter environment. 
To do this, the following notation is adopted: 
##EQU124## 
It can be seen that CardC!=2.sup.I and that s.sub.i,k corresponds to the 
cell k of the estimation window for the event C.sub.i. 
Then, it will be assumed that: 
EQU X=P.sub.p 
P.sub.p being the total forward and backward prediction error energy for 
the cells of the majority population in the estimation window 
EQU X.sub.p-1.sup.2 =(X.sub.1,p-1.sup.2, . . . ,X.sub.I,p-1.sup.2) 
with 
##EQU125## 
Applying the formula (32) to P.sub.p on condition that X.sub.p-1.sup.2 is 
known, we get: 
##EQU126## 
for 
##EQU127## 
with 
##EQU128## 
Taking into account the fact that: 
EQU P(X.andgate.C.sub.i /Z)=P(C.sub.i .andgate.Z). 
and replacing X by P.sub.p and Z by X.sub.p-1.sup.2 the operation is 9 
reduced, as with (32), to the expression (33). 
For the rest of development, Bayes' formula will be applied to P(C.sub.i 
/X.sub.p-1.sup.2): 
##EQU129## 
now 
##EQU130## 
so that 
##EQU131## 
To simplify the reading hereinafter, it will be taken that: 
##EQU132## 
It is possible to express as a synthetis: E(P.sub.p /C.sub.i 
.andgate.X.sub.p-1.sup.2). 
##EQU133## 
with 
##EQU134## 
Namely, reducing the operation to the expression (33): 
##EQU135## 
now 
##EQU136## 
for 
##EQU137## 
now as 
EQU P(G.orgate.H/K)=P(G/K)+P(H/K)-P(G.andgate.H/K) 
and 
##EQU138## 
then 
##EQU139## 
that is 
##EQU140## 
now 
##EQU141## 
for 
##EQU142## 
making it possible to write: 
##EQU143## 
Hereinafter, it will be preferred to write: 
##EQU144## 
On the basis of the expression (34), we come to the new expression: 
##EQU145## 
now 
##EQU146## 
Here, s.sub.j is the reference given to the cell j of the estimation 
window. 
It has been seen here above that: 
##EQU147## 
This results in a new writing for the expression (35): 
##EQU148## 
with 
##EQU149## 
By derivation of (36) with respect to .mu..sub.I,p, the expression of the 
reflection coefficient is found: 
##EQU150## 
with 
##EQU151## 
In the above expression K(X.sub.i,p-1.sup.2) is given by: 
##EQU152## 
P(s.sub.i =1)=p.sub.1 et P(s.sub.i =0)=p.sub.0 are a priori probabilities 
such that p.sub.0 +p.sub.1 =1. 
For the case p=0, we have: 
##EQU153## 
with 
##EQU154## 
For the cases (s.sub.i =1) and (p=0), namely the case of clutter in the 
cell, for example chaff clutter, there is a log-normal relationship for 
the amplitude: 
##EQU155## 
with .nu.=N-2: relationship of X.sup.2 having .nu. degrees of liberty 
##EQU156## 
for m.sub.p and .rho..sub.p defined by: 
##EQU157## 
and a Gaussian law in phase 
##EQU158## 
with .nu.=N-3: relationship of X.sup.2 having u degrees of liberty 
##EQU159## 
m.sub..PHI. and .sigma..sub..PHI. being defined by: 
##EQU160## 
For the case (s.sub.i =1) and (p.noteq.0), namely the case of clutter in 
the cell, the statistical weighting of the cells is given by: 
##EQU161## 
with .nu.=2.(N-p): relationship of X.sup.2 with .nu. degees of liberty 
##EQU162## 
For the case (s.sub.i =0) and (p=0), namely the case of the target in the 
cell, the assumption taken is that the energy borne by each component of 
the signal follows a centered normal relationship of variance 
.sigma..sub.c.sup.2 (mean theoretical energy re-emitted by the target). 
The law of distribution of amplitude is then: 
##EQU163## 
with .nu.=N: relationship of x.sup.2 with .nu. degrees of liberty 
##EQU164## 
while the law of distribution in phase is Gaussian: 
##EQU165## 
with .nu.=N-3: relationship of X.sup.2 having u degrees of liberty 
##EQU166## 
m.sub..PHI. and .sigma..sub..PHI. being defined by: 
##EQU167## 
For the case (s.sub.i =1) and (p.noteq.0), namely the case of the target in 
the cell, the assumption taken is that the target is not the same, in 
terms of its Doppler spectrum, as with the clutter (if the contrary were 
the case, the cell would be likened to a cluttered cell which would then 
no longer involve penalties). In this case, the outputs of the filter at 
each p order possess the characteristics of the target that has not been 
filtered. P(X.sub.i,p-1.sup.2 /s.sub.i =0) may then be broken down as for 
the case (p=0), into two relationships: a law of distribution in amplitude 
and a law of distribution in phase: 
EQU P(X.sub.i,p-1.sup.2 /s.sub.i =0)=.rho.(X.sub.i,p-1,.rho..sup.2 /s.sub.i 
=0)..PHI.(X.sub.i,p-1,.PHI..sup.2 /s.sub.i =0) 
The law of distribution in amplitude is equal to: 
##EQU168## 
with .nu.=2.(N-p): relationship of X.sup.2 with .nu. degrees of liberty 
##EQU169## 
The law of distribution in phase is Gaussian and is equal to: 
##EQU170## 
with: .nu.=2.(N-p-1)-2: relationship of the x.sup.2 with .nu. degrees of 
liberty. 
##EQU171## 
.epsilon..sub.f,.PHI. and .epsilon..sub.b,.PHI. being defined by: 
##EQU172## 
with 
##EQU173## 
and .sigma..sub..PHI. being defined by: 
##EQU174## 
with: 
##EQU175## 
The expression of the weighted multisegment reflection coefficient (37) can 
be regularized: 
##EQU176## 
with 
##EQU177## 
By derivation of the expression (38) with respect to .mu..sub.I,p, the 
statistically weighted and regularized multisegment reflection coefficient 
is found: 
##EQU178## 
with 
##EQU179## 
If this method of regularization of Burg's method is to be compared with 
Kitagawa's and Gersch's regularization method applied to the standard 
least error squares method, it can be observed, by resuming the reading of 
the relationship (4): 
##EQU180## 
that we have: 
##EQU181## 
for: 
##EQU182## 
so that the reflection coefficient of Burg's regularized method defined in 
the relationship (29) can also be written as follows: 
##EQU183## 
or, by using Kronecker's symbols: 
##EQU184## 
This can also be written as: 
##EQU185## 
with: 
##EQU186## 
The terms of an implicit autocorrelation matrix are recognized in the terms 
c'.sub.k-l : 
EQU .OMEGA.'.sub.p =.OMEGA..sub.p +D.sup.T.D: 
where .OMEGA..sub.p is the usual autocorrelation matrix of the samples of 
signals: 
EQU .OMEGA..sub.p ={c.sub.k-l }.sub.k,l.epsilon.(1,p) 
and D is the diagonal matrix: 
##EQU187## 
This implicit autocorrelation matrix corresponds precisely to the explicit 
autocorrelation matrix of the method of autoregressive predictive model 
spectral analysis obtained by the regularized least errors square method 
developed by Kitagawa and Gersch. 
This shows that the proposed new method for the regularization of Burg's 
maximum entropy method, provides a smoothing of the spectrum similar to 
that obtained with the Kitagawa's and Gersch's method of regularized least 
error squares while at the same time retaining the advantages related to 
Burg's maximum entropy method, namely the fact of having low computation 
costs and of being capable of being implanted in real time with a lattice 
structure, being extendable to the multisegment case introduced by Haykin 
in frequency filtering or temporal smoothing of spectral analysis and 
being robust with respect to computation noise, quantization errors and 
rounding-off operations. 
Compared with the method of regularization of Burg's maximum entropy method 
by Silverstein and Pimbley's free energy minimization method, it can be 
shown that the regularization obtained with the new proposed method of 
regularization is better. Indeed, the method of regularization of Burg's 
maximum entropy method by the minimization of free energy brings into play 
only the residual energy of the prediction errors and therefore 
regularizes the solution comprehensively by preventing the denominator of 
the reflection coefficient from getting cancelled while the new method of 
regularization proposed places a direct constraint on the flattening of 
the spectrum. 
To be certain of this, it may be recalled that the solution based on 
minimum free energy consists in making a search for the root, with a 
modulus smaller than 1, of the third degree equation with real 
coefficients in .gamma..sub.p of the relationship (10): 
EQU (1-.gamma..sub.p.sup.2).(.gamma..sub.p.G.sub.p +.vertline.D.sub.p 
.vertline.)=-2.alpha..gamma..sub.p 
D.sub.p and G.sub.p having their usual definitions of the relationships 
(10). This root can be put in the following form: 
##EQU188## 
which gives a reflection coefficient having the form: 
##EQU189## 
When p is close to the optimum order or is overestimated and when the 
prediction model becomes low or noisy, we have: 
##EQU190## 
where .sigma..sup.2 is a variance of white noise. 
Now, for .rho..sub..alpha. =0, namely without regularization, the 
computation noises b.sub.c modify the ratio of D.sub.p * and G.sub.p and 
therefore the reflection coefficient .mu..sub.p,: 
##EQU191## 
and introduce'spurious poles into the transfer function H.sub.p (z) of the 
associated predictive filter. 
By contrast, in the case of the minimum free energy (.rho..sub..alpha. 
.noteq.0), O.sub.2 (b.sub.c) becomes negligible in relation to 
.rho..sub..alpha. and thus the reflection coefficient .mu..sub.p 
approaches zero in overcoming the pollution introduced by the computation 
noise. This prevents problems of singularity which induce spurious poles 
in the transfer function H.sub.p (z) of the associated predictive filter 
(spurious peaks in the spectrum). 
It can then be noted that it is not necessary to resolve the third degree 
equation with real coefficients in .gamma..sub.p : 
EQU (1-.gamma..sub.p.sup.2).(.gamma..sub.p. G.sub.p +.vertline.D.sub.p 
.vertline.)=-2.alpha..gamma..sub.p 
but that it is enough to find a .rho..sub..alpha. that is great enough such 
that: 
##EQU192## 
is close to zero when p tends towards p.sub.opt. 
FIG. 1 shows a diagram of a computation circuit 1 for the computation of 
the complex coefficients {a.sub.p,0. . . ,a.sub.p,p } of a p order 
prediction model on the basis of the reflection coefficients {.mu..sub.1, 
. . . ,.mu..sub.p }, according to Burg's maximum entropy method which may 
or may not be regularized by the new method proposed. This computation 
circuit implements the relationship (12) of recurrence on the p order of 
the model: 
EQU a.sub.p,k =a.sub.p-1,k +.mu..sub.p.a*.sub.p-1,p-k with k .epsilon.1, . . . 
,p-1! (40) 
* designating the conjugate operator with: 
##EQU193## 
The first line which determines all the coefficients a.sub.k,o having an 
index zero of all the 0 to p order predictive models contains no operator 
since all the coefficients have a unit value. 
The following lines which determine the coefficients a.sub.p,k having an 
index k varying from 1 to p-1 enclose a succession of conjugate operators 
(*), multipliers (.times.) and summators (.SIGMA.) implementing the 
recursive relationship (30). 
For example, for the term a.sub.2,1 the recursive relationship (40) gives 
the definition: 
EQU a.sub.2,1 =a.sub.1,1 +.mu..sub.2.a.sub.1,1 * 
so that it is obtained by means of a two-input summator (.SIGMA.) 
receiving, at a first input, the term a.sub.1,1, this term a.sub.1,1 being 
none other than the reflection coefficient .mu..sub.1 and, at a second 
input, the term .mu..sub.2.a.sub.1,1 * delivered by a multiplier 
(.times.). This multiplier (.times.) receives firstly the term .mu..sub.2 
and secondly the conjugate term a.sub.1,1 or the conjugate term .mu..sub.1 
delivered by means of a conjugate operator (*) interposed between the 
input of the term .mu..sub.1 and the input of the multiplier (.times.). 
For the term a.sub.3,1 the recursive relationship (29) gives the definition 
: 
EQU a.sub.3,1 =a.sub.2,1 +.mu..sub.3.a.sub.2,2 * 
which is obtained by means of a summator (.SIGMA.), a multiplier (.times.) 
and a conjugate operator (*) that are appropriately connected. 
This is so up to the term a.sub.p,p which is taken to be equal to the 
reflection coefficient .mu..sub.p. 
FIG. 2 shows the drawing of a circuit for the computation of a reflection 
coefficient .mu..sub.m, with m varying from 1 to p, according to the new 
proposed method of regularization of Burg's maximum entropy method. 
At the bottom of FIG. 2, there is a lattice filter 2 giving the forward and 
backward prediction errors .epsilon..sub.f,n,m and .epsilon..sub.b,n,m, m 
varying from 1 to p, as a function of the value of the sample of the 
signal x.sub.n. This lattice filter 2 implements the recursive 
relationships on the order which, in Burg's maximum entropy method, link 
the forward and backward prediction errors of a p order model to the 
forward and backward prediction errors of the immediately lower p-1 order 
model: 
##EQU194## 
with: 
EQU .epsilon..sub.f,n,0 =.epsilon..sub.b,n,0 =X.sub.n 
This lattice filter 2 has a sequence of p identical stages defining two 
parallel channels, one on which there are available the forward prediction 
errors .epsilon..sub.f,n,m, with m varying from 1 to p, and the other on 
which there are available the backward prediction errors 
.epsilon..sub.b,n-1,m, with m varying from 1 to p. Each stage has a 
two-input summator (.SIGMA.) on the forward prediction errors channel and 
a delay cell (Z.sup.-1) on the backward prediction error channel, This 
delay cell delays the value of the error of a signal sample and is 
followed by a summator (.SIGMA.). Between these two channels, there are 
two multipliers (.times.) taking the crossed products. One is the 
multiplier of the forward prediction error available at input of the stage 
with the conjugate of the reflection coefficient .mu.*.sub.m, m being the 
order of the stage, to apply it to the summator (.SIGMA.) of the backward 
prediction error channel. The other is the multiplier of the backward 
prediction error available at output of the delay cell of the stage with 
the reflection coefficient .mu..sub.m, m being the order of the stage, to 
apply it to the summator (.SIGMA.) of the forward prediction error 
channel. 
In addition to the lattice filter 2, there are circuits 3 and 4 for the 
computation of the terms D*.sub.m and G.sub.m, with m varying from 1 to p, 
of the numerator and the denominator of the fraction for the definition of 
the reflection coefficient .mu..sub.m of Burg's maximum entropy method: 
##EQU195## 
The circuit 3 for the computation of the term D*.sub.m has a summator 30 
with N-m inputs connected to the outputs of the N-m last stages of two 
N-stage shift registers 31 and 32 by means of multiplier circuits 
(.times.) and conjugate operators (*). The shift register 31 is connected 
to the inputs of the summator 30 by conjugate operators (*) and multiplier 
circuits (.times.). It is connected by its series input, at the lattice 
filter 2, to the output of the m-1.sup.th stage on the forward prediction 
error channel. The shift register 32 is connected to the inputs of the 
summator 30 solely by the multiplier circuits (.times.) and is connected 
by its series input, at the lattice filter 2, between the delay circuit 
(Z.sup.-1) and the summator (.SIGMA.) of the m.sup.th stage, on the 
backward prediction error channel. A multiplier circuit (.times.) 
connected to the output of the summator 30 enables the application of a 
2/(N-m) weighting after the summation to obtain the term D*.sub.m. 
The circuit 4 for the computation of the term G.sub.m has a summator 40 
with 2(N-m) inputs connected to the outputs of the N-m last stages of two 
N-stage shift registers 41 and 42 by means of individual circuits for the 
computation of the square formed by a multiplier circuit (.times.) with 
two inputs combined with a conjugate operator (*) interposed in one of the 
inputs. One of the shift registers 41 is connected by its series input, at 
the lattice filter 2, to the output of the m-1.sup.the stage on the 
channel of the forward prediction error while the other shift register 42 
is connected by its series input, at the lattice filter 2, between the 
delay circuit (Z.sup.-1) and the summator (.SIGMA.) of the m.sup.th stage 
on the backward prediction error channel. A multiplier circuit (.times.) 
connected at output of the summator 40 enables the application of a 
1/(N-m) weighting after the summation to obtain the term G.sub.m. 
In addition to the computation circuits 1, 2, 3 and 4 which are in the 
devices for the real time implementation of the spectral analysis by the 
non-regularized Burg's maximum entropy method, there are specific 
computation circuits 5, 6 and 7 enabling the implementation of the new 
regularization method proposed. 
The computation circuits 5 enable the computation of the intermediate 
parameters .beta..sub.m-1,k while the computation circuits 6 and 7 enable 
the computation of the corrective regularization terms to be applied to 
the numerator and denominator of the fraction defining the reflection 
coefficient .mu..sub.m. 
The computation circuits 5 enable the computation of the intermediate 
parameters: 
EQU .beta..sub.m-1,k =.lambda..sub.0 +.lambda..sub.1 (2 .pi.).sup.2.(m-k).sup.2 
with k .epsilon.1, . . . ,m-1! 
They are each formed by a summator (.SIGMA.) with two inputs directly 
receiving the constant .lambda..sub.0 at a first input and the constant 
.lambda..sub.1 at a second input by means of a multiplier (.times.) taking 
the product of the constant .lambda..sub.1 by the coefficient (2 
.pi.).sup.2.(m-k).sup.2. 
The computation circuit 6 enables the computation of the corrective term: 
##EQU196## 
applied to the numerator of the quotient defining the reflection 
coefficient .mu..sub.m. It consists of a summator 60 with m-1 inputs, the 
k.sup.th input receiving the coefficient a.sub.m-1,k by means of two 
multipliers (.times.) arranged in sequence, one taking a product by the 
coefficient a.sub.m-1,m-k, the other taking the product by 
.beta..sub.m-1,k. A multiplier (.times.) placed at output of the summator 
60 enables the application of a weighting by 2 to arrive at the full 
expression of the corrective term of the numerator. 
The computation circuit 7 enables the computation of the corrective term: 
##EQU197## 
applied to the denominator of the quotient defining the reflection 
coefficient .mu..sub.m. It consists of of a summator 70 with m-1 inputs, 
the k.sup.th input receiving the coefficient a.sub.m-1,k by means of an 
individual circuit for the computation of the square constituted by a 
multiplier circuit (.times.) with two inputs combined with a conjugate 
operator (*) interposed in one of its inputs. A multiplier (.times.) 
placed at output of the summator 70 enables the application of a weighting 
by 2 to arrive at the complete expression of the corrective term of the 
denominator. 
The outputs of the computation circuits 3 and 6 at which there are 
available the term D*.sub.m of the numerator of the quotient for defining 
the reflection quotient .mu..sub.m in Burg's maximum entropy method and 
the corrective term of this numerator for the regularization are connected 
to two inputs of a summator 8 fitted out at output with a multiplier 
circuit (.times.) 9 operating a multiplication by -1 to take account of 
the minus sign in the quotient for defining the reflection coefficient. 
The outputs of the computation circuits 4 and 7 at which there are 
available the term G.sub.m of the denominator of the coefficient for 
defining the reflection coefficient .mu..sub.m in Burg's maximum entropy 
method and the corrective term of this denominator for the regularization 
are connected to the two inputs of a summator 10. 
To obtain the reflection coefficient .mu..sub.m, the regularized term of 
the numerator of its definition quotient, available at output of the 
multiplier circuit 9, is divided by the regularized term of the 
denominator of its definition quotient, available at output of the 
summator 10, by means of a two-input multiplier circuit 11, one of which 
is connected directly to the output of the multiplier circuit 9 while the 
other is connected to the output of the summator 10 by means of an 
inverter circuit 12. 
FIGS. 3, 4 and 5 show the different frequency spectra obtained by Burg's 
method respectively non-regularized and regularized by the free energy 
minimum method and regularized by the new method proposed, for a complex 
signal having two frequency lines from which, in the course of ten 
noise-infested draws, ten evenly distributed samples have been drawn. 
FIG. 3, which depicts the spectra obtained for the ten noise-infested draws 
by implementing the non-regularized Burg's method up to the maximum 
possible order of nine (10-1), shows that the different versions of the 
spectrum obtained for the different drawings do not overlap except for the 
two frequency lines that actually exist and are affected by high amplitude 
spurious lines that differ according to the version. The result thereof is 
a mean spectrum that is far too noise-infested to be capable of being 
exploited. 
FIG. 4, which shows the spectra obtained for the ten noise-infested drawing 
operations by implementing Burg's method up to the maximum possible order 
of nine, with a regularization by the minimum free energy method, taking 
account of a coefficient .alpha. equal to 1, shows that the different 
versions of the spectrum obtained for the different drawings always 
overlap for the two frequency lines that effectively exist and have 
spurious lines with far smaller amplitudes that tend to be superimposed. 
The resultant mean spectrum has become exploitable. However it is observed 
that the regularization has caused the loss of no small degree of 
amplitude in both of the really existing frequency lines. 
FIG. 5 shows the spectra obtained for the ten noise-infested draws 
implementing the Burg method up to the maximum possible order of nine, 
with regularization according to the new method proposed, taking account 
of a coefficient .lambda..sub.0 equal to 1 and a coefficient 
.lambda..sub.1 equal to 0.001. This figure shows that the different 
spectra henceforth have only the two actually existing lines with minor 
ripples. With respect to the spectra of FIG. 4, it is noted that this new 
method of regularization has prompted a lower attenuation of the two 
actually existing frequency lines and a better attenuation of the spurious 
frequency lines. It is therefore more efficient.