Improved ultrasonic detection using the analytic signal magnitude

Improved ultrasonic detection using the analyticsignal magnitude P.M. GAMME

LL

An alternative to rectification is proposed for detection of an ultrasonic signal for medical or non-destructive evaluation applications. With this technique the envelope of the ultrasonic waveform is obtained by calculation of the magnitude of the complex analytic signal. An implementation, particularly suited for digital data processing, is presented which utilizes a fast Fourier transform for computation of the quadrature component of the received signal using the Hilbert transform. An ordinary transducer, pulser, and receiver are used together with a digital data acquisition system. The results of this processing are compared with conventional rectification, both with and without smoothing. The superior resolvability of two closely spaced interfaces on the A-mode presentation is clearly evident. Similar improvement in B-mode quality will result when an image is required. (where a and b are real constants) complex form

Introduction Each line of an ultrasonic B-mode image is formed by brightening the trace according to the magnitude of the corresponding pulse-echo return. Although this A-mode signal is often not displayed, its characteristics nevertheless limit the quality that can be obtained in the B-mode display. The properties of the pulser circuit, transducer, amplifier, and detector circuit all are factors determining the characteristics of the A-mode signal. These factors affect the longitudinal (range) resolution and the accuracy of boundary delineation. Most commercial pulse-echo systems use full-wave rectitication as part of the signal processing chain. The envelope of the output from the rectifier is determined by the successive monopolar peaks. Since there are only a few peaks in a single reflected pulse, the envelope is poorly defined. It is customary to smooth the signal by an R-C filter. This improves interpretability at the expense of resolution of closely spaced interfaces. It will be shown in this paper that there is a way to provide a signature that is free of carrier modulation effects. Operator recognition of closely spaced interfaces will be easier with such a system. This method uses a comparatively little-known concept introduced by Gabor in 1947’ known as the analytic signal. Background In his original paper Gabor’ defined the analytic signal. When treating an electronic signal as an analytic signal, the real signal is replaced by a complex form. For example, with a simple harmonic function, the real signal f(t)

= a cos wt + b sin at

(1)

The author is at the Jet Propulsion Laboratory, California Institute of Technology, Pasadena, California 91103, and the Radiation Physics Section, Department of Radiology, University of Southern California Medical School, Los Angeles, California 90033, USA. Paper received 8 September 1980.

0041-624X/81/020073-04/$02.00 ULTRASONICS.

MARCH

1981

h(t) = f(t)

is replaced by the

t ig(t) = (a - ib) exp(iot).

(2)

The functiong(t) is produced fromf(t) by replacing cos wt by sin ot and sin wt by -cos wt. The function g(t) thus represents a signal in quadrature (retarded by 7r/2 rad) with the signal represented by f(t). If s(t) is not a simple harmonic function, the complex (analytic) signal can be obtained by similar treatment to each component of its Fourier expansion’. Heyser3 has investigated the relationship between the analytic signal and the rate-of-arrival of energy. He has shown that the square of the magnitude of the analytic signal is proportional to the instantaneous rate-of-arrival of the total energy. The square of the real signal, on the other hand, is proportional to the rate-of-arrival of one of the components of the energy. Such a component may be the kinetic energy, the potential energy, or some linear combination of the two. The square of the real signal will thus be zero at any instant when one of the component energies is zero, whereas the square of the analytic signal magnitude will only be zero when the total (kinetic plus potential) instantaneous energy is zero. Since the analytic signal magnitude is directly related to the rate of energy arrival, it is the optimal estimator of interface location for echo signals of the type commonly used in ultrasound analysis. A swept-frequency ultrasound system using the magnitude and phase of the analytic signal has been in use at the Jet Propulsion Laboratory since 1970. Images produced by this system were presented by Heyser and Le Croissette4 in 1974. A more detailed discussion of the mathematical relationships in the time and frequency domain, showing how the frequency spectrum can be compressed to obtain precise timeof-arrival information, was published in 1976’. Theoretical The technique 0

1981

approach proposed here is to replace conventional

IPC BusinessPress 73

rectification of the ultrasound signal by a means of computing the magnitude of the analytic signal. The resulting A-mode signal can then either be displayed directly or can be used to produce an improved B-mode image. The detection scheme advocated here, which uses the magnitude of the analytic signal instead of rectification, produces a signal which is proportional to the square root of the rate of arrival of energy at the transducer. The transducer will have the same directivity pattern as when measured by conventional means since the pressure (as opposed to the square of the pressure) is still integrated over the face of the transducer. Some characteristics will differ from those of the energy-sensitive transducer,6 which performs a spatial integration of the square of the pressure across the transducer surface. Clearly, each of these techniques will have merit in different applications. If the real part of a signal, YJ~), has time that it is not of zero amplitude, part, vi(t) can be obtained from the is the Cauchy principal value of the

