METHOD AND SYSTEM FOR ESTIMATING A QUANTITY REPRESENTATIVE OF SOUND ENERGY

Abstract

A method and associated system for estimating a quantity representative of the sound energy at at least one point of a three-dimensional space where a plurality of antennas are situated, each including at least K acoustic sensors, K being higher than or equal to 2, includes for each antenna of the plurality of antennas, production of a plurality of signals representative of the sound field at the antenna in question, for each antenna of the plurality of antennas, determination of a raw value of the quantity at the point based on at least K+1 elements of a matrix that are based respectively on pairwise combinations of representative signals produced by the antenna in question, and determination of an estimated value of the quantity at the point by combining the raw values of the quantity at the point determined respectively for the various antennas of the plurality of antennas.

Claims

1. A method for estimating a quantity representative of the sound energy at at least one point of a three-dimensional space where a plurality of arrays are located, each of the plurality of arrays including at least K acoustic sensors, K being greater than or equal to 2, the method comprising: for each array of the plurality of arrays, producing a plurality of signals representative of the sound field at the respective array; for each array of the plurality of arrays, determining a raw value of said quantity at said point based on at least K+1 elements of a matrix that are based respectively on pairwise combinations of representative signals produced by the respective array; and determining an estimated value of said quantity at said point by combining the respective raw values of said quantity at said point, determined for the respective arrays of the plurality of arrays.

2. The estimation method according to claim 1, wherein the determining the raw value for the respective array comprises: determining, based on said matrix, a directional value of the quantity representative of the sound energy received at the respective array from a direction connecting said point and the respective array, and determining the raw value for the respective array based on the determined directional value determined.

3. The estimation method according to claim 1, further comprising, for each array of the plurality of arrays, determining, based on said matrix, a plurality of directional values of the quantity representative of the sound energy received at the respective array, from a plurality of directions respectively.

4. The method according to claim 3, further comprising, for at least one array of the plurality of arrays, refining the directional values by a beamforming technique.

5. The method according to claim 3, further comprising, for each array of the plurality of arrays, determining raw values of said quantity at a plurality of points based on the determined directional values determined for the respective array.

6. The estimation method according to claim 5, further comprising, for each point of said plurality of points, determining an estimated value of said quantity at the respective point by combining the respective raw values determined for the respective arrays of the plurality of arrays at the respective point.

7. The method according to claim 6, further comprising refining the raw values by a beamforming technique using the respective estimated values estimated for the respective points of the plurality of points.

8. The estimation method according to claim 1, wherein the estimated value of said quantity is determined by applying to the raw values, a multi-variable function whose having an image that is zero for any antecedent comprising at least one zero variable.

9. The estimation method according to claim 8, wherein the estimated value of said quantity is equal to the inverse of a sum of the inverses of the raw values.

10. The estimation method according to claim 8, wherein the estimated value of said quantity is equal to the M-th root of the product of the raw values, where M is the number of arrays of the plurality of arrays.

11. The estimation method according to claim 1, wherein each of said pairwise combinations of representative signals are an estimation of the mathematical expectation of the product of the respective representative signals.

12. The estimation method according to claim 1, wherein said representative signals are produced by processing measurements respectively acquired by the acoustic sensors of the respective array.

13. The estimation method according to claim 1, wherein said quantity is an acoustic power.

14. A system for estimating a quantity representative of the sound energy at at least one point of a three-dimensional space, the system comprising: a plurality of arrays, each of the arrays comprising at least K acoustic sensors and configured to produce a plurality of signals representative of the sound field at the respective array and to determine a raw value of said quantity at said point based on at least K+1 elements of a matrix that are based respectively on pairwise combinations of representative signals produced by the respective array, K being greater than or equal to 2; and a processor configured to determine an estimated value of said quantity at said point by combining the raw values of said quantity at said point, determined for the respective arrays of the plurality of arrays.

15. The estimation method according to claim 2, further comprising, for each array of the plurality of arrays, determining, based on said matrix, a plurality of directional values of the quantity representative of the sound energy received at the respective array, from a plurality of directions respectively.

16. The estimation method according to claim 2, wherein the estimated value of said quantity is determined by applying, to the raw values, a multi-variable function having an image that is zero for any antecedent comprising at least one zero variable.

17. The estimation method according to claim 2, wherein each of said pairwise combinations of representative signals are an estimation of the mathematical expectation of the product of the respective representative signals.

