Abstract:
The invention relates to a method of processing a noisy sound signal and to a device for implementing said method. The inventive method includes the following steps which are performed for each frame of a noisy sound signal, consisting in: applying a frequency-domain transform; estimating a power spectral density (PSD) of the noise for the frame; calculating a first noise reduction filter from the estimated PSD of the noise and from an estimation of the PSD of a wanted signal corresponding to the frame; filtering the frame with the aid of the calculated noise reduction filter, in order to obtain a first noise-corrected estimation for the frame; and obtaining a frame of a second signal from the first noise-corrected estimation of the aforementioned frame, comprising harmonics which are essentially at the same positions as the wanted signal corresponding to the frame of the noisy sound signal.

Description:
This application claims priority from PCT/FR2005/002284 filed Sep. 14, 2005, which claims priority from French Application No. 04 09819, filed Sep. 16, 2004, both of which are hereby incorporated by reference in their entirety. 
   BACKGROUND OF THE INVENTION 
   The present invention relates to the processing of noisy sound signals. It relates in particular to the reduction of the noise present in such signals. 
   Techniques for reducing noise, that is to say a disturbing signal, within a sound signal are known. They are aimed at taking account of the acoustic environment in which the sound signal appears so as to improve the quality and the intelligibility of this signal. These techniques consist in extracting the useful information from the sound signal considered by performing a processing on this noisy signal. Such techniques apply for example to spoken communications, in applications such as telephony, teleconferencing, videoconferencing where the sound signal is then transmitted between several talkers. They apply moreover in respect of applications of sound pick-up in noisy surroundings, or else of voice recognition, the performance of which is greatly altered when the voice signal is pronounced in a noise-filled environment. 
   These techniques usually consist in estimating a transfer function of a noise reduction filter, then in carrying out a filtering processing on the basis of a multiplication in the spectral domain. They come within approaches termed “noise reduction by short-term spectral attenuation”. 
   According to these techniques, the sound signal considered x(n) comprises a useful signal component s(n) and a noise component b(n), n representing a temporal index in discrete time. It will however be noted that a representation of the signal in continuous time could also be adopted. The signal x(n) is organized as successive frames x(n,k) of constant length and of index k. Each of these frames is firstly multiplied by a weighting window making it possible to improve the later estimation of the spectral quantities necessary for the calculation of the noise reduction filter. Each frame thus windowed is then analyzed in the spectral domain, for example with the aid of a discrete or fast Fourier transformation. This operation is called short-term Fourier transformation (STFT). 
   The frequency representation X(k,f) thus obtained of the signal observed, where f is a frequency index, makes it possible at one and the same time to estimate the transfer function H(k,f) of the noise reduction filter, and to apply this filter in the spectral domain by simple multiplication between this transfer function and the short-term spectrum of the noisy signal. The result of the filtering may thus be written:
 
 Ŝ ( k,f )= H ( k,f ) X ( k,f ).
 
   A return to the time domain of the signal obtained is then performed by an inverse spectral transform. The corresponding temporal signal is finally synthesized by a block overlap and add technique (OLA standing for “overlap add”) or else by a block save technique (OLS standing for “overlap save”). This operation of reconstructing the signal in the time domain is called inverse short-term Fourier transformation (ISTFT). 
   A detailed description of the methods of short-term spectral attenuation will be found in the references: J. S. Lim, A. V. Oppenheim, “Enhancement and bandwidth compression of noisy speech”, Proceedings of the IEEE, Vol. 67, pp. 1586-1604, 1979; and R. E. Crochiere, L. R. Rabiner, “Multirate digital signal processing”, Prentice Hall, 1983. 
   The short-term spectral attenuation H(k,f) applied to the observation signal X(k,f) over the temporal segment of index k and with the frequency component f, is generally determined on the basis of the estimate of the local signal-to-noise ratio SNR(k,f). A characteristic common to all suppression rules resides in their asymptotic behavior, given by:
 
 H ( k,f )≈1 for SNR( k,f )&gt;&gt;1
 
 H ( k,f )≈0 for SNR( k,f )&lt;&lt;1.
 
   In most techniques, the following assumptions are made: the noise and the useful signal are statistically uncorrelated, the useful signal is intermittent (presence of periods of silence) and the human ear is not sensitive to the signal phase (which is not in general modified by the processing). 
   Among the suppression rules commonly employed may be cited by way of example: power spectral subtraction, amplitude spectral subtraction and direct implementation of the Wiener filter. For these rules, the short-term estimate of the frequency component f of the useful speech signal may be written respectively: 
                       S   ^     SSP     ⁡     (     k   ,   f     )       =             γ   ss     ⁡     (     k   ,   f     )             γ   ss     ⁡     (     k   ,   f     )       +       γ   bb     ⁡     (     k   ,   f     )             ⁢     X   ⁡     (     k   ,   f     )                 (   1   )               
for the power spectral subtraction (see the aforesaid article by J. S. Lim and A. V. Oppenheim);
 
                       S   ^     SSA     ⁡     (     k   ,   f     )       =       [     1   -           γ   bb     ⁡     (     k   ,   f     )             γ   ss     ⁡     (     k   ,   f     )       +       γ   bb     ⁡     (     k   ,   f     )               ]     ⁢     X   ⁡     (     k   ,   f     )                 (   2   )               
for amplitude spectral subtraction (see S. F. Boll, “Suppression of acoustic noise in speech using spectral subtraction”, IEEE Trans. on Audio, Speech and Signal Processing, Vol. 27, No. 2, pp. 113-120, April 1979); and
 
                       S   ^     Wiener     ⁡     (     k   ,   f     )       =           γ   ss     ⁡     (     k   ,   f     )             γ   ss     ⁡     (     k   ,   f     )       +       γ   bb     ⁡     (     k   ,   f     )           ⁢     X   ⁡     (     k   ,   f     )                 (   3   )               