been recorded for all then the imaginary Hilbert transform, which integral:

+Vi(t)

=

i L

v,(i)

(n)-’ (t-t’)-’

dt’

The uniqueness of this relationship is related to minimum phase considerations, pole locations, and causality2*3. Equation (3) represents a convolution of the signal with the kernel I/(nt). By the convolution theorem, convolution of two functions can be implemented in the frequent? domain by multiplication of their Fourier transforms. Generally, direct computation of the convolution is more time consuming for digital computations than the process of transformation to the frequency domain using a fast Fourier transform (fft) program, multiplication of the spectra, and then performing the inverse fft. The kernel l/(nt) is not of class L2 (that is, it is not squareintegrable). It is, however, a well-behaved functiona17, whose Fourier transform is -i sgnm, (where sgn(j) = -1 for fO). Thus, the convolution of (3) can be performed in the frequency domain by multiplying the Fourier transform of vr(t)by -i sgn(j) and then applying an inverse Fourier transform to the product. If a real fft program is used, multiplication of the fft of the real signal by -i sgnm to obtain the fft of the imaginary part, then performing the inverse fft, is the most appropriate approach. If a complex fft program is used, however, a shortcut can be applied, which takes advantage of certain symmetry properties of the Fourier transform of analytic signals. The Fourier transform of the complex analytic signal can be obtained from the Fourier transform of one of its components by suppressing the negative frequency contributions’. This extremely simple means of obtaining the full (complex) analytic signal consists of: 1 - Fourier transforming the real data using a complex fft with the imaginary part set equal to zero. 2 - Setting all negative frequency components to zero. 3 - Taking an inverse fft. The real and imaginary parts of the complex analytic signal are then in the appropriate data arrays. With the full analytic signal available, its magnitude may be calculated as the square root of the sum of the squares of the real and the imaginary parts at each point in time. The

74

analytic signal magnitude can then be used either to produce a better A-mode display or can be further processed by grey scale modification (logarithmic amplification, thresholding, etc) and combined with information on the position and orientation of the transducer (or of the array steering) to produce an improved B-mode image. The analytic signal magnitude would be expected to be superior to rectification for automated signal identification in medicinegl” and in non-destructive evaluation”. Experimental

verification

Experiments were performed to compare this processing scheme with conventional rectification and filtering. A digital data acquisition system was used which digitized the unprocessed (raw rf) ultrasonic signal. For comparison of the effect of the signal processing schemes on the same input signal, the same digital record that was used for analytic signal magnitude processing was also digitally rectified and smoothed The amplified echo was digitized using a Biomation 8 100 transient recorder. This instrument was interfaced with a Digital Equipment Corporation LSI-11 microprocessor, which stored the data on diskettes and performed all of the later processing. The plots were produced using the graphics mode of the Diablo Hytype 1641 terminal. All programming was done in I_abforth12 , which is an RT-11 resident version of Forth.13,14 This language was chosen because assembly language routines can be written directly in the main program and because it offers the flexibility and immediate accessibility of an interpreter. The floating point error traps were modified to allow processing to continue when underflow errors were encountered. The A-scans were produced using a Panametrics” 5052-PR pulser-receiver. The damping control, which varies the shunt resistance across the transducer, was set for maximum damping. A moderately damped 2.25 MHz, 12.7 mm diameter unfocused transducer was used. The range of the target from the transducer was approximately 10 cm. The target consisted of a large sheet of acrylic plastic (trade names: Lucite, Plexiglass, or Perspex) which was carefully aligned to obtain the maximum specular echo. In some experiments a 1 cm section of formalin-fured hog liver was placed between the transducer and the target to simulate the general type of signal distortion (attenuation, refraction, and scattering) that is encountered when imaging through biological media. This tissue was present in all of the data presented in Figs 14. The data were digitized to eight-bit resolution at 10 ns intervals for records of 2048 points. Although this is considerably oversampled by the Nyquist criterion, economy of data was not an important consideration and this slight extravagance provides a more recognizable signal than does sampling at a rate that gives only two samples per cycle of the highest frequency. For economy of data storage, a 2048 point real-valued fft algorithm was used. To minimize leakage effects, 205 points at the beginning of the record and the same number at the end were weighted (apodized) by multiplication by a raised cosine function. The reflected signal from a single interface is shown in Fig. la. Fig. lb shows the computed Hilbert transform of that signal. The vertical arrows show that the peak of each of these waveforms occurs at the same epoch as the inflection of the other, as required by the dispersion relations16.