18. The estimation method according to claim 2, wherein said representative signals are produced by processing measurements respectively acquired by the acoustic sensors of the respective array.

19. The estimation method according to claim 2, wherein said quantity is an acoustic power.

Description

DISCLOSURE OF THE INVENTION

[0005] In this context, the present invention proposes a method for estimating a quantity representative of the sound energy at at least one point of a three-dimensional space where a plurality of arrays are located, each comprising at least K acoustic sensors, K being higher than or equal to 2, comprising the following steps: [0006] for each array of the plurality of arrays, producing a plurality of signals representative of the sound field at the array in question; [0007] for each array of the plurality of arrays, determining a raw value of said quantity at said point based on at least K+1 elements of a matrix that are based respectively on pairwise combinations of representative signals produced by the array in question; [0008] determining an estimated value of said quantity at said point by combining the raw values of said quantity at said point, determined for the various arrays of the plurality of arrays respectively.

[0009] The use of arrays each comprising at least 2 acoustic sensors (and preferably at least 4 acoustic sensors) allows a fine analysis of the sound field at the array. The various signals resulting from this analysis allow generating a matrix that renders accurately the sound field present at the array. The sound field analysis is thus both rich and made in a compact way, so that it is possible to map correctly the sound field at the array.

[0010] In the case where 2nd-order ambisonic signals are used, for example, the number K of acoustic sensors per array is higher than or equal to 9. In the case where 3rd-order ambisonic signals are used, the number K of acoustic sensors per array is higher than 16.

[0011] The step of determining a raw value for a given array may comprise the following sub-steps: [0012] determining, based on said matrix, a directional value of the quantity representative of the sound energy received at the given array from a direction connecting said point and the given array; [0013] determining the raw value for the given array based on the directional value determined.

[0014] In order to cover a set of directions about each array, the estimation method may comprise, for each array of the plurality of arrays, a step of determining, based on said matrix, a plurality of directional values of the quantity representative of the sound energy received at the array in question, from a plurality of directions respectively.

[0015] Moreover, it may be provided, for at least one array of the plurality of arrays, a step of refining the directional values by means of a beamforming technique.

[0016] The estimation method may further comprise, in this case, for each array of the plurality of arrays, a step of determining raw values of said quantity at a plurality of points based on the directional values determined for the array in question.

[0017] The method may then comprise, for each point of said plurality of points, a step of determining an estimated value of said quantity at the point in question by combining the raw values determined for the various arrays of the plurality of arrays at the point in question.

[0018] A mapping of the quantity representative of the sound energy is thus carried out.

[0019] The method may further comprise a step of refining the raw values by means of a beamforming technique using the values estimated for the various points of the plurality of points.

[0020] In practice, the estimated value of said quantity may be determined by applying to the raw values a multi-variable function whose image is zero for any antecedent comprising at least one zero variable, which makes it possible to determine relatively simply the estimated value of said quantity on the basis of the raw values.

[0021] The estimated value of said quantity may for example be equal to the inverse of the sum of the inverses of the raw values.

[0022] According to another conceivable possibility, the estimated value of said quantity may be equal to the M-th root of the product of the raw values, where M is the number of arrays of the plurality of arrays.

[0023] The pairwise combination of representative signals are for example each an estimation of the mathematical expectation of the product of the representative signals in question.

[0024] The above-mentioned representative signals may be produced by processing measurements respectively acquired by the acoustic sensors of the array in question.

[0025] The above-mentioned quantity is for example the acoustic power. As an alternative, it could be the acoustic pressure (defined on the basis of the square root of the acoustic power).

[0026] The present invention also relates to a system for estimating a quantity representative of the sound energy at at least one point of a three-dimensional space comprising: [0027] a plurality of arrays each comprising at least K acoustic sensors and each adapted to produce a plurality of signals representative of the sound field at the array in question and to determine a raw value of said quantity at said point on the basis of at least K+1 elements of a matrix that are based respectively on pairwise combinations of representative signals produced by the array in question, K being higher than or equal to 2; and [0028] a processor adapted to determine an estimated value of said quantity at said point by combining the raw values of said quantity at said point, determined for the various arrays of the plurality of arrays respectively.

[0029] Of course, the different features, alternatives and embodiments of the invention can be associated with each other according to various combinations, insofar as they are not mutually incompatible or exclusive.

DETAILED DESCRIPTION OF THE INVENTION

[0030] Moreover, various other features of the invention will be apparent from the appended description made with reference to the drawings that illustrate non-limiting embodiments of the invention, and wherein:

[0031] FIG. 1 schematically shows a system comprising a processor and a plurality of arrays;

[0032] FIG. 2 is a flowchart showing the main steps of a method for estimating a quantity representative of the sound energy according to the invention;

[0033] FIG. 3 schematically shows a particular direction in space relative to an array of the plurality of arrays of FIG. 1;

[0034] FIG. 4 schematically shows a meshing of the space about the array of FIG. 3;

[0035] FIG. 5 is a flowchart showing a method for refining directional values by means of a beamforming technique; and

[0036] FIG. 6 is a flowchart showing a method for refining raw and estimated values by means of a beamforming technique.

[0037] The system shown in FIG. 1 comprises a processor P and a plurality of arrays (here M arrays) A.sub.1, A.sub.m, A.sub.M.

[0038] As schematically shown in FIG. 1, the various arrays A.sub.1, A.sub.m, A.sub.M are respectively located at various points of a three-dimensional space E.

[0039] Each array A.sub.m comprises several acoustic sensors S.sub.i each capable of making a measurement of a sound field present at the acoustic sensor S.sub.i in question. In FIG. 1 is schematically shown an acoustic wave transmitted by a sound source G, but the invention applies regardless the number of sound sources.

[0040] In the example described herein, each array A.sub.m comprises exactly K acoustic sensors (K being higher than or equal to 2, preferably K being higher than or equal to 4), for example 35 acoustic sensors. As an alternative, however, certain acoustic arrays could comprise more than K acoustic sensors.

[0041] Each array A.sub.m also comprises a processing unit U adapted to process the signals measured by the acoustic sensors S.sub.i of the array in question, as explained hereinafter.

[0042] Each array A.sub.m can moreover communicate with processor P (for example by means of a wireless link or, as an alternative, a wire link) in order to allow data exchanges between the processing unit U of this array A.sub.m and processor P.

[0043] FIG. 2 shows the main steps of a method for estimating a quantity representative of the sound energy according to the invention. In the example described herein, the quantity used to represent the sound energy is the acoustic power.

[0044] Steps E2 to E8 that will now be described are implemented in each of the arrays A.sub.1, A.sub.m, A.sub.M. However, for the sake of brevity, a single array reference is given below: A.sub.m.

[0045] The method starts by a step E2 of acquiring respective measurements by the K acoustic sensors S.sub.i of each array A.sub.m of the plurality of arrays.

[0046] In the example described hereinafter, step E2 further comprises a processing (by the processing unit U of each array A.sub.m) of the measurements acquired by the K acoustic sensors S.sub.i of the array A.sub.m in question in order to produce signals s.sub.k(t) representative of the sound field at the array A.sub.m in question. According to the representation used, these signals s.sub.k(t) may be complex signals (i.e. represented as complex number in order to define a module, or amplitude, and a phase) or real signals.

[0047] These signals s.sub.k(t) are for example L-order ambisonic signals. The L-order ambisonic representation indeed allows representing the sound field at the array A.sub.m in question by means of N signals s.sub.k(t) with N=(L+1).sup.2. Generally, the number K of acoustic sensors is higher than or equal to the number N of signals s.sub.k(t) produced.

[0048] The method then continues, at each array A.sub.m (and by means of the processing unit U of the array A.sub.m in question), by a step E4 of determining directional values p.sup.(m)(Ω) of the acoustic power received at the array A.sub.m from a plurality of directions Ω.

[0049] As schematically shown in FIG. 3, processing unit U determines for example the directional value p.sup.(m)(Ω) for a set of directions Ω forming an angular meshing about the array A.sub.m in question, each direction Ω corresponding to particular spherical angular coordinates (θ, φ), where θ is the elevation (between 0 and π) and φ is the azimuth (between 0 and 2π). For example, in practice, a meshing (such as a Lebedev meshing) is used, which comprises a number of directions between a few tens and several thousands (i.e. in practice between 50 and 5,000).

[0050] For each direction Ω (and at each time for which the estimation is made), each processing unit U determines for that purpose the elements of a covariance matrix C.sub.ss in which: [0051] each diagonal element is an estimation of the mathematical expectation of the square of the module of one of the signals s.sub.k(t) representative of the sound field (the covariance matrix C.sub.ss here comprises N diagonal elements); [0052] each non-diagonal element is an estimation of the mathematical expectation of the product of one s.sub.i(t) of the signals representative of the sound field by the conjugate of another s.sub.i(t) of the signals representative of the sound field (the covariance matrix C.sub.ss here comprises N(N−1) non-diagonal elements).

