EP0422158A1 - Procede et dispositif d'extraction d'estimateurs des caracteres d'un ecoulement instationnaire a partir d'un signal doppler a valeurs complexes - Google Patents

Procede et dispositif d'extraction d'estimateurs des caracteres d'un ecoulement instationnaire a partir d'un signal doppler a valeurs complexes

Info

Publication number
EP0422158A1
EP0422158A1 EP90905554A EP90905554A EP0422158A1 EP 0422158 A1 EP0422158 A1 EP 0422158A1 EP 90905554 A EP90905554 A EP 90905554A EP 90905554 A EP90905554 A EP 90905554A EP 0422158 A1 EP0422158 A1 EP 0422158A1
Authority
EP
European Patent Office
Prior art keywords
frequency
coefficients
average
signal
correlation
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Withdrawn
Application number
EP90905554A
Other languages
German (de)
English (en)
Inventor
Alain Herment
Guy 3 résidence de Chevreuse DEMOMENT
Jean-Paul Guglielmi
Philippe Dumee
Pierre Peronneau
Alain Barbet
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Institut National de la Sante et de la Recherche Medicale INSERM
Original Assignee
Institut National de la Sante et de la Recherche Medicale INSERM
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Institut National de la Sante et de la Recherche Medicale INSERM filed Critical Institut National de la Sante et de la Recherche Medicale INSERM
Publication of EP0422158A1 publication Critical patent/EP0422158A1/fr
Withdrawn legal-status Critical Current

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01P—MEASURING LINEAR OR ANGULAR SPEED, ACCELERATION, DECELERATION, OR SHOCK; INDICATING PRESENCE, ABSENCE, OR DIRECTION, OF MOVEMENT
    • G01P5/00—Measuring speed of fluids, e.g. of air stream; Measuring speed of bodies relative to fluids, e.g. of ship, of aircraft
    • G01P5/24—Measuring speed of fluids, e.g. of air stream; Measuring speed of bodies relative to fluids, e.g. of ship, of aircraft by measuring the direct influence of the streaming fluid on the properties of a detecting acoustical wave
    • G01P5/241—Measuring speed of fluids, e.g. of air stream; Measuring speed of bodies relative to fluids, e.g. of ship, of aircraft by measuring the direct influence of the streaming fluid on the properties of a detecting acoustical wave by using reflection of acoustical waves, i.e. Doppler-effect
    • G01P5/244—Measuring speed of fluids, e.g. of air stream; Measuring speed of bodies relative to fluids, e.g. of ship, of aircraft by measuring the direct influence of the streaming fluid on the properties of a detecting acoustical wave by using reflection of acoustical waves, i.e. Doppler-effect involving pulsed waves