ULTRASONICS.

MARCH 1981

.a-. .

The phase information contains the sign of the reflection coefficient and the Doppler shift.

. . . .

a

,-. __

The results of three types of detection (rectification alone, rectification with smoothing, and computation of the analytic signal amplitude) have been compared in Fig. 3 and Fig. 4 for the case of one or two planar interfaces with the tissue between the transducer and the reflector. The echoes from the tissue, which would be far to the left of these traces, are not shown. Similar results were obtained when only water was in the path between the transducer and the reflector. The results with tissue in the path are shown as they represent a more realistic case which includes distortions induced by the propagation media.

. .

, ,-; : .. .. : . . .

In Fig. 3, the results of these three types of detection are shown for a single planar reflector. The absolute value of the received signal, which would be the result of full wave rectification with little or no smoothing, is shown in Fig. 3a. It is certainly difficult to identify the precise location of the peak signal on this A-mode presentation.

.. . .

b . . .

.

.

.z: .

. . .:

Fig. 1 Quadrature components of the analytic signal of the reflection from a single interface: a - real part; b - imaginary part

Fig. 2 is an isometric plot of the real and imaginary parts of this complex function against time. Note that this complex function is a spiral of slowly varying amplitude. This circular nature is due to the sine-cosine duality of the analytic signal. This complex function may be considered as a complex exponential, of approximately 2.25 MHz frequency, modulated by the more slowly varying amplitude function. A detection scheme that recovers this slowly varying amplitude function is equivalent to shifting this signal to baseband and contains exactly the same positional information as the original signal; only the phase information has been suppressed.

Fig. 3b is the result of applying an R-C filter, simulated by digital processing, to Fig. 3a. The tune constant was chosen so that exp(-t/KC) was equal to 0.5 for t equal to approximately two half-cycles of the 2.25 MHz transducer centre frequency. This choice of a time constant produces a smoothing that represents a typical compromise: it is not quite enough smoothing to give a smooth rising edge to the signal, although it is high enough to produce a long trailing edge, which sacrifices resolution. A longer or shorter time constant would result in improvement of one of these aspects but degradation of the other. Multi-pole filters are an alternative which may produce better results for some signals. They are not considered here since the main purpose of this study is

Fig. 3 Processed signatures of the signal from a single interface: a - rectified; b - smoothed; c - magnitude of the analytic signal

;:,T-.. ;- :... ;A? /

c Complex representation Fig. 2 from a single interface

ULTRASONICS.

MARCH

of the analytic

1981

signal of the reflection

/

,ts

\. :

.’ i *,

i

‘\:,,j

b

J?,;*,

i

’

2

Time r@l

\

‘\-

Fig. 4 Processed signatures of a signal from two closely spaced interfaces: a - rectified; b -smoothed; c -magnitude of the analytic signal

75

to

provide a simple comparison of this new technique with the typical limitations of conventional processing. Note, incidentally, that smoothing delays the peak by approximately one half-cycle, as compared to Fig. 3a.

The analytic signal magnitude has been successfully used in our laboratory for the production of images of biological specimens that include the soft internal echoes.

The analytic magnitude of the same signal, as shown in Fig. 3c, gives an A-scan that is easier to read than the rectified signal and has better resolution than the rectified and smoothed signal. This signal is much smoother than the signal that is only rectified and yet its peak occurs at approximately the same epoch. The analytic signal magnitude is expected to be the best measure of the time at which the peak occurs. Since the maximum amplitude of the true analytic signal of Fig. 2 could occur at any phase angle, including along the imaginary axis, the peak of the rectified real part of Fig. 3a could differ from the true energy maximum by as much as a quarter of a cycle. Notice that the weak tail that follows the main part of the echo can now be identified as not being due to the simple ring-down of the transducer; it is suspected as being an echo from the transducer backing.

Acknowledgements