[0053] The covariance matrix C.sub.ss provides a set of statistical information about the spatial properties of the sound field, in particular about the position of the sound sources and the more or less strong correlation of the signals transmitted by them. From this point of view, each element of the matrix enriches the information and thus allows refining the analysis performed.

[0054] In the case where ambisonic signals s.sub.k(t) with real values are used as described here, the covariance matrix C.sub.ss is written:

[00001]$C_{\mathrm{ss}}=\left[\begin{array}{ccc}E\left\{{\left[s_1\left(t\right)\right]}^2\right\}& .Math.& E\left\{s_1\left(t\right).Math.s_N\left(t\right)\right\}\\ .Math.& \ddots & .Math.\\ E\left\{s_N\left(t\right).Math.s_1\left(t\right)\right\}& .Math.& E\left\{{\left[s_N\left(t\right)\right]}^2\right\}\end{array}\right],$

[0055] where E is a function estimating the mathematical expectation of the signal in question.

[0056] In practice, the function E may be an indicator of central tendency of the signal in question over a predetermined number of samples of this signal (the samples used in the calculation of the central tendency indicator being generally the last samples produced). The function E is for example the (sliding) average of the signal over this predetermined number of (last) samples.

[0057] The directional value p.sup.(m)(Ω) of the acoustic power received from a direction Ω is then written:

[00002]$p^{\left(m\right)}\left(\Omega \right)=\frac{1}{a^H\left(\Omega \right)C_{\mathrm{ss}}^{-1}a\left(\Omega \right)},$

[0058] where (.).sup.H is the transpose-conjugate operator and a(Ω) is a steering vector of the direction Ω defined as follows in the case of ambisonic signals s.sub.k(t):

[00003]$a\left(\Omega \right)=\left(\begin{array}{c}{Y}_0^0\left(\theta ,\phi \right)\\ Y_l^q\left(\theta ,\phi \right)\\ Y_L^L\left(\theta ,\phi \right)\end{array}\right),$

[0059] where Y.sub.l.sup.q is the spherical harmonic function of real value of order l and degree q and the variables θ and φ represent the direction Ω in spherical coordinates. (The number of spherical harmonic functions of order lower or equal to L being equal to (L+1).sup.2, the vector a(Ω) is of dimension N=(L+1).sup.2 and the above-mentioned covariance matrix C.sub.ss of size N×N.)

[0060] The directional values p.sup.(m)(Ω) obtained for a given array A.sub.m can potentially be refined by means of a beamforming technique, as described hereinafter with reference to FIG. 5.

[0061] The processing unit U of each array A.sub.m then performs a step E6 of determining raw values p.sup.(m)(r) of the acoustic power at a plurality of points of the three-dimensional space E, the position of a point being given by a coordinate vector r.

[0062] The points where the raw values p.sup.(m)(r) of acoustic power are determined are for example predefined and are the same for all the arrays A.sub.m. These points form for example a meshing of the area of interest of the three-dimensional space (area that thus comprises all the arrays A.sub.m).

[0063] A portion of this meshing around a particular array A.sub.m is schematically shown in FIG. 4. (For the sake of clarity, the meshing shown in FIG. 4 is two-dimensional, but this meshing may in practice be three-dimensional.) In practice, for example, a number of points is used, which is included (depending on the size of the meshing) between a few tens and tens of thousands (i.e. in practice between 50 and 50,000). As an alternative, it could be possible to consider only points located in a given plane (as shown in FIG. 4).

[0064] For each point indicated by its coordinates r (comprising for example three coordinates (x,y,z) as shown in FIG. 4), the raw value p.sup.(m)(r) of acoustic power at that point is determined based on the directional value p.sup.(m)(w) determined at step E4 for the direction w connecting this point and the array A.sub.m in question, for example by interpolation of the directional values p.sup.(m)(ω) determined at step E4. In the case described herein, in which the possible directions vary in azimuth and elevation, it may be for example an interpolation by spherical spline or using spherical harmonic functions. In the case where the meshing is two-dimensional and the considered directions therefore extend in the plane of the meshing, a linear or quadratic interpolation can be used.

[0065] For each array A.sub.m, the raw values p.sup.(m)(r) determined by this array A.sub.m (precisely by the processing unit U of this array A.sub.m) are then transmitted at step E8 to processor P.

[0066] Processor P thus receives at step E10 all the raw values p.sup.(m)(r) determined by all the arrays A.sub.m of the plurality of arrays.

[0067] Processor P can hence determine, at step E12, for all the points considered, an so estimated value p.sup.(all)(r) of the acoustic power at the point in question by combining the raw values p.sup.(m)(r) for this point received from the various arrays A.sub.m.

[0068] The estimated value p.sup.(all)(r) for a given point (of coordinates r) is for example determined by applying, to the raw values p.sup.(m)(r) for this given point, a multi-variable function f (the number of variables x.sub.1, x.sub.2, . . . , x.sub.M being equal to the number of arrays) and whose image f(x.sub.1, x.sub.2, . . . , x.sub.M) is equal to zero for any antecedent (x.sub.1, x.sub.2, . . . , x.sub.M) comprising at least one zero variable x.sub.i.

[0069] In other words, the function f verifies: f(x.sub.1, x.sub.2, . . . , x.sub.M)=0 if (at least) one index i exists between 1 and M such that x.sub.i=0.

[0070] The estimated value p.sup.(all)(r) for a given point is then equal to:

p.sup.(all)(r)=f(p.sup.(1)(r),p.sup.(2)(r), . . . ,p.sup.(M)(r)).

[0071] The use of such a function f is interesting in that it allows obtaining an estimated value p.sup.(all)(r) that is zero (or very low in practice) whenever one of the raw values p.sup.(m)(r) is zero (or very low in practice), which tends to reduce the occurrence of noise in the estimation process.

[0072] In practice, the estimated value p.sup.(all)(r) for a given point may be determined as follows:

[00004]$p^{\left(\mathrm{all}\right)}\left(r\right)=\frac{1}{{{\Sigma }_{m=1}^M\left[p^{\left(m\right)}\left(r\right)\right]}^{-1}}.$

[0073] In other words, in this case, the estimated value p.sup.(all)(r) is equal to the inverse of the sum of the inverses of the M raw values p.sup.(m)(r).

[0074] This solution, based on the hypothesis of absence of interaction between the various arrays, is simple to implement and gives good results in practice.

[0075] It can be noticed that this possible embodiment corresponds to the case of a function f as proposed hereinabove due to the fact that the above expression of p.sup.(all)(r) may also be written:

[00005]$p^{\left(\mathrm{all}\right)}\left(r\right)=\frac{p^{\left(1\right)}\left(r\right).Math.p^{\left(2\right)}\left(r\right).Math..Math.p^{\left(M\right)}\left(r\right)}{p^{\left(2\right)}\left(r\right).Math.p^{\left(3\right)}\left(r\right).Math..Math.p^{\left(M\right)}\left(r\right)+.Math.+p^{\left(1\right)}\left(r\right).Math.p^{\left(2\right)}\left(r\right).Math..Math.p^{\left(M-1\right)}\left(r\right)}.$

[0076] According to a conceivable alternative, the processor determines as follows the estimated value p.sup.(all)(r) for a given point:

p.sup.(all)(r)=(Π.sub.m=1.sup.Mp.sup.(m)(r)).sup.1/M.

[0077] In other words, in this case, the estimated value p.sup.(all)(r) is equal to the M-th root of the product of the raw values (p.sup.(m)(r)).

[0078] Once the estimated values p.sup.(all)(r) determined for the plurality of considered points, it is possible to refine the corresponding raw vales p.sup.(m)(r) associated with the various arrays A.sub.m, and thus the estimated values p.sup.(all)(r) themselves, by means of a beamforming technique, as described hereinafter with reference to FIG. 6.

[0079] FIG. 5 shows a method for refining the directional values p.sup.(m)(Ω) by means of a beamforming technique. The implementation of this refining method is described for a particular array A.sub.m, but the method may be implemented (for example, by the processing unit U of the array in question) for several of said arrays (or even for all said arrays).

[0080] As indicated hereinabove, this refining method may take place when a set of D directional values p.sup.(m)(Ω.sub.i) has been determined (as indicated hereinabove in relation with step E4) for D directions Ω.sub.1, . . . , Ω.sub.D, respectively.

[0081] This refining method starts by a step E20 of determining a matrix V.sup.(m) defined as followed:

V.sup.(m)=W.sup.(m)A.sup.H(AW.sup.(m)A.sup.H+R).sup.−1, with

[0082] A a matrix obtained by concatenating the pointing vectors a(Ω.sub.i) defined as hereinabove, each for a direction Ω.sub.i, and associated with the D directions Ω.sub.i, . . . , Ω.sub.D, respectively,

[0083] W.sup.(m) the diagonal matrix comprising (in diagonal) the directional values p.sup.(m)(Ω.sub.i) previously determined for the D directions Ω.sub.1, . . . , Ω.sub.D,

[0084] R a regularization matrix that allows taking account of the presence of diffuse noise in the measured signals.

[0085] Reference may be made to the book “Geophysical Data Analysis: Diverse Inverse Theory, 4th Edition”, Academic Press, 2008, p. 62 for more details about this technique indicated as the solution to the “weighted damped least-square problem”.

[0086] The refining method continues with a step E22 in which the processing unit U in question determines a refined version of the matrix W.sup.(m) (and therefore of the directional values p.sup.(m)(Ω.sub.i) present on the diagonal of this matrix) as follows:

Z.sup.(m)=V.sup.(m)C.sub.ss.sup.(m)V.sup.(m)H

W.sup.(m)=diag(Z.sup.(m))

[0087] where C.sub.ss.sup.(m) is the covariance matrix determined (as indicated hereinabove at step E4) for the array A.sub.m in question (and at the moment in question),

[0088] where diag is the operator that, with matrix Z.sup.(m), associates the diagonal matrix W.sup.(m), whose diagonal elements are identical to those of matrix Z.sup.(m) (and whose other elements are zero).

[0089] The new directional values p.sup.(m)(Ω.sub.i) present on the diagonal of the so-obtained matrix W.sup.(m) may be used for the following of the method.

[0090] Steps E20 and E22 may in practice be repeated several times to further refine the directional values p.sup.(m)(Ω.sub.i).

[0091] FIG. 6 shows a method for refining the directional values p.sup.(m)(r) and the estimated values p.sup.(all)(r) by means of a beamforming technique.

[0092] This refining method starts with a step E30 in which processor P determines, for each array A.sub.m, a matrix V.sup.(m) as follows:

V.sup.(m)=W.sup.(all)A.sup.(m)H(A.sup.(m)W.sup.(all)A.sup.(m)H+R).sup.−1, with

[0093] A.sup.(m) a matrix obtained by concatenating the pointing vectors a(ω.sub.i) defined, for a set of T points (of the area of interest) indicated by vectors r.sub.1, r.sub.2, . . . , r.sub.T, by the direction Wi connecting the array A.sub.m to the point r.sub.i in question (the pointing vector associated with a particular direction being defined hereinabove),

[0094] W.sup.(m) the diagonal matrix comprising (in diagonal) the estimated values p.sup.(all)(r.sub.i) previously determined for the T points of coordinates r.sub.1, r.sub.2, . . . , r.sub.T,

[0095] R a regularization matrix that allows taking account of the presence of diffuse noise in the measured signals and of sound sources present out of the area of interest.

[0096] This solution is of the same type as that proposed hereinabove for refining the directional values and reference can therefore also be made to the above-mentioned book for more details on this subject.

[0097] The refining method continues with a step E32 of determining, for each array A.sub.m, refined raw values p.sup.(m)(r.sub.i). For that purpose, processor P determines the matrix V.sup.(m)C.sub.ss.sup.(m)V.sup.(m)H, the refined raw values p.sup.(m)(r.sub.1), p.sup.(m)(r.sub.2), . . . , p.sup.(m)(r.sub.T) then being the diagonal elements of this matrix V.sup.(m)C.sub.ss.sup.(m)V.sup.(m)H (matrix C.sub.ss.sup.(m) being as hereinabove the covariance matrix determined at step E4 for the array A.sub.m in question).

[0098] Processor P can then obtain at step E34, for each point r.sub.i of the plurality of T points of coordinates r.sub.1, r.sub.2, . . . , r.sub.T, a refined estimated value p.sup.(all)(r.sub.i) by combining the M refined raw values p.sup.(m)(r.sub.i) obtained for this point r.sub.i for the various arrays A.sub.m, respectively, for example by the combination method described hereinabove at step E12:

p.sup.(all)(r.sub.i)=f(p.sup.(1)(r.sub.i),p.sup.(2)(r.sub.i), . . . ,p.sup.(M)(r.sub.i)).

[0099] Steps E30 to E34 may in practice be repeated several times to further refine the raw values p.sup.(m)(r.sub.i) and the estimated values p.sup.(all)(r.sub.i).