Alias Suppression in Pulsed Doppler Systems

The invention relates to methods and apparatus for measuring velocity by analysis of Doppler shift induced in broadband pulsed waves. The invention is particularly useful in conjunction with echo ultrasound systems which measure blood flow within a human or animal body.

Pulsed Doppler ultrasound systems are commonly used to measure and map the velocity of blood flow within human and animal bodies. Pulses of ultrasound energy are directed into the body along a path which intersects blood vessels or a coronary chamber. Ultrasound energy from the pulses is backscattered from blood within the vessel or chamber and returns to a transducer where it is converted into an electrical signal. If the blood has a velocity component along the direction of propagation of the ultrasound waves, the frequency of the scattered echoes will be shifted, in relation to the frequency of the incident ultrasound energy. The Doppler shift which is thus induced in the echoes can be analyzed to yield a numeric estimate of the blood velocity and/or to produce a map of blood velocity as a function of position within the body.

U.S. Pat. 4,930,513 describes a system wherein the amplitude of echoes along an ultrasound A-line are sampled at discrete times at a rate which is above the Nyquist frequency. RF sample vectors at a selected range from successive A-lines, taken in the same direction through a region of the body, form a two-dimensional matrix with element positions described by a first, fast-time, variable which specifies the range of a data sample along its A-line, and a second, slow-time, variable which specifies the position of the A-line within a group of collected A-lines. The data matrix is processed with a two-dimensional discrete Fourier transform, taken with respect to the fast-time variable and the slow-time variable, which maps the data set into a discrete two-dimensional Fourier frequency space wherein constant velocity Doppler shifts appear on radial lines whose angle with respect to the vertical axis is a function of velocity. The slow-time axis corresponds to a slow-frequency (Doppler frequency) axis in the two-dimensional Fourier frequency space and the fast-time axis corresponds to a fast frequency (radio frequency) axis in the Fourier frequency space. Echoes of wideband pulses which are scattered from a moving target are mapped as generally elliptical shapes in the two-dimensional Fourier space. The angle between the major axis of an ellipse and the coordinate axes is a high quality measure of the velocity of the scattering medium. The angular distribution of the spectrum in two-dimensional Fourier frequency space (and hence velocity components in the region of the object) are estimated by computing a radial projection of the transformed data array.

The discrete Fourier transform is periodic in nature. As a result, the spectral components produced by high velocity scatterers will tend to "wrap around" unit cells in the Fourier frequency space. This is the familiar "aliasing" problem associated with a Fourier transform of any sampled signal. Thus, a typical radial projection of the two dimensional spectrum generated by pulse echoes from a moving medium will have a number of peaks at different velocities. The peak with the greatest maximum and least width generally corresponds to the true velocity spectrum of the scattering medium. The other peaks, which are alias terms, generally have smaller maxima and greater widths. Because of the random nature of ultrasound scattering from blood, an alias term may occasionally have a higher maximum or a narrower width than the true term. However, it is very unusual for an alias term to have both a higher maximum and a narrower width than the true term. Thus, a human observer of a plotted radial projection can almost always identify the term in the radial projection that corresponds to the true velocity spectrum of the scattering medium.

The invention is an apparatus and method for processing two-dimensional Fourier transform data from pulsed Doppler systems to automatically identify the true velocity components of a radial projection and suppress alias terms.

The method is based on an algorithm which uses a weighted combination of the maximum value and a width measure (such as variance) of each term in the radial projection to identify a term that corresponds to the true velocity. To do this, it is necessary to appropriately isolate each term in the radial projection so that the maximum value and width can be computed.

To compute the width of a term in the radial projection it is necessary to only use data samples in a "window" that is centered on each term and extends half way to adjacent terms on each side. If this data window is not properly centered, contributions from adjacent terms can cause significant error in the computed width.

One way to locate the windows is to center them around the peak values; unfortunately, because of the random nature of the scattering, the terms tend to be ragged and unsymmetrical. Therefore the mean of each term in the radial projection is used as a more appropriate point on which to center the window. In order to compute the mean for a term, however, it is necessary to have the window appropriately centered, otherwise the computation of the mean will be contaminated by contributions from adjacent terms.

An alternate approach is to locate the troughs between the terms for use as the window edges. This can be done by taking advantage of the periodic nature of the terms as well as the fact that the energy in each term, although distributed differently, is approximately the same. This can be accomplished by an iterative calculation that successively approximates the location of the mean of the spectrum such that the terms are evenly divided by this mean.

A preferred embodiment of the invention uses an algorithm which assumes only that there are an even number of terms in the spectrum. The locations of the terms, even when a portion of a term is located on the other side of the spectrum, does not matter. In practice, the algorithm has been found to converge very rapidly, yielding an accurate value of the mean after only a few iterations.

Initially, the mean is computed starting at a first data sample in the spectrum. If the computed mean is at the center of the spectrum, then the starting point has been found; the width and the maximum value of each term can then be computed. However, the first computed mean will often not be at the center of the spectrum because of wrap around of terms in the radial projection. In this case, the data are rotated so that the first computed mean is at the center and the mean of the rotated spectrum is recomputed. The process is repeated --compute mean, rotate data, compute mean, rotate data, etc.--until the computed mean is at the center of the data. When the computed mean is at the center of the data, then the mean has been accurately located.

Once the mean has been located, the method divides the radial projection into adjacent non-overlapping windows of the same size, one for each term, with the first window starting at the mean. Within each window the maximum and width are computed, and in turn, they are used to compute a score for each window term in the radial projection according to a function that is maximized by a large maximum and a small width, for example, the formula:

where .beta. is a parameter that controls the relative weight given to max and variance and maxmax is the maximum of the entire spectrum. The window with the lowest score is taken to be the window containing the true velocity spectrum; the spectra in the other windows are assumed to be aliases and are suppressed.

The invention is described with reference to the attached drawings in which:

FIG. 1 is an ultrasound examination system in accordance with the invention;

FIG. 2 illustrates the radial projection of an aliased signal in two-dimensional Fourier frequency space;

FIG. 3 is an example of a radial projection of a two-dimensional Doppler spectrum indicating blood flow velocity;

FIG. 4 is a radial projection of a two-dimensional Doppler spectrum in which part of one term has wrapped around; and

FIG. 5 is a radial projection with correctly located windows for computation of the mean and width of each spectral component.

U.S. Pat. No. 4,390,513 is incorporated herein by reference.

FIG. 1 illustrates an ultrasound blood velocity measurement system which incorporates the invention. A transmitter 100 excites an ultrasound transducer 110 via a TR switch 120 to transmit broadband pulses of ultrasound energy into a body 130 long a narrow beam 140. The ultrasound energy in the beam is backscattered from structures in the body, for example from blood flowing in an artery 150, to produce echoes which return to and are detected by transducer 110. The frequency of the ultrasound signals scattered from structures in the body, which are moving in the direction of the axis of the beam 140 will be shifted with respect to the frequency of the ultrasound energy produced by the transducer in accordance with the Doppler equation. Echo signals produced by the transducer are coupled to a receiver 160 via the TR switch 120 and are amplified therein. The output of the receiver is connected to a radio frequency digitizer 170 which extracts and digitizes samples of the echo signals. The output of the digitizer is connected to the input of a range and line gate 180. Each transmitted RF pulse produces an A-line signal which maps range along the axis of beam 140 into signal arrival time. Control circuit 190 operates range gate 180 to select signal samples which originate within the region of the body defined by range, R, and beam 140. The digitizer 170 periodically samples the RF signal of the A-line at a rate which is at least twice the frequency of the highest component in the signal. Successive samples which lie within the range, R, along each A-line are stored as a column vector in a matrix array in memory 200. Successive A-lines are stored at successive column in the matrix so that each row of the matrix represents samples taken at corresponding ranges along successive A-lines. Individual elements stored in the matrix in memory 200 may thus be identified and addressed by a "fast-time" index, (which identifies a sample along each A-line and thus a row in the matrix) in conjunction with a "slow-time" index (which identifies the specific A-line and thus the column in the matrix).

The array stored in memory 200 is processed in a two-dimensional discrete Fourier transform processor 210 and the results are stored in a second processor 220. Depending on speed and data constraints of the system, the processor 210 may either be accomplished as a software routine in a general purpose digital computer or a microprocessor or, alternately, may be effected by dedicated fast Fourier transform processor chips.

The output of the processor 210 is an array of data elements which represent the transform of the sampled A-line data with respect to a first (slow-time) variable which identifies individual A-line vectors in a data set and a second (fast-time) variable which represents the discrete sample position along each A-line. The array in memory 220 may be interpreted as a two-dimensional frequency spectrum representation with respect to a fast frequency variable which corresponds to the RF frequency spectrum of the received echo signal and a slow-frequency variable which corresponds to the Doppler frequency components.

The two dimensional discrete Fourier transform data provides information about the velocity of the scattering medium in a region corresponding to the range gate R. If the medium in the range gate is moving towards the transmitter with a velocity, v, the two-dimensional discrete Fourier transform will have a generally elliptical form, T, as shown in FIG. 2. For typical parameter values the major axis of the ellipse passes through the origin of the two-dimensional frequency plane making an angle, .theta., with the vertical axis. The spectral data stored in memory 220 are processed in an analyzer 230 to extract velocity information. In a preferred embodiment, the velocity spectrum is determined by radial projection in the two-dimensional Fourier frequency space. The angle, .theta., and hence the velocity, v, is estimated by computing the radial projection of the magnitude of the two-dimensional Fourier transform as taught by U.S. Pat. No. 4,930,513.

If components of the velocity are sufficiently high, the two-dimensional Fourier transform data in memory 220 will include alias spectral components, A, as well as true spectral components, T (FIG. 2). In general, the radial projections, a, of the alias spectra, A, will have smaller magnitude and be wider than the radial projection, t, of the true spectrum, T. This is further illustrated in FIG. 3 which is a cartesian plot of a radial projection velocity spectrum. In a practical system, all of the alias components are not necessarily contiguous. As illustrated in FIG. 4 one of the alias components, a*, may wrap around the end of the spectrum.

In order to suppress alias spectral components, the analyzer 230 (FIG. 1) includes a program which divides the radial projection (FIG. 5) spectrum into windows, W.sub.1 ---- W.sub.n, using the troughs between terms as the window edges. This is accomplished by an iterative calculation that successively approximates the location of the mean, M, of the spectrum such that the terms are evenly divided. Initially the mean is computed starting at the first data sample in the spectrum. If the computed mean is at the center of the spectrum, then the starting point has been found. If the mean is not at the center, the spectrum is rotated so that the initially computed mean is at the center and the new mean of the rotated spectrum is calculated. The process is iterated until the mean of the rotated spectrum is at its center.

The maximum value and width of each term can then be computed and tested to identify the true spectrum term, t, and the alias terms are suppressed.

Appendix A contains an example of a C language computer program with a module alias.sub.13 sup which implements the alias suppression techniques of the invention. ##SPC1##