The results of applying the same three processing schemes to overlapping echoes from the front and back surfaces of a thin piece of Lucite are shown in Fig. 4. If the general transducer response were known, the existence of two discrete interfaces might be surmised by close inspection of Fig. 4a. The smoothing of Fig. 4b actually serves to make interpretation more difficult. It would be difficult to image these interfaces separately on a B-mode. To separate these echoes, the decision threshold would have to be between the lower peak and the valley, which represents a tolerance of 28% of the signal amplitude. These two interfaces are clearly resolved by the analytic signal magnitude processing of Fig. 4c. For this signal to be imaged as two separate points the decision threshold could be anywhere between the lower peak and the deep valley between the peaks, a range which represents 77% of the signal amplitude. Clearly the analytic signal magnitude processing makes identification of the two interfaces easier than either rectification alone or than rectification with smoothing.

The work described in this paper was carried out at the Jet Propulsion Laboratory, California Institute of Technology, jointly sponsored by the National Aeronautics and Space Administration under Contract NAS7-100, and by the National Institutes of Health, Biomedical Research Support Grant Program Division of Research Resources, under BRSG Grant RR07003, by agreement with the National Aeronautics and Space Administration. The analysis programs and the microprocessor used in this work were originally funded by the NASA Office of Life Sciences. The author is grateful to Richard C. Heyser for his enlightening discussions on the significance of the analytic signal. References 1 2

3 4

5

6

Conclusions This study clearly demonstrates that the analytic signal magnitude provides better resolution than does either rectification alone or rectification with smoothing, at least for the case of simple boundaries that may be closely spaced. In these experiments the analytic signal magnitude was calculated by computer-processing of sampled data. These calculations serve to demonstrate the utility of the method. Although digital implementation of such processing is currently prohibitively expensive, in a few years production of a clinical or industrial instrument incorporating this technique may well become commercially viable. Alternative implementations, using analogue circuits, are currently being investigated.

76

10 11

12 13

14 15 16

Gabor, D., Theory of Communication, J. Znst ElectricalErzgrs (London), 93 (1946) 429-457 Blake, W.K., Waterhouse, R.V., The use of cross-spectral Density Measurements in Partially Reverberant Sound Fields, J. Sound Vib 54 (1977) 589-599 Heyser, R.C., Determination of Loudspeaker Signal Arrival Times: Part 111, J. Audio Eng. Sot., 19 (1971) 902-905 Heyser, R.C., Le Croissette, D.H., A New Ultrasonic Imaging System Using Time Delay Spectrometry, Ultrasound Med Biol l(1974) 119-131 Le Croissette, D.H., Heyser, R.C., Attenuation and Velocity Measurements Using Time Delay Spectrometry, Proc. Conf. Ultrasonic Tissue Characterization, National Bureau of Standards Spec. Pub. 453 (1976) 81-95 Busse, L.J., Miller, J.G., Yuhas, D.E., Mimbs, J.W., Weiss, A.N., Sobel, B.E., Phase Cancellation Effects: A Source of Attenuation Artifact Eliminated by a CdS Acoustoelectric Receiver, Ultrasound in Medicine, 3, Denis White, Ed., Plenum Press, New York, 1519-1535 Lighthill, M.J., Introduction to Fourier Analysis and Generalized Functions, Cambridge, The University Press (1964) 16-21 Bracewell, R., The Fourier Transform and its Applications, New York, McGraw-Hill (1965) 267-272 Hudson, A.C., Muller, H.R., A Pocket-Size Echoencephalograph: Description and First Clinical Results, Surgical Neurology, 2 (1974) 201-206 Galicich, J.H., Williams, J.B., A Computerized Echoencephalograph, J. Neurosurg. 35 (1971) 453459 Mucciardi, A.N. Adaptive Learning Network (ALN) Hardware, PIOC. ARPA/AFML Review of Progress in Quantitative NDE (July 17-21’1978, Scripps Institute), DARPA Report No. AFML-TR-78-205 (1979) available NTIS. 88-91 Laboratory Software Systems, Inc., 3634 Mandeville Canyon Rd., Los Angeles, CA 90049, USA Ewing, M.S., The Caltech Forth Manual, Second Edition, Owens Valley Radio Observatory, California Institute of Technology, Pasadena, California 91125 (1978) Moore, C.H., Forth: A New Way to Program a Mini-computer, Astronomy and Astrophysics Supplement, 15 (1974) 497-511 Panametrics, Inc., Waltham, MA 02154, USA Heyser, R.C. Some Useful Graphical Relationships, J. Audio Eng. Sot. 23 (1975) 562-565

ULTRASONICS.

MARCH 1981