Definitions

  • the invention relates to the extraction of character estimators of an unsteady flow from a Doppler signal with complex values.
  • the invention relates to the extraction of character estimators of an unsteady flow, from a signal.
  • Sampled complex value Doppler designates a parameter providing an indication representative of a characteristic of the flow at a determined time, such as in particular the average Doppler frequency (representing the average speed of the flow) and the average bandwidth of the spectrum. Doppler (representing the dispersion of velocities in the flow).
  • the invention finds a particularly important application in medical ultrasound imaging with Doppler effect: the average frequency then gives an indication of the speed in a vessel or a cavity.
  • the bandwidth makes it possible to determine the flow regime and, in particular, to detect the jet regimes.
  • the average frequency is given by
  • patent application FR n ° 89 03761 for "Method and device for real-time spectral analysis of complex unsteady signals” makes it possible to represent the signals concerned by the present invention using an autoregressive parametric model, whose coefficients remain stationary over short time intervals. We could therefore calculate the estimators mentioned above from the coefficients of the autoregressive model. But this approach results in a large volume of computation and difficulties of computation of integral.
  • the invention aims to provide a method for extracting estimators making it possible to arrive at high precision by calculations on a small number of samples, therefore with a relatively reduced calculation volume for a significant estimator.
  • the invention uses the covariance function of the sampled signal as a calculation aid.
  • the latter being locally stationary, its covariance is stationary in each block of samples and identifies with an autocorrelation function whose coefficients C Y (n), for whole delays, are the Fourier coefficients of the dsp g y (f) which is itself a periodic function of the frequency. So we have :
  • the invention proposes a method for extracting estimators according to which: the Doppler spectrum of a flow is measured and the Doppler signal is sampled, characterized in that: the coefficients of the autocorrelation function are calculated for each block d 'samples; the sum of the correlation coefficients corresponding to the different available delays is calculated, each weighted by the integral of the product of the frequency by a negative exponential as a function of the delay to obtain the average frequency.
  • formula (8) shows a real average frequency, which corresponds to physical reality and that, in the case of a signal with real values, whose correlation function is real and even, the frequency mean is zero.
  • the invention consequently also proposes a method according to which the sum of the imaginary parts of the autocorrelation coefficients of the signal is each assigned a weighting coefficient.
  • the number p of coefficients is limited to a low value which leads to deviations for the frequencies close to the sampling frequency. Since it is often sufficient to obtain a satisfactory representation in a limited frequency range, where the dsp is maximum, the factors 1 / n can be replaced by factors a n chosen to improve the representation in the zone judged essential, at the cost of degradation in other areas.
  • the selection of a n will be carried out by a process of optimization of a polynomial approximation of the part considered essential.
  • the average frequency is thus obtained directly.
  • the invention also provides a device for implementing the above-defined method, which can be entirely wired and therefore capable of operating in real time at high speed, comprising only elements of a simple nature (multipliers and accumulators).
  • FIG. 1 is a block diagram of a circuit for calculating the average frequency and the average bandwidth of an unsteady Doppler signal with complex value, representing the Doppler spectrum of a flow from the real parts and imaginary of its correlation function;
  • FIG. 2 is a block diagram showing a circuit for calculating the imaginary part of the coefficients of the signal autocorrelation function, coefficients acting as inputs for the calculation of the average value in the circuit of Figure 1;
  • Figure 3 similar to Figure 2, shows a circuit for calculating the real part of the coefficients of the autocorrelation function, involved in the calculation of the average bandwidth.
  • the average frequency is calculated by bringing into play all the terms of formula (8). But it will often be possible to truncate the sum by omitting the terms corresponding to a value of n greater than a determined threshold. Indeed, the weight of successive correlation coefficients varies in 1 / n. These coefficients themselves cannot have a module greater than unity and are often of the same order of magnitude. Consequently, their contribution often decreases rapidly with n.
  • the notation Im (Cxx) designates the coefficients Im which are involved in formula (8).
  • the successive coefficients pass through a cascade constituted by a memory S and a network of flip-flops or "latch" L.
  • the coefficients corresponding to the delays successive are applied one after the other to a multiplier 10 whose second input receives the multiplier coefficients 0, 1 / ⁇ , ... (-1) n + 1 / ⁇ n stored in a register 11 or an addressable memory for form the successive terms of the sum of formula (8).
  • the average frequency appears in a memory
  • this value is normalized as a function of the signal strength of which the auto-correlation moment for a zero delay constitutes a representation.
  • this zero delay autocorrelation moment is also used for the calculation of the average bandwidth.
  • the autocorrelation coefficients of the signal are taken directly from the input of the real parts Re. Normalization can be carried out using a divider or, better, as indicated in FIG. 1, using a table 18 constituted by a logical network programmable.
  • This fundamental process can be refined to further improve the quality of the estimate of the average frequency, by acting on the estimator and / or the iterations of the estimator.
  • the estimator can be modified to correct various effects and improve the compromise between the robustness with respect to the approximations introduced by the windowing of the signal and the frequency dynamics.
  • Correcting the windowing effect on the correlation means amounts to correcting the bias introduced into the calculation of the correlation coefficients by the windowing of the signal.
  • the calculation performed is not that of Cy (n) but that of ((pn -) / p) Cy (n), where: p is the number of points in the window,
  • the summation is performed for n varying from 1 to p.
  • the correction of the windowing effect on the weighted sum of the correlation coefficients is carried out by using more correlation moments in the calculation of the mean frequency.
  • f ⁇ [(1 / ⁇ n) (-1) n + 1 Im ⁇ Cy (n)>)] with the summation performed for n varying from 1 to p by the sum: f 2 - ⁇ [(1 / ⁇ ) (-1) 3m-1 + 1 a ⁇ (3m-1) im ⁇ C'y (3m-1) >].
  • the summation is carried out for m varying from 1 to p and 1 varying from 0 to 2.
  • R (n) [(Re ⁇ Cy (n)>) 2 + (Im ⁇ Cy (n)>) 2 ] 1/2
  • R (n-1) [(Re ⁇ Cy (n-1)>) 2 + (Im ⁇ Cy (n-1)>) 2 ] 1/2
  • R (3m-1) [2R (n) + R (n-1)] / 3
  • ⁇ (n-1) tg -1 [(Im ⁇ Cy (n-1)>) / (Re ⁇ Cy (n-1)>)]
  • ⁇ (3m-1) [2 ⁇ (n) + ⁇ (n-1)] / 3
  • ⁇ (3M-2) [ ⁇ (n) + 2 ⁇ (n-1)] / 3
  • C'y (3m-2) R (3m-2) [cos ⁇ (3m-2) + i.sin ⁇ (3m-2)].
  • f u k + 1 f k + 1 , f 1 k + 1 , f 2 k + 1 or f 3 k + 1
  • the average frequency already estimated at the previous instant, at the previous pixel of the same line or finally at the previous line of the same image is used to calculate the average frequency.
  • f u k represents one of the estimators f, f 1 , f 2 or fg at time k or at pixel k or at line k.
  • the estimated average frequency is not to vary too much between k and k + 1.
  • ⁇ ⁇ + ⁇ 0 ⁇ 0 , ⁇ 1 and ⁇ 2 are set by the user.
  • ⁇ 0 is the increment used to modify ⁇ .
  • ⁇ 1 prohibits rapid variations of the average frequency around zero.
  • the average frequency is calculated by the estimator f 1 .
  • the average frequency at time k + 1 is calculated as the sum of the average frequency at time k and a change in average frequency at time k + 1 according to l '' following algorithm: Initialization:
  • f u k + 1 f k + 1 f 1 k + 1 f 2 k + 1 , or f 3 k + 1
  • the average bandwidth used to calculate the circuit illustrated in Figure 1 is defined by formula (2) for a standardized dsp.
  • the average frequency being known, it remains to calculate (or the variance ⁇ ) which constitutes the second order moment of the dsp:
  • the part of the circuit of Figure 1 providing the mean of the square of the frequency has a constitution similar to that of the part of calculation of The only difference is to add 1/12 to each input value from the looped register 12a, before normalization.
  • the average value can be calculated from only the imaginary parts of the autocorrelation coefficients which can be calculated directly from the samples of the signal, which has a real part Re and an imaginary part Im.
  • the circuit shown in Figure 2 has two channels with the same constitution.
  • the samples Re and Im of the Doppler signal are applied to two respective inputs 22 and 24.
  • each channel comprises, in FIG. 2, two sets of shift registers on each track, operating alternately.
  • the samples Re of the real part are alternately directed, by the switch 30, to two registers 26 and 26a and two registers 28 and 28a.
  • the samples Im are alternately directed, by a switch 32 controlled in synchronism with the switch 30, to registers 34 and 34a and registers 36 and 36a.
  • Multipliers 38 and 40 receive the signals from the registers through switches controlled in synchronism with switches 30 and 32 and carry out the partial sums of the samples of the same block, designated by v1 to vn in Figure 2, are first in phase, then with successive shifts.
  • the circuit shown in Figure 1 can be completed to go beyond the second order moment, in particular in the case of a highly disturbed signal or with a very asymmetrical spectrum.
  • An alternative embodiment of the invention uses, as input data, not the samples of the real and imaginary parts of the complex signal, but the parameters of an autoregressive model obtained for example by implementing the method according to the request âe Patent "Method and device for real time spectral analysis of complex unsteady signals" FR n ° 89 03761.
  • the inverse Levinson algorithm is a fast algorithm requiring order of p 2 elementary arithmetic operations to pass from a prediction vector a p to the corresponding correlation matrix.
  • the matrix obtained is a Toeplitz matrix satisfying the optimality conditions and guaranteeing high precision, given the absence of approximation.
  • the application of the method described above to medical Doppler imaging has shown satisfactory results on short blocks, of only eight samples, representative of current working conditions in color Doppler imaging.
  • the method according to the invention revealed a lower variance than that of the conventional methods for determining the average frequency by correlation methods, on the estimation of the average frequency in diastole.

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Acoustics & Sound (AREA)
  • Multimedia (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • General Physics & Mathematics (AREA)
  • Ultra Sonic Daignosis Equipment (AREA)

Abstract

Le procédé et le dispositif sont notamment utilisables pour extraire les estimateurs du signal Doppler en imagerie médicale ultrasonore. Le dispositif utilise comme paramètre d'entrée, pour calculer la fréquence moyenne, les coefficients de corrélation Im(Cxx) correspondant à la partie imaginaire du signal. Ces coefficients sont appliqués en séquence à une entrée d'un multiplieur (10). L'autre entrée reçoit des coefficients de pondération de la forme 0, 1/0, 1/pi,...,(-1)p+1/pip, où p est le nombre de coefficients de pondération, stockés dans un registre (11). Les produits partiels sont accumulés et fournissent une somme représentative de la fréquence moyenne.

Description

Procédé et dispositif d'extraction d'estimateurs des caractères d'un écoulement instationnaire à partir d'un signal Doppler à valeurs complexes L'invention concerne l'extraction d'estimateurs des caractères d'un écoulement instationnaire, à partir d'un signal Doppler à valeurs complexes échantillonné. Le terme "estimateur" désigne un paramètre fournissant une indication représentative d'une caractéristique de l'écoulement à un instant déterminé, telle notamment que la fréquence Doppler moyenne (représentant la vitesse moyenne de l'écoulement) et la largeur de bande moyenne du spectre Doppler (représentant la dispersion des vitesses dans l'écoulement).
L'invention trouve une application particulièrement importante en imagerie médicale ultrasonore à effet Doppler : la fréquence moyenne donne alors une indication sur la vitesse dans un vaisseau ou une cavité. La largeur de bande permet de déterminer le régime d'écoulement et, en particulier, de déceler les régimes de jet.
II existe des procédés de calcul permettant de déterminer la fréquence moyenne et la largeur de bande moyenne d'un signal Doppler, pour un bloc déterminé de données sur lesquelles le signal est quasi stationnaire, à partir de la densité spectrale de puissance ou d.s.p., gY(f) pour le bloc de données.
La fréquence moyenne est donnée par
(1) La largeur de bande moyenne Δf est donnée par (2)
' J
où est la fréquence moyenne et
la moyenne des carrés de donnée par
/
( 2bis )
Les formules (2bis) et (1) supposent que la d.s.p. est préalablement normalisée, c'est-à-dire que l'on a
On connaît une méthode de calcul de la fréquence moyenne par une technique d'auto-corrélation (Kasai et al, "Real time two-dimensional blood flow imaging using an auto-correlation technique", IEEE TRANSACTIONS on sonics and ultrasonics, Vol. 32, pp. 458-464, mai 1985). Elle a l'inconvénient de partir de l'hypothèse que le module de la fonction de corrélation du signal considéré a une dérivée nulle au point zéro, alors que cette hypothèse est loin d'être toujours vérifiée.
Par ailleurs, la demande de brevet FR n° 89 03761 pour "Procédé et dispositif d'analyse spectrale en temps réel de signaux instationnaires complexes" permet de représenter les signaux concernés par la présente invention à l'aide d'un modèle paramétrique autorégressif, dont les coefficients restent stationnaires sur les intervalles de temps courts. On pourrait en conséquence calculer les estimateurs mentionnés plus haut à partir des coefficients du modèle autorégressif. Mais cette approche se traduit par un volume de calcul important et des difficultés de calcul d'intégrale.
L'invention vise à fournir un procédé d'extraction d'estimateurs permettant d'arriver à une précision élevée par des calculs sur un petit nombre d'échantillons, donc avec un volume de calcul relativement réduit pour un estimateur significatif.
Dans ce but, l'invention utilise comme auxiliaire de calcul la fonction de covariance du signal échantillonné. Ce dernier étant localement stationnaire, sa covariance est stationnaire dans chaque bloc d'échantillons et s'identifie à une fonction d'autocorrélation dont les coefficients CY(n), pour des retards entiers, sont les coefficients de Fourier de la d.s.p. gy(f) qui est elle-même une fonction périodique de la fréquence. On a donc :
dans laquelle :
On voit qu'on remplace les intégrales (1) et (2bis) par des sommes pondérées des coefficients de la fonction d'autocorrélation du signal. Dans la formule (4), le poids affecté aux valeurs successives CY(n) de la fonction de corrélation décroit très rapidement avec n, de sorte qu'un petit nombre de termes peut suffire, lorsque les coefficients CY ont le même ordre de grandeur.
On arrive ainsi à une valeur de la fréquence moyenne donnée par :
On voit que la fréquence moyenne peut être calculée de façon simple à partir des coefficients de la fonction de corrélation.
Ces coefficients peuvent être calculés :
- soit directement à partir des échantillons,
- soit à partir des coefficients d'un modèle paramétrique autorégressif, par un algorithme rapide, par exemple du type "Levinson inverse" qui sera rappelé plus loin.
Dans ce dernier cas, on peut utiliser, sans les inconvénients mentionnés plus haut, les résultats fournis dans la demande de brevet du même jour déjà mentionnée.
L'invention propose un procédé d'extraction d'estimateurs suivant lequel : on mesure le spectre Doppler d'un écoulement et on échantillonne le signal Doppler , caractérisé en ce que : on calcule les coefficients de la fonction d'autocorrélation pour chaque bloc d'échantillons ; on calcule la somme de coefficients de corrélation correspondant aux différents retards disponibles, pondérés chacun par l'intégrale du produit de la fréquence par une exponentielle négative fonction du retard pour obtenir la fréquence moyenne.
On peut considérer la fréquence moyenne comme le moment d'ordre 1 du spectre du signal Doppler complexe, la largeur de bande du spectre comme le moment d'ordre 2 du spectre et d'autres caractéristiques (dissymétrie du spectre ou aplatissement) comme représentées par des moments d'ordre supérieur.
Il est possible de simplifier notablement les calculs d'extraction en analysant et en adaptant la formule (5). Cette formule fait apparaître que le calcul de la fréquence moyenne se ramène au calcul d'une somme de coefficients de corrélation pondérés par des intégrales de la forme :
(6)
Or, on a, pour n = 0, une intégrale de la forme f.df qui est nulle et on peut écrire, de façon générale, les intégrales pour n 0 sous la forme :
( 7 ) On en déduit, en tenant compte des propriétés de symétrie hermitienne des coefficients d'autocorrélation d'un signal à valeurs complexes, la formule (8) :
( 8 ) On voit qu'il suffit, pour calculer la fréquence moyenne ee disposer de la partie imaginaire Im(Cxx) des coefficients d'autocorrélation de la fonction d ' entrée échantillonnée.
On peut encore remarquer que la formule (8) fait apparaître une fréquence moyenne réelle, ce qui correspond à la réalité physique et que, dans le cas d'un signal à valeurs réelles, dont la fonction de corrélation est réelle et paire, la fréquence moyenne est nulle.
L'invention propose en conséquence également un procédé suivant lequel on calcule la somme des parties imaginaires des coefficients d'autocorrélation du signal affectées chacune d'un coefficient de pondération.
Si le nombre de coefficients disponibles était infini, la valeur du coefficient de pondération serait fonction uniquement de n pour l'approximation la meilleure possible.
Mais, dans la pratique, le nombre p de coefficients est limité à une valeur faible qui conduit à des écarts pour les fréquences proches de la fréquence d'échantillonnage. Etant donné qu'il est souvent suffisant d'obtenir une représentation satisfaisante dans une plage de fréquence limitée, là où la d.s.p. est maximale, les facteurs 1/n peuvent être remplacés par des facteurs an choisis pour améliorer la représentation dans la zone jugée essentielle, au prix d'une dégradation dans d'autres zones. La sélection de an s'effectuera par un processus d'optimisation d'une approximation polynomiale de la partie jugée essentielle.
On obtient ainsi directement la fréquence moyenne.
L'invention propose également un dispositif permettant de mettre en oeuvre le procédé ci-dessus défini, pouvant être entièrement câblé et donc pouvant fonctionner en temps réel à vitesse élevée, ne comportant que des éléments de nature simple (multiplieurs et accumulateurs). L'invention sera mieux comprise à la lecture de la description qui suit de modes particuliers de mise en oeuvre, donnés à titre d'exemples non limitatifs. La description se réfère aux dessins qui l'accompagnent, dans lesquels :
- la Figure 1 est un synoptique d'ensemble d'un circuit de calcul de la fréquence moyenne et de la largeur de bande moyenne d'un signal Doppler instationnaire à valeur complexe, représentant le spectre Doppler d'un écoulement à partir des parties réelles et imaginaires de sa fonction de corrélation ;
- la Figure 2 est un synoptique montrant un circuit de calcul de la partie imaginaire des coefficients de la fonction d'autocorrélation du signal, coefficients intervant en tant qu'entrées pour le calcul de la valeur moyenne dans le circuit de la Figure 1 ;
- la Figure 3, similaire à la Figure 2, montre un circuit de calcul de la partie réelle des coefficients de la fonction d'autocorrélation, intervenant dans le calcul de la largeur de bande moyenne.
On décrira successivement le calcul de la fréquence moyenne et celui de la largeur de bande moyenne Δf à l'aide du dispositif selon l'invention montré en Figure 1. Il faut cependant noter que, souvent, il n'est pas nécessaire pour l'utilisateur du dispositif de disposer à la fois de la fréquence moyenne et de la largeur de bande et que, dans ce cas, le dispositif peut être simplifié par suppression des composants inutiles.
Calcul de la fréquence moyenne
Dans le mode de réalisation de la Figure 1, la fréquence moyenne est calculée en mettant en jeu tous les termes de la formule (8). Mais il sera souvent possible de tronquer la somme en omettant les termes correspondant à une valeur de n supérieure à un seuil déterminé. En effet, le poids des coefficients de corrélation successifs varie en 1/n. Ces coefficients eux-mêmes ne peuvent avoir un module supérieur à l'unité et sont souvent du même ordre de grandeur. En conséquence, leur contribution décroît souvent rapidement avec n.
Sur la Figure 1, la notation Im (Cxx) désigne les coefficients Im qui interviennent dans la formule (8). Pour permettre la synchronisation de l'ensemble du circuit (à partir d'une horloge centrale non représentée), les coefficients successifs transitent par une cascade constituée par une mémoire S et un réseau de bascules ou "latch" L. Les coefficients correspondant aux retards successifs sont appliqués l'un après l'autre à un multiplieur 10 dont la seconde entrée reçoit les coefficients multiplicateurs 0, 1/π,... (-1)n+1/πn stockés dans un registre 11 ou une mémoire adressable pour former les termes successifs de la somme de la formule (8).
Les produits successifs sont accumulés pour toutes les valeurs de n = 1 à n = +p dans un sommateur constitué par exemple par un registre à décalage 12 bouclé par un réseau de bascules 14 formant "latch".
La fréquence moyenne apparaît dans une mémoire
16 à la sortie du registre. Suivant la formule (8), on normalise cette valeur en fonction de la puissance du signal dont le moment d'auto-corrélation pour un retard nul constitue une représentation. Sur la Figure 1, ce moment d'auto-corrélation à retard nul est également utilisé pour le calcul de la largeur de bande moyenne. Il est directement prélevé sur l'entrée des parties réelles Re des coefficients d'auto-corrélation du signal. La normalisation peut s'effectuer à l'aide d 'un diviseur ou, mieux, comme indiqué sur la Figure 1, à l'aide d'une table 18 constituée par un réseau logique programmable.
Ce processus fondamental peut être affiné pour améliorer encore la qualité de l'estimation de la fréquence moyenne, en agissant sur l'estimateur et/ou les itérations de l'estimateur.
I - L'estimateur peut être modifié pour corriger divers effets et améliorer le compromis entre la robustesse vis-à-vis des approximations introduites par le fenêtrage du signal et la dynamique fréquentielle.
1. La correction de l'effet de fenêtrage sur les moyens de corrélation revient à corriger le biais introduit dans le calcul des coefficients de corrélation par le fenêtrage du signal. Le calcul effectué n'est pas celui de Cy(n) mais celui de ((pn-)/p) Cy(n), où : p est le nombre de points dans la fenêtre,
n le retard du coefficient de corrélation. L'effet de Cy(1) étant de loin le plus important dans l'estimation de f, on propose comme estimateur : f1 = (P/(P-1)) f
fχ = (P/(P-1)) ∑ [(1/πn) (-1)n+1 Im<Cy(n)>]
La sommation est effectuée pour n variant de 1 à p.
2. La. correction de l'effet de fenêtrage sur la somme pondérée des coefficients de corrélation est effectuée en utilisant davantage de moments de corrélation dans le calcul de la fréquence moyenne.
Pour une meilleure estimation, on remplace la somme : f = Σ [(1/πn) (-1)n+1 Im<Cy(n)>)] avec la sommation effectuée pour n variant de 1 à p par la somme : f2 - ∑ [(1/π) (-1)3m-1+1 aα(3m-1) im<C'y(3m-1)>].
La sommation est effectuée pour m variant de 1 à p et 1 variant de 0 à 2.
Les valeurs C'y(m-1) sont obtenues par extrapolation non linéaire à partir des valeurs Cy(n).
Dans ce qui suit, la fonction de corrélation
Cy(n) sera considérée comme une fonction complexe de la forme :
R(n) exp(iθ(n)) de module R(n) et d'argument θ.
R(n) = [(Re<Cy(n)>)2 + (Im<Cy(n)> )2]1/2
R(n-1) = [(Re<Cy(n-1)>)2 + (Im<Cy(n-1)>)2]1/2
R(3m) = R(n)
R(3m-1) = [2R(n)+R(n-1)]/3
R(3m-2) = [R(n)+2R(n-1)]/3 θ(n) = tg-1[(Im<Cy(n)>) / (Re<Cy(n)>)]
θ(n-1)= tg-1[(Im<Cy(n-1)>) / (Re<Cy(n-1)>)]
θ(3m) = θ(n)
θ(3m-1) = [2θ(n)+θ(n-1)]/3
θ(3M-2) = [θ(n)+2θ(n-1)]/3
C'y(3m) = R(3m) [cosθ(3m)+i.sinθ(3m)]
C'y(3m-1) = R(3m-1) [cosθ(3m-1)+i.sinθ(3m-1)]
C'y(3m-2) = R(3m-2) [cosθ(3m-2)+i.sinθ(3m-2)].
Les coefficients aα(m-1) sont calculés comme suit, pour α réel compris entre [0 et 1] et fixé par l'utilisateur. aα ( 3m ) = 3 / 3m = 1 / m
aα( 3m-1 ) = 3α / ( 3m-1 )
aα( 3m-2 ) = 3α / ( 3m-2 ) . II faut noter que, si α = 0 : f2 = Σ (1/π) [(-1)3m+1 (1/3m) Im<C'y(3m)>]
f2 = f car n - 3m. 3. Enfin, la correction de l'effet de fenêtrage sur les moments de corrélation et sur la somme pondérée des coefficients de corrélation s'effectue en cumulant les deux améliorations précédentes 1 et 2 : f3 = [(1+2α)p/((1+2α)p-1)] f2
f3 = [(1+2α)p/((1+2α)p-1)] x Σ [(1/π)(-1)3m-1+1 aα(3m-1) Im <C'y(3m-1)>]. II faut noter que, si α = 0 : f3 = f1
Si a = 1 f3 _ (1/2π) tg-1 [Im<Cy(1)>/Re<Cy(1)>].
II - Les itérations de l'estimateur peuvent être modifiées pour dépasser le compromis Robustesse/Dynamique fréquentielle des estimateurs précédents.
1. Au niveau de l'estimation locale, on utilise un décalage itératif du domaine d'intégration dans le calcul de la fréquence moyenne de manière à s'affranchir du phénomène de Gibbs lié à la troncature de la somme sur les moments de corrélation. On dispose de p échantillons du signal isolés comme précédemment.
On utilise un des quatre estimateurs précédents f, f1 , f2 ou f3 appelé fu dans la suite.
On boucle l'estimateur sur lui-même selon l'enchaînement suivant :
Initialisation :
k = 0
fu0 = 0
Incrémentation de l'indice :
k = k+1 Décalage fréquentiel du signal Doppler :
yk+1(t) = yk(t) exp (-2πfu kt)
Calcul des moments de corrélation du signal Doppler : Cy(n) ou C'y(m-1)
Calcul de la fréquence moyenne :
fu k+1 = fk+1, f1 k+1, f2 k+1 ou f3 k+1
Test de fin d'itérations :
Signe [fu k+1] 4 Signe [fu k] ou :
k > k0
k0 défini à l'avance par l'utilisateur, ou :
les deux tests combinés.
Si le test est validé, alors :
f = fu k+1
Sinon, retour à l'incrémentation de l'indice. 2. Au niveau de l'estimation dynamique, on utilise la fréquence moyenne déjà estimée à l'instant précédent, au pixel précédent de la même ligne ou enfin à la ligne précédente de la même image pour calculer la fréquence moyenne.
fu k représente l'un des estimateurs f, f1, f2 ou fg à l'instant k ou au pixel k ou à la ligne k.
(a) Pour améliorer la "douceur" fréquentielle, on impose à la fréquence moyenne estimée de ne pas trop varier entre k et k+1.
On utilise pour cela le synoptique suivant :
Initialisation :
k = 0
fu0 = 0
α0, α1, a2
α = α + α0 α0, α1 et α2 sont fixés par l'utilisateur.
α0 est l'incrément utilisé pour modifier α.
α1 interdit les variations rapides de la fréquence moyenne autour de zéro.
α2 règle la contrainte de douceur. Incrémentation de l'indice de fenêtre, de pixel ou de ligne : k = k+1
Calcul de la valeur de o : α = α-α0. Calcul des coefficients de corrélation du signal Doppler : C'y(m-1)
Calcul de la fréquence moyenne : fu k+1 = f2 k+1 ou f3 k+1 Test de fin d'itération : abs [(fu k+1-fu k)/(fu k+1l)] <α2 ou α < 0 Si le test est validé : f = fu k+1
Sinon, retour à l'incrémentation de l'indice. II faut noter que, si α0 = 1, l'enchaînement des opérations se réduit comme indiqué plus haut en 1.3.
Estimation de la fréquence moyenne par l'estimateur : f4 = (1/2π) tg-1 [Im<Cy(1)>/Re<Cy(1)>]
Test : abs[(fu k+1-fu k)/(fu k+α1)] <α2. Si le test est validé, on garde f.
Sinon, on calcule la fréquence moyenne par l'estimateur f1.
(b) Pour améliorer le suivi fréquentiel, la fréquence moyenne à l'instant k+1 est calculée comme la somme de la fréquence moyenne à l'instant k et d'une variation de fréquence moyenne à l'instant k+1 selon l'algorithme suivant : Initialisation :
k = 0
fu 0 = 0
Incrémentation de l'indice de fenêtre, de pixel ou de ligne :
k = k+1
Décalage fréquentiel du signal Doppler :
yk+1(t) = yk(t) exp (-2πfu kt)
Calcul des moments de corrélation du signal Doppler : Cy( n) ou C ' y( m-1 )
Calcul de la fréquence moyenne incrémentale :
fu k+1 = fk+1 f1 k+1 f2 k+1, ou f3 k+1
Calcul de la fréquence moyenne
f = fu k+1+fu k
Calcul de la largeur de bande moyenne
La largeur de bande .moyenne que permet de calculer le circuit illustré en Figure 1 est définie par la formule (2) pour une d.s.p. normalisée. La fréquence moyenne étant connue, il reste à calculer (ou la variance σ) qui constitue le moment d'ordre 2 de la d.s.p. :
(9) qui peut s'écrire
(10)
Le calcul est très similaire à celui de : il consiste à faire la somme de coefficients d'autocorrélation CY(n) pondérés par les intégrales :
( 11 ) L'intégrale est égale à 1/12 pour n = 0 ; pour n 0, l'intégrale d'ordre n est égale à :
(12)
Toujours en utilisant les propriétés de symétrie hermitienne des fonctions d'auto-corrélation, on arrive à l'expression (13) similaire à la formule (8) :
(13) Les termes 1/n2 peuvent être remplacés par des termes bn 2 pour assurer une approximation optimale dans une zone particulièrement intéressante (à d.s.p. maximale).
On voit que la moyenne du carré de la fréquence se réduit au terme constant 1/12 si le signal est un bruit blanc, ce qui est normal puisque la d.s.p. reste égale à l'unité sur l'intervalle [-1/2, + 1/2] ; la fréquence moyenne est alors nulle, Δf prend sa valeur maximale qui est égale à 1/
La partie du circuit de la Figure 1 fournissant la moyenne du carré de la fréquence a une constitution similaire à celle de la partie de calcul de La seule différence consiste à ajouter 1/12 à chaque valeur d'entrée provenant du registre bouclé 12a, avant normalisation.
Le calcul de Δf à partir des valeurs disponibles de et 2 peut s'effectuer de façon simple dans un circuit 20.
Comme on l'a indiqué plus haut, la valeur moyenne peut être calculée à partir des seules parties imaginaires des coefficients d'autocorrélation qui peuvent être calculés directement à partir des échantillons du signal, qui a une partie réelle Re et une partie imaginaire Im.
Le circuit montré en Figure 2 comporte deux voies ayant la même constitution. Les échantillons Re et Im du signal Doppler sont appliqués sur deux entrées respectives 22 et 24. Pour permettre au circuit de traiter des blocs d'échantillons successifs et adjacents, chaque voie comporte, sur la Figure 2, deux jeux de registres à décalage sur chaque voie, fonctionnant en alternance. Les échantillons Re de la partie réelle, sont alternativement dirigés, par le commutateur 30, vers deux registres 26 et 26a et deux registres 28 et 28a. De façon similaire, les échantillons Im sont alternativement dirigés, par un commutateur 32 commandé en synchronisme avec le commutateur 30, vers des registres 34 et 34a et des registres 36 et 36a.
Des multiplieurs 38 et 40 reçoivent les signaux des registres à travers des commutateurs commandés en synchronisme avec les commutateurs 30 et 32 et effectuent les sommes partielles des échantillons d'un même bloc, désignés par v1 à vn sur la Figure 2, sont d'abord en phase, puis avec des décalages successifs.
Le schéma montré en Figure 2 ne comporte pas de moyens arrêtant le calcul lorsque les termes deviennent négligeables ; le calcul est donc poursuivi avec n décalages successifs, en nombre égal à celui des échantillons disponibles.
Le calcul des termes successifs Im(Cxx) de la fonction de corrélation s'effectue en chaque point à l'aide d'un soustracteur 42 ; les termes sont accumulés dans un registre 44 bouclé par un réseau de bascules 46 jusqu'à ce que la totalisation ait été effectuée sur les n points disponibles.
On trouve ainsi successivement, sur la sortie Im
(Cxx) la partie imaginaire du premier moment (après la première séquence de calcul) à l'issue de la première séquence, puis les parties imaginaires des autres moments, correspondant chaque fois à un décalage déterminé, au cours des séquences suivantes. Ces moments sont successivement appliqués, dans le schéma de la Figure 1, au multiplieur 10.
S'il est nécessaire de calculer les coefficients d'auto-corrélation réels, on peut utiliser un circuit du genre montré en Figure 3, qui ne sera pas décrit étant donné sa similitude avec celui de la Figure 2.
Le circuit montré en Figure 1 peut être complété pour aller au-delà du moment d'ordre 2, notamment dans le cas d'un signal fortement perturbé ou à spectre très dissymétrique.
Une variante de réalisation de l'invention utilise, comme données d'entrée, non pas les échantillons des parties réelle et imaginaire du signal complexe, mais les paramètres d'un modèle autorégressif obtenu par exemple par mise en oeuvre du procédé suivant la demande âe brevet "Procédé et dispositif d'analyse spectrale en temps réel de signaux instationnaires complexes" FR n° 89 03761.
Pour cela, elle met en oeuvre un algorithme de Levinson qui permet de calculer les coefficients constituant la matrice de corrélation à partir de la matrice des paramètres du modèle autorégressif du signal Doppler. On pourra trouver une description de l'algorithme de Levinson dans G. Demoment, Algorithmes Rapides, Cours de l'Ecole Supérieure d'Electricité, n° 3152, 1987. L'algorithme de Levinson inverse est un algorithme rapide nécessitant de l'ordre de p2 opérations arithmétiques élémentaires pour passer d'un vecteur de prédiction ap à la matrice de corrélation correspondante. La matrice obtenue est une matrice de Toeplitz satisfaisant les conditions d'optimalite et garantissant une précision élevée, étant donné l'absence d'approximation. L'application du procédé décrit ci-dessus à l'imagerie Doppler médicale a montré des résultats satisfaisants sur des blocs courts, de huit échantillons seulement, représentatifs de conditions de travail courantes en imagerie Doppler couleurs. En particulier, le procédé suivant l'invention a révélé une variance inférieure à celle des méthodes classiques de détermination de la fréquence moyenne par des méthodes de corrélation, sur l'estimation de la fréquence moyenne en diastole.

Claims

REVENDICATIONS
1. Procédé d'extraction d'estimateurs des caractères, et notamment de la vitesse d'un écoulement instationnaire, suivant lequel on mesure le spectre Doppler d'un écoulement et on échantillonne le signal Doppler, caractérisé en ce que : on calcule les coefficients de la fonction d'autocorrélation pour chaque bloc d'échantillons ; on calcule la somme de coefficients de corrélation correspondant aux différents retards disponibles, pondérés chacun par l'intégrale du produit de la fréquence par une exponentielle négative fonction du retard pour obtenir la fréquence moyenne.
2. Procédé selon la revendication 1, caractérisé en ce qu'on calcule la somme des parties imaginaires des coefficients de corrélation du signal affectées, pour chaque coefficient de corrélation, d'un coefficient de pondération dont la valeur est fonction de l'ordre du coefficient d'autocorrélation.
3. Procédé selon la revendication 1 ou 2, caractérisé en ce qu'on détermine la largeur de bande par calcul de la moyenne des carrés de la fréquence et soustraction.
4. Dispositif d'extraction d'estimateurs des caractères d'un écoulement instationnaire, à partir d'un signal Doppler échantillonné à valeurs complexes, caractérisé en ce qu'il comprend : des moyens pour appliquer en séquence les parties imaginaires des coefficients de corrélation du signal pendant chaque fois ledit intervalle, à une entrée d'un multiplieur (10) ; des moyens (11) pour appliquer à l'autre entrée du multiplieur (10), en synchronisme, des coefficients de pondération de la forme 0, 1/π,..., (-1)P+1/πp, où p est le nombre de coefficients de pondération ; et des moyens (12, 14) pour accumuler les produits partiels et fournir une somme représentative de la fréquence moyenne.
5. Dispositif selon la revendication 4, caractérisé en ce qu'il comporte de plus des moyens (18) pour normaliser la fréquence moyenne, comprenant une table de transcodage (18) qui reçoit la sortie des moyens d'accumulation et une valeur représentative du coefficient d'autocorrélation de la fonction pour un retard nul.
6. Dispositif selon la revendication 4 ou 5, caractérisé en ce qu'il comprend de plus des moyens de détermination de la moyenne des carrés de la fréquence, comprenant des moyens pour appliquer en séquence les parties réelles des coefficients de corrélation du signal pendant chaque fois ledit intervalle de temps à une entrée d'un multiplieur ; des moyens pour appliquer à l'autre entrée du multiplieur, en synchronisme, des coefficients de pondération de la forme 0, -1/π2,..., (-1)P/p2π2 ; et des moyens pour accumuler les produits partiels.
7. Dispositif selon la revendication 6, caractérisé en ce qu'il comprend de plus des moyens de normalisation de la moyenne des carrés de la fréquence, comprenant une table de transcodage (18a) recevant le coefficient d'autocorrélation pour un retard nul et prévus pour ajouter une valeur constante à la sortie des moyens d'accumulation et diviser le résultat par ledit coefficient d'autocorrélation à retard nul.
8. Dispositif selon l'une quelconque des revendications 4 à 7, caractérisé en ce qu'il comporte de plus des moyens pour déterminer directement les parties imaginaires des coefficients de corrélation, et éventuellement les parties réelles des coefficients de corrélation, à partir des échantillons de la partie réelle et de la partie imaginaire du signal, chaque fois sur ledit intervalle de temps.
9. Dispositif selon l'une quelconque des revendications 4 à 7, caractérisé en ce qu'il comprend, de plus, des moyens pour déterminer ledit coefficient de corrélation à partir des paramètres d'un modèle autorégressif du signal par application de l'algorithme de Levinson inverse.
10. Dispositif selon la revendication 4 ou 6, caractérisé en ce que les coefficients de pondération sont ajustés pour optimiser l'approximation de la fréquence moyenne et/ou de la largeur de bande dans une plage déterminée et limitée de fréquence.
EP90905554A 1989-03-22 1990-03-22 Procede et dispositif d'extraction d'estimateurs des caracteres d'un ecoulement instationnaire a partir d'un signal doppler a valeurs complexes Withdrawn EP0422158A1 (fr)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
FR8903763A FR2644895B1 (fr) 1989-03-22 1989-03-22 Procede et dispositif d'extraction d'estimateurs des caracteres d'un ecoulement instationnaire a partir d'un signal doppler a valeurs complexes
FR8903763 1989-03-22

Publications (1)

Publication Number Publication Date
EP0422158A1 true EP0422158A1 (fr) 1991-04-17

Family

ID=9379957

Family Applications (1)

Application Number Title Priority Date Filing Date
EP90905554A Withdrawn EP0422158A1 (fr) 1989-03-22 1990-03-22 Procede et dispositif d'extraction d'estimateurs des caracteres d'un ecoulement instationnaire a partir d'un signal doppler a valeurs complexes

Country Status (3)

Country Link
EP (1) EP0422158A1 (fr)
FR (1) FR2644895B1 (fr)
WO (1) WO1990011528A1 (fr)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107137285A (zh) * 2017-06-06 2017-09-08 广东康王日化有限公司 一种牙膏用三七提取物的制备方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
See references of WO9011528A1 *

Also Published As

Publication number Publication date
FR2644895A1 (fr) 1990-09-28
WO1990011528A1 (fr) 1990-10-04
FR2644895B1 (fr) 1991-07-05

Similar Documents

Publication Publication Date Title
EP0749083B1 (fr) Procédé et dispositif de détermination du spectre de fréquence d&#39;un signal
EP3238210B1 (fr) Procede de traitement et d&#39;analyse d&#39;un signal, et dispositif mettant en oeuvre ledit procede
EP1615049B1 (fr) Traitement cohérent rapide pour codes à spectre de raies périodiques
EP1908401A1 (fr) Méthode et dispositif de mesure d&#39;une pulsation cardiaque lors de la pratique d&#39;un sport rythmique
EP0223657A1 (fr) Dispositif de calcul d&#39;une transformée de fourier discrète, et son application à la compression d&#39;impulsion dans un système radar
FR2583883A1 (fr) Dispositif et procede pour produire un facteur de merite de rapport signal a bruit pour des donnees en codage numerique
EP0166836B1 (fr) Procédé de caractérisation de la structure d&#39;un milieu, et dispositif mettant en oeuvre ce procédé
EP3671250B1 (fr) Interféromètre numérique à sous-échantillonnage
EP0140726A1 (fr) Procédé de mesure des paramètres d&#39;écoulement d&#39;un fluide et dispositif mettant en oeuvre ce procédé
WO2005004002A2 (fr) Procede de traitement d’une sequence sonore, telle qu’un morceau musical
WO1990011494A1 (fr) Procede et dispositif d&#39;analyse spectrale en temps reel de signaux instationnaires complexes
EP0309336B1 (fr) Procédé de mesure de la fréquence d&#39;un signal périodique et fréquencemètre pour la mise en oeuvre du procédé
CN112652290B (zh) 产生混响音频信号的方法及音频处理模型的训练方法
EP1172980B1 (fr) Dispositif pour la classification de signaux complexes contenant une modulation numérique linéaire
EP0422158A1 (fr) Procede et dispositif d&#39;extraction d&#39;estimateurs des caracteres d&#39;un ecoulement instationnaire a partir d&#39;un signal doppler a valeurs complexes
JP2016500847A (ja) デジタルプロセッサベースの複素音響共鳴デジタル音声分析システム
Kusuma et al. Multichannel sampling of parametric signals with a successive approximation property
EP1895433A1 (fr) Procédé d&#39;estimation de phase pour la modélisation sinusoidale d&#39;un signal numérique
FR2854247A1 (fr) Procede de traitement de cubes sismiques correspondant pour une meme zone au sol, a differentes valeurs de deports source/recepteur et/ou d&#39;angles d&#39;incidence
FR2737578A1 (fr) Dispositif radar doppler a impulsions avec determination complete du vecteur vitesse de la cible
EP0740165B1 (fr) Procédé et dispositif de traitement de signel pour lever l&#39;ambiguité d&#39;un radar Doppler
WO2005019841A2 (fr) Procede de traitement du signal pour le calcul de spectres de signaux sous-echantillonnes
FR2678738A1 (fr) Procede de detection autoregressive d&#39;un signal sinusouidal complexe dans du bruit et d&#39;estimation de sa frequence pour un radar a impulsions.
EP1233521A1 (fr) Procédé de synchronisation d&#39;un signal d&#39;horloge avec un signal de référence
EP0298786A1 (fr) Dispositif d&#39;analyse spectrale

Legal Events

Date Code Title Description
PUAI Public reference made under article 153(3) epc to a published international application that has entered the european phase

Free format text: ORIGINAL CODE: 0009012

17P Request for examination filed

Effective date: 19901116

AK Designated contracting states

Kind code of ref document: A1

Designated state(s): AT BE CH DE DK ES FR GB IT LI LU NL SE

17Q First examination report despatched

Effective date: 19911213

STAA Information on the status of an ep patent application or granted ep patent

Free format text: STATUS: THE APPLICATION IS DEEMED TO BE WITHDRAWN

18D Application deemed to be withdrawn

Effective date: 19931001