for Wiener filtering (cf. aforesaid article by J. S. Lim and A. V. Oppenheim).
 
   In these expressions, γ ss (k,f) and γ bb (k,f) respectively represent the power spectral densities of the useful signal and of the noise that are present within the frequency component f of the observation signal X(k,f) over the time window of index k. 
   On the basis of the above expressions, it is possible to study, as a function of the local signal-to-noise ratio measured on a given frequency component f, the behavior of the spectral attenuation applied to the noisy signal. These curves are plotted in  FIG. 1  for the abovementioned three short-term suppression rules. It may be noted that the set of rules gives a substantially identical attenuation when the local signal-to-noise ratio is significant (right-hand part of  FIG. 1 ). The optimal power subtraction rule, optimal in the sense of the maximum likelihood for Gaussian models (see O. Cappé, “Elimination of the musical noise phenomenon with the Ephraim and Malah noise suppressor”, IEEE Trans. on Speech and Audio Processing, Vol. 2, No. 2, pp. 345-349, April 1994) is still the one for which the power of the noise remains the most significant at the output of the processing. For the three suppression rules, we may note that a small variation in the local signal-to-noise ratio about a cutoff value suffices for switching from the case of total attenuation (H(k,f)≈0) to the case of negligible spectral modification (H(k,f)≈1). 
   This latter property constitutes one of the causes of the phenomenon dubbed “musical noise”. Specifically, the ambient noise, comprising at one and the same time deterministic and random components, can be characterized only during the periods of vocal non-activity. On account of the presence of random components, there are very strong variations between the real contribution of a frequency component f of the noise during the periods of vocal activity and its average estimate made over several frames during the instants of vocal non-activity. On account of this difference, the estimate of the local signal-to-noise ratio may fluctuate about the cutoff level and hence give rise at the output of the processing to spectral components which appear and then disappear and whose average lifetime does not statistically exceed the order of magnitude of the analysis window considered. The generalization of this behavior over the whole of the passband introduces audible and annoying residual noise. 
   Several studies have endeavored to reduce the influence of this residual noise. The solutions advocated follow several avenues: averaging of the short-term estimates (cf. aforesaid article by S. F. Boll), overestimation of the noise power spectrum (see M. Berouti et al., “Enhancement of speech corrupted by acoustic noise”, Int. Conf. on Speech, Signal Processing, pp. 208-211, 1979; and P. Lockwood, J. Boudy, “Experiments with a non-linear spectral subtractor, hidden Markov models and the projection for robust speech recognition in cars”, Proc. of EUSIPCO&#39;91, pp. 79-82, 1991), or else tracking of the minima of the noise spectral density (see R. Martin, “Spectral subtraction based on minimum statistics”, in Signal Processing VII: Theories and Applications, EUSIPCO&#39;94, pp. 1182-1185, September 1994). 
   A relatively effective solution for suppressing musical noise consists of an estimator of the power spectral density of the useful signal termed “directed-decision” (see Y. Ephraim, and D. Malah, “Speech enhancement using a minimum mean square error short-time spectral amplitude estimator”, IEEE Trans. on Audio, Speech and Signal Processing, Vol. 32, No. 6, pp. 1109-1121, 1984 and the aforesaid article by O. Cappé). This estimator effects a compromise between the instantaneous and long-term power spectral density of the useful signal, thereby making it possible to effectively eliminate the musical noise. It is moreover known to improve this solution by making up the delay inherent in this estimator (see FR2820227 and C. Plapous, C. Marro, L. Mauuary, P. Scalart, “A Two-Step Noise Reduction Technique”, ICASSP, May 2004). 
   Several studies have also pertained to the establishing of new suppression rules based on statistical models of the speech and additive noise signals. These studies have made it possible to introduce new algorithms dubbed “soft-decision” algorithms since they possess an additional degree of freedom with respect to the conventional methods (see R. J. Mac Aulay, M. L. Malpass, “Speech enhancement using a soft-decision noise suppression filter”, IEEE trans. on Audio, Speech and Signal Processing, Vol. 28, No. 2, pp. 137-145, April 1980, Y. Ephraim, D. Malah, “Speech enhancement using optimal non-linear spectral amplitude estimation”, Int. Conf. on Speech, Signal Processing, pp. 1118-1121, 1983, and Y. Ephraim, D. Malah article, “Speech enhancement using a minimum mean square error short-time spectral amplitude estimator”, stated above). 
   As was mentioned above, the calculation of the short-term spectral attenuation relies on the estimation of the signal-to-noise ratio on each of the spectral components. By way of example, the equations given above each involve the following quantity: 
   
     
       
         
           
             SNR 
             ⁡ 
             
               ( 
               
                 k 
                 , 
                 f 
               
               ) 
             
           
           = 
           
             
               
                 
                   γ 
                   ss 
                 
                 ⁡ 
                 
                   ( 
                   
                     k 
                     , 
                     f 
                   
                   ) 
                 
               
               
                 
                   γ 
                   bb 
                 
                 ⁡ 
                 
                   ( 
                   
                     k 
                     , 
                     f 
                   
                   ) 
                 
               
             
             . 
           
         
       
     
   
   Thus, the performance of the noise reduction technique, especially in terms of distortions and of effective reduction of the noise level, are governed by the relevance of this estimator of the signal-to-noise ratio. 
   This defect constitutes the major limitation of the known speech denoising systems. Specifically, the current denoising systems are incapable of denoising the harmonics characterized by too low a signal-to-noise ratio. In practice, the denoising algorithms use the SNR to detect the presence or the absence of a speech component for each frequency. If the estimated SNR is too unfavorable, then the algorithm considers that there is no signal component and suppresses it. Thus, harmonics may be destroyed by the known denoising systems, although it is known a priori that such harmonics must exist. Now, it should be noted that, in the majority of languages, the voiced sounds (harmonics) represent a very large part of the sounds uttered. 
   An object of the present invention is to overcome the limitation of the known denoising systems. 
   Another object of the invention is to improve the performance of noise reduction methods. 
   Another object of the invention is to propose a sound signal processing which does not distort the signal excessively. In particular, the processing of the signal performed makes it possible to preserve all or part of the harmonics included in this signal. 
   Another object of the invention is to limit the appearance of musical noise on completion of the sound signal processing. 
   Another object of the invention is to obtain a good estimate of the harmonic comb of a useful signal. 
   SUMMARY OF THE INVENTION 
   The invention thus proposes a method of processing a noisy sound signal organized as successive frames. The process comprises the following steps relating to at least one of said frames:
         applying a transform to the frequency domain to said frame of the noisy sound signal;   estimating a power spectral density of the noise for said frame;   calculating a first noise reduction filter on the basis of the estimated power spectral density of the noise and of an estimate of the power spectral density of a useful signal corresponding to said frame;   filtering said frame of the noisy sound signal with the aid of the first noise reduction filter calculated, to obtain a first denoised estimate of said frame; and   obtaining a frame of a second signal on the basis of the first denoised estimate of said frame of the noisy sound signal, said frame of the second signal comprising harmonics substantially at the same positions as the useful signal corresponding to said frame of the noisy sound signal.       

   Such a processing of the noisy sound signal thus makes it possible to obtain a regeneration of harmonicity at the output of a first noise reduction filter applied to the signal. The thus-obtained frame of the second signal is constructed in such a way as to eliminate the distortions of the frame of the denoised signal which may have appeared during the first filtering. 
   This processing makes it possible to obtain a good estimate of the harmonic comb of the useful signal, which may be the subject of a subsequent utilization. 
   Advantageously, in the domain of the denoising, the method furthermore comprises the following steps:
         calculating a second noise reduction filter on the basis of the estimated power spectral density of the noise and of a combination of the powers of the first denoised estimate of said frame and of the frame of the second signal obtained;   filtering said frame of the noisy sound signal, with the aid of the second noise reduction filter calculated, to obtain a second denoised estimate of said frame; and   synthesizing the second denoised estimate of said frame.       

   In this embodiment, the second noise reduction filter is calculated in such a way as to preserve the harmonics, since it is driven by the signal emanating from the harmonicity regeneration. The second denoised estimate of said frame obtained on completion of the second filtering is thus more efficacious than those obtained in the conventional noise reduction systems, where the harmonics of the input signal are destroyed or at least impaired. 
   The calculation of the first and, as appropriate, of the second noise reduction filter may comprise a first pass implementing a technique of short-term spectral attenuation. By way of illustration, the following techniques may be cited for example: power spectral subtraction, amplitude spectral subtraction and an open-loop Wiener filter, etc. Advantageously, it may also comprise a second pass implementing a technique of short-term spectral attenuation, such as for example a power spectral subtraction, an amplitude spectral subtraction and an open-loop Wiener filter, etc., and use an estimate of the power spectral density of the useful signal corresponding to said frame taking account of the calculation performed during the first pass. 
   Advantageously, the obtaining of a frame of a second signal comprising harmonics substantially at the same positions as the useful signal corresponding to said frame of the noisy sound signal comprises the application of a nonlinear function to the first denoised estimate of said frame of the noisy sound signal when said first estimate is in the time domain, and the application of a circular convolution between the first denoised estimate of said frame of the noisy sound signal and of a transform into the frequency domain of a nonlinear function when said first estimate is in the frequency domain. 
   Said nonlinear function may be, for example, one among: a single-wave rectification function, an absolute value, a maximum between said first denoised estimate of said frame of the noisy sound signal and a threshold, and a minimum between said first denoised estimate of said frame of the noisy sound signal and a threshold. 
   The invention furthermore proposes a device for processing noisy sound signals, comprising means designed to implement the abovementioned method. 
   The invention also proposes a computer program on an information medium, this program comprising instructions adapted for the implementation of the abovementioned method, when the program is loaded and executed by computing means. 

   
     BRIEF DESCRIPTION OF THE DRAWINGS 
       FIG. 1 , already commented on, is a diagram comparing the behavior of three known short-term suppression rules; 
       FIG. 2  is a diagram showing various functional blocks able to implement an embodiment of the invention; 
       FIG. 3  is a diagram showing an example of restoration of the lost harmonics in a denoised signal; 
       FIG. 4  is a diagram showing the effect of the restoration of harmonics in a signal according to an embodiment of the invention; and 
       FIG. 5  is a diagram showing an example of estimation of a second noise reduction filter according to an embodiment of the invention. 
   

   DESCRIPTION OF PREFERRED EMBODIMENTS 
     FIG. 2  shows a device  1  according to an embodiment of the invention. A current frame  2  of a noisy sound signal enters the device  1  so as to undergo a processing therein. This frame is designated by x(k,n) according to the notation introduced above. 
   The device  1  comprises an analysis unit  3  whose role is to switch the current frame of the noisy signal x(k,n) into the frequency domain X(k,f). Firstly, the frame is multiplied by a weighting window w(n) which delivers the signal x w (k,n):
 
 x   w ( k,n )= w ( n ) x ( k,n ).
 
   The switch to the frequency domain is then carried out using a discrete Fourier transform (DFT), i.e.:
 
 X ( k,f )=DFT(( x   w ( k,n )).
 
   The DFT in question may advantageously be implemented by a fast Fourier transform (or FFT). However, other transformations to the frequency domain, such as the wavelet transform are also possible. The same holds for the corresponding inverse operations, that is to say the inverse discrete Fourier transforms (IDFT) which will be mentioned later, for returning to the time domain. 
   A function  4  for detecting vocal activity (DVA) may moreover advantageously be implemented on the current frame  2 . The DVA makes it possible to know when to update the estimate of the power spectral density (PSD) of the noise. Thus, for each “noise only” frame k b  detected by the DVA, the noise power spectral density {circumflex over (γ)} bb (k b ,f) is estimated by a functional block  5  according to the recursive expression:
 
{circumflex over (γ)} bb ( k   b   ,f )=α( k   b ){circumflex over (γ)} bb ( k   b −1 ,f )+((1−α( k   b )| X ( k   b   ,f )| 2 .
 
   The parameter α(k b ) controls the smoothing factor. It can vary over time. 
   If the current frame  2  is not detected as a noise frame, then the estimate of the power spectral density of the noise is pegged, that is to say we retain the last value of PSD obtained for a previous noise frame. 
   It will be noted that the estimate of {circumflex over (γ)} bb (k b ,f) is not limited to this exponential-smoothing estimator, any other power spectral density estimator may be used. 
   Any type of DVA may be used, whether the latter operates in the time domain or in the frequency domain. It is even possible to dispense with such detection of vocal activity. 
   An estimation of a first noise reduction filter is then performed by the functional block  7  of  FIG. 2 . The transfer function of this first noise reduction filter is advantageously estimated in the frequency domain by virtue of a two-pass technique (see FR2820227 and C. Plapous, C. Marro, L. Mauuary, P. Scalart, “A Two-Step Noise Reduction Technique”, ICASSP, May 2004 stated above). In the first pass, the filter having the following transfer function is calculated:
 
 Ĥ   step1 ( k,f )= f   step1 (({circumflex over (γ)} ss ( k,f ),{circumflex over (γ)} bb ( k,f )).
 
   This expression for the transfer function of the filter is a mathematical expression which depends on the PSD of the useful signal {circumflex over (γ)} ss (k,f) and that of the noise {circumflex over (γ)} bb( k,f). In this sense, the function f step1  may be chosen according to a short-term spectral attenuation technique, so that Ĥ step1 (k,f) implements for example the power spectral subtraction rule (equation (1)), the amplitude spectral subtraction rule (equation (2)), the open-loop Wiener filter rule (equation (3)), etc. Any other frequency-based noise suppression rule may also be implemented to estimate Ĥ step1 (k,f). 
   The calculation of the noise power spectral density {circumflex over (γ)} bb (k b ,f) was detailed above. As far as the spectral quantity {circumflex over (γ)} ss (k,f) is concerned, it cannot be obtained directly on account of the mixing of the signal and of the noise during the periods of vocal activity. To calculate it, a directed-decision estimator is used (see Y. Ephraim, D. Malah, “Speech enhancement using a minimum mean square error short-time spectral amplitude estimator”, stated above), in accordance with the following expression:
 
{circumflex over (γ)} ss ( k,f )=β( k )| Ŝ ( k− 1 ,f )| 2 +(1−β( k )) P[|X ( k,f )| 2 −{circumflex over (γ)} bb ( k,f )],
 
where β(k) is a barycentric parameter that can vary over time and Ŝ(k−1,f) is the spectrum of the useful signal estimated in relation to the previous frame of index k−1. The function P, which ensures the thresholding of the quantity (|X(k,f)| 2 −{circumflex over (γ)} bb (k,f)) which may possibly be negative in case of an estimation error, is given by:
 
   
     
       
         
           
             P 
             ⁡ 
             
               [ 
               
                 z 
                 ⁡ 
                 
                   ( 
                   
                     k 
                     , 
                     f 
                   
                   ) 
                 
               
               ] 
             
           
           = 
           
             { 
             
               
                 
                   
                     z 
                     ⁡ 
                     
                       ( 
                       
                         k 
                         , 
                         f 
                       
                       ) 
                     
                   
                 
                 
                   
                     
                       if 
                       ⁢ 
                       
                           
                       
                       ⁢ 
                       
                         z 
                         ⁡ 
                         
                           ( 
                           
                             k 
                             , 
                             f 
                           
                           ) 
                         
                       
                     
                     &gt; 
                     0 
                   
                 
               
               
                 
                   0 
                 
                 
                   
                     otherwise 
                     . 
                   
                 
               
             
           
         
       
     
   
   It will be noted that the estimation of {circumflex over (γ)} ss (k,f) is not limited to this directed-decision estimator. Specifically, an exponential-smoothing estimator or any other power spectral density estimator may be used. 
   The transfer function Ĥ step1 (k,f) may thereafter be reused to refine the estimate of the PSD of the useful signal {circumflex over (γ)} ss (k,f). We then obtain the quantity {circumflex over (γ)} ss2 (k,f) given by:
 
{circumflex over (γ)} ss2 ( k,f )=| Ĥ   step1 ( k,f ) X ( k,f )| 2 .
 
   The second pass then consists in calculating the estimator Ĥ 1 (k,f) of the transfer function of the first noise reduction filter on the basis of {circumflex over (γ)} ss2 (k,f), i.e.:
 
 Ĥ   1 ( k,f )= f   step2 ({circumflex over (γ)} ss2 ( k,f ),{circumflex over (γ)} bb ( k,f )).
 
   This two-pass calculation, the particular feature of which resides in a “faster” updating of the PSD of the useful signal {circumflex over (γ)} ss (k,f), confers two advantages on the first noise reduction filter Ĥ 1 (k,f). On the one hand, one obtains faster tracking of the non-stationarities of the useful signal, in particular during fast variations of its temporal envelope (e.g. attacks or extinctions of the speech signal during a silence/speech transition). On the other hand, the noise reduction filter is better estimated, this being manifested through an enhancement of the performance of the method. 
   Ĥ step1 (k,f) and Ĥ 1 (k,f) may implement a short-term spectral attenuation technique, such as for example the power spectral subtraction rule (equation (1)), the amplitude spectral subtraction rule (equation (2)), or else the open-loop Wiener filter rule (equation (3)). Any other frequency-based noise suppression rule may also be implemented to estimate Ĥ step1 (k,f) and H 1 (k,f). Of course, it is also possible to limit oneself to the first pass, without implementing the second. 
   A filtering is then carried out by the functional block  6  of the device  1 , according to the first filter calculated H 1 (k,f). It may be performed in the frequency domain, this corresponding to multiplying two spectra. As things stand, this is equivalent to an operation of circular convolution in the time domain. It is thus necessary to take certain precautions to avoid distortions due to temporal aliasing which are manifested when listening by clicks in tempo with the frames. Thus, to satisfy the linear convolution constraint, it is necessary both to add a certain number of null samples to each input frame (so-called “zero padding” technique) and to limit the temporal support of the impulse response of the noise reduction filter (this may be performed in the time domain or frequency domain). 
   It is noted that, to limit the temporal support of the impulse response of the noise reduction filter, it is possible to introduce a constraint in the time domain, this requiring:
     i) a first “inverse” spectral transformation making it possible to obtain the impulse response h 1 (k,n) on the basis of the knowledge of the transfer function of the filter H 1 (k,f).   ii) a limitation of the number of points of this impulse response, to obtain a truncated temporal filter h 1 ′(k,n),   iii) a second “direct” spectral transformation making it possible to obtain the modified transfer function of the filter H 1 ′(k,f) on the basis of the constrained impulse response h 1 ′(k,n).   

   The (frequency) transfer function of the noise reduction filter Ĥ 1 (k,f) being available, the (temporal) impulse response ĥ 1 (k,n) is obtained by means of an inverse discrete Fourier transform (IDFT), i.e.:
 
 ĥ   1 ( k,n )=IDFT( Ĥ   1  ( k,f )).
 
   This impulse response is then limited temporally by selecting and weighting the most significant L filt1  coefficients by a window w filt1 :
 
 ĥ   1 ( k,n )= w   filt1 ( n )ĥ 1 ( k,n ).
 
   This limitation of the temporal support of the noise reduction filter presents a dual advantage. On the one hand it makes it possible to avoid the problems of temporal aliasing (compliance with the linear convolution). On the other hand, it ensures a smoothing making it possible to avoid the effects of too aggressive a filter. 
   The transfer function of the noise reduction filter H 1 ′(k,f) is therefore obtained by discrete Fourier transform (DFT):
 
 Ĥ   1 ′( k,f )=DFT( ĥ   1 ′( k,n )).
 
   Finally, the first denoised version of the noisy signal is obtained by frequency filtering, that is to say by multiplying the spectrum of the noisy signal and the transfer function of the noise reduction filter, i.e.:
 
 Ŝ   1 ( k,f )= Ĥ   1 ′( k,f ) X ( k,f ).
 
   This frequency filtering step may also be carried out in an equivalent manner by a temporal filtering, that is to say by directly filtering x w (k,n) by ĥ 1 ′(k,n), then doing a discrete Fourier transform (DFT) of the result. 
   This first estimate of the useful signal Ŝ 1 (k,f) possesses audible distortions since certain harmonics have been suppressed by the noise reduction filter H 1 ′(k,f), as explained in the introduction. Of course, the level of distortion is directly related to the SNR. The more powerful the noise in a frequency zone, the more prone to degradation are the harmonics of this zone. 
   To alleviate this drawback, a functional block  8  of the device  1  represented in  FIG. 2  creates an artificial signal Ŝ harmo (k,f) which possesses harmonics just where they had been destroyed or very greatly attenuated. 
   The signal Ŝ harmo (k,f) may be obtained in the time domain or in the frequency domain. In the time domain it is possible to restore harmonicity to a signal by applying a nonlinearity to it, for example a function from among: a single-wave rectification function, an absolute value, a minimum/maximum with respect to a threshold, etc. If NL denotes the nonlinear function used, as illustrated in step  14  of  FIG. 3 , we have:
 
 ŝ   harmo ( k,n )= NL ( ŝ   1 ( k,n )),
 
   where ŝ 1 (k,n) is obtained, according to step  13  of  FIG. 3 , by inverse discrete Fourier transform (IDFT):
 
 ŝ   1 ( k,n )=IDFT( Ŝ   1 ( k,f )).
 
   We then obtain Ŝ harmo (k,f), as indicated in step  15  of  FIG. 3 , by discrete Fourier transform (DFT):
 
 Ŝ   harmo ( k,f )=DFT( ŝ   harmo ( k,n )).
 
   It will be noted that the harmonics are recreated at the same positions as those of the useful signal. This is implicitly ensured by the fact that a temporal nonlinearity is used to restore them. 
   It will be noted moreover that, in the case where the nonlinearity is applied in the time domain, it is possible to avoid the switch to the frequency domain by applying the nonlinearity to the signal ŝ 1 (k,n) obtained directly through the convolution of x w (k,n) and ĥ 1 ′(k,n). 
     FIG. 4  illustrates the effect and the benefit of the nonlinearity. The first spectrum (curve  16 ) corresponds to a frame of the non-noisy useful signal and serves as reference. The second spectrum (curve  17 ) corresponds to the signal Ŝ 1 (k,f). It is clearly apparent in this second spectrum that certain harmonics have been completely destroyed and others degraded. The last spectrum (curve  18 ) corresponds to the signal Ŝ harmo (k,f) mentioned above. It is apparent that the nonlinearity applied to the signal ŝ 1 (k,n) has effectively made it possible to recreate the missing harmonics in the signal ŝ harmo (k,n). The spectrum of this signal, Ŝ harmo (k,f), therefore possesses harmonicity information that is very useful in particular for calculating a new noise reduction filter which will be capable of preserving the harmonics that the conventional noise reduction techniques destroy. 
   The temporal nonlinearity possesses a frequency equivalent which consists in carrying out a circular convolution between the spectrum of the signal Ŝ 1 (k,f) and the spectrum of a signal which possesses the same harmonic structure (at the level of the positions of the harmonics) as the useful signal. Such a spectrum may for example be obtained by virtue of the temporal signal d(k,n) which corresponds to a function of ŝ 1 (k,n). It is therefore also possible to obtain the signal Ŝ harmo (k,f) in the frequency domain, in the following manner:
 
 Ŝ   harmo ( k,f )= Ŝ   1 ( k,f )           DFT( d ( k,n )),
 
where the sign           corresponds to the circular convolution operator. If the temporal nonlinearity chosen is the absolute value function, for example, then the signal d(k,n) may be expressed thus, the function sgn designating the sign of the value to which it is applied:
 
 d ( k,n )= sgn ( ŝ   1 ( k,n )).

   Thus, the signal Ŝ harmo (k,f), where it be obtained through a temporal nonlinearity or through a frequency circular convolution, possesses harmonics just where those of Ŝ 1 (k,f) were destroyed or degraded. It constitutes per se a good estimate of the harmonic comb of the denoised signal, in which the harmonics are preserved at the same positions as in the useful signal. This estimate of the harmonic comb is particularly beneficial and may be utilized, in particular, but not exclusively, for denoising purposes. 
   Advantageously, this signal can moreover make it possible to calculate a noise reduction filter capable of preserving the harmonics normally destroyed by the conventional algorithms, as is described hereinbelow. 
   It will be noted that the present invention is not limited to the examples described hereinabove to obtain the signal Ŝ harmo (k,f) It extends in fact to any other method aimed at recreating the harmonics suppressed by the noise reduction filter H 1 ′(k,f). 
   Advantageously, a second filter is thereafter estimated by the functional block  9  of  FIG. 2 . The transfer function of this second noise reduction filter Ĥ 2 (k,f) is estimated in the frequency domain, as illustrated in  FIG. 5 , in the following manner, g representing a function:
 
 Ĥ   2 ( k,f )= g ({circumflex over (γ)} ss     —     harmo ( k,f ),{circumflex over (γ)} bb ( k,f )).
 
   This expression for the transfer function of the filter is a mathematical expression which depends on the PSD of the noise {circumflex over (γ)} bb (k,f) and that of the useful signal estimated by virtue of the harmonic restoration step, i.e.:
 
{circumflex over (γ)} ss     —     harmo ( k,f )=ρ( k,f )| Ŝ   1 ( k,f )| 2 +(1−ρ( k,f ))| Ŝ   harmo ( k,f )| 2 .
 
The parameter ρ(k,f) is used to control the level of reinjection of the signal Ŝ harmo (k,f) in the calculation of the PSD γ ss     —     harmo (k,f) as a function of the nonlinearity NL which was chosen to create the signal Ŝ harmo (k,f). It will be noted that this parameter may be fixed or else depend on frequency and/or time.
 
   The function g may be chosen so that Ĥ 2 (k,f) implements for example the power spectral subtraction rule (equation (1)), the amplitude spectral subtraction rule (equation (2)), or else the open-loop Wiener filter rule (equation (3)). Any other frequency-based noise suppression rule may of course be implemented to estimate Ĥ 2 (k,f). 
   This noise reduction transfer function Ĥ 2 (k,f) possesses the advantage of preserving the harmonics of the signal which are normally destroyed by the conventional noise reduction techniques. This therefore limits the distortions of the denoised signal. It should also be noted that this filter retains the good properties of the filter Ĥ 1 (k,f), that is to say good tracking of the non-stationarities and very little residual musical noise. 
   A second filtering is thereafter carried out by the functional block  10  of the device  1 , on the basis of the second estimated noise reduction filter. Thus, as for the first transfer function, the (temporal) impulse response ĥ 2 (k,n) equivalent to the (frequency) transfer function of the noise reduction filter Ĥ 2 (k,f) is obtained by means of an inverse discrete Fourier transform (IDFT), that is to say:
 
 ĥ   2 ( k,n )=IDFT( Ĥ   2 ( k,f )).
 
   This impulse response is thereafter advantageously limited temporally by selecting and weighting the most significant L filt2  coefficients by a window w filt2 , i.e.:
 
 ĥ   2 ′( k,n )= w   filt2 ( n ) ĥ   2 ( k,n ).
 
   The transfer function of the noise reduction filter Ĥ 2 ′(k,f) is thereafter obtained by discrete Fourier transform (DFT):
 
 Ĥ   2 ′( k,f )=DFT( ĥ   2 ′( k,n )).
 
   Finally, the denoised frequency frame is obtained by frequency filtering, that is to say by multiplying the spectrum of the noisy signal and the transfer function of the noise reduction filter, i.e.:
 
 Ŝ   2 ( k,f )= Ĥ   2 ′( k,f ) X ( k,f ).
 
   This frequency filtering step may also be carried out in an equivalent manner by a temporal filtering, that is to say by directly filtering x w (k,n) by ĥ 2 ′(k,n). In this case, we obtain directly the signal ŝ(k,n) in the time domain. 
   When the filtering is performed in the frequency domain, it is then advisable to return to the time domain. This is carried out by the functional block  11  of the device  1 , by applying an inverse spectral transform (IDFT), i.e.:
 
 ŝ ( k,n )=IDFT( Ŝ   2 ( k,f )).
 
   The denoised output signal is thereafter advantageously synthesized by a block overlap and add technique (OLA standing for “overlap add”) or else a block save technique (OLS standing for “overlap save”). This operation of reconstructing the signal in the time domain is called inverse short-term Fourier transformation (ISTFT). 
   On completion of these steps, the frame  12  is obtained at the output of the device  1 . This output frame  12  is a denoised estimate of the input frame  2 , in which the distortion is limited. 
   In the rest of the description, advantageous choices of implementation of the invention are described. These choices are exemplary embodiments and should not be regarded as having a limiting nature. 
   For the implementation of the analysis and the synthesis that are carried out by the functional blocks  3  and  11  of the device  1 , the OLA technique is advantageously chosen. There exists moreover an overlap between the frames analyzed. For example, an overlap rate of 50% may be chosen. The current frame  2  of the noisy signal x(k,n) {n=0, . . . , L−1}, of length L, is firstly multiplied by a weighting window w(n) which delivers the signal x w (k,n):
 
 x   w ( k,n )= w ( n ) x ( k,n ).
 
   In the example, the window used is a Hanning window of size L: 
               w   ⁡     (   n   )       =       0.5   -     0.5   ⁢     cos   ⁡     (       π   ⁢           ⁢   n     L     )       ⁢           ⁢   n       =   0       ,   …   ⁢           ,     L   -   1     ,         
where L=256, this corresponding to frames of length 32 ms at the sampling frequency Fe=8 kHz.
 
   The switch to the frequency domain is carried out advantageously by using the fast Fourier transform (FFT), numerically optimized version of the discrete Fourier transform (DFT). To get a better frequency resolution and to satisfy the linear convolution constraints, an FFT of length L FFT =512 is used. It is therefore necessary beforehand to lengthen the windowed block x w (k,n) by 256 null samples (“zero-padding”) 
   
     
       
         
           
             
               x 
               FFT 
             
             ⁡ 
             
               ( 
               
                 k 
                 , 
                 n 
               
               ) 
             
           
           = 
           
             { 
             
               
                 
                   
                     
                       
                         
                           x 
                           w 
                         
                         ⁡ 
                         
                           ( 
                           
                             k 
                             , 
                             n 
                           
                           ) 
                         
                       
                       , 
                     
                   
                   
                     
                       
                         n 
                         = 
                         0 
                       
                       , 
                       … 
                       ⁢ 
                       
                           
                       
                       , 
                       
                         L 
                         - 
                         1 
                       
                     
                   
                 
                 
                   
                     
                       0 
                       , 
                     
                   
                   
                     
                       
                         n 
                         = 
                         L 
                       
                       , 
                       … 
                       ⁢ 
                       
                           
                       
                       , 
                       
                         
                           L 
                           FFT 
                         
                         - 
                         1 
                       
                     
                   
                 
               
               . 
             
           
         
       
     
   
   The FFT X(k,f) of the input frame is then obtained by:
 
 X ( k,f )=FFT(( x   FFT ( k,n )).
 
   A detection of vocal activity may be carried out, optionally, by the functional block  4 . 
   As far as the estimation of the noise power spectral density {circumflex over (γ)} bb (k b ,f) is concerned, carried out by the functional block  5 , it is then possible to proceed advantageously in the following manner. For each “noise only” frame k b  detected by the DVA, the noise power spectral density {circumflex over (γ)} bb (k b ,f) is estimated by the following recursive expression:
 
{circumflex over (γ)} bb ( k   b   ,f )=α{circumflex over (γ)} bb ( k   b −1 ,f )+(1−α)| X ( k   b   ,f )| 2 .
 
   If the current frame is not detected as a noise frame, then the estimate of the power spectral density of the noise is pegged. The smoothing quantity α is chosen constant and equal to α=0.8825, this corresponding to a time constant of 128 ms (having regard to the 50% overlap and the sampling frequency Fe=8 kHz), judged sufficient to ensure a compromise between a reliable estimate and a tracking of the temporal variations of the statistics of the noise. 
   The transfer function Ĥ step1 (k,f) of the first noise reduction filter estimated by the functional block  7  is advantageously that of an open-loop Wiener filter such that: 
                   H   ^       step   ⁢           ⁢   1       ⁡     (     k   ,   f     )       =       η   ⁡     (     k   ,   f     )         1   +     η   ⁡     (     k   ,   f     )             ,         
where η(k,f) represents the a priori signal-to-noise ratio (SNR) defined theoretically by the ratio of the estimates of the PSDs of the useful signal {circumflex over (γ)} ss (k,f) and of the noise {circumflex over (γ)} bb (k,f), that is to say:
 
   
     
       
         
           
             η 
             ⁡ 
             
               ( 
               
                 k 
                 , 
                 f 
               
               ) 
             
           
           = 
           
             
               
                 
                   
                     γ 
                     ^ 
                   
                   ss 
                 
                 ⁡ 
                 
                   ( 
                   
                     k 
                     , 
                     f 
                   
                   ) 
                 
               
               
                 
                   
                     γ 
                     ^ 
                   
                   bb 
                 
                 ⁡ 
                 
                   ( 
                   
                     k 
                     , 
                     f 
                   
                   ) 
                 
               
             
             . 
           
         
       
     
   
   The PSD of the useful signal {circumflex over (γ)} ss (k,f) is advantageously obtained using a directed-decision estimator, i.e.:
 
{circumflex over (γ)} ss ( k,f )=β| Ŝ ( k− 1 ,f )| 2 +(1−β) P[|   X ( k,f )| 2 −{circumflex over (γ)} bb ( k,f )],
 
where the barycentric parameter β is chosen constant and equal to β=0.98. As indicated above, the function P ensures the thresholding of the quantity weighted by (1−β):
 
   
     
       
         
           
             P 
             ⁡ 
             
               [ 
               
                 z 
                 ⁡ 
                 
                   ( 
                   
                     k 
                     , 
                     f 
                   
                   ) 
                 
               
               ] 
             
           
           = 
           
             { 
             
               
                 
                   
                     
                       z 
                       ⁡ 
                       
                         ( 
                         
                           k 
                           , 
                           f 
                         
                         ) 
                       
                     
                   
                   
                     
                       
                         if 
                         ⁢ 
                         
                             
                         
                         ⁢ 
                         
                           z 
                           ⁡ 
                           
                             ( 
                             
                               k 
                               , 
                               f 
                             
                             ) 
                           
                         
                       
                       &gt; 
                       0 
                     
                   
                 
                 
                   
                     0 
                   
                   
                     otherwise 
                   
                 
               
               . 
             
           
         
       
     
   
   The transfer function Ĥ step1 (k,f) is thereafter advantageously reused to refine the estimate of the PSD of the useful signal {circumflex over (γ)} ss (k,f). One thus obtains the quantity {circumflex over (γ)} ss2 (k,f) given by:
 
{circumflex over (γ)} ss2 ( k,f )=| Ĥ   step1 ( k,f ) X ( k,f )| 2 ,
 
which is used to do a second estimate of the a priori signal-to-noise ratio η 2 (k,f), given by:
 
   
     
       
         
           
             
               η 
               2 
             
             ⁡ 
             
               ( 
               
                 k 
                 , 
                 f 
               
               ) 
             
           
           = 
           
             
               
                 
                   
                     γ 
                     ^ 
                   
                   
                     ss 
                     ⁢ 
                     
                         
                     
                     ⁢ 
                     2 
                   
                 
                 ⁡ 
                 
                   ( 
                   
                     k 
                     , 
                     f 
                   
                   ) 
                 
               
               
                 
                   
                     γ 
                     ^ 
                   
                   bb 
                 
                 ⁡ 
                 
                   ( 
                   
                     k 
                     , 
                     f 
                   
                   ) 
                 
               
             
             . 
           
         
       
     
   
   A second pass then advantageously consists in the calculation of the estimator of the transfer function Ĥ 1 (k,f) of the first noise reduction filter on the basis of η 2 (k,f), i.e.: 
   
     
       
         
           
             
               
                 H 
                 ^ 
               
               1 
             
             ⁡ 
             
               ( 
               
                 k 
                 , 
                 f 
               
               ) 
             
           
           = 
           
             
               
                 
                   η 
                   2 
                 
                 ⁡ 
                 
                   ( 
                   
                     k 
                     , 
                     f 
                   
                   ) 
                 
               
               
                 1 
                 + 
                 
                   
                     η 
                     2 
                   
                   ⁡ 
                   
                     ( 
                     
                       k 
                       , 
                       f 
                     
                     ) 
                   
                 
               
             
             . 
           
         
       
     
   
   The filtering carried out by the functional block  6  may then be carried out as indicated hereinbelow. The (frequency) transfer function of the noise reduction filter Ĥ 1 (k,f) being available, the corresponding temporal response ĥ 1 (k,n) is obtained by advantageously using the inverse fast Fourier transform (IFFT), which is a numerically optimized version of the inverse discrete Fourier transform (IDFT), i.e.:
 
 ĥ   1 ( k,n )=IFFT( Ĥ   1 ( k,f )).
 
   Because the signal ĥ 1 (k,n) is real, the corresponding temporal filter is firstly rendered causal. The L filt1 =256 coefficients of this filter corresponding to the samples that are significant for the application envisaged are selected thereafter. The impulse response thus obtained is thereafter weighted by a Hanning window of length L filt1 . Finally, the weighted impulse response is supplemented with L FFT −L filt1  zeros, to give the impulse response ĥ 1 ′(k,n) which complies with the aliasing-free linear convolution constraint. 
   The transfer function of the noise reduction filter Ĥ 1 ′(k,f) is thereafter obtained by advantageously using the fast Fourier transform (FFT), numerically optimized version of the discrete Fourier transform (DFT), i.e.:
 
 Ĥ   1 ′( k,f )=FFT( ĥ   1 ′( k,n )).
 
   Finally, the first denoised version of the noisy signal is obtained by frequency filtering, that is to say by multiplying the spectrum of the noisy signal and the transfer function of the noise reduction filter:
 
 Ŝ   1 ( k,f )= Ĥ   1 ′( k,f ) X ( k,f ).
 
   In the step of restoring the lost harmonics that is performed by the functional block  8 , the signal Ŝ 1 (k,f) is thereafter used to create the artificial signal Ŝ harmo (k,f). Firstly, the temporal signal ŝ 1 (k,n) is obtained by advantageously using the inverse fast Fourier transform (FFTI):
 
 ŝ   1 ( k,n )=IFFT( ŝ   1 ( k,f )).
 
   Thereafter, the single-wave nonlinear rectification function is applied to this signal, this having the aim of restoring harmonicity to it, i.e.:
 
 ŝ   harmo ( k,n )=max( ŝ   1 ( k,n ),0).
 
   Of course, any other nonlinearity than the single-wave rectification may be used as indicated above (absolute value, maximum/minimum with respect to a threshold, etc.). Finally, the signal Ŝ harmo (k,f) is obtained by advantageously using the fast Fourier transform (FFT), i.e.:
 
 Ŝ   harmo ( k,f )=FFT( ŝ   harmo ( k,n )).
 
   This signal is thereafter exploited to advantageously calculate the transfer function of a second noise reduction filter. 
   In the present exemplary embodiment, the estimate of such a second filter is carried out by the functional block  9  of the device  1  in the following manner. The transfer function of the second noise reduction filter Ĥ 2 (k,f) is that of an open-loop Wiener filter such that: 
                   H   ^     2     ⁡     (     k   ,   f     )       =         η   harmo     ⁡     (     k   ,   f     )         1   +       η   harmo     ⁡     (     k   ,   f     )             ,         
where ηharmo(k,f) represents the a priori signal-to-noise ratio (SNR) defined theoretically by the ratio of the estimates of the PSD of the useful signal {circumflex over (γ)} ss     —     harm (k,f) and of the PSD of the noise {circumflex over (γ)} bb (k,f), i.e.:
 
                 η   harmo     ⁡     (     k   ,   f     )       =           γ   ^       ss_harmo   ⁢               ⁡     (     k   ,   f     )             γ   ^     bb     ⁡     (     k   ,   f     )           ,         
where the PSD of the useful signal {circumflex over (γ)} ss     —     harmo (k,f) is obtained thus:
 {circumflex over (γ)} ss     —     harmo ( k,f )=ρ( k,f )| Ŝ   1 ( k,f )| 2 +(1−ρ( k,f ))| Ŝ   harmo ( k,f )| 2 . 
   The parameter ρ(k,f) for controlling the level of reinjection of the signal Ŝ harmo (k,f) is chosen fixed and equal to 0.5 for the application envisaged. 
   The (frequency) transfer function of the noise reduction filter Ĥ 2 (k,f) being available, the corresponding temporal response ĥ 2 (k,n) is obtained by advantageously using the inverse fast Fourier transform (FFTI), numerically optimized version of the inverse discrete Fourier transform (IDFT):
 
 ĥ   2 ( k,n )=IFFT( Ĥ   2 ( k,f )).
 
   Because the signal ĥ 2 (k,n) is real, the corresponding temporal filter is firstly rendered causal. The L filt2 =256 coefficients of this filter corresponding to the samples that are significant for this application are selected thereafter. The impulse response thus obtained is thereafter weighted by a Hanning window of length L filt2 . Finally, the weighted impulse response is supplemented with L FFT =L filt2  zeros to give the impulse response ĥ 2 (k,n) which complies with the aliasing-free linear convolution constraint. 
   The transfer function of the noise reduction filter Ĥ 2 ′(k,f) is thereafter obtained by advantageously using the fast Fourier transform (FFT), numerically optimized version of the discrete Fourier transform (DFT):
 
 Ĥ   2 ′( k,f )=FFT( ĥ   2 ′( k,n )).
 
   The values of the modulus of the noise reduction filter Ĥ 2 ′(k,f) may optionally be thresheld to avoid too aggressive a noise reduction. 
   Finally, the denoised frequency frame is obtained by frequency filtering carried out by the functional block  10 , that is to say by multiplying the spectrum of the noisy signal and the transfer function of the noise reduction filter, i.e.:
 
 Ŝ   2 ( k,f )= Ĥ   2 ′( k,f ) X ( k,f ).
 
   The return to the time domain of the signal thus obtained is thereafter performed, at the level of the functional block  11 , by inverse spectral transform, by advantageously using the inverse fast Fourier transform (FFTI), that is to say:
 
 ŝ ( k,n )=IDFT( Ŝ   2 ( k,f )).
 
   The denoised output signal is thereafter synthesized by the block overlap and add technique OLA. 
   It will be noted moreover that the device  1  illustrated in  FIG. 2  may be placed in a particular piece of equipment, as a function of the application envisaged, such as a piece of sound pick-up equipment, a piece of communication equipment or else a piece of voice recognition equipment. 
   It will be noted moreover that the processing implemented by the device may be performed in a local manner, as in the case of a communication terminal for example, or else in a manner centralized in a network, for example in the case of a voice recognition server. 
   Additionally, in a practical manner, a device  1  for processing noisy sound signals, according to the invention, consists of (electronic) hardware means and/or software means suitable for implementing a method of processing a signal, according to the invention. 
   According to a preferred implementation, the steps of the method of processing a noisy sound signal, according to the invention, are determined by the instructions of a computer program used in such a piece of equipment according to the invention. 
   The method according to the invention is then implemented when the aforesaid program is loaded into computing means incorporated in the equipment, and the operation of which is then controlled by the execution of the program. 
   Here, the expression “computer program” is understood to mean one or more computer programs forming a (software) set, the purpose of which is the implementation of the invention, when it is executed by appropriate computing means. 
   Consequently, the subject matter of the invention is also such a computer program, in particular in the form of a piece of software stored on an information medium. Such an information medium may consist of any entity or device capable of storing a program according to the invention. 
   For example, the medium in question may comprise a hardware storage means, such as a ROM, for example a CD ROM or a microelectronic circuit ROM, or else a magnetic recording means, for example a hard disk. As a variant, the information medium may be an integrated circuit in which the program is incorporated, the circuit being suitable for executing or for being used in the execution of the method in question. 
   Moreover, the information medium may also be a transmissible immaterial medium, such as an electrical or optical signal that can be conveyed via an electrical or optical cable, by radio or by other means. A program according to the invention may in particular be downloaded from an Internet-type network. 
   From a design point of view, a computer program according to the invention can use any programming language and be in the form of source code, object code, or of code intermediate between source code and object code (for example a partially compiled form), or in any other form desirable for implementing a method according to the